Iterative Hessian Mixing

What the solver mixes

Write the data as \(\mathbf X\in\mathbb R^{n\times d}\), with one feature vector \(\mathbf x_i^\top\) per person, and responses \(\mathbf y\in\mathbb R^n\). Least squares minimizes

\[ L(\boldsymbol\theta)=\tfrac12\|\mathbf y-\mathbf X\boldsymbol\theta\|_2^2. \]

Its curvature is the Gram matrix \(\mathbf X^\top\mathbf X\). An iterative Hessian sketch approximates that matrix while retaining the exact gradient [PW, JMLR ’16]. Draw a matrix \(\mathbf S\in\mathbb R^{k\times n}\) of independent standard Gaussians, with \(k<n\), and form the sketch \(\mathbf S\mathbf X\). The approximation has the correct mean:

\[ \mathbb E\!\left[\frac1k(\mathbf S\mathbf X)^\top(\mathbf S\mathbf X)\right] =\mathbf X^\top\mathbf X, \]

because \(\mathbb E[\mathbf S^\top\mathbf S]=k\mathbf I_n\). The iterative solver uses fresh sketches to compute its successive updates.

For privacy, imagine releasing the sketch itself. The mixing weights \(\mathbf S\) stay internal; the released object is only \(\mathbf S\mathbf X\). Each of its rows is a Gaussian vector:

\[ \mathbf X^\top\mathbf s\sim\mathcal N(\mathbf0,\boldsymbol\Sigma), \qquad \boldsymbol\Sigma=\mathbf X^\top\mathbf X =\sum_{i=1}^n\mathbf x_i\mathbf x_i^\top. \]

So a \(k\)-row sketch is \(k\) independent samples whose covariance contains the data. Blanking out one person’s row \(\mathbf x^\top\) changes that covariance to

