Multiplicative Updates as Preconditioning

I started thinking about multiplicative updates after reading about ScaledGD. Both algorithms scale their gradients using the current factors. Rewriting MU in this form reveals how its diagonal scaling responds to coordinate units and component correlation.

Seeing the gradient inside the multiplicative update

Scaled gradient descent (ScaledGD) makes a familiar promise in low-rank recovery. By adapting the gradient to the current factors, it can achieve convergence rates independent of the condition number of the target matrix under appropriate observation and initialization assumptions Tong, Ma, and Chi (2021). This invites a question about an older factorization algorithm. Multiplicative updates also adapt each entry to the current factors. Could they already provide a similar form of preconditioning? We can investigate this by rewriting the update itself, then comparing the matrix that multiplies its gradient with the one used by ScaledGD.

Consider a nonnegative data matrix \(X\in\mathbb{R}_+^{m\times n}\) and factors \(W\in\mathbb{R}_+^{m\times r}\), \(H\in\mathbb{R}_+^{r\times n}\). We minimize \(f(W,H)=\tfrac{1}{2}\|WH-X\|_F^2\). Its gradients are

\[ G_W=\nabla_W f=WHH^\top-XH^\top, \qquad G_H=\nabla_H f=W^\top WH-W^\top X. \]

Ordinary gradient descent subtracts \(\eta G_W\) from \(W\). MU is usually written as \(W^+=W\odot XH^\top/(WHH^\top)\), where \(\odot\) and the quotient act entrywise. To expose its gradient, subtract the current iterate and combine the numerator. For positive entries and denominators, the exact algebra is

\[ \begin{aligned} W^+-W &=W\odot\left(\frac{XH^\top}{WHH^\top}-1\right)\\ &=W\odot\frac{XH^\top-WHH^\top}{WHH^\top}\\ &=-\underbrace{\frac{W}{WHH^\top}}_{\text{entrywise gradient scaling}}\odot G_W. \end{aligned} \]

Thus each entry of the gradient receives its own step size, \(W_{ij}/(WHH^\top)_{ij}\). The same calculation gives an \(H\)-block displacement \(-H/(W^\top WH)\odot G_H\). Classical MU performs these two block updates in sequence, evaluating the second expression with the updated \(W^+\). This rescaled-gradient interpretation already appears in Lee and Seung’s original analysis Lee and Seung (2001).

We can now put MU beside the ScaledGD directions for this objective. The following expressions compare directions evaluated at the same factor pair.

\[ \begin{array}{c|cc} &\Delta W&\Delta H\\ \hline \mathrm{MU}&-\dfrac{W}{WHH^\top}\odot G_W &-\dfrac{H}{W^\top WH}\odot G_H\\[6pt] \mathrm{ScaledGD}&-G_W(HH^\top)^{-1} &-(W^\top W)^{-1}G_H \end{array} \]

Both rules adapt the gradient using the factors. MU rescales individual entries. ScaledGD multiplies by an inverse Gram matrix, allowing the update of one component to use gradient information from other components. To understand how much conditioning each operation removes, we need to identify the curvature these two scalings act on. Holding one factor fixed makes that curvature explicit.

What the two preconditioners do to a block

For the squared reconstruction loss, the Gram scaling used by ScaledGD has a direct interpretation through the curvature of each factor update.

Hold \(H\) fixed. The squared Frobenius loss separates across rows of \(W\), so each row can be analyzed independently. Write \(w=W_{k,:}^\top\) and \(x=X_{k,:}^\top\) for one row, and abbreviate \(C=HH^\top\) and \(b=Hx\). Its objective is

\[ \phi(w)=\frac{1}{2}\|H^\top w-x\|_2^2 =\frac{1}{2}w^\top Cw-b^\top w+\frac{1}{2}x^\top x. \]

Because the residual is linear in \(w\), this is a quadratic problem with gradient \(g(w)=Cw-b\) and Hessian \(\nabla^2\phi(w)=C\). The Gram-scaled matrix update is \(W^+=W-\eta G_WC^{-1}\). Right multiplication acts independently on each row. Transposing row \(k\), using the symmetry of \(C^{-1}\), gives

