Differentially Private Coordinate Descent for Composite Empirical Risk Minimization

Paul Mangold, Aurélien Bellet, Joseph Salmon, Marc Tommasi

Introduction

Machine learning fundamentally relies on the availability of data, which can be sensitive or confidential. It is now well-known that preventing learned models from leaking information about individual training points requires particular attention . A standard approach for training models while provably controlling the amount of leakage is to solve an empirical risk minimization (ERM) problem under a differential privacy (DP) constraint . In this work, we aim to design a differentially private algorithm which approximates the solution to a composite ERM problem of the form:

Differential privacy constraints induce a trade-off between the privacy and the utility (i.e., optimization error) of the solution of (1). This trade-off was made explicit by , who derived lower bounds on the achievable error given a fixed privacy budget. To solve the DP-ERM problem in practice, the most popular approaches are based on Differentially Private variants of Stochastic Gradient Descent (DP-SGD) , in which random perturbations are added to the (stochastic) gradients. analyzed DP-SGD in the non-smooth DP-ERM setting, and then proposed an efficient DP-SVRG algorithm for composite DP-ERM. Both algorithms match known lower bounds. SGD-style algorithms perform well in a wide variety of settings, but also have some flaws: they either require small (or decreasing) step sizes or variance reduction schemes to guarantee convergence, and they can be slow when gradients’ coordinates are imbalanced. These flaws propagate to the private counterparts of these algorithms. Despite a few attempts at designing other differentially private solvers for ERM under different setups , the differentially private optimization toolbox remains limited, which undoubtedly restricts the resolution of practical problems.

In this paper, we propose and analyze a Differentially Private proximal Coordinate Descent algorithm (DP-CD), which performs updates based on perturbed coordinate-wise gradients (i.e., partial derivatives). Coordinate Descent (CD) methods have encountered a large success in non-private machine learning due to their simplicity and effectiveness , and have seen a surge of practical and theoretical interest in the last decade . In contrast to SGD, they converge with constant step sizes that adapt to the coordinate-wise smoothness of the objective. Additionally, CD updates naturally tend to have a lower sensitivity. Operating with partial gradients thus enables our private algorithm to reduce the perturbation required to guarantee privacy without resorting to amplification by subsampling .

We propose a novel analysis of proximal CD with perturbed gradients to derive optimal upper bounds on the privacy-utility trade-off achieved by DP-CD. We prove a recursion on distances of CD iterates to an optimal point that keeps track of coordinate-wise regularity constants in a tight manner and allows to use large, constant step sizes that yield high utility. Our results highlight the fact that DP-CD can exploit imbalanced gradient coordinates to outperform DP-SGD. They also improve upon known convergence rates for inexact CD in the non-private setting . We assess the optimality of DP-CD by deriving lower bounds that capture coordinate-wise Lipschitz regularity measures, and show that DP-CD matches those bounds up to logarithmic factors. Our lower bounds also suggest interesting perspectives for future work on DP-CD algorithms.

Our theoretical results have important consequences for practical implementations, which heavily rely on gradient clipping to achieve good utility. In contrast to DP-SGD, DP-CD requires to set coordinate-wise clipping thresholds, which can lead to impractical coordinate-wise hyperparameter tuning. We instead propose a simple rule for adapting these thresholds from a single hyperparameter. We also show how the coordinate-wise smoothness constants used by DP-CD can be estimated privately. We validate our theory with numerical experiments on real and synthetic datasets. These experiments further show that even in balanced problems, DP-CD can still improve over DP-SGD, confirming the relevance of DP-CD for DP-ERM.

Our main contributions can be summarized as follows:

We propose the first proximal CD algorithm for composite DP-ERM, formally prove its utility, and highlight regimes where it outperforms DP-SGD.

We show matching lower bounds under coordinate-wise regularity assumptions.

We give practical guidelines to use DP-CD, and show its relevance through numerical experiments.

Preliminaries

In this section, we introduce important technical notions that will be used throughout the paper.

We start by defining two conjugate norms that will be crucial in our analysis, for they allow to keep track of coordinate-wise quantities. Let ⟨u,v⟩=∑j=1puivi\left\langle u,v\right\rangle=\sum_{j=1}^{p}u_{i}v_{i} be the Euclidean dot product, let M=diag⁡(M1,…,Mp)M=\operatorname{diag}(M_{1},\dots,M_{p}) with M1,…,Mp>0M_{1},\dots,M_{p}>0, and

The above component-wise regularity hypotheses are not restrictive: Λ\Lambda-Lipschitzness implies (Λ,…,Λ)(\Lambda,\dots,\Lambda)-component-Lipschitzness and β\beta-smoothness implies (β,…,β)(\beta,\dots,\beta)-component-smoothness. Yet, the actual component-wise constants of a function can be much lower than what can be deduced from their global counterparts. This will be crucial for our analysis and in the performance of DP-CD.

When ψ\psi is the characteristic function of a convex set (with separable components), the regularity assumptions only need to hold on this set. This allows considering problem (1) with a smooth objective under box-constraints.

Let D\mathcal{D} be a set of datasets and F\mathcal{F} a set of possible outcomes. Two datasets D,D′∈DD,D^{\prime}\in\mathcal{D} are said neighboring (denoted by D∼D′D\sim D^{\prime}) if they differ on at most one element.

A randomized algorithm A:D→F\mathcal{A}:\mathcal{D}\rightarrow\mathcal{F} is (ϵ,δ)(\epsilon,\delta)-differentially private if, for all neighboring datasets D,D′∈DD,D^{\prime}\in\mathcal{D} and all S⊆FS\subseteq\mathcal{F} in the range of A\mathcal{A}:

In this paper, we consider the classic central model of DP, where a trusted curator has access to the raw dataset and releases a model trained on this datasetIn fact, our privacy guarantees hold even if all intermediate iterates are released (not just the final model)..

Differentially Private Coordinate Descent

In this section, we introduce the Differentially Private proximal Coordinate Descent (DP-CD) algorithm to solve problem (1) under (ϵ,δ)(\epsilon,\delta)-DP constraints. We first describe our algorithm, show how to parameterize it to satisfy the desired privacy constraint, and prove corresponding utility results. Finally, we compare these utility guarantees with DP-SGD.

Update (3) only requires the computation of the jj-th entry of the gradient. To satisfy differential privacy, we perturb this gradient entry with additive Gaussian noise of variance σj2\sigma_{j}^{2}. The complete DP-CD procedure is shown in Algorithm 1. At each iteration, we pick a coordinate uniformly at random and update according to (3), albeit with noise addition (see line 7). For technical reasons related to our analysis, we use a periodic averaging scheme (line 8). This scheme is similar to DP-SVRG , although no variance reduction is required since DP-CD computes coordinate gradients over the whole dataset.

2 Privacy Guarantees

For Algorithm 1 to satisfy (ϵ,δ)(\epsilon,\delta)-DP, the noise scales σ=(σ1,…,σp)\sigma=(\sigma_{1},\dots,\sigma_{p}) can be calibrated as given in Theorem 3.1.

The dependence of the noise scales on ϵ\epsilon, δ\delta, nn and TKTK (the number of updates) in Theorem 3.1 is standard in DP-ERM. However, the noise is calibrated to the loss function’s component-Lipschitz constants. These can be much lower their global counterpart, the latter being used to calibrate the noise in DP-SGD algorithms. This will be crucial for DP-CD to achieve better utility than DP-SGD in some regimes. We also note that, unlike DP-SGD, DP-CD does not rely on privacy amplification by subsampling , and thereby avoids the approximations required by these schemes to bound the privacy loss.

Theorem 3.1 assumes ϵ∈(0,1]\epsilon\in(0,1] to give a simple closed form for the noise scales. In practice we compute tighter values numerically using Rényi DP formulas directly (see Eq. 18 in Appendix B), removing the need for this assumption.

3 Utility Guarantees

We now state our central result on the utility of DP-CD for the composite DP-ERM problem. As done in previous work, we use the asymptotic notation O~\widetilde{O} to hide non-significant logarithmic factors. Non-asymptotic utility bounds can be found in Appendix C.

For FF convex, K=O(RMpnϵ∥L∥M−1)K=O\left(\frac{R_{M}\sqrt{p}n\epsilon}{\left\lVert L\right\rVert_{M^{-1}}}\right), and T=1T=1, then:

