## 3 Method for Estimating the NMF Rank  
### Rank Estimation Using [PERSON]'s Unbiased Risk Estimator  
If we define \(\mathbf{V}_{0}\triangleq\mathbf{W}\mathbf{H}\), NMF Model (2) can be rewritten as follows:  
\[\mathbf{V}=\mathbf{V}_{0}+\mathbf{E} \tag{8}\]  
According to the work of [PERSON] and [PERSON] [9] and the example shown in Section 2, a well-decomposed basis corresponds to a musical note event. Because different musical notes have different fundamental frequencies (\(f_{0}\)) and because the spectrum of a note event forms a harmonic line spectrum comprising \(nf_{0}\) (\(n\) is an integer), the columns of the matrix \(\mathbf{W}\) are linearly independent (Figure 5). Hence, the nonnegative ranks of all of the matrices \(\mathbf{W}\), \(\mathbf{H}\), and \(\mathbf{V}_{0}=\mathbf{W}\mathbf{H}\) are equal to \(R\). Therefore, the estimation of the number of basis vectors is the same as the estimation of a nonnegative rank in a noisy environment. Unfortunately, the estimating of the nonnegative rank is NP-hard problem [24]. Because the rank and the nonnegative rank are related (the matrix rank is the lower bound of the nonnegative rank [25]), we decide to estimate the rank instead of the nonnegative rank. Figure 4: Example of the NMF analysis results for (**a**) the column vectors of \(\mathbf{W}\) and (**b**) the row vectors of \(\mathbf{H}\) when the NMF is excessively small. Figure 3: Example of NMF analysis for (**a**) the column vectors of \(\mathbf{W}\) and (**b**) the row vectors of \(\mathbf{H}\), when the number of basis vectors is excessively large. If we assume that the rank of \(\mathbf{V_{0}}\) is \(R\), \(\mathbf{R_{V_{0}}}=(1/N)\mathbf{V_{0}}\mathbf{V_{0}^{T}}\) has \(R\) non-zero eigenvalues but \(\mathbf{R_{V}}=(1/N)\mathbf{V}\mathbf{V}^{T}\) has more than \(R\) non-zero eigenvalues because of \(\mathbf{E}\). Therefore, the rank estimation problem is the same as the problem of identifying \(R\) signal eigenvalues in a noisy environment, which can be handled using noisy principal component analysis (PCA). The SURE algorithm is a well-known method for model order estimation in the noisy PCA problem; it selects the model order as the value that minimizes a risk function. To apply the SURE method to the NMF problem, Model (2) can be rewritten as follows:  
\[\mathbf{v}(n)=\mathbf{Wh}(n)+\mathbf{e}(n)=\mathbf{v_{0}}(n)+\mathbf{e}(n) \tag{9}\]  
where \(\mathbf{v}(n)\), \(\mathbf{v_{0}}(n)\), \(\mathbf{h}(n)\), and \(\mathbf{e}(n)\) are the \(n\)th columns of \(\mathbf{V}\), \(\mathbf{V_{0}}\), \(\mathbf{H}\), and \(\mathbf{E}\), respectively. \(n\) denotes both the column index and the discrete-time index. The SURE method selects the model order \(R\), which is same as the number of column vectors of the matrix \(\mathbf{W}\) that minimizes the risk function [16; 18]  
\[\mathfrak{R}_{R} = E\|\mathbf{v_{0}}(n)-\hat{\mathbf{v}_{0}}(n)\|^{2} \tag{10}\] \[= E\|\mathbf{e}(n)-\hat{\mathbf{e}}(n)\|^{2}\] \[= E\|\hat{\mathbf{e}}(n)\|^{2}-2E(\hat{\mathbf{e}}^{T}(n)\mathbf{ e}(n))+K\sigma_{e}^{2}\]  
where \(\hat{\mathbf{e}}(n)=\mathbf{v}(n)-\hat{\mathbf{v}}_{0}(n)\), \(\hat{\mathbf{v}}_{0}(n)\) provides the estimation of \(\mathbf{v}_{0}(n)\), and \(\sigma_{e}^{2}\) is the variance of \(\mathbf{e}(n)\). \(E(\cdot)\) and \(\|\cdot\|\) represent the expectation and the L2-norm, respectively, \(K\) is the length of the vector \(\mathbf{v}(n)\), and \(\hat{\mathbf{e}}^{T}(n)\) indicates the transpose of \(\hat{\mathbf{e}}(n)\). According to [PERSON]'s lemma [26], the second term has the relationship [18]  
\[E\left(\hat{\mathbf{e}}^{T}(n)\mathbf{e}(n)\right)=\sigma_{e}^{2}E\left\{ \mathrm{tr}\left(\frac{\partial\hat{\mathbf{e}}(n)}{\partial\mathbf{e}^{T}(n) }\right)\right\} \tag{11}\]  
where \(\mathrm{tr}(\cdot)\) represents the trace operator of a matrix, if we assume that \(\mathbf{e}(n)\) has a Gaussian distribution.