\[ w^+=(W_{k,:}-\eta G_{W,k,:}C^{-1})^\top =w-\eta C^{-1}g(w). \]

For the block comparison set \(\eta=1\) and assume \(C\succ0\). Introduce the two gradient multipliers with matching notation, \(P_{\rm MU}(w)=\operatorname{diag}(w_i/(Cw)_i)\) and \(P_{\rm Gram}=C^{-1}\). The complete updates are

\[ \begin{aligned} w^+_{\rm MU}&=w-P_{\rm MU}(w)g(w),\\ w^+_{\rm Gram}&=w-P_{\rm Gram}g(w). \end{aligned} \]

The Gram multiplier is the inverse Hessian of this row problem. Its unit step therefore coincides with a Newton step and reaches the unconstrained least-squares solution,

\[ w^+_{\rm Gram}=w-[\nabla^2\phi(w)]^{-1}g(w) =w-C^{-1}(Cw-b)=C^{-1}b. \]

When this solution is nonnegative it also solves the constrained row problem. A simultaneous ScaledGD iteration computes both factor updates from the current pair.

What curvature approximation does MU invert? Define \(D(w)=P_{\rm MU}(w)^{-1}=\operatorname{diag}((Cw)_i/w_i)\). This puts MU in the same inverse-curvature form as the Gram update,

\[ \begin{aligned} w^+_{\rm MU}&=w-D(w)^{-1}g(w),\\ w^+_{\rm Gram}&=w-C^{-1}g(w). \end{aligned} \]

For positive \(w\), the entrywise nonnegativity of \(C\) gives

\[ v^\top(D(w)-C)v =\frac{1}{2}\sum_{i,j}C_{ij}w_iw_j \left(\frac{v_i}{w_i}-\frac{v_j}{w_j}\right)^2\geq0. \]

Thus \(D(w)\succeq C\). MU inverts a diagonal upper bound on the true Hessian. This also connects the update directly to majorization-minimization. Holding \(D(w)\) fixed, replace \(C\) in the exact quadratic expansion by its upper bound to obtain

\[ \phi(w+v)\leq Q(v\mid w) :=\phi(w)+g(w)^\top v+\frac{1}{2}v^\top D(w)v. \]

The minimizer satisfies \(g(w)+D(w)v_*=0\). Substituting it gives the exact equivalence

\[ w+\arg\min_v Q(v\mid w) =w-D(w)^{-1}g(w) =w\odot\frac{b}{Cw}=w^+_{\rm MU}. \]

The last expression is nonnegative, and the upper model guarantees descent through \(\phi(w^+_{\rm MU})\leq Q(v_*\mid w)\leq Q(0\mid w)=\phi(w)\). The preconditioner in the gradient update is precisely the inverse curvature of the MM upper model.

To measure what preconditioning achieves, track the distance to the unconstrained block solution \(w_*=C^{-1}b\). Write the current error as \(e=w-w_*\) and the error after one update as \(e^+=w^+-w_*\). The superscript + means “after the update.” Since \(Cw_*=b\), the gradient is \(g(w)=C(w-w_*)=Ce\). Substituting into the unit-step update gives

\[ \begin{aligned} e^+&=w^+-w_*\\ &=w-Pg(w)-w_*\\ &=e-PCe=(I-PC)e. \end{aligned} \]

Thus the product PC determines how much error remains after a step. To read its effect in orthogonal directions, express the error in coordinates scaled by the fixed positive definite preconditioner P. Set \(z=P^{-1/2}e\) and \(z^+=P^{-1/2}e^+\). The same update becomes

\[ z^+=\left(I-\underbrace{P^{1/2}CP^{1/2}}_{B}\right)z. \]

The symmetric matrix B is the curvature in these scaled coordinates. Along an eigenvector with eigenvalue \(\lambda\), the corresponding error component is multiplied by \(1-\lambda\). An eigenvalue of one removes that component in one step; an eigenvalue close to zero leaves most of it. Substituting the two preconditioners gives

\[ \begin{aligned} B_{\rm Gram}&=C^{-1/2}CC^{-1/2}=I,\\ B_{\rm MU}(w)&=D(w)^{-1/2}CD(w)^{-1/2}. \end{aligned} \]