where RM=max⁡(F(w0)−F(w∗),∥w0−w∗∥M)R_{M}=\max(\sqrt{F(w^{0})-F(w^{*})},\left\lVert w^{0}-w^{*}\right\rVert_{M}) and more simply RM=∥w0−w∗∥MR_{M}=\left\lVert w^{0}-w^{*}\right\rVert_{M} when ψ=0\psi=0.

For FF μM\mu_{M}-strongly convex w.r.t. ∥⋅∥M\smash{\left\lVert\cdot\right\rVert_{M}}, K=O(p/μM)K=O\left(p/\mu_{M}\right), and T=O(log⁡(nϵμM/p∥L∥M−1))T=O\left(\log(n\epsilon\mu_{M}/p\left\lVert L\right\rVert_{M^{-1}})\right), then:

Expectations are over the randomness of the algorithm.

The above inequality shows that coordinate-wise updates leave a fraction p−1p\frac{p-1}{p} of the function “unchanged”, while the remaining part decreases (up to additive noise). Importantly, all quantities are measured in MM-norm. When summing (4) for k=0,…,K−1k=0,\dots,K-1, its left hand side simplifies and its right hand side is simplified as a telescoping sum:

Our novel convergence proof of CD is also useful in the non-private setting. In particular, we improve upon known convergence rates for inexact CD methods with additive error , under the hypothesis that gradients are noisy and unbiased. In their formalism, we have α=0\alpha=0 and β=∥σ∥M−12/p\beta=\left\lVert\sigma\right\rVert_{M^{-1}}^{2}/p. With our analysis, the algorithm requires 2pRM2/(ξ−pβ)2pR_{M}^{2}/(\xi-p\beta) (resp. 4p/μMlog⁡((F(w0)−F∗)/(ξ−pβ))4p/\mu_{M}\log((F(w^{0})-F^{*})/(\xi-p\beta))) iterations to achieve expected precision ξ>pβ\xi>p\beta when FF is convex (resp. μM\mu_{M}-strongly-convex w.r.t. ∥⋅∥M\left\lVert\cdot\right\rVert_{M}), improving upon ’s results by a factor pβ/2RM2\sqrt{p\beta/2R_{M}^{2}} (resp. μM/2\mu_{M}/2). See Section C.3 for details. Moreover, unlike this prior work, our analysis does not require the objective to decrease at each iteration, which is essential to guarantee DP.

Our utility guarantees stated in Theorem 3.3 directly depend on precise coordinate-wise regularity measures of the objective function. In particular, the initial distance to optimal, the strong convexity parameter and the overall sensitivity of the loss function are measured in the norms ∥⋅∥M\smash{\left\lVert\cdot\right\rVert_{M}} and ∥⋅∥M−1\smash{\left\lVert\cdot\right\rVert_{M^{-1}}} (i.e., weighted by coordinate-wise smoothness constants or their inverse). In the remainder of this section, we thoroughly compare our utility results with existing ones for DP-SGD. We will show the optimality of our utility guarantees in Section 4.

4 Comparison with DP-SGD and DP-SVRG

When the smoothness constants MM are all equal, ∥L∥M−1RM=∥L∥2RI\left\lVert L\right\rVert_{M^{-1}}R_{M}=\left\lVert L\right\rVert_{2}R_{I} and ∥L∥M−12/μM=∥L∥22/μI{\left\lVert L\right\rVert_{M^{-1}}^{2}}/{\mu_{M}}={\left\lVert L\right\rVert_{2}^{2}}/{\mu_{I}}. This boils down to comparing ∥L∥2\left\lVert L\right\rVert_{2} to Λ\Lambda. As Λ≤∥L∥2≤pΛ\Lambda\leq\left\lVert L\right\rVert_{2}\leq\sqrt{p}\Lambda, DP-CD can be up to pp times worse than DP-SGD. This can only happen when features are extremely correlated, which is generally not the case in machine learning. We show empirically in Section 6.2 that, even in balanced regimes, DP-CD can still significantly outperform DP-SGD.

Lower Bounds

If FF is μI\mu_{I}-strongly-convex w.r.t. ∥⋅∥2\left\lVert\cdot\right\rVert_{2}:

We recover the lower bounds of for Λ\Lambda-Lipschitz losses as a special case of ours by setting L1=⋯=Lp=Λ/pL_{1}=\cdots=L_{p}={\Lambda}/{\sqrt{p}}. In this case, the loss function used in our proof is indeed (∑j=1pLj2)1/2=Λ(\sum_{j=1}^{p}L_{j}^{2})^{1/2}=\Lambda-Lipschitz. To relate these lower bounds to the performance of DP-CD, consider a suboptimal version of our algorithm where the step sizes are set to γ1=⋯=γp=(max⁡jMj)−1\gamma_{1}=\cdots=\gamma_{p}=({\max_{j}M_{j}})^{-1}. In this setting, results from Theorem 3.3 still hold, and match the lower bounds from Theorem 4.1 up to logarithmic factors. We leave open the question of the optimality of DP-CD under the additional hypothesis of smoothness.

We note that the assumption on the sum of the LjL_{j}’s over a set of indices J\mathcal{J} in Theorem 4.1 can be eliminated at the cost of an additional factor of Lmin⁡/Lmax⁡{L_{\min}}/{L_{\max}} for convex losses and (Lmin⁡/Lmax⁡)2({L_{\min}}/{L_{\max}})^{2} for strongly-convex losses, making the bound looser. Although the aforementioned assumption may seem solely technical, we conjecture that better utility is possible when a few coordinate-wise Lipschitz constants dominate the others. We discuss this further in Section 8.

DP-CD in Practice

We now discuss practical questions related to DP-CD. First, we show how to implement coordinate-wise gradient clipping using a single hyperparameter. Second, we explain how to privately estimate the smoothness constants. Finally, we discuss the possibility of standardizing the features and how this relates to estimating smoothness constants for the important problem of fitting generalized linear models.

In DP-CD, gradients are released one coordinate at a time and should thus be clipped in a coordinate-wise fashion. Using the same threshold for each coordinate would ruin the ability of DP-CD to account for imbalance across gradient coordinates, whereas tuning coordinate-wise thresholds as pp individual hyperparameters {Cj}j=1p\{C_{j}\}_{j=1}^{p} is impractical.

2 Private Smoothness Constants

3 Feature Standardization

Numerical Experiments

For DP-SGD, we use constant step sizes and standard gradient clipping. For DP-CD, we adapt the coordinate-wise clipping thresholds from one hyperparameter, as described in Section 5.1. Similarly, coordinate-wise step sizes are set to γj=γ/Mj\gamma_{j}=\gamma/M_{j}, where γ\gamma is a hyperparameter. When the coordinate-wise smoothness constants are not all equal, we also consider DP-CD with privately computed MjM_{j}’s, as described in Section 5.2. For each dataset and each algorithm, we simultaneously tune the clipping threshold, the number of passes over the dataset and, for DP-CD and DP-SGD, the step sizes. After tuning these parameters, we report the relative error to the (non-private) optimal objective value. The complete tuning procedure is described in Section G.1, where we also give the best error for various numbers of passes for each algorithm and dataset. The code used to obtain all our results is available in a public repository https://gitlab.inria.fr/pmangold1/private-coordinate-descent/ and in the supplementary material.

In the Electricity and California datasets, features are naturally imbalanced. DP-CD can exploit this through the use of coordinate-wise smoothness constants. We also consider a variant of DP-CD (DP-CD-P) which dedicates 10%10\% of the privacy budget ϵ\epsilon to estimate these constants (see Section 5.2) from a crude upper bound on each feature (twice their maximal absolute value). It then uses the resulting private smoothness constants in step sizes and clipping thresholds. Figure 1 shows that DP-CD outperforms DP-SGD and DP-SCD by an order of magnitude on both datasets, even when the smoothness constants are estimated privately.

2 Balanced Datasets

To assess the performance of DP-CD when coordinate-wise smoothness constants are balanced, we standardize the Electricity and California datasets (see Section 5.3). As standardization is done for all algorithms, we do not account for it in the privacy budget. On standardized datasets, coordinate-wise smoothness constants are all equal, removing the need of estimating them privately. We report the results in Figure 2. Although our theory suggests that DP-CD may do worse than DP-SGD in balanced regimes, we observe that it still improves over DP-SGD (and DP-SCD) in practice. Similar observations hold in our challenging Sparse LASSO problem, where DP-SGD is barely able to make any progress. We believe these results are in part due to the beneficial effect of clipping in DP-CD, and the fact that DP-SGD relies on amplification by subsampling, for which privacy accounting is not perfectly tight. Additionally, CD methods are known to perform well on fitting linear models: our results show that this transfers well to private optimization.

