# Adaptive Clustering via Symmetric Nonnegative Matrix Factorization of the Similarity Matrix
## 2 The Alternating Nonnegative Least Squares for NMF  
Problem (3) is nonconvex and finding its global minimum is NP-hard [15]. In general, nonconvex optimization algorithms guarantee only the stationarity of the limit points, and only local minima can be detected. Most iterative algorithms for solving (3) are based on a simple block nonlinear Gauss-Seidel scheme [4; 5]. An initial \(W_{0}\in\mathbb{R}_{+}^{m\times k}\) is chosen, and the sequences  
\[\left\{ \begin{array}{l}H_{\
u}=\underset{H\geq 0}{\text{argmin}}\ \frac{1}{2}\ \|A-W_{\
u-1}\,H^{T}\|_{F}^{2},\\ W_{\
u}=\underset{W\geq 0}{\text{argmin}}\ \frac{1}{2}\ \|A^{T}-H_{\
u}\,W^{T}\|_{F}^{2}, \end{array} \right. \tag{4}\]  
are computed for \(\
u=1,2,\ldots\), until a suitable stopping condition is satisfied, such as by checking whether a local minimum of the objective function \(f(W,H)\) has been sufficiently well-approximated. The convergence of this scheme, called Alternating Nonnegative Least Squares (ANLS), follows from [16]. Although the original problem (3) is nonconvex, the two subproblems of (4) are convex. Let us denote by \(\varphi(X)\) the objective function to be minimized in both cases,  
\[\varphi(X)=f(Y,X)=\tfrac{1}{2}\ \|B-Y\,X^{T}\|_{F}^{2}, \tag{5}\]  
where \(B=A\) and \(Y=W_{\
u-1}\) for the first subproblem in (4), and \(B=A^{T}\) and \(Y=H_{\
u}\) for the second subproblem in (4). Problem (5) can be dealt with by applying aprocedure for nonnegatively constrained least squares optimization, such as an Active-Set-like method [13; 17; 18; 19] or an iterative inexact method as a gradient descent method modified to satisfy the nonnegativity constraints. We have performed in [14] a preliminary ad hoc experimentation, comparing an Active-Set-like method (the one called BPP and coded as Algorithm 2 in [13]) with the Greedy Coordinate Descent method, called GCD in [7]. The results have shown that in our context, GCD is faster and equally reliable. GCD, specifically designed for nonnegative constrained least squares problems, implements a selection policy of the elements updated during the iteration, based on the largest decrease of the objective function. We give here a brief description of GCD (the corresponding code can be found as Algorithm 1 in [7]). The gradient of \(\varphi(X)\) is  
\[G(X)=\big{(}X\,Y^{T}-B^{T}\big{)}Y=XQ-B^{T}Y,\quad\text{where}\quad Q=Y^{T}Y.\]  
GCD computes a sequence of matrices \(X_{j}\), \(j=1,2,\ldots\), until termination. Initially, \(X_{0}\) is set equal to \(H_{\
u-1}\) for the first subproblem in (4) and to \(W_{\
u-1}\) for the second subproblem in (4) (with \(H_{0}=O\), which ensures that \(G(X_{0})\leq O\)). At the \(j\)th iteration, matrix \(X_{j}\) is obtained by applying a single coordinate correction to the previous matrix \(X_{j-1}\) according to the rule \(X_{j}=X_{j-1}+s\,\mathbf{e}_{r}\mathbf{e}_{i}^{T}\), where \(r\) and \(i\) are suitably selected indices and \(\mathbf{e}_{r}\) and \(\mathbf{e}_{i}\) are the \(r\)th and \(i\)th canonical vectors of compatible lengths. The scalar \(s\) is determined by imposing that \(\varphi(X_{j})\), as a function of \(s\), is the minimum on the set  
\[S=\{s\text{ such that }(X_{j-1})_{r,j}+s\geq 0\}.\]Denoting \(g_{r,i}=\big{(}G(X_{j-1})\big{)}_{r,i}\) and \(q_{i,i}\) the \(i\)th principal element of \(Q\), the value which realizes the minimum of \(\varphi(X_{j})\) on \(S\) is  
\[\widehat{s}=-\,\frac{g_{r,i}}{q_{i,i}}\quad\text{if}\quad\frac{g_{r,i}}{q_{i,i} }\leq(X_{j-1})_{r,i}\quad\text{and}\quad\widehat{s}=-\,(X_{j-1})_{r,i}\quad \text{otherwise}.\]  
In correspondence, the objective function is decreased by  
\[\theta_{r,i}^{(j)}=\varphi(X_{j-1})-\varphi(X_{j})=-\,g_{r,i}\,\widehat{s}- \frac{1}{2}\,q_{i,i}\,\widehat{s}^{2}.\]  
A natural choice for indices \((r,i)\) would be the one that maximizes \(\theta_{r,i}^{(j)}\) on all the pairs \((r,i)\), but it would be too expensive. So it is suggested to proceed by rows.