The Gram inverse makes every eigenvalue one. The bound \(D(w)\succeq C\) gives \(B_{\rm MU}(w)\preceq I\), while its smallest eigenvalue can still approach zero. MU changes its multiplier at every iteration, so the next section evaluates it at a positive optimum to derive the local contraction factor. The square identity explains how additional curvature enters MU’s model when correlated components move by different relative amounts. It vanishes for a common relative change \(v=aw\). We next examine component rescaling and correlation to determine which conditioning effects remain.

Which conditioning does MU remove?

A factorization has a built-in freedom of scale. For any positive diagonal matrix \(S\), the factors \(\widetilde{W}=WS\) and \(\widetilde{H}=S^{-1}H\) give exactly the same product. Scaling one component can make a Gram matrix badly conditioned even though the reconstruction is unchanged. MU absorbs this change of representation. Under the transformation, both \(XH^\top\) and \(WHH^\top\) acquire the right factor \(S^{-1}\). Their entrywise ratio is unchanged, so the updated factor satisfies \(\widetilde{W}^+=W^+S\). The following \(H\) update transforms in the same way. Starting from rescaled factors therefore produces identical reconstruction trajectories in exact arithmetic, even when the component scales differ greatly.

The corresponding property of the Gram direction extends to invertible changes of basis. Substituting \(W\mapsto WT\) and \(H\mapsto T^{-1}H\) into the Gram inverses gives \(\Delta W\mapsto\Delta WT\) and \(\Delta H\mapsto T^{-1}\Delta H\). Such transformations can mix components as well as rescale them. For NMF, feasible transformed factors must remain nonnegative.

To see what remains after MU’s scaling, return to the fixed-row problem and a positive optimum \(w_*\). There \(g(w_*)=0\). Differentiating \(w-D(w)^{-1}g(w)\) therefore leaves

\[ \begin{aligned} D_*&=D(w_*),\\ J_{\rm MU}&=I-D_*^{-1}C,\\ B&=D_*^{-1/2}CD_*^{-1/2}. \end{aligned} \]

The derivative of the changing preconditioner multiplies the zero gradient. This is why its value at the optimum determines the local linearization. The upper-bound property gives \(0\prec B\preceq I\). Also, \(D_*w_*=Cw_*\) implies that \(D_*^{1/2}w_*\) is an eigenvector of \(B\) with eigenvalue one. Hence the largest eigenvalue is exactly one, and the local spectral factor becomes

\[ q_{\rm MU}=1-\lambda_{\min}(B) =1-\frac{1}{\kappa(B)}. \]

This formula identifies the condition number relevant to MU’s fixed-block dynamics. MU removes a family of ill-conditioning when \(\kappa(B)\) stays bounded as the original problem becomes ill-conditioned. If \(C\) is diagonal, then \(D(w)=C\) and \(B=I\), however different the diagonal entries are. In that case MU and the unit Gram step coincide. More generally, a bound \(C\succeq\alpha D_*\) with a uniform \(\alpha>0\) yields a uniform contraction factor at most \(1-\alpha\).

Correlated components give a contrasting family. Choose

\[ C=\begin{pmatrix}1&c\\c&1\end{pmatrix},\qquad w_*=\begin{pmatrix}1\\1\end{pmatrix}. \]

Here \(0\leq c<1\) and \(b=Cw_*\). This quadratic comes from a nonnegative least-squares problem. We can take the rows of \(H\) to be \((1,0)\) and \((c,\sqrt{1-c^2})\), and set \(x=H^\top w_*\). The diagonal entries of \(C\) already agree, while its off-diagonal correlation approaches one. At the optimum,

\[ \begin{aligned} D_*&=(1+c)I,\\ \operatorname{spec}(B)&=\left\{1,\frac{1-c}{1+c}\right\},\\ \kappa(B)&=\kappa(C)=\frac{1+c}{1-c}. \end{aligned} \]

MU leaves the entire block condition number in this family. The full Gram inverse maps both curvature eigenvalues to one. The off-diagonal information determines their difference here.

The scale and correlation controls have a general matrix interpretation. Let \(S=\operatorname{diag}(\sqrt{C_{ii}})\) and \(R=S^{-1}CS^{-1}\). Then \(C=SRS\) and \(R_{ii}=1\). The diagonal matrix \(S\) describes component magnitudes, while \(R\) contains their normalized inner products. Their effects interact through a congruence transformation. In particular,