3 Running Time

The results above showed that DP-CD yields better utility than DP-SGD. We also observe that DP-CD tends to reach these results in up to 1010 times fewer passes on the data than DP-SGD (see Section G.1 for detailed results). Additionally, when accounting for running time, DP-CD significantly outperforms DP-SGD: we refer to Section G.2 for the counterparts of Figure 1 and 2 as a function of the running time instead of the number of passes.

Related Work

Differentially Private Empirical Risk Minimization was first studied by , using output perturbation (adding noise to the solution of the non-private ERM problem) and objective perturbation (adding noise to the ERM objective itself). then proposed DP-SGD and proved its near-optimality. obtained faster convergence rates using a DP version of the SVRG algorithm . DP-SGD has become the standard approach to DP-ERM. In our work, we show that coordinate-wise updates can have lower sensitivity than DP-SGD updates and propose a DP-CD algorithm achieving competitive results. A private variant of the Frank-Wolfe algorithm (DP-FW) was also proposed to solve constrained DP-ERM problems . Although these algorithms achieve a good privacy-utility trade-off in theory, we are not aware of any empirical evaluation. DP-FW algorithms access gradients indirectly through a linear optimization oracle over a constrained set. Restricting to a constrained set is not necessary in DP-CD, allowing its use for a different family of problems.

Recent work has also studied algorithms and utility guarantees for stochastic convex optimization under differential privacy constraints, a problem very similar to DP-ERM. [27, 4, following work from] extended results known for DP-ERM to this setting, showing that the population risk of DP-SCO is asymptotically equivalent to the one of non-private SCO. Efficient algorithms for DP-SCO were proposed by , and studied stochastic variants of DP-FW. As detailed by results from DP-ERM can be converted to DP-SCO.

Coordinate descent (CD) algorithms have a long history in optimization. have shown convergence results for (block) CD algorithms for nonsmooth optimization. later proved a global non-asymptotic 1/k1/k convergence rate for CD with random choice of coordinates for a convex, smooth objective. Parallel, proximal variants were developed by , while further considered non-separable non-smooth parts. introduced Dual CD algorithms for smooth ERM, showing performance similar to SVRG. We refer to and for detailed reviews on CD. Inexact CD was studied by , but their analysis requires updates not to increase the objective, which is hardly compatible with DP. We obtain tighter results for inexact CD with noisy gradients (see Remark 3.4).

Conclusion and Discussion

We presented the first differentially private proximal coordinate descent algorithm for composite DP-ERM. Using an original approach to analyze proximal CD with perturbed gradients, we derived optimal upper bounds on the privacy-utility trade-off achieved by DP-CD. We also prove new lower bounds under a component-Lipschitzness assumption, and showed that DP-CD matches these bounds. Our results demonstrate that DP-CD strongly outperforms DP-SGD when gradients’ coordinates are imbalanced. Numerical experiments show that DP-CD also performs very well in balanced regimes. The choice of coordinate-wise clipping thresholds is crucial for DP-CD to achieve good utility in practice, and we provided a simple rule to set them.

Although DP-CD already achieves good utility when most coordinates have small sensitivity, our lower bounds suggest that even better utility could be achieved by dynamically allocating more privacy budget to coordinates with largest sensitivities. A promising direction is to design DP-CD algorithms that leverage active set methods , which could provide practical alternatives to recent DP-SGD approaches that use a subspace assumption . Finally, we believe that adaptive clipping techniques may help to further improve the practical performance of DP-CD when coordinate-wise smoothness constants are more balanced.

Acknowledgments

The authors would like to thank the anonymous reviewers who provided useful feedback on previous versions of this work, which helped to improve the paper.

This work was supported in part by the Inria Exploratory Action FLAMED and by the French National Research Agency (ANR) through grant ANR-20-CE23-0015 (Project PRIDE) and ANR-20-CHIA-0001-01 (Chaire IA CaMeLOt).

References

Appendix A Lemmas on Sensitivity

In this section, we let X\mathcal{X} be the universe where the data is drawn from. To upper bound the sensitivities of a function’s gradient, we start by recalling in Lemma A.1 that (coordinate) gradients are bounded by (coordinate-wise-)Lipschitz constants. We then link this upper bound with gradients’ sensitivities in Lemma A.2.

which is the claim of the first statement. To prove the second statement, we proceed similarly: the triangle inequality and Lemma A.1 give the following upper bounds:

We obtain the inequality (2) stated in Section 2 as a corollary.

Appendix B Proof of Theorem 3.1

To track the privacy loss of an adaptive composition of KK Gaussian mechanisms, we use Rényi Differential Privacy [38, RDP]. We note that similar results are obtained with zero Concentrated Differential Privacy . This flavor of differential privacy, gives tighter privacy guarantees in that setting, as it reduces the noise variance by a multiplicative factor of log⁡(K/δ)\log(K/\delta) in comparison to the usual advanced composition theorem of differential privacy . Importantly, RDP can be translated back to differential privacy.

In this section, we recall the definition and main properties of zCDP. We denote by D\mathcal{D} the set of all datasets over a universe X\mathcal{X} and by F\mathcal{F} the set of possible outcomes of the randomized algorithms we consider.

We will use the Rényi divergence (Definition B.1), which gives a distribution-oriented vision of privacy.

For two random variables YY and ZZ with values in the same domain C\mathcal{C}, the Rényi divergence is, for α>1\alpha>1,

We now define RDP in Definition B.2. RDP provides a strong privacy guarantee that can be converted to classical differential privacy (Lemma B.3 and Corollary B.8).

A randomized algorithm A:D→F\mathcal{A}:\mathcal{D}\rightarrow\mathcal{F} is (α,ϵ)(\alpha,\epsilon)-Rényi-differentially private (RDP) if, for all all datasets D,D′∈DD,D^{\prime}\in\mathcal{D} differing on at most one element,

If a randomized algorithm A:D→F\mathcal{A}:\mathcal{D}\rightarrow\mathcal{F} is (α,ϵ)(\alpha,\epsilon)-RDP, then it is (ϵ+log⁡(1/δ)α−1,δ)(\epsilon+\frac{\log(1/\delta)}{\alpha-1},\delta)-differentially private for all 0<δ<10<\delta<1.

The above (α,ϵ)(\alpha,\epsilon)-RDP guarantees hold for multiple values of α,ϵ\alpha,\epsilon. As such, ϵ=ϵ(α)\epsilon=\epsilon(\alpha) can be seen as a function of α\alpha, and Lemma B.3 ensures that the algorithm is (ϵ′,δ)(\epsilon^{\prime},\delta)-DP for

We can now restate in Theorem B.5 the composition theorem of RDP, which is key in designing private iterative algorithms.

Let A1,…,AK:D→F\mathcal{A}_{1},\dots,\mathcal{A}_{K}:\mathcal{D}\rightarrow\mathcal{F} be K>0K>0 randomized algorithms, such that for 1≤k≤K1\leq k\leq K, Ak\mathcal{A}_{k} is (α,ϵk(α))(\alpha,\epsilon_{k}(\alpha))-RDP, where these algorithms can be chosen adaptively (i.e., Ak\mathcal{A}_{k} can use to the output of Ak′\mathcal{A}_{k^{\prime}} for all k′<kk^{\prime}<k). Let A:D→FK\mathcal{A}:\mathcal{D}\rightarrow\mathcal{F}^{K} such that for D∈DD\in\mathcal{D}, A(D)=(A1(D),…,AK(D))\mathcal{A}(D)=(\mathcal{A}_{1}(D),\dots,\mathcal{A}_{K}(D)). Then A\mathcal{A} is (α,∑k=1Kϵk(α))\left(\alpha,\sum_{k=1}^{K}\epsilon_{k}(\alpha)\right)-RDP.

Finally, we define the Gaussian mechanism (Definition B.6), as used in Algorithm 1, and restate in Lemma B.7 the privacy guarantees that it satisfies in terms of RDP.

