$\boldsymbol{V} \in \mathbb{R}^{m \times n}$を入力とし,これを$\boldsymbol{V} \approx \boldsymbol{W} \boldsymbol{H}^T$と分解することを目的とする．ただし，$\boldsymbol{W} \in \mathbb{R}^{m \times k}$, $\boldsymbol{H} \in \mathbb{R}^{n \times k}$とする．なお，今回は単純のため，正則化は考えない．拡張は難しくない（はず）．
さて，解きたい最適化問題は以下で与えられる．
\begin{align*}
\min_{\boldsymbol{W}, \boldsymbol{H}} &\ \left\|\boldsymbol{V} - \boldsymbol{W}\boldsymbol{H}^T \right\|^2_{F} \\
s.t. &\ \boldsymbol{W} \succeq \boldsymbol{0} \\
&\ \boldsymbol{H} \succeq \boldsymbol{0} \\
\end{align*}
なお，今回はFrobeniusで行列間の距離を考えるが，他の指標でも構わない（KL Divergenceとか）．Lagrangianは以下のようになる．
\begin{align*}
L(\boldsymbol{W}, \boldsymbol{H}, \boldsymbol{A}, \boldsymbol{B})
&= \left\|\boldsymbol{V} - \boldsymbol{W}\boldsymbol{H}^T \right\|^2_{F} - \sum_{i=1}^m \sum_{j=1}^k A_{ij} W_{ij} - \sum_{i=1}^n \sum_{j=1}^k B_{ij} H_{ij} \\
&= tr\left(\boldsymbol{H} \boldsymbol{W}^T \boldsymbol{W} \boldsymbol{H}^T\right) - 2 tr \left(\boldsymbol{V}^T \boldsymbol{W} \boldsymbol{H}^T\right)
- tr\left(\boldsymbol{W} \boldsymbol{A}^T\right) - tr\left(\boldsymbol{H} \boldsymbol{B}^T\right) + const.
\end{align*}
以下文章が適当になる．
KKT condition (Stationarity):
\begin{align}
\frac{\partial L}{\partial \boldsymbol{W}} = 0 &\Rightarrow 2 \boldsymbol{W} \boldsymbol{H}^T \boldsymbol{H} - 2 \boldsymbol{V} \boldsymbol{H} - \boldsymbol{A} = 0 \cdots (1) \\
\frac{\partial L}{\partial \boldsymbol{H}} = 0 &\Rightarrow 2 \boldsymbol{H} \boldsymbol{W}^T \boldsymbol{W} - 2 \boldsymbol{V}^T \boldsymbol{W} - \boldsymbol{B} = 0 \cdots (2) \\
\end{align}
KKT condition (Complementary slackness):
\begin{align*}
-A_{ij}W_{ij} &= 0 \ (\forall i=1, \ldots, m \ \forall j=1, \ldots, k) \cdots (3)\\
A_{ij} &\geq 0 \ (\forall i=1, \ldots, m \ \forall j=1, \ldots, k) \cdots (4)\\
-B_{ij}H_{ij} &= 0 \ (\forall i=1, \ldots, n \ \forall j=1, \ldots, k) \cdots (5)\\
B_{ij} &\geq 0 \ (\forall i=1, \ldots, n \ \forall j=1, \ldots, k) \cdots (6)\\
\end{align*}
\begin{align*}
(1), (3) &\Rightarrow W_{ij} \left(\left[\boldsymbol{W} \boldsymbol{H}^T \boldsymbol{H}\right]_{ij} - \left[\boldsymbol{V} \boldsymbol{H}\right]_{ij} \right) = 0
\Rightarrow W_{ij} = W_{ij} \frac{\left[\boldsymbol{V} \boldsymbol{H}\right]_{ij}}{\left[\boldsymbol{W} \boldsymbol{H}^T \boldsymbol{H}\right]_{ij}} \cdots (7) \\
(2), (5) &\Rightarrow H_{ij} \left( \left[\boldsymbol{H} \boldsymbol{W}^T \boldsymbol{W}\right]_{ij} - \left[\boldsymbol{V}^T \boldsymbol{W}\right]_{ij} \right)= 0
\Rightarrow H_{ij} = H_{ij} \frac{\left[\boldsymbol{V}^T \boldsymbol{W}\right]_{ij}}{\left[\boldsymbol{H} \boldsymbol{W}^T \boldsymbol{W}\right]_{ij}} \cdots (8)
\end{align*}
(7), (8)から解は右辺の関数のfixed pointであることがわかる．なので，以下の更新式（fixed point iteration）で解を求めることができる．
\begin{align*}
W_{ij} \leftarrow W_{ij} \frac{\left[\boldsymbol{V} \boldsymbol{H}\right]_{ij}}{\left[\boldsymbol{W} \boldsymbol{H}^T \boldsymbol{H}\right]_{ij}} \\
H_{ij} \leftarrow H_{ij} \frac{\left[\boldsymbol{V}^T \boldsymbol{W}\right]_{ij}}{\left[\boldsymbol{H} \boldsymbol{W}^T \boldsymbol{W}\right]_{ij}}
\end{align*}
ちなみに，
W_{ij} \leftarrow W_{ij} - \epsilon_{ij} \left(\left[\boldsymbol{W} \boldsymbol{H}^T \boldsymbol{H}\right]_{ij} - \left[\boldsymbol{V} \boldsymbol{H}\right]_{ij}\right)
がGradient descentなわけですが，
\epsilon_{ij} = \frac{W_{ij}}{\left[\boldsymbol{W} \boldsymbol{H}^T \boldsymbol{H}\right]_{ij}}
とおくと，fixed point iterationは，adaptive stepsizeを用いた，Gradient descentだと見ることもできる．$\boldsymbol{H}$も同様．