\[ \kappa(C)\leq\kappa(S)^2\kappa(R). \]

This bound relates coordinate scales and normalized correlations to the Gram condition number. MU is equivariant to changes of \(S\), whereas the Gram inverse also normalizes the correlations within this fixed block. MU’s remaining conditioning additionally depends on the relative coefficients at the optimum.

This is a finite-iteration guarantee throughout an explicit neighborhood. For any family with \(\kappa_{\rm MU}\leq K\), it gives a common geometric rate even if \(\kappa(C)\) diverges through coordinate rescaling. Both the initialization condition and \(\kappa_{\rm MU}\) are invariant under positive diagonal changes of coordinates.

Proof of Theorem 1

Write the error as \(e=w-w_*\). The coordinate bound is

\[ |e_i|\leq\sqrt{(C^{-1})_{ii}}\,\|e\|_C. \]

The initialization condition therefore places the initial point in a relative neighborhood of the optimum,

\[ (1-\delta)w_*\leq w_0\leq(1+\delta)w_*. \]

MU decreases the quadratic objective, whose gap equals half the squared error in the C norm. The same energy bound and relative neighborhood consequently hold at every iterate. Entrywise nonnegativity of C bounds both the numerator and denominator of each diagonal entry of D. Writing \(\beta=(1+\delta)/(1-\delta)\) gives

\[ D(w)\preceq\beta D_*. \]

The normalized curvature at the optimum has largest eigenvalue one and smallest eigenvalue equal to the reciprocal of the MU condition number. Hence

\[ C\succeq\frac{D_*}{\kappa_{\rm MU}}, \qquad D(w)^{-1}\succeq\frac{C^{-1}}{\beta\kappa_{\rm MU}}. \]

Set the objective gap to \(E(w)=\phi(w)-\phi(w_*)\). The MM decrease and the preceding inverse bound imply

\[ \begin{aligned} E(w^+)&\leq E(w)-\tfrac{1}{2}g(w)^\top D(w)^{-1}g(w)\\ &\leq E(w)-\frac{g(w)^\top C^{-1}g(w)}{2\beta\kappa_{\rm MU}}\\ &=\left(1-\frac1{\beta\kappa_{\rm MU}}\right)E(w). \end{aligned} \]

The last equality uses the quadratic identity \(g(w)^\top C^{-1}g(w)=2E(w)\). Iterating this bound proves the theorem. ∎

The explorer extends this comparison to three components. Let \(T=\operatorname{diag}(s,1,1)\) and let R have diagonal entries one and off-diagonal entries c. We use \(C=T^{-1}RT^{-1}\) and \(w_*=T\mathbf{1}\). Then \(D_*=(1+2c)T^{-2}\), so the MU-normalized curvature is \(R/(1+2c)\). Its condition number is \((1+2c)/(1-c)\), independent of the coordinate scale s. Gram scaling gives the identity.

Each shape is a three-dimensional equal-energy ellipsoid for the corresponding curvature matrix. Increasing scale stretches the original geometry. Increasing correlation also distorts the MU geometry, while the Gram geometry remains spherical. The three 3D views rotate together. The section controls rotate a central plane through all three ellipsoids. Its highlighted intersection forms the outer contour in the 2D view above. Nested equal-energy contours show the curvature within the plane, with faint lines marking the long and short axes. The displayed axis ratio measures conditioning within that section; its square equals the condition number of the curvature restricted to the plane. The larger condition-number display describes the full three-dimensional curvature. Each ellipsoid is uniformly resized to fit its panel, preserving its axis ratios. These shapes use MU’s scaling at the optimum.

This comparison also clarifies which condition number is being discussed. The Gram condition number \(\kappa(C)\) depends on the factors’ coordinates. The remaining condition number \(\kappa(B)\) describes MU’s local block update. ScaledGD’s recovery results concern the target matrix condition number \(\kappa_X=\sigma_1(W_*H_*)/\sigma_r(W_*H_*)\), with observation and initialization assumptions that connect the factor geometry to recovery Tong, Ma, and Chi (2021).