# Sparse Envelope Spectra for Feature Extraction of Bearing Faults Based on NMF
## 2 Materials and Methods  
### Nmf  
With the non-negative constraints, NMF is also widely used in image recognition, data mining etc. [23; 24; 25]. Given a non-negative matrix \(\mathbf{V}=\left[\begin{array}{cccc}\mathbf{v_{1}}&\mathbf{v_{2}}&\cdots&\mathbf{v_{n}}\\ \end{array}\right]\in\mathbf{R}^{m\times n}\), where \(m\) and \(n\) are the numbers of samples and features respectively, NMF is designed to decompose the matrix \(\mathbf{V}\) into the product of two non-negative matrices \(\mathbf{W}=\left[\begin{array}{cccc}\mathbf{w_{1}}&\mathbf{w_{2}}&\cdots&\mathbf{w_{r}}\\ \end{array}\right]\in\mathbf{R}^{m\times r}\) and \(\mathbf{H}=\left[\begin{array}{cccc}\mathbf{h_{1}}&\mathbf{h_{2}}&\cdots&\mathbf{h_{r}}\\ \end{array}\right]^{T}\in\mathbf{R}^{r\times n}\), such that Equation (1) holds:  
\[\mathbf{V_{n\times m}}\approx\mathbf{W_{n\times r}}\mathbf{H_{r\times m}}\] \[\mathbf{s.t.}\ \mathbf{W}\geq 0,\ \mathbf{H}\geq 0 \tag{1}\]where \(W\) and \(H\) are the basis matrix and the weight matrix, respectively, and \(r\) is the low rank of the embedding space. Normally, the choice of \(r\) is to satisfy \((m+n)<mn\). \(W\) can be regarded as a set of bases that linearly approximate the observed array \(V\), and \(H\) is the non-negative projection weight of the observed array \(V\) on the basis \(W\).  
Due to the purely additive expression, NMF provides a good explanation of the local characteristics of data. In addition, from the perspective of clustering, the decomposition result of NMF is clustering for the non-negative data [26]. [PERSON] [27] proved the consistency among NMF, K-means and spectral clustering methods. We can interpret the columns of \(W\) as the posterior probabilities for row clustering and the rows of \(H\) as the posterior probabilities for column clustering.  
For the solution of Equation (1), the objective function based on Euclidean distance is written as:  
\[\begin{split}&\mathbf{min}:f_{NMF}(\mathbf{V}||\mathbf{W}\mathbf{ H})=\frac{1}{2}||\mathbf{V}-\mathbf{W}\mathbf{H}||_{F}^{2}=\frac{1}{2}\sum \limits_{ij}\left(\mathbf{V}_{ij}-(\mathbf{W}\mathbf{H})_{ij}\right)^{2}\\ &\textbf{s.t.}\ \ \mathbf{W}\geq 0,\mathbf{H}\geq 0\end{split} \tag{2}\]  
where \(||\bullet||_{F}^{2}\) denotes the Frobenius norm. Using the gradient descent method for this optimization problem, the classical multiplicative iterative rules are given by:  
\[\mathbf{H}\leftarrow\mathbf{H}\frac{\mathbf{W}^{T}\mathbf{V}}{\mathbf{W}^{T} \mathbf{W}\mathbf{H}}, \tag{3}\]  
\[\mathbf{W}\leftarrow\mathbf{W}\frac{\mathbf{V}\mathbf{H}^{T}}{\mathbf{W} \mathbf{H}\mathbf{H}^{T}}. \tag{4}\]  
The multiplicative iterative algorithm is simple and easy to use, but the poor performance and convergence are also pointed out in paper [28; 29]. Therefore, improved algorithms are proposed, such as the alternating least squares algorithm [30], the projected quasi-Newton algorithm [31], the projected gradient algorithm [28], and the block principal pivoting algorithm [32]. The block principal pivoting algorithm is developed with the non-negative least squares (NNLS) framework.  
\[\left\{\begin{array}{l}\min\limits_{\mathbf{W}\geq 0}\left\|\mathbf{H}^{T} \mathbf{W}^{T}-\mathbf{V}^{T}\right\|_{F}^{2}\text{, by fixing }\mathbf{H}\\ \min\limits_{\mathbf{H}\geq 0}\left\|\mathbf{W}\mathbf{H}-\mathbf{V} \right\|_{F}^{2}\text{, by fixing }\mathbf{W}\end{array}\right. \tag{5}\]  
To obtain good convergence, the main iterative rules are as follows:  
\[\mathbf{H}\leftarrow\left(\mathbf{W}^{T}\mathbf{W}\right)^{-1}\mathbf{W}^{T} \mathbf{V}, \tag{6}\]  
\[\mathbf{W}\leftarrow\mathbf{V}\mathbf{H}^{T}\left(\mathbf{H}\mathbf{H}^{T} \right)^{-1}. \tag{7}\]  
Meanwhile, block principal rotation [33] is used to improve the convergence of the algorithm. Hence, the block principal pivoting algorithm is chosen for analyzing complex vibration signals in this paper.