🤖 machine learning

An Efficient Newton Algorithm for Nonnegative Matrix Factorization with the Kullback-Leibler Divergence

This paper proposes a new efficient Newton-type algorithm for Kullback-Leibler Nonnegative Matrix Factorization that utilizes a second-order Taylor expansion and a generalized HALS approach to overcome the limitations of existing separable majorant methods, achieving provable convergence and competitive performance across diverse datasets.

Original authors: Damien Lesens, Jérémy E. Cohen, Bora Uçar

Published 2026-07-16

📖 4 min read☕ Coffee break read

CC BY 4.0

Original authors: Damien Lesens, Jérémy E. Cohen, Bora Uçar

Original paper licensed under CC BY 4.0 (http://creativecommons.org/licenses/by/4.0/). ✨ This is an AI-generated explanation of the paper below. It is not written or endorsed by the authors. For technical accuracy, refer to the original paper. Read full disclaimer

Imagine you are trying to solve a giant jigsaw puzzle, but with a twist: you don't have the picture on the box, and you can't see the pieces clearly. All you have is a blurry, messy pile of data. In the world of computer science, this is called Nonnegative Matrix Factorization (NMF). It's a tool used to take a big, complicated table of numbers (like a spreadsheet of song lyrics or a photo made of light pixels) and break it down into two smaller, simpler tables that, when multiplied together, recreate the original picture. The "nonnegative" part just means all the numbers have to be zero or positive—no negative numbers allowed, because you can't have "minus three" apples or "minus five" words in a sentence.

But here's the tricky part: how do you know your simplified tables are a good fit? If the data you are looking at comes from counting things—like how many times a word appears in a book, or how many photons hit a camera sensor—the math gets a little weird. The errors aren't like the smooth, bell-shaped curves of a standard math class; they are more like the jittery, unpredictable nature of raindrops hitting a roof. To measure the fit in these cases, scientists use a special ruler called the Kullback-Leibler (KL) divergence. Think of it as a "surprise meter." If your model predicts a word will appear 10 times, but it actually appears 100 times, the surprise meter goes off the charts. The goal is to find the two small tables that make this surprise meter read as low as possible.

For a long time, the best way to solve this puzzle was to take tiny, cautious steps, checking the surprise meter after every single move. This method, known as "Multiplicative Updates," has been the champion for years. But what if there was a way to take a giant leap, looking ahead to see where the path leads, rather than just shuffling your feet? That is exactly what this paper explores.

The authors, Damien Lesens, Jérémy E. Cohen, and Bora Uçar, argue that the old "tiny step" method has hit a wall. They propose a new, bolder strategy: a Newton-type algorithm. In the world of math, a Newton method is like a hiker who doesn't just look at the ground under their feet, but looks at the shape of the entire hill to decide the best direction to run. Instead of just looking at the slope (the first derivative), this new method looks at the curvature (the second derivative) to predict exactly where the bottom of the valley is.

However, there's a catch. The math for this "big leap" is incredibly complex and doesn't play nicely with the rule that all numbers must be positive. Most attempts to use this powerful tool in the past have been too slow or too messy to be useful. The authors' main breakthrough is showing how to tame this complex math. They invented a new way to solve the problem efficiently by adapting an existing technique called HALS (Hierarchical Alternating Least Squares). They essentially created a "generalized" version of this tool that can handle the heavy lifting of the second-order math without getting bogged down.

The result is an algorithm they call KL-HALS. In their tests, this new method proved to be a powerhouse on audio recordings and synthetic data, often finding better solutions faster than the current state-of-the-art methods. However, the results were more nuanced on other data types. On image datasets, the new method was actually the second-best, trailing behind a simpler algorithm that uses a different type of math (Frobenius norm), and on large document datasets with high complexity, it sometimes converged more slowly than older methods. This suggests that while the "big leap" strategy is powerful, the terrain of the data matters; sometimes the old "tiny steps" are still the most efficient path.

Interestingly, the authors also proved mathematically that the old "tiny step" method (Multiplicative Updates) is actually the best possible version of that specific type of cautious approach. This means that to get any faster, you must stop being cautious and start using the "big leap" strategy they developed, even if it requires more computing power per step. They also found that starting the process with a smart "warm-up" (scaling the initial numbers correctly) helps the algorithm find its footing much faster. In short, this paper doesn't just offer a slightly better tool; it suggests a fundamental shift in how we should approach this specific type of data puzzle, proving that sometimes, taking a calculated giant leap is better than a million tiny shuffles—provided you are on the right kind of terrain.

Technical Summary: An Efficient Newton Algorithm for Nonnegative Matrix Factorization with the Kullback-Leibler Divergence

Problem Statement
Nonnegative Matrix Factorization (NMF) approximates a nonnegative matrix $V$ as the product of two low-rank nonnegative factors, $W$ and $H$ ( $V \approx WH$ ). When data follows a Poisson distribution (common in count datasets like term-document matrices or photon-counting images), the Kullback-Leibler (KL) divergence is the appropriate loss function. The resulting optimization problem, (NMF-KL), is non-convex jointly in $W$ and $H$ but convex in each separately.

Existing state-of-the-art algorithms for (NMF-KL) predominantly rely on the Majorization-Minimization (MM) framework. These methods minimize a separable upper-bound (majorant) of the loss function to derive update rules, such as the Multiplicative Updates (MU) algorithm. While efficient ($O(MNR)$ per iteration), the authors argue that this approach has reached its theoretical limits regarding the tightness of the surrogate function.

Methodology
The proposed approach, KL-HALS, departs from the separable surrogate paradigm by utilizing the second-order Taylor expansion of the loss function. This results in a non-separable quadratic surrogate. The core challenges addressed are:

Non-negativity Constraints: Handling constraints without resorting to expensive projected Newton steps or active-set methods that lack general convergence guarantees.
Efficient Minimization: Minimizing the non-separable quadratic surrogate efficiently.
Step Size Selection: Choosing step sizes that ensure convergence without relying on Lipschitz smoothness assumptions (which do not hold for KL divergence).

The authors leverage the framework of self-concordant functions (Tran-Dinh et al. [22]) to address step size selection and convergence guarantees. They prove that the (NNKL) subproblem is self-concordant, allowing for the use of damped Newton steps with specific step sizes ( $\alpha_k = 1/(1+\lambda_k)$ ) that guarantee a decrease in the loss.

To minimize the resulting non-separable quadratic surrogate under non-negativity constraints, the paper proposes a generalization of the Hierarchical Alternating Least Squares (HALS) algorithm, termed GHALS.

GHALS: Instead of solving for a single scalar entry or a separable sum, GHALS solves a sum of multivariate quadratics with different Hessian matrices for each column of $H$ . It updates rows of $H$ cyclically using a closed-form expression derived from the specific structure of the quadratic terms.
Algorithm Structure: The full algorithm alternates between updating $H$ and $W$ . For each factor update, it computes first-order ( $A$ ) and second-order ( $T$ ) coefficients, solves the surrogate problem using GHALS, and applies a step size.
Complexity: The per-iteration complexity is $O(MNR^2)$ , dominated by the computation of the tensor $T$ (the Hessian approximation). This is higher than the $O(MNR)$ complexity of separable MM methods, but the authors argue the trade-off is favorable due to faster convergence.
Sparse Implementation: For large sparse datasets, the authors implement a sparse version that avoids forming the dense $WH$ matrix, reducing complexity to depend on the number of non-zero entries ($nnz(V)$).

Key Contributions

Theoretical Limit of Separable Majorants: The authors prove that the separable upper-bound used by the MU algorithm is the tightest possible global separable majorant for the NNKL loss. This theoretically justifies moving away from separable surrogates to achieve better performance.
Sub-optimality of MU Burg: They demonstrate theoretically that the MU Burg algorithm (based on Bauschke et al. [17]) uses a looser upper-bound than MU, explaining its empirically observed inferior performance.
Efficient Non-Separable Minimization: The paper introduces GHALS, a generalization of HALS that efficiently minimizes non-separable quadratic surrogates under non-negativity constraints, a task previously not solved efficiently in the NMF-KL literature.
Convergence Guarantees: By utilizing the self-concordance properties of the KL loss, the authors provide a version of the algorithm (KL-HALS-descent) with provable convergence guarantees using damped Newton steps.
Initialization Strategy: The paper proposes using Sinkhorn's algorithm to optimally scale the rows and columns of the initial factors, ensuring the starting point is close to the basin of quadratic convergence.

Experimental Results
The authors evaluated KL-HALS and its damped variant against state-of-the-art algorithms (MU, SN, CCD, AmSOM, AMUSOM, FPA, HALS) on a diverse set of datasets:

Audio: On the MAPS dataset, KL-HALS converged to the best loss but was slower than MU for high ranks due to $R^2$ complexity. On a second audio excerpt, KL-HALS outperformed all others.
Images: On face and hyperspectral datasets, HALS (optimized for Frobenius loss) surprisingly performed best, suggesting Gaussian and Poisson models may be equivalent for these specific data types. KL-HALS was the second-best performer, consistently outperforming MU and CCD.
Documents: On the CLUTO and Verb datasets (large, sparse count data), KL-HALS achieved the best final loss values across various relative errors, significantly outperforming MU and CCD. However, its convergence speed degraded as the rank $R$ increased, consistent with its higher computational complexity.
Synthetic Data: On both low-rank and full-rank synthetic datasets, KL-HALS consistently converged to the best loss values, often outperforming other methods, particularly in noisy settings.

Significance and Claims
The paper claims that exploring non-separable surrogates is essential for further progress in NMF-KL, as the space for improvement on separable upper-bounds is exhausted. The proposed method demonstrates that minimizing a second-order Taylor expansion under non-negativity constraints is feasible and effective.

The authors modestly note that while their $O(MNR^2)$ method is competitive and often superior in final accuracy, the computational cost remains a barrier for very high-rank problems compared to $O(MNR)$ first-order methods. They suggest that future work could focus on approximating the Hessian (quasi-Newton approaches) or exploring extrapolated updates to accelerate convergence, potentially bridging the gap between the accuracy of second-order methods and the speed of first-order algorithms. The work establishes that second-order information, when handled correctly via self-concordance and specialized solvers like GHALS, yields high-quality iterates that justify the additional computational investment.

Technical Summary: An Efficient Newton Algorithm for Nonnegative Matrix Factorization with the Kullback-Leibler Divergence

More like this