The Gaussian mechanism with noise σ2\sigma^{2} is (α,Δ(f)2α2σ2)(\alpha,\frac{\Delta(f)^{2}\alpha}{2\sigma^{2}})-RDP, where Δ(f)=sup⁡D,D′∥f(D)−f(D′)∥2\Delta(f)=\sup_{D,D^{\prime}}\left\lVert f(D)-f(D^{\prime})\right\rVert_{2} (for neighboring D,D′D,D^{\prime}) is the sensitivity of ff.

The function h=fΔ(f)h=\frac{f}{\Delta(f)} has sensitivity 11, thus for any s>0s>0, the Gaussian mechanism MhGauss(⋅;s)\mathcal{M}_{h}^{Gauss}(\cdot;s) is (α,α2σ2)(\alpha,\frac{\alpha}{2\sigma^{2}})-RDP [38, Corollary 1]. As f=Δ(f)×hf=\Delta(f)\times h, we have MfGauss(⋅;σ)=Δ(f)×MhGauss(⋅;σΔ(f))\mathcal{M}_{f}^{Gauss}(\cdot;\sigma)=\Delta(f)\times\mathcal{M}_{h}^{Gauss}(\cdot;\frac{\sigma}{\Delta(f)}). This mechanism is thus (α,Δ(f)2α2σ2)(\alpha,\frac{\Delta(f)^{2}\alpha}{2\sigma^{2}})-RDP. ∎

Let 0<ϵ≤1,0<δ<130<\epsilon\leq 1,0<\delta<\tfrac{1}{3}. If a randomized algorithm A:D→F\mathcal{A}:\mathcal{D}\rightarrow\mathcal{F} is (α,γα2σ2)(\alpha,\frac{\gamma\alpha}{2\sigma^{2}})-RDP with γ>0\gamma>0 and σ=3γlog⁡(1/δ)ϵ\sigma=\frac{\sqrt{3\gamma\log(1/\delta)}}{\epsilon} for all α>1\alpha>1, it is also (ϵ,δ)(\epsilon,\delta)-DP.

From Remark B.4 it holds that A\mathcal{A} is (ϵ′,δ)(\epsilon^{\prime},\delta)-DP with ϵ′=min⁡α>1{γα2σ2+log⁡(1/δ)α−1}.\epsilon^{\prime}=\min_{\alpha>1}\left\{\frac{\gamma\alpha}{2\sigma^{2}}+\frac{\log(1/\delta)}{\alpha-1}\right\}. This minimum is attained when the derivative of the objective is zero, which is the case when γ2σ2=log⁡(1/δ)(α−1)2\frac{\gamma}{2\sigma^{2}}=\frac{\log(1/\delta)}{(\alpha-1)^{2}}, resulting in α=1+2log⁡(1/δ)σ2γ\alpha=1+\sqrt{\frac{2\log(1/\delta)\sigma^{2}}{\gamma}}. A\mathcal{A} is thus (ϵ′,δ)(\epsilon^{\prime},\delta)-DP with

Choosing σ=3γlog⁡(1/δ)ϵ\sigma=\frac{\sqrt{3\gamma\log(1/\delta)}}{\epsilon} now gives

where the first inequality comes from ϵ≤1\epsilon\leq 1, thus ϵ2≤ϵ\epsilon^{2}\leq\epsilon and δ<1/3\delta<1/3 thus 1log⁡(1/δ)≤1\frac{1}{\log(1/\delta)}\leq 1. The second inequality follows from 1/6+2/3≈0.983<11/6+\sqrt{2/3}\approx 0.983<1. ∎

B.2 Proof of Theorem 3.1

We are now ready to prove Theorem 3.1. From the privacy perspective, Algorithm 1 adaptively releases and post-processes a series of gradient coordinates protected by the Gaussian mechanism. We thus start by proving Lemma B.9, which gives an (ϵ,δ)(\epsilon,\delta)-differential privacy guarantee for the adaptive composition of KK Gaussian mechanisms.

Let σ>0\sigma>0. Lemma B.7 guarantees that the kk-th Gaussian mechanism with noise scale σk=Δ(fk)σ>0\sigma_{k}=\Delta(f_{k})\sigma>0 is (α,α2σ2)(\alpha,\frac{\alpha}{2\sigma^{2}})-RDP. Then, the composition of these KK mechanisms is, according to Theorem B.5, (α,kα2σ2)(\alpha,\frac{k\alpha}{2\sigma^{2}})-RDP. This can be converted to (ϵ,δ)(\epsilon,\delta)-DP via Corollary B.8 with γ=K\gamma=K, which gives σk=Δ(fk)3klog⁡(1/δ)ϵ\sigma_{k}=\frac{\Delta(f_{k})\sqrt{3k\log(1/\delta)}}{\epsilon} for k∈[K]k\in[K]. ∎

thus by Lemma B.9 and the post-processing property of DP, Algorithm 1 is (ϵ,δ)(\epsilon,\delta)-differentially private. ∎

Appendix C Proof of Utility (Theorem 3.3)

Let D∈XnD\in\mathcal{X}^{n} be a dataset of nn elements drawn from a universe X\mathcal{X}. Recall that we consider the following composite empirical risk minimization problem:

C.2 Proof of Theorem 3.3

In this section, we prove our central theorem that guarantees the utility of the DP-CD algorithm. To this end, we start by proving a lemma that upper bounds the expected value of F(θk+1)F(\theta^{k+1}) in Algorithm 1. Using this lemma, we prove sub-linear convergence for the inner loop of DP-CD. This gives the sub-linear convergence of our algorithm for convex losses. Under the additional hypothesis that FF is strongly convex, we show that iterates of the outer loop of DP-CD converge linearly towards the (unique) minimum of FF.

To avoid notational clutter, we will write γjgj\gamma_{j}g_{j} instead of γjgjej\gamma_{j}g_{j}e_{j} throughout this section.

The regularization terms can now be reorganized using the separability of ψ\psi, as done by . Indeed, we notice that

Plugging (28) in (26) results in the following:

which gives the lemma since F=f+ψF=f+\psi. ∎

To exploit this result, we need to upper bound the right hand side of (C.1) for the realizations of θk\theta^{k} in Algorithm 1. This is where our proof differs from classical convergence proofs for coordinate descent methods. Namely, we rewrite the right hand side of (C.1) so as to obtain telescopic terms plus a bias term resulting from the addition of noise, as shown in Lemma C.3.

where ∥σ∥Γ2=∑j=1pγjσj2\left\lVert\sigma\right\rVert_{\Gamma}^{2}=\sum_{j=1}^{p}\gamma_{j}\sigma_{j}^{2} and the expectations are taken over the random choice of jj and η\eta, conditioned upon the realization of θk\theta^{k}.

From Lemma C.1 with θ=θk\theta=\theta^{k}, w=ww=w and g=gg=g as defined above we obtain

We can upper bound the right hand term of (C.2.1) using the convexity of ff and ψ\psi:

where we use the slight abuse of notation ∂ψ(θk−Γg)\partial\psi(\theta^{k}-\Gamma g) to denote any vector in the subdifferential of ψ\psi at the point θk−Γg\theta^{k}-\Gamma g. We now rewrite the dot product:

where the second equality follows from ⟨g,−Γg⟩=−∥g∥Γ2\left\langle g,-\Gamma g\right\rangle=-\left\lVert g\right\rVert_{\Gamma}^{2} and ∥Γg∥M2=∥g∥Γ2M2\left\lVert\Gamma g\right\rVert_{M}^{2}=\left\lVert g\right\rVert_{\Gamma^{2}M}^{2}. We split (38) into two terms: a “descent” term and a “noise” term.

Rewriting the “descent” term. We first focus on the “descent” term. As γj=1Mj\gamma_{j}=\frac{1}{M_{j}} for all j∈[p]j\in[p], it holds that γj2Mj=γj\gamma_{j}^{2}M_{j}=\gamma_{j} which gives −∥g∥Γ2+12∥g∥Γ2M2=−∥g∥Γ2+12∥g∥Γ2=−12∥g∥Γ2-\left\lVert g\right\rVert^{2}_{\Gamma}+\frac{1}{2}\left\lVert g\right\rVert_{\Gamma^{2}M}^{2}=-\left\lVert g\right\rVert^{2}_{\Gamma}+\frac{1}{2}\left\lVert g\right\rVert_{\Gamma}^{2}=-\frac{1}{2}\left\lVert g\right\rVert^{2}_{\Gamma}. We can now rewrite the “descent” term as a difference of two norms, materializing the distance to ww, weighted by the inverse of the step sizes Γ−1\Gamma^{-1}:

where we factorized the norm to obtain the last inequality. We can rewrite (42) as an expectation over the random choice of the coordinate jj (drawn uniformly in [p][p]), given the realizations of θk\theta^{k} and of the noise η\eta (which determines gg):

Finally, we remark that γj−1∣θjk−wj∣2−γj−1∣θjk−γjgj−wj∣2=∥θk−w∥Γ−12−∥θk−γjgj−w∥Γ−12\gamma_{j}^{-1}\left|\theta^{k}_{j}-w_{j}\right|^{2}-\gamma_{j}^{-1}\left|\theta_{j}^{k}-\gamma_{j}g_{j}-w_{j}\right|^{2}=\left\lVert\theta^{k}-w\right\rVert_{\Gamma^{-1}}^{2}-\left\lVert\theta^{k}-\gamma_{j}g_{j}-w\right\rVert_{\Gamma^{-1}}^{2}, as only one coordinate changes between the two vectors, and the squared norm ∥⋅∥Γ−12\left\lVert\cdot\right\rVert_{\Gamma^{-1}}^{2} is separable. We thus obtain

For an update of the coordinate j∈[p]j\in[p], the optimality condition of the proximal operator gives, for ηj\eta_{j} the realization of the noise drawn at the current iteration when coordinate jj is chosen:

and we now separate this term in two using g~\widetilde{g}:

It is now time to consider the expectation with respect to the noise of these terms. First, as g~j\widetilde{g}_{j} is not dependent on the noise anymore, it simply holds that

The last step of our proof now takes care of the following term:

where each inequality comes from the triangle inequality. The non-expansiveness property of the proximal operator (see , Section 2.3) is now key to our result, as it yields

We now have everything to prove the lemma by plugging (56) and (53) into expected value of (52), and then (52) and (42) back into (38) to obtain, after using the Tower property of conditional expectations:

C.2.2 Convergence Lemma

Lemma C.3 allows us to prove a result on the mean of KK consecutive noisy coordinate-wise gradient updates, by simply summing it and rewriting the terms. This gives Lemma C.4, which is the key lemma of our proof.

The term F(wˉt)−F(w∗)F(\bar{w}^{t})-F(w^{*}) essentially remains in the inequality due to the composite nature of FF. When ψ=0\psi=0, MM-component-smoothness of f(⋅;d)f(\cdot;d) (for d∈Xd\in\mathcal{X}) gives

and the result of Lemma C.4 further simplifies as:

Summing Lemma C.3 for k=0k=0 to k=Kk=K and w=w∗w=w^{*}, taking expectation with respect to all choices of coordinate and random noise and using the tower property gives:

As wˉt+1=1K∑k=1Kθk\bar{w}^{t+1}=\frac{1}{K}\sum_{k=1}^{K}\theta^{k}, the convexity of FF gives F(wˉt+1)≤1K∑k=1KF(θk)−F(w∗)F(\bar{w}^{t+1})\leq\frac{1}{K}\sum_{k=1}^{K}F(\theta^{k})-F(w^{*}). Plugging this inequality into (65) and combining the result with (64) gives

We conclude the proof by using the fact that Γj=Mj−1\Gamma_{j}=M_{j}^{-1} for all j∈[p]j\in[p], thus ∥⋅∥Γ=∥⋅∥M−1\left\lVert\cdot\right\rVert_{\Gamma}=\left\lVert\cdot\right\rVert_{M^{-1}} and ∥⋅∥Γ−1=∥⋅∥M\left\lVert\cdot\right\rVert_{\Gamma^{-1}}=\left\lVert\cdot\right\rVert_{M}. ∎

C.2.3 Convex Case

Setting K=RMpnϵ∥L∥M−18log⁡(1/δ)K=\frac{R_{M}\sqrt{p}n\epsilon}{\left\lVert L\right\rVert_{M^{-1}}\sqrt{8\log(1/\delta)}} yields:

In the convex case, we iterate only once in the inner loop (since T=1T=1). As such, wpriv=wˉ1w^{priv}=\bar{w}^{1}, and applying Lemma C.4 with wˉt+1=wˉ1\bar{w}^{t+1}=\bar{w}^{1}, wt=wˉ0w^{t}=\bar{w}^{0} and σj\sigma_{j} chosen as in Theorem 3.1 gives the result. Taking K=RMpnϵ∥L∥M−18log⁡(1/δ)K=\frac{R_{M}\sqrt{p}n\epsilon}{\left\lVert L\right\rVert_{M^{-1}}\sqrt{8\log(1/\delta)}} then gives

and the result follows from 28+128≈8.48<92\sqrt{8}+\frac{12}{\sqrt{8}}\approx 8.48<9. ∎

C.2.4 Strongly Convex Case

Setting T=log⁡2(32n2ϵ2(F(wˉ0)−F(w∗))p(1+1/μM)∥L∥M−12log⁡(1/δ))T=\log_{2}\left(\frac{32n^{2}\epsilon^{2}(F(\bar{w}^{0})-F(w^{*}))}{p(1+1/\mu_{M})\left\lVert L\right\rVert_{M^{-1}}^{2}\log(1/\delta)}\right) yields:

It remains to set K=2p(1+1/μM)K=2p(1+1/\mu_{M}) to obtain

Recursive application of this inequality gives

where we upper bound the sum by the value of the complete series. It remains to replace ∥σ∥M2\left\lVert\sigma\right\rVert_{M}^{2} by its value to obtain the result. Taking T=log⁡2((F(wˉ0)−F(w∗))n2ϵ224p(1+1/μM)∥L∥M−12log⁡(1/δ))T=\log_{2}\left(\frac{(F(\bar{w}^{0})-F(w^{*}))n^{2}\epsilon^{2}}{24p(1+1/\mu_{M})\left\lVert L\right\rVert_{M^{-1}}^{2}\log(1/\delta)}\right) then gives

C.3 Proof of Remark 1

where the expectation is taken over the random noise ηj\eta_{j}, and −γjgj=prox⁡γjψj(θjk−γj(∇jf(θk)+ηj))−θjk-\gamma_{j}g_{j}=\operatorname{prox}_{\gamma_{j}\psi_{j}}(\theta^{k}_{j}-\gamma_{j}(\nabla_{j}f(\theta^{k})+\eta_{j}))-\theta^{k}_{j} as defined in the analysis of Algorithm 1. We need to link the proximal operator we use in DP-CD with the quantity VjηjV_{j}^{\eta_{j}} that we just defined:

which can be rewritten as Vj(θk,−γjgj)≤min⁡tVj(θk,t)+⟨ηj,γj(gj−gj∗)⟩V_{j}(\theta^{k},-\gamma_{j}g_{j})\leq\min_{t}V_{j}(\theta^{k},t)+\left\langle\eta_{j},\gamma_{j}(g_{j}-g_{j}^{*})\right\rangle. Taking the expectation yields

Finally, we remark that ∣gj−gj∗∣≤∣γjηj∣\left|g_{j}-g_{j}^{*}\right|\leq\left|\gamma_{j}\eta_{j}\right| and the non-expansiveness of the proximal operator gives

When the objective function FF is convex, we use Lemma C.4 to obtain, since ∥σ∥M−12=βp\left\lVert\sigma\right\rVert_{M^{-1}}^{2}=\beta p,

Therefore, when FF is convex, we get F(w1)−F(w∗)≤ξF(w^{1})-F(w^{*})\leq\xi, for ξ>βp\xi>\beta p, as long as 2pRM2K≤ξ−βp\frac{2pR_{M}^{2}}{K}\leq\xi-\beta p, that is K≥2pRM2ξ−βpK\geq\frac{2pR_{M}^{2}}{\xi-\beta p}.

In comparison, [51, Theorem 5.1 therein ] gives convergence to ξ>2pRM2β\xi>\sqrt{2pR_{M}^{2}\beta} when K≥2pRM2ξ−2pRM2βK\geq\frac{2pR_{M}^{2}}{\xi-\sqrt{2pR_{M}^{2}\beta}}. We thus gain a factor βp/2RM2\sqrt{\beta p/2R_{M}^{2}} in utility. Importantly, our utility upper bound does not depend on initialization in that setting, whereas the one of does.

