I am working with Nonnegative Matrix Factorization, where I factor a nonnegative matrix $\mathbf{X}$ into a basis matrix, $\mathbf{W}$ and a coefficient matrix, $\mathbf{H}$.
$$\mathbf{X}\approx \mathbf{W H}$$
Furthermore, I aim for the columns of $\mathbf{H}$ to be sparse and to calculate it I therefore use an iterative update scheme, Nonnegative Matching Pursuit (NMP), as described in http://doi.org/10.1109/MLSP.2010.5589219.
My problem arises when I try to update $\mathbf{W}$. In the algorithm, I firstly initialize $\mathbf{W}$ with random values (between 0 and 1). $\mathbf{H}$ is as mentioned calculated from NMP. Then the multiplicative update rules as described by D. Lee and H. Seung here is used to update $\mathbf{W}$ and $\mathbf{H}$ sequentially. The update rules are as follows:
$$ \mathbf{W} \leftarrow \mathbf{W} \otimes \frac{\mathbf{X}\mathbf{H}^T}{\mathbf{W}\mathbf{H}\mathbf{H}^T} $$
$$ \mathbf{H} \leftarrow \mathbf{H} \otimes \frac{\mathbf{W}^T\mathbf{X}}{\mathbf{W}^T\mathbf{W}\mathbf{H}} $$ where $\otimes$ and $\frac{\cdots}{\cdots}$ in this case are notations for elementwise multiplcation and division respectively.
The issue is that if I have rows of $\mathbf{H}$ or columns of $\mathbf{W}$ which are entirely zeros, I get division by 0. A row of zeros in $\mathbf{H}$ indicate that there are 'atoms' (columns) of $\mathbf{W}$ which are not used in the reconstruction of any of the signals (columns) in $\mathbf{X}$.
So far, my solution has been to firstly identify the problematic rows of $\mathbf{H}$ and deleting them outright along with the corresponding atoms in $\mathbf{W}$. But is there a more "correct" or elegant solution, where I end up with matrices with the same dimensions as I specified before running the algorithm? Columns of exclusively $0$'s in $\mathbf{W}$ I can simply reinitialize by setting all elements to $\sqrt{D}$ ($D$ being the number of rows in $\mathbf{W}$), but I struggle to find a good solution for the $0$-rows of $\mathbf{H}$.
I hope you can help me with your suggestions.