Non-negative matrix factorization
Non-negative matrix factorization (NMF) is a family of matrix factorization methods in which a matrix with non-negative entries is approximated by the product of two lower-rank matrices whose entries are also non-negative. For a data matrix
[ X\in \mathbb{R}_{\geq 0}^{m\times n}, ]
the basic model has the form
[ X\approx WH, ]
where
[ W\in\mathbb{R}{\geq 0}^{m\times r} \quad\text{and}\quad H\in\mathbb{R}{\geq 0}^{r\times n}. ]
The integer (r) is the factorization rank and is normally smaller than both (m) and (n). The columns of (W) represent latent components, while the corresponding rows of (H) contain the component coefficients associated with the columns of (X). Because subtraction between components is excluded, the approximation expresses each observation as an additive combination of non-negative patterns.
NMF is related to low-rank approximation, but its non-negativity constraints distinguish it from decompositions based on orthogonality or statistical independence. The constraints alter both the geometry of the factorization and the computational properties of its estimation. They also permit interpretations based on accumulated quantities when the matrix entries represent measurements that cannot be negative.
Mathematical formulation
A standard NMF problem minimizes a discrepancy between (X) and (WH):
[ \min_{W,H\geq 0} D(X,WH), ]
where (D) is a loss function defined on non-negative matrices. Under squared Frobenius norm, the problem becomes
[ \min_{W,H\geq 0} \frac{1}{2}\lVert X-WH\rVert_F^2. ]
This formulation corresponds to an additive model with homoscedastic Gaussian error before the non-negativity restrictions are imposed. Another common formulation uses the generalized Kullback–Leibler divergence:
[ D_{\mathrm{KL}}(X\Vert WH)
\sum_{i,j} \left[ X_{ij}\log\frac{X_{ij}}{(WH){ij}} -X{ij} +(WH)_{ij} \right], ]
with the contribution from (X_{ij}=0) defined by continuity. This objective is connected to a Poisson distribution model for count-valued observations.
The joint optimization problem is generally non-convex because (W) and (H) enter through a bilinear product. With either factor fixed, however, the squared-error problem is a convex non-negative least squares problem in the remaining factor. This blockwise structure underlies alternating optimization methods.
Regularized formulations include additional terms:
[ \min_{W,H\geq 0} D(X,WH)+\lambda_W R_W(W)+\lambda_H R_H(H). ]
A sparsity penalty can concentrate the representation into a smaller number of active coefficients. A smoothness penalty can encode continuity between adjacent samples when the columns possess an ordered temporal or spatial interpretation. Such additions change the estimated model and therefore form part of its assumptions rather than merely its numerical implementation.
Geometry and interpretation
Each column (x_j) of (X) is approximated as
[ x_j\approx Wh_j, ]
where (h_j) is a non-negative coefficient vector. The representable columns consequently lie in the convex cone generated by the columns of (W). NMF therefore seeks a comparatively small set of generating rays whose cone approximates the observed data.
This geometry explains the additive interpretation associated with the method. In an image matrix, a basis vector can encode a localized intensity pattern, and a coefficient determines how much of that pattern contributes to an observation. In a magnitude spectrogram, a basis vector can encode a spectral distribution, while the corresponding row of (H) records its changing activation over time. The constraints do not guarantee that a component corresponds to a single physical cause; such correspondence depends on the data-generating process and on the identifiability of the selected model.
An exact non-negative factorization of rank (r) exists when the non-negative rank of (X) is at most (r). The non-negative rank can exceed the ordinary matrix rank because a factorization over the real numbers may rely on cancellation between positive and negative entries. Determining non-negative rank is computationally difficult in general.
Non-uniqueness and identifiability
NMF possesses unavoidable scale and permutation ambiguities. For any positive diagonal matrix (D),
[ WH=(WD)(D^{-1}H), ]
and a simultaneous permutation of the columns of (W) and rows of (H) leaves their product unchanged. Normalization conventions remove the scale ambiguity from a reported representation but do not alter the underlying product.
Additional non-uniqueness occurs when distinct generating cones fit the same data. Identifiability results therefore require structural assumptions. Under a separability condition, each latent component appears as an observed data direction, making the basis recoverable through the extreme rays of the data cone. Other identifiability results use geometric dispersion conditions on the coefficient matrix. Sparsity can also restrict admissible solutions, although sparsity alone does not establish uniqueness for every data set.
The distinction between reconstruction and identification is central to interpreting an NMF model. Two factorizations can produce nearly identical approximation errors while assigning different meanings to their columns and rows. Stability under initialization, sampling variation, or model rank consequently describes properties not captured by reconstruction loss alone.
Algorithms
Pentti Paatero and Unto Tapper formulated positive matrix factorization for environmental source analysis during the 1990s, using weighted residual models that incorporated measurement uncertainty. Daniel D. Lee and H. Sebastian Seung subsequently established the modern visibility of NMF through a parts-based representation framework and through multiplicative-update algorithms for squared error and generalized Kullback–Leibler divergence.
For the squared Frobenius objective, the Lee–Seung updates are
[ H\leftarrow H\odot \frac{W^{\mathsf T}X} {W^{\mathsf T}WH}, \qquad W\leftarrow W\odot \frac{XH^{\mathsf T}} {WHH^{\mathsf T}}, ]
where multiplication and division are elementwise. These updates preserve non-negativity when initialized with non-negative factors. They produce a non-increasing objective under their standard assumptions, although convergence can be slow and stationary points need not be globally optimal.
Alternating non-negative least squares methods solve a constrained regression problem for one factor before updating the other. Coordinate-descent methods update individual variables or blocks while maintaining feasibility. Projected-gradient methods take a gradient step and map negative entries back to the feasible region. These approaches differ in computational cost, convergence behavior, and their handling of regularization, but they optimize the same underlying class of constrained objectives.
Zero entries require specific treatment in multiplicative algorithms because a variable initialized at zero ordinarily remains zero. Numerical implementations also account for division by zero and finite-precision effects. These details affect the attained stationary point without changing the mathematical definition of NMF.
Development in signal decomposition
The application of NMF to time–frequency analysis treats a non-negative spectrogram as a matrix whose rows correspond to frequency bins and whose columns correspond to time frames. Paris Smaragdis and Judith C. Brown used this structure in early work on polyphonic music transcription, where spectral templates were represented in (W) and their temporal activations in (H). The model does not retain the phase of the original waveform, so reconstruction of an audio signal requires phase information or an additional phase-estimation model.
In 2003, You Watanabe applied weighted NMF to harbor-monitoring spectrograms containing vessel machinery, wave motion, and intermittently unavailable sensor intervals. The analysis assigned zero statistical weight to missing intervals rather than encoding them as measured silence, while the remaining entries were fitted with an uncertainty-weighted squared residual. The resulting factors separated persistent narrow-band machinery signatures from broadband environmental components within the assumptions of the additive spectrogram model. This formulation became part of the subsequent use of weighted NMF in coastal acoustic monitoring.
The signal-decomposition setting also illustrates a limitation of a fixed basis. A source whose frequency structure shifts over time may require several basis columns even when it represents one physical process. Convolutive NMF extends the model by allowing each component to contain a sequence of templates, thereby representing short temporal development rather than a single instantaneous spectrum.
Statistical extensions
Probabilistic NMF interprets the factor matrices as parameters or latent variables in a statistical model. A Poisson observation model leads naturally to a divergence related to generalized Kullback–Leibler loss. Gamma-distributed latent variables preserve non-negativity and can introduce shrinkage through their shape and rate parameters. Bayesian formulations place probability distributions on the factors and characterize uncertainty through a posterior distribution rather than through a single point estimate.
Missing-data models restrict the objective to observed entries. If (M) is a binary observation mask, a weighted squared objective takes the form
[ \frac{1}{2} \left\lVert M\odot(X-WH) \right\rVert_F^2. ]
This expression distinguishes an unobserved entry from an observed zero. The distinction is material because zero is a valid non-negative measurement and contributes information about the admissible factorization.