When the objective function FF is μM\mu_{M}-strongly-convex w.r.t. to ∥⋅∥M\left\lVert\cdot\right\rVert_{M}, then from (75) we obtain, as long as K≥4/μMK\geq 4/\mu_{M}, that

Appendix D Comparison with DP-SGD

We start by the scenario where coordinate-wise smoothness constants are balanced and all equal to M=M1=⋯=MpM=M_{1}=\cdots=M_{p}. We observe that

We then consider the convex and strongly-convex functions separately:

Convex functions: it holds that RM=MRIR_{M}=\sqrt{M}R_{I}, which yields the equality ∥L∥M−1RM=∥L∥2RI\left\lVert L\right\rVert_{M^{-1}}R_{M}=\left\lVert L\right\rVert_{2}R_{I}.

which means that ff is MμMM\mu_{M}-strongly-convex with respect to ∥⋅∥2\left\lVert\cdot\right\rVert_{2}. This gives ∥L∥M−12μM=∥L∥22/MμI/M=∥L∥22μI\frac{\left\lVert L\right\rVert_{M^{-1}}^{2}}{\mu_{M}}=\frac{\left\lVert L\right\rVert_{2}^{2}/M}{\mu_{I}/M}=\frac{\left\lVert L\right\rVert_{2}^{2}}{\mu_{I}}.

In light of the results summarized in Table 1, it remains to compare ∥L∥2=∑j=1pLj2\left\lVert L\right\rVert_{2}=\sqrt{\sum_{j=1}^{p}L_{j}^{2}} with Λ\Lambda, for which it holds that Λ≤∑j=1pLj2≤pΛ\Lambda\leq\sqrt{\sum_{j=1}^{p}L_{j}^{2}}\leq\sqrt{p}\Lambda, which is our result.

When smoothness constants are disparate, we discuss the case where

the first coordinate of wˉ0\bar{w}^{0} is already very close to its optimal value so that M1∣wˉ10−w1∗∣≪∑j≠1Mj∣wˉj0−wj∗∣M_{1}\left|\bar{w}^{0}_{1}-w^{*}_{1}\right|\ll\sum_{j\neq 1}M_{j}\left|\bar{w}^{0}_{j}-w^{*}_{j}\right|. Under this hypothesis,

which implies that when ff is μM\mu_{M} strongly-convex with respect to ∥⋅∥M\left\lVert\cdot\right\rVert_{M}, it is Mmin⁡μMM_{\min}\mu_{M} strongly-convex with respect to ∥⋅∥2\left\lVert\cdot\right\rVert_{2}. This yields, under our hypotheses, ∥L∥M−12μM≈Λ2/Mmax⁡μI/Mmin⁡=Mmin⁡Mmax⁡Λ2μI\frac{\left\lVert L\right\rVert_{M^{-1}}^{2}}{\mu_{M}}\approx\frac{\Lambda^{2}/M_{\max}}{\mu_{I}/M_{\min}}=\frac{M_{\min}}{M_{\max}}\frac{\Lambda^{2}}{\mu_{I}}. In both cases, DP-CD can get arbitrarily better than DP-SGD, and gets better as the ratio Mmax⁡/Mmin⁡{M_{\max}}/{M_{\min}} increases.

The two hypotheses we describe above are of course very restrictive. However, it gives some insight about when and why DP-CD can outperform DP-SGD. Our numerical experiments in Section 6 confirm this analysis, even in less favorable cases.

Appendix E Proof of Lower Bounds

To prove lower bounds on the utility of LL-component-Lipschitz functions, we extend the proof of to our setting (that is, LL-component-Lipschitz functions and unconstrained composite optimization). There are three main difficulties in adapting their proof:

Second, Lemma 5.1 of must be extended to our LL-component-Lipschitz setting. To do so, we consider datasets with points in ∏j=1p{−Lj,Lj}\prod_{j=1}^{p}\{-L_{j},L_{j}\} rather than {−1/p,1/p}p\{-1/\sqrt{p},1/\sqrt{p}\}^{p}, and carefully adapt the construction of the dataset DD so that ∥∑i=1ndi∥2=Ω(min⁡(n∥L∥2,p∥L∥2/ϵ))\left\lVert\sum_{i=1}^{n}d_{i}\right\rVert_{2}=\Omega(\min(n\left\lVert L\right\rVert_{2},{\sqrt{p}\left\lVert L\right\rVert_{2}}/{\epsilon})), which is essential to prove our lower bounds.

Third, the lower bounds of rely on fingerprinting codes, and in particular on the result of which uses such codes to prove that (when nn is smaller than some n∗n^{*} we describe later) differential privacy is incompatible with precisely and simultaneously estimating all pp counting queries defined over the columns of the dataset DD. In our construction, since all columns of DD now have different scales, we need an additional hypothesis on the repartition of the LjL_{j}’s, (i.e., that ∑j∈JLj2=Ω(∥L∥2)\sum_{j\in\mathcal{J}}L_{j}^{2}=\Omega(\left\lVert L\right\rVert_{2}) for all J⊆[p]\mathcal{J}\subseteq[p] of a given size), which is not required in existing lower bounds (where all columns have equal scale).

We start our proof by recalling and extending to our setting the notions of counting queries (Definition E.1) and accuracy (Definition E.2), as described by . The main feature of our definitions is that we allow the set X\mathcal{X} to have different scales for each of its coordinates, and that we account for this scale in the definition of accuracy. We denote by conv⁡(X)\operatorname{conv}(\mathcal{X}) the convex hull of a set X\mathcal{X}.

Let n>0n>0. A counting query on X\mathcal{X} is a function q:Xn→conv⁡(X)q:\mathcal{X}^{n}\rightarrow\operatorname{conv}(\mathcal{X}) defined using a predicate q:X→Xq:\mathcal{X}\rightarrow\mathcal{X}. The evaluation of the query qq over a dataset D∈Xn\mathcal{D}\in\mathcal{X}^{n} is defined as the arithmetic mean of qq on D\mathcal{D}:

In our proof, we will use a specific class of queries: one-way marginals (Definition E.3), that compute the arithmetic mean of a dataset along one of its column.

Let X=∏j=1p{−Lj;Lj}\mathcal{X}=\prod_{j=1}^{p}\{-L_{j};L_{j}\} or X={0,Lj}p\mathcal{X}=\{0,L_{j}\}^{p}. The family of one-way marginals on X\mathcal{X} is defined by queries with predicates qj(x)=xjq_{j}(x)=x_{j} for x∈Xx\in\mathcal{X}. For a dataset D∈XnD\in\mathcal{X}^{n} of size nn, we thus have qj(D)=1n∑i=1ndi,jq_{j}(D)=\frac{1}{n}\sum_{i=1}^{n}d_{i,j}.

E.2 Lower Bound for One-Way Marginals

We can now restate a key result from , which shows that there exists a minimal number n∗n^{*} of records needed in a dataset to allow achieving both accuracy and privacy on the estimation of one-way marginals on X=({0,1}p)n\mathcal{X}=(\{0,1\}^{p})^{n}. This lemma relies on the construction of re-identifiable distribution (see [11, Definition 2.10]). One can then use this distribution to find a dataset on which a private algorithm can not be accurate (see [11, Lemma 2.11]).

For ϵ>0\epsilon>0 and p>0p>0, there exists a number n∗=Ω(pϵ)n^{*}=\Omega(\frac{\sqrt{p}}{\epsilon}) such that for all n≤n∗n\leq n^{*}, there exists no algorithm that is both (1/3,1/75)(1/3,1/75)-accurate and (ϵ,o(1n))(\epsilon,o\left(\frac{1}{n}\right))-differentially private for the estimation of one-way marginals on ({0,1}p)n(\{0,1\}^{p})^{n}.