\[ \boldsymbol\Sigma'=\boldsymbol\Sigma-\mathbf x\mathbf x^\top. \]

The privacy question has become a question about distinguishing two Gaussians.

One person’s share of the covariance

Assume for now that \(\boldsymbol\Sigma\) is positive definite. Rescale the coordinates by \(\boldsymbol\Sigma^{-1/2}\), which makes the reference covariance the identity. Writing \(\mathbf w=\boldsymbol\Sigma^{-1/2}\mathbf x\), the two covariances become

\[ \mathbf I_d \quad\text{and}\quad \mathbf I_d-\mathbf w\mathbf w^\top. \]

Every direction perpendicular to \(\mathbf w\) is unchanged. Along \(\mathbf w\), the variance falls from \(1\) to \(1-u\), where

\[ \boxed{u=\|\mathbf w\|_2^2 =\mathbf x^\top\boldsymbol\Sigma^{-1}\mathbf x.} \]

This is the row’s leverage score. It measures the person’s contribution relative to the variation already present in their direction. Ten identical one-dimensional rows each have share \(1/10\); a row that supplies the only variation in some direction has share \(1\). In the latter case, removing the row eliminates a direction entirely, so the sketch can reveal its presence.

For \(0\le u<1\), the remaining calculation is one-dimensional. Let \(p\) be the density of \(\mathcal N(0,1)\) and \(q\) that of \(\mathcal N(0,1-u)\). The KL divergence averages the log likelihood ratio under \(q\):

\[ D_{\mathrm{KL}}(q\|p) =\mathbb E_{z\sim q}\!\left[\log\frac{q(z)}{p(z)}\right]. \]

Substituting the two Gaussian densities gives

\[ \log\frac{q(z)}{p(z)} =-\frac12\log(1-u) -\frac{z^2}{2}\left(\frac1{1-u}-1\right). \]

Since \(\mathbb E_q[z^2]=1-u\), taking the expectation leaves

\[ \boxed{ D_{\mathrm{KL}}(q\|p) =\frac{-u-\log(1-u)}2. } \]

The invertible rescaling preserves likelihood ratios, and the unchanged Gaussian coordinates contribute zero to the divergence. So this is also the divergence between the original Gaussian vectors with covariances \(\boldsymbol\Sigma-\mathbf x\mathbf x^\top\) and \(\boldsymbol\Sigma\).

Expanding the logarithm gives:

\[ \begin{aligned} -\log(1-u)&=u+\frac{u^2}{2}+O(u^3),\\ D_{\mathrm{KL}}&=\frac{u^2}{4}+O(u^3). \end{aligned} \]

The linear term cancels. Halving a small share therefore cuts this average evidence by about four. For \(k\) independent sketch rows, log likelihood ratios add, so the divergence is \(k\) times as large: releasing more mixtures also provides more evidence.

Noise as extra rows

A dataset may contain a person with a large share. We can reduce it by adding independent Gaussian noise to the sketch:

\[ \widetilde{\mathbf X}=\mathbf S\mathbf X+\eta\boldsymbol\xi, \]

where \(\boldsymbol\xi\in\mathbb R^{k\times d}\) has independent standard Gaussian entries and \(\eta\ge0\) is the noise scale. Each released row now has covariance \(\mathbf X^\top\mathbf X+\eta^2\mathbf I_d\). There is another way to generate exactly the same distribution:

\[ \begin{bmatrix}\mathbf S&\boldsymbol\xi\end{bmatrix} \begin{bmatrix}\mathbf X\\\eta\mathbf I_d\end{bmatrix} =\mathbf S\mathbf X+\eta\boldsymbol\xi. \]

Adding noise is equivalent to appending \(d\) artificial rows before mixing. These phantom rows contribute variance in every direction, so no person needs to supply it alone. If \(\|\mathbf x_i\|_2\le1\), the new share satisfies

\[ u_i=\mathbf x_i^\top(\mathbf X^\top\mathbf X+\eta^2\mathbf I_d)^{-1}\mathbf x_i \le\frac1{\lambda_{\min}(\mathbf X^\top\mathbf X)+\eta^2}. \]

A higher covariance floor makes every share smaller. This is the basis of the Gaussian mixing privacy analysis [LSSLSW, 2025].

In the example below, the covariance before adding noise is the identity and one person’s share is \(r\). Adding variance \(\eta^2\) changes their share to \(u=r/(1+\eta^2)\). Start with no added noise and vary \(r\); then add noise while keeping the person fixed.

Person's share before noise, r: 0.60

Added variance, η²: 0.00

━ Person present ┄ Person removed

━ Exact KL ┄ Quadratic approximation, u²/4

One sketch row. Both densities are rescaled by the standard deviation with the person present, so their variances are 1 and 1 − u. Axes stay fixed as the controls move. KL is average evidence for this pair of distributions, not a differential privacy budget.

Small KL divergence is an explanation of the mechanism, not yet a differential privacy guarantee. KL controls an average; privacy also has to control rare outputs with unusually large likelihood ratios, in both directions between neighboring datasets.

From the covariance calculation to a privacy guarantee

A mechanism is \((\varepsilon,\delta)\)-differentially private if, for every pair of neighboring datasets \(D,D'\) and every output event \(A\),

\[ \Pr(M(D)\in A)\le e^\varepsilon\Pr(M(D')\in A)+\delta. \]

Here neighbors differ by blanking one person’s record, and the requirement includes both orders of the pair.

First fix \(\eta\) independently of the data, and suppose every dataset under consideration has covariance floor \(\lambda_{\min}(\mathbf X^\top\mathbf X)+\eta^2\ge\lambda_0>1\), with row norms at most one. The Gaussian mixing bound [LSSLSW, Lemma 1] controls the Rényi divergence between neighboring output laws \(P,Q\) in both directions:

\[ D_\alpha(P\|Q)\le R_\alpha, \qquad R_\alpha= \frac{k}{2(\alpha-1)} \left[ \alpha\log\left(1-\frac1{\lambda_0}\right) -\log\left(1-\frac\alpha{\lambda_0}\right) \right], \qquad 1<\alpha<\lambda_0. \]

Rényi divergence controls an exponential moment of the privacy loss \(\ell(z)=\log(p(z)/q(z))\):

\[ \mathbb E_{z\sim P}[e^{(\alpha-1)\ell(z)}] =e^{(\alpha-1)D_\alpha(P\|Q)}. \]

Markov’s inequality therefore bounds \(\Pr_P(\ell>\varepsilon)\) by \(\delta\) when \(\varepsilon=R_\alpha+\log(1/\delta)/(\alpha-1)\). Outside that exceptional set, \(p\le e^\varepsilon q\); integrating over any event gives the required privacy inequality. Thus a valid bound is [Mironov, 2017]

\[ \varepsilon= \inf_{1<\alpha<\lambda_0} \left\{R_\alpha+\frac{\log(1/\delta)}{\alpha-1}\right\}. \]

The range of orders matters: the logarithm in \(R_\alpha\) ceases to be finite at \(\alpha=\lambda_0\). An optimization that ignores this boundary can claim a privacy bound the argument does not support.

This is a guarantee for the sketch at a fixed noise scale and a uniform covariance floor. IHM also accounts for selecting the scale privately and for the gradient releases across iterations.

For arbitrary bounded data, choosing \(\eta^2\ge\lambda_0\) supplies the floor without inspecting the data at all. IHM can add less noise when the data already supplies enough variance, using a private lower estimate of its smallest Gram eigenvalue. That estimate has its own privacy cost; substituting the exact, unreleased eigenvalue into a data-dependent noise rule would not establish the same guarantee.

The gradient still needs protection

The sketch handles curvature, but an iterative solver also uses the responses through the negative gradient

\[ \mathbf g(\boldsymbol\theta) =\mathbf X^\top(\mathbf y-\mathbf X\boldsymbol\theta) =\sum_i\mathbf x_i(y_i-\mathbf x_i^\top\boldsymbol\theta). \]

This sum is exact, so the sketch’s randomness does not protect it. Condition on the previous private outputs, which fix the current iterate \(\boldsymbol\theta\). Clip each residual to \([-c,c]\) and add fresh noise:

\[ \widetilde{\mathbf g} =\sum_i\mathbf x_i\operatorname{clip}_c(y_i-\mathbf x_i^\top\boldsymbol\theta) +\sigma\mathbf z, \qquad \mathbf z\sim\mathcal N(\mathbf0,\mathbf I_d). \]

Blanking a person changes the clipped sum by at most \(c\) in Euclidean norm, because \(\|\mathbf x_i\|_2\le1\). That bound lets the Gaussian mechanism calibrate \(\sigma\) [Mironov, 2017]. IHM combines this private gradient with the private curvature sketch and accounts for both costs over its iterations.

How the two pieces enter an update

With fresh independent randomness each round, the paper’s update is

\[ \widehat{\boldsymbol\theta}_{t+1} =\widehat{\boldsymbol\theta}_t +\left(\frac1k\widetilde{\mathbf X}_t^\top\widetilde{\mathbf X}_t\right)^{-1} \left(\widetilde{\mathbf g}_t-\eta^2\widehat{\boldsymbol\theta}_t\right). \]

Take \(k\ge d\) and a positive definite row covariance so the sketched Gram matrix is invertible almost surely. The extra \(-\eta^2\widehat{\boldsymbol\theta}_t\) term matches the ridge penalty introduced by the phantom rows. If the gradient noise is zero and clipping is inactive, its stationary point solves

\[ (\mathbf X^\top\mathbf X+\eta^2\mathbf I_d)\boldsymbol\theta =\mathbf X^\top\mathbf y. \]

With active clipping, the deterministic gradient instead belongs to a clipped-residual loss, so that ridge least-squares conclusion need not hold. The paper’s accuracy analysis chooses parameters to control clipping and optimization error; its privacy argument accounts for every release [LSSLW, 2026, Algorithm 1 and Theorem 1].

For the full regression algorithm and its accuracy bounds, see Lev et al. (2026). Earlier routes include Sheffet’s Gaussian-projection regression and Wang’s AdaSSP, which perturbs the sufficient statistics directly.