Differentially Private Eigenvalues

3.1. Allocation of pair rows

The PMWM procedure consists of two stages: quantile estimation and winsorized mean estimation. The \(m\) pair rows of \(W\) can either be reused in both stages or divided into disjoint subsets.

Let \(I_q\) denote the indices used for quantile estimation and \(I_M\) those used for mean estimation, with \(m_q=|I_q|\) and \(m_M=|I_M|\). The two stages may reuse all \(m\) pair rows, in which case \(I_q=I_M=\{1,\ldots,m\}\), or use disjoint subsets, with \(m_q=\lfloor m/2\rfloor\) and \(m_M=m-m_q\).

3.2. Private quantiles and cutoff bounds

For each component, define the trimming proportion as

\[ \zeta = \min\left\{ \max\left( \frac{c_{\mathrm{trim}}}{m_Q}, \eta_{\mathrm{trim}} \right), 0.49 \right\}, \]

where \(c_{\mathrm{trim}}>0\) is a baseline trimming constant and \(0\le\eta_{\mathrm{trim}}<0.5\) is a lower-bound tuning value.

The lower and upper cutoffs are estimated using the private quantiles at levels \(\zeta\) and \(1-\zeta\), respectively:

\[ L_\ell^{\mathrm{priv}} = \mathcal Q_{\varepsilon_Q,\beta} (W_{I_Q,\ell},\zeta), \qquad U_\ell^{\mathrm{priv}} = \mathcal Q_{\varepsilon_Q,\beta} (W_{I_Q,\ell},1-\zeta). \]

Here, \(\mathcal Q_{\varepsilon_Q,\beta}\) denotes the one-sided unbounded private quantile mechanism of Durfee (2023), with privacy parameter \(\varepsilon_Q\) and geometric grid parameter \(\beta>1\). Since the squared pair differences are nonnegative, the quantile mechanism uses zero as a public lower bound for its search. The parameter \(\beta\) controls the grid resolution, with values closer to one producing a finer grid.

The implementation searches the public grid \(t_i=\beta^i-1\) for \(i=1,\ldots,5000\). If no noisy count reaches the noisy target, it issues a warning and returns \(t_{5000}\) before restricting the cutoffs to \([a,b]\). With the default \(\beta=1.001\), the largest grid value is \(1.001^{5000}-1\approx147.04\). At this fixed step limit, choosing \(\beta\) closer to one refines the grid but also reduces its search range, which can limit cutoff accuracy for larger-scale squared pair differences. The public pmwm_control() interface does not expose the step limit.

After quantile estimation, the private cutoffs are restricted to a prespecified public interval \([a,b]\), where \(a<b\):

\[ L_\ell = \min\{b,\max\{a,L_\ell^{\mathrm{priv}}\}\}, \qquad U_\ell = \min\{b,\max\{a,U_\ell^{\mathrm{priv}}\}\}. \]

If the resulting upper cutoff is below the lower cutoff, the interval is reset to \([a,b]\). These are bounds on the estimated cutoffs; the original paired values do not have to lie inside \([a,b]\).

3.3. Winsorized mean and Gaussian perturbation

Define the winsorization rule and component means by

\[ \phi_{L_\ell,U_\ell}(w)=\min\{U_\ell,\max\{L_\ell,w\}\}, \qquad \bar w_\ell^{\mathrm{PMWM}}= \frac1{m_M}\sum_{r\in I_M}\phi_{L_\ell,U_\ell}(w_{r\ell}). \]

Conditional on the privately selected cutoffs, the joint mean sensitivity is bounded by

\[ \Delta_2\le \frac1{m_M}\sqrt{\sum_{\ell=1}^k (U_\ell-L_\ell)^2}\le\frac{\sqrt{k}(b-a)}{m_M}. \]

To ensure DP, Gaussian noise is added to the vector of winsorized means:

\[ \widetilde{\lambda}^{\mathrm{PMWM,raw}}= \bar w^{\mathrm{PMWM}}+Z,\qquad Z_\ell\sim N(0,\sigma_\ell^2),\qquad \sigma_\ell=\frac{\sqrt{k}(U_\ell-L_\ell)}{m_M\mu_M}, \]

Here, \(\mu_M\) is the Gaussian-DP parameter calibrated from the mean-estimation privacy allocation \((\varepsilon_M,\delta_M)\).

3.4. Privacy allocation and tuning choices

The scree epsilon budget is divided equally between private quantile estimation and noisy winsorized mean estimation. For each of the \(k\) components, two private quantiles are released to estimate the lower and upper cutoffs. Thus, each of the \(2k\) quantile releases uses \(\varepsilon_Q =\frac{\varepsilon_{\mathrm{scree}}}{4k}\), with no delta allocation.

The \(k\)-dimensional winsorized mean vector is then released using the Gaussian mechanism with privacy budget

\[ (\varepsilon_M,\delta_M) = \left( \frac{\varepsilon_{\mathrm{scree}}}{2}, \delta_{\mathrm{scree}} \right). \]

By sequential composition, the two stages together satisfy \((\varepsilon_{\mathrm{scree}},\delta_{\mathrm{scree}})\)-DP, provided that the respective mechanism guarantees hold.

The PMWM tuning parameters are:

4.1. Huber loss and score

For a threshold \(\tau>0\), the Huber loss function and its derivative are defined as

\[ \rho_\tau(u)= \begin{cases} \frac12u^2,& |u|\le\tau,\\ \tau|u|-\frac12\tau^2,& |u|>\tau, \end{cases} \qquad \psi_\tau(u)=\rho_\tau'(u) =\min\{\tau,\max\{-\tau,u\}\}. \]

For each component \(\ell=1,\ldots,k\), the corresponding non-private Huber mean estimator is defined by

\[ \widehat\mu_{\tau_\ell,\ell} \in \underset{\theta\in\mathbb R}{\arg\min} \frac1m \sum_{r=1}^m \rho_{\tau_\ell}(w_{r\ell}-\theta). \]

Since \(|\psi_{\tau_\ell}(u)|\le\tau_\ell\), extreme residuals have bounded influence on the iterative updates.

4.2. Private scale and Huber thresholds

Before estimating the Huber mean, we privately estimate the variability of the paired squared differences for each component. Then, the resulting private scale proxy \(\widetilde{m}_{2,\ell}\) is used to determine the Huber threshold \(\tau_\ell\).

To obtain this scale proxy, the \(m\) paired squared differences are paired again to form additional half squared differences. Their block medians are then summarized using a private histogram over dyadic scales \(2^r\). The public range \(r_{\min}\le r\le r_{\max}\) determines the candidate scales, while \(M\) controls the number of blocks. The detailed procedure is provided in the private robust estimator for \(m_2\).

Let \(\rho\in(0,1)\) denote the proportion of the scree epsilon budget allocated to private scale estimation. For \(k\) components, the privacy allocations are

\[ \varepsilon_{\mathrm{scale},\ell} = \frac{\rho\varepsilon_{\mathrm{scree}}}{k}, \qquad \varepsilon_{\mathrm{GD}} = (1-\rho)\varepsilon_{\mathrm{scree}}. \]

Using the private scale proxy and the privacy budget allocated to gradient descent, the Huber threshold for each component is defined as

\[ \tau_\ell = \sqrt{\widetilde{m}_{2,\ell}} \sqrt{ \frac{ (\varepsilon_{\mathrm{GD}}/k)m }{ \sqrt{(1+\log m)\log m} } }. \]

4.3. Noisy gradient descent

Let \(\theta_\ell^{(0)}=\mu_0\) be a public initial value, \(\eta_0>0\) the step size, and \(T\) the number of iterations. The private iterates are projected to the nonnegative half-line:

\[ \theta_\ell^{(t+1)}= \max\left\{0,\,\theta_\ell^{(t)}+ \frac{\eta_0}{m}\sum_{r=1}^m \psi_{\tau_\ell}(w_{r\ell}-\theta_\ell^{(t)})+ Z_\ell^{(t)}\right\}. \]

The Gaussian perturbations are drawn independently across coordinates and iterations with

\[ Z_\ell^{(t)}\sim N(0,\sigma_{\ell,H}^2),\qquad \sigma_{\ell,H}= \frac{2\eta_0\tau_\ell\sqrt{kT}}{m\mu_{\mathrm{GD}}}, \]

where \(\mu_{\mathrm{GD}}\) is the Gaussian-DP parameter calibrated from the gradient-descent privacy allocation \((\varepsilon_{\mathrm{GD}},\delta_{\mathrm{scree}})\). After \(T\) steps, \(\widetilde{\lambda}_\ell^{\mathrm{Huber,raw}}=\theta_\ell^{(T)}\).

The scale step consumes \(\rho\varepsilon_{\mathrm{scree}}\) and no delta; the gradient step consumes \(((1-\rho)\varepsilon_{\mathrm{scree}},\delta_{\mathrm{scree}})\). Their releases must be accounted for together. The method requires at least four initial pair rows, corresponding to at least eight observations.

The Huber tuning parameters are:

The iteration and block counts may be chosen as functions of public sample size; the corresponding automatic choices used by the package are described in huber_control()