To leverage this result in our setting of private empirical risk minimization, we start by extending it to queries on X=∏j=1p{−Lj;Lj}\mathcal{X}=\prod_{j=1}^{p}\{-L_{j};L_{j}\}. Before stating the main theorem of this section (Theorem E.5), we describe a procedure χL:({0,1}p)n→X3n\chi_{L}:(\{0,1\}^{p})^{n}\rightarrow\mathcal{X}^{3n} (with L1,…,Lp>0L_{1},\dots,L_{p}>0), that takes as input a dataset D∈({0,1}p)nD\in(\{0,1\}^{p})^{n} and outputs an augmented and rescaled version. This procedure is crucial to our proof and is defined as follows. First, it adds 2n2n rows filled with 11’s to DD, which ensures that the sum of each column of DD is Θ(n)\Theta(n) (which gives the lower bound on MM in Theorem E.5). Then it rescales each of these columns by subtracting 1/21/2 to each coefficient and multiplying the jj-th column of DD (j∈[p]j\in[p]) by 2Lj2L_{j}. The resulting dataset DLaug=χL(D)D_{L}^{aug}=\chi_{L}(D) is a set of 3n3n points with values in X=∏j=1p{−Lj,Lj}\mathcal{X}=\prod_{j=1}^{p}\{-L_{j},L_{j}\}, with the property that, for all j∈[p]j\in[p], 3nLj≥∑i=1n(DLaug)i,j≥nLj3nL_{j}\geq\sum_{i=1}^{n}(D_{L}^{aug})_{i,j}\geq nL_{j}. For D∈({0,1}p)nD\in(\{0,1\}^{p})^{n}, we show how to reconstruct qj(χL(D))q_{j}(\chi_{L}(D)) from qj(D)q_{j}(D) in 1.

where we use the slight abuse of notation by denoting the one-way marginals qj:X3n→conv⁡(X)q_{j}:\mathcal{X}^{3n}\rightarrow\operatorname{conv}(\mathcal{X}) and qj:({0,1}p)n→pq_{j}:(\{0,1\}^{p})^{n}\rightarrow^{p} in the same way.

Let D∈({0,1}p)nD\in(\{0,1\}^{p})^{n}, and let Daug∈({0,1}p)3nD^{aug}\in(\{0,1\}^{p})^{3n} constructed by adding 2n2n rows of 11’s at the end of DD. Let DLaug=χL(D)D_{L}^{aug}=\chi_{L}(D). We remark that

Then, we link qj(Daug)q_{j}(D^{aug}) with qj(DLaug)q_{j}(D^{aug}_{L}):

combining (97) and (98) gives the result. ∎

Let M=Ω(min⁡(n∥L∥2,p∥L∥2ϵ))M=\Omega\left(\min\left(n\left\lVert L\right\rVert_{2},\frac{\sqrt{p}\left\lVert L\right\rVert_{2}}{\epsilon}\right)\right), and define the set of queries Q\mathcal{Q} composed of pp queries qj(D)=1n∑i=1ndi,jq_{j}(D)=\frac{1}{n}\sum_{i=1}^{n}d_{i,j} for j∈[p]j\in[p]. Let A\mathcal{A} be a (ϵ,δ)(\epsilon,\delta)-differentially-private randomized algorithm. Let α,β∈\alpha,\beta\in. We will show that there exists a dataset DD such that ∥∑i=1ndi∥2∈[M−1,M+1]\left\lVert\sum_{i=1}^{n}d_{i}\right\rVert_{2}\in[M-1,M+1] for which A(D)\mathcal{A}(D) is not (α,β)(\alpha,\beta)-accurate.

Assume, for the sake of contradiction, that A:X3n→conv⁡(X)\mathcal{A}:\mathcal{X}^{3n}\rightarrow\operatorname{conv}(\mathcal{X}) is (13α,β)(\tfrac{1}{3}\alpha,\beta)-accurate for Q\mathcal{Q}. Then, for each dataset D′∈X3nD^{\prime}\in\mathcal{X}^{3n}, we have

Importantly, for all D∈({0,1})p)nD\in(\{0,1\})^{p})^{n}, the randomized algorithm A\mathcal{A} satisfies (100) for the dataset DLaug=χL(D)∈X3nD_{L}^{aug}=\chi_{L}(D)\in\mathcal{X}^{3n}. We now construct the mechanism A~:({0,1}p)n→p\widetilde{\mathcal{A}}:(\{0,1\}^{p})^{n}\rightarrow^{p} that takes a dataset D∈({0,1}p)nD\in(\{0,1\}^{p})^{n}, constructs DLaug=χL(D)D_{L}^{aug}=\chi_{L}(D) and runs A\mathcal{A} on it. It then outputs A~(D)\widetilde{\mathcal{A}}(D) such that, for j∈[p]j\in[p], A~j(D)=32LjAj(DLaug)−Lj3\widetilde{\mathcal{A}}_{j}(D)=\frac{3}{2L_{j}}\mathcal{A}_{j}(D_{L}^{aug})-\frac{L_{j}}{3}. Using 1, the results of A~\widetilde{\mathcal{A}} and be linked to the ones of A\mathcal{A}, as

Therefore, if A\mathcal{A} satisfies (100) and (101), then A~:({0,1}p)n→p\widetilde{\mathcal{A}}:(\{0,1\}^{p})^{n}\rightarrow^{p} satisfies, for all D∈({0,1}p)nD\in(\{0,1\}^{p})^{n},

which is exactly the definition of (α,β)(\alpha,\beta)-accuracy for A~\widetilde{\mathcal{A}}. Remark that since A~\widetilde{\mathcal{A}} is only a post-processing of A\mathcal{A}, without additional access to the dataset itself, A~\widetilde{\mathcal{A}} is itself (ϵ,δ)(\epsilon,\delta)-differentially-private. We have thus constructed an algorithm that is both accurate and private for n≤n∗n\leq n^{*}, which contradicts the result of Lemma E.4 when β=175\beta=\frac{1}{75}. This proves the existence of a dataset D∈({0,1}p)nD\in(\{0,1\}^{p})^{n} such that for DLaug=χL(D)D_{L}^{aug}=\chi_{L}(D), A(DLaug)\mathcal{A}(D_{L}^{aug}) is not (13α,β)(\tfrac{1}{3}\alpha,\beta)-accurate on Q\mathcal{Q}, which means that with probability at least 1/31/3, there exists a subset J⊆[p]\mathcal{J}\subseteq[p] of cardinal ∣J∣≥⌈βp⌉\left|\mathcal{J}\right|\geq\lceil\beta p\rceil such that

where the second inequality comes from the fact that ∣J∣≥⌈βp⌉=⌈p75⌉\left|\mathcal{J}\right|\geq\lceil\beta p\rceil=\lceil\frac{p}{75}\rceil and our hypothesis on ∑j∈JLj2\sum_{j\in\mathcal{J}}L_{j}^{2}. Notice that when L1=⋯=Lp=1pL_{1}=\cdots=L_{p}=\frac{1}{\sqrt{p}}, we recover the result of , since ∥L∥2=1\left\lVert L\right\rVert_{2}=1 it holds with probability at least 1/31/3 that

and in that case, since all LjL_{j}’s are equal, it indeed holds that ∑j∈JLj2=Ω(∥L∥2)\sqrt{\sum_{j\in\mathcal{J}}L_{j}^{2}}=\Omega(\left\lVert L\right\rVert_{2}). Finally, we remark that the sum of each column of DLaugD_{L}^{aug} is ∑i=1ndi,j≥nLj\sum_{i=1}^{n}d_{i,j}\geq nL_{j}, and as such, we have ∥∑i=1ndi∥2=∑j=1p(∑i=1ndi,j)2≥∑j=1pn2Lj2=n∥L∥2\left\lVert\sum_{i=1}^{n}d_{i}\right\rVert_{2}=\sqrt{\sum_{j=1}^{p}(\sum_{i=1}^{n}d_{i,j})^{2}}\geq\sqrt{\sum_{j=1}^{p}n^{2}L_{j}^{2}}=n\left\lVert L\right\rVert_{2}.

We get the result in that case by augmenting the dataset D∗D^{*} that we constructed in the first part of this proof. To do so, we follow the steps described by in the proof of their Lemma 5.1. The construction consists in choosing a vector c∈Xc\in\mathcal{X}, and adding ⌈n−n∗2⌉\lceil\frac{n-n^{*}}{2}\rceil rows with cc, and ⌊n−n∗2⌋\lfloor\frac{n-n^{*}}{2}\rfloor rows with −c-c to the dataset D∗D^{*}. This results in a dataset D′D^{\prime} such that ∥∑i=1ndi∥=Ω(n∗∥L∥2)=Ω(p∥L∥2ϵ)\left\lVert\sum_{i=1}^{n}d_{i}\right\rVert=\Omega(n^{*}\left\lVert L\right\rVert_{2})=\Omega(\frac{\sqrt{p}\left\lVert L\right\rVert_{2}}{\epsilon}), since the contributions of rows −c-c and cc (almost) cancel out. The theorem follows from observing that (n∗nα,β)(\frac{n^{*}}{n}\alpha,\beta)-accuracy on this augmented dataset implies (α,β)(\alpha,\beta)-accuracy on the original dataset. As such, if an algorithm is both private and (n∗nα,β)(\frac{n^{*}}{n}\alpha,\beta)-accurate on the dataset D′D^{\prime}, we get a contradiction, which gives the theorem as n∗n=pnϵ\frac{n^{*}}{n}=\frac{\sqrt{p}}{n\epsilon}. ∎

Without the assumption on the distribution of the LjL_{j}’s, we can still get an inequality that resembles (103): ∥A(DLaug)−q(DLaug)∥2≥\eqrefthm:lower−bound−one−way−marginals:accuracy−A∑j∈J4Lj29≥227Lmin⁡Lmax⁡∥L∥2\left\lVert\mathcal{A}(D_{L}^{aug})-q(D_{L}^{aug})\right\rVert_{2}\overset{\eqref{thm:lower-bound-one-way-marginals:accuracy-A}}{\geq}\sqrt{\sum_{j\in\mathcal{J}}\frac{4L_{j}^{2}}{9}}\geq\frac{2}{27}\frac{L_{\min}}{L_{\max}}\left\lVert L\right\rVert_{2}, with probability at least 1/31/3, and we get a result similar to Theorem E.5, except with an additional multiplicative factor Lmin⁡/Lmax⁡L_{\min}/L_{\max}.

E.3 Lower Bound for Convex Functions

To find the solution of (105), we look for w∗w^{*} so that the objective’s gradient is zero, that is

so that ∥w∗∥2=β∥∑i=1ndi∥2∥∑i=1ndi∥2=β\left\lVert w^{*}\right\rVert_{2}=\frac{\beta}{\left\lVert\sum_{i=1}^{n}d_{i}\right\rVert_{2}}\left\lVert\sum_{i=1}^{n}d_{i}\right\rVert_{2}=\beta. To prove the lower bound, we remark that

At this point, we can proceed similarly to to relate this quantity to private estimation of one-way marginals. We let M=Ω(min⁡(n∥L∥2,∥L∥2p/ϵ))M=\Omega(\min(n\left\lVert L\right\rVert_{2},\left\lVert L\right\rVert_{2}\sqrt{p}/\epsilon)) and A\mathcal{A} be an (ϵ,δ)(\epsilon,\delta)-differentially private mechanism that outputs a private solution wprivw^{priv} to (105). Suppose, for the sake of contradiction, that for every dataset DD with ∥∑i=1ndi∥2∈[M−1;M+1]\left\lVert\textstyle{\sum_{i=1}^{n}d_{i}}\right\rVert_{2}\in[M-1;M+1], it holds with probability at least 2/32/3 that

We now derive from A\mathcal{A} a mechanism A~\widetilde{\mathcal{A}} to estimate one-way marginals. To do this, A~\widetilde{\mathcal{A}} runs A\mathcal{A} to obtain wprivw^{priv} and outputs Mnβwpriv\frac{M}{n\beta}w^{priv}. We obtain that with probability at least 2/32/3,

where q(D)=1n∑i=1ndiq(D)=\frac{1}{n}\sum_{i=1}^{n}d_{i}. This is in contradiction with Theorem E.5. We thus proved that ∥wpriv−w∗∥=Ω(β)\left\lVert w^{priv}-w^{*}\right\rVert=\Omega(\beta), with probability at least 1/31/3. As a consequence, we now obtain that with probability at least 1/31/3,

which gives the desired result on the expectation of F(wpriv;D)−F(w∗;D)F(w^{priv};D)-F(w^{*};D).

with probability at least 1/31/3, which is in contradiction with Remark E.6. We thus get an additional factor of Lmin⁡/Lmax⁡L_{\min}/L_{\max} in the lower bound:

E.4 Lower Bound for Strongly-Convex Functions

To prove a lower bound for strongly-convex functions, we let μI>0\mu_{I}>0, L1,…,Lp>0L_{1},\dots,L_{p}>0, W=∏j=1p[−Lj2μI,+Lj2μI]\mathcal{W}=\prod_{j=1}^{p}[-\frac{L_{j}}{2\mu_{I}},+\frac{L_{j}}{2\mu_{I}}] and D={d1,…,dn}∈∏j=1p{±Lj2μI}D=\{d_{1},\dots,d_{n}\}\in\prod_{j=1}^{p}\{\pm\frac{L_{j}}{2\mu_{I}}\}. We consider the following problem, which fits in our setting:

It remains to apply Theorem E.5 to obtain that, with probability at least 1/31/3,

which gives the lower bound on the expected value of F(wpriv;D)−F(w∗)F(w^{priv};D)-F(w^{*}). Note that without the additional assumption on the distribution of the LjL_{j}’s, Remark E.6 directly gives the result with an additional multiplicative factor (Lmin⁡/Lmax⁡)2(L_{\min}/L_{\max})^{2}:

Appendix F Private Estimation of Smoothness Constants

Assuming that the practitioner knows an approximate upper bound bjb_{j} over the Mj(i)M_{j}^{(i)}’s, they can enforce it by clipping Mj(i)M_{j}^{(i)} to bjb_{j} for each i∈[n]i\in[n]. The sensitivity of the average of the clipped Mj(i)M_{j}^{(i)}’s is thus 2bj/n2b_{j}/n. One can then compute an estimate of M1,…,MpM_{1},\dots,M_{p} under ϵ\epsilon-DP using the Laplace mechanism as follows:

where the factor pp in noise scale comes from using the simple composition theorem , and Lap(λ)\text{Lap}(\lambda) is a sample drawn in a Laplace distribution of mean zero and scale λ\lambda. The computed constant can then directly be used in DP-CD, allocating the remaining budget ϵ−ϵ′\epsilon-\epsilon^{\prime} to the optimization procedure.

Appendix G Additional Experimental Details and Results

We simultaneously tune these three hyperparameters for each algorithm across the following grid:

step size: 10 logarithmically-spaced values between 10−610^{-6} and 11 for DP-SGD, and between 10−210^{-2} and 1010 for DP-CD.Recall that step sizes for CD algorithms are coordinate-wise, and thus larger than in SGD algorithms. We empirically verify that the best step size always lies strictly inside the considered interval for both DP-CD and DP-SGD.

clipping threshold: 100 logarithmically-spaced values, between 10−310^{-3} and 10610^{6}.

number of passes: 5 values (2, 5, 10, 20 and 50).

We run each algorithm on each dataset 5 times on each combination of hyperparameter values. We then keep the set of hyperparameters that yield the lowest value of the objective at the last iterate, averaged across the 55 runs.

In Table 2, we report the best relative error (in comparison to optimal objective value) at the last iterate, averaged over five runs, for each dataset, algorithm, and total number of passes on the data. As such, each cell of this table corresponds to the best value obtained after tuning the step size and clipping hyperparameters for a given number of passes.

G.2 Running Time

In this section, we report the running times of DP-CD and DP-SGD. We implemented DP-CD and DP-SGD in C++, with Python bindingsThe code is available at https://gitlab.inria.fr/pmangold1/private-coordinate-descent/.. The design matrix and the labels are kept in memory as dense matrices of the Eigen library. No special code optimization nor tricks is applied to the algorithms, except for the update of residuals at each iteration of DP-CD, which prevents from accessing the complete dataset at each step. All experiments were run on a laptop with 16GB of RAM and an Intel(R) Core(TM) i7-10610U CPU @ 1.80GHz.

Figure 3 shows the same experiments as in Figure 1 and Figure 2, but as a function of the running time. In our implementation, DP-CD runs about 44 times as fast as DP-SGD for a given number of iterations (see Figure 3(a) and Figure 3(b) for 5050 iterations). On the three other plots, Figure 3(c), Figure 3(d) and Figure 3(e), DP-CD yields better results in less iterations. DP-CD is thus particularly valuable in these scenarios: combined with its faster running time, it provides accurate results extremely fast. For completeness, we provide in Table 3 the full table of running time, corresponding to Table 2 and Figure 3. These results show that, for a given number of passes on the data, DP-CD consistently runs about 55 times faster than DP-SGD.