The Discrete Gaussian for Differential Privacy

Introduction

Differential Privacy [DMNS06] provides a rigorous standard for ensuring that the output of an algorithm does not leak the private details of individuals contained in its input. A standard technique for ensuring differential privacy is to evaluate a function on the input and then add a small amount of random noise to the result before releasing it. Specifically, it is common to add noise drawn from a Laplace or Gaussian distribution, which is scaled according to the sensitivity of the function – i.e., how much one person’s data can change the function value. These are two of the most fundamental algorithms in differential privacy, which are used as subroutines in almost all differentially private systems. For example, differentially private algorithms for convex empirical risk minimization and deep learning are based on adding noise to gradients [BST14, ACGMMTZ16].

However, the Laplace and Gaussian distributions are both continuous over the real numbers. As such, it is not possible to even represent a sample from them on a finite computer, much less produce such a sample. One might suppose that such issues are purely of theoretical interest, and that they can be resolved in practice by simply using standard floating-point arithmetic and representations. Unfortunately, this is not the case: Mironov [Mir12] demonstrated that the naïve use of finite-precision approximations can result in catastrophic failures of privacy. In particular, by examining the low-order bits of the noisy output, the noiseless value can often be determined. Mironov demonstrated that this information allows the entire input dataset can be rapidly reconstructed, while only a negligible privacy loss is recorded by the system. Despite this demonstration, the flawed methods continue to appear in open source implementations of differentially private mechanisms. This demonstrates a real need for us to provide safe and practical solutions to enable the deployment of differentially private systems in real-world privacy-critical settings. In this work, we carefully consider how to securely implement these basic differentially private methods on finite computers that cannot faithfully represent real numbers.

One solution to this problem involves instead sampling from a discrete distribution that can be sampled on a finite computer. For many natural queries, the output of the function to be computed is naturally discrete – e.g., counting how many records in a dataset satisfy some predicate – and hence there is no loss in accuracy when adding discrete noise to it. Otherwise, the function value must be rounded before adding noise.

The (continuous) Gaussian distribution has many advantages over the (continuous) Laplace distribution (and also some disadvantages), making it better suited for many applications. For example, the Gaussian distribution has lighter tails than the Laplace distribution. In settings with a high degree of composition – i.e., answering many queries with independent noise, rather than a single query – the scale (e.g., variance) of Gaussian noise is also lower than the scale of Laplace noise required for a comparable privacy guarantee. The privacy analysis under composition of Gaussian noise addition is typically simpler and sharper; in particular, these privacy guarantees can be cleanly expressed in terms of concentrated differential privacy (CDP) [DR16, BS16] and related variants of differential privacy [Mir17, BDRS18, DRS19]. (See Section 4 for further discussion.)

Thus, it is natural to wonder whether a discretization of the Gaussian distribution retains the privacy and utility properties of the continuous Gaussian distribution, as is the case for the Laplace distribution. In this paper, we show that this is indeed the case.

Our investigations focus on three aspects of the discrete Gaussian: privacy, utility, and sampling. In summary, we demonstrate that the discrete Gaussian provides the same level of privacy and utility as the continuous Gaussian. We also show that it can be efficiently sampled on a finite computer, thus addressing the shortcomings of continuous distributions discussed earlier. Along the way, we both prove and empirically demonstrate a number of additional properties of the discrete Gaussian of interest and useful to those deploying it for privacy purposes (and even otherwise). We proceed to elaborate on our contributions and findings.

Utility.

Sampling.

We can practically sample a discrete Gaussian on a finite computer. We present a simple and efficient exact sampling procedure that only requires access to uniformly random bits and does not involve any real-arithmetic operations or non-trivial function evaluations (Algorithm 3). As there are previous methods (see, e.g., Karney’s algorithm [Kar16], which was an inspiration for our work, and the more recent work of Du, Fan, and Wei [DFW20]), we do not consider this to be one of our primary contributions. Nonetheless, we include these results as we consider our methods to be simpler, and in order to make our paper self-contained; we also provide open source code implementing our algorithm [Dga]. Our results on how to sample are provided in Section 5.

We provide a thorough comparison between the discrete Gaussian and the discrete Laplace distribution in Section 4. This includes statements of the privacy and utility guarantees for the discrete Laplace, and discussing its performance under composition in depth.

On a technical note, while the takeaway of many of our conclusions is that the discrete and continuous Gaussian are qualitatively similar, we comment that such statements are non-trivial to prove, in particular relying upon methods such as the Poisson summation formula and Fourier analysis. For instance, even basic statements on the stability property of Gaussians under linear combinations do not hold for the discrete counterpart, with their approximate versions being highly involved to establish (see, e.g., [AR16]).

2 Related Work

As originally observed and demonstrated by Mironov [Mir12], naïve implementations of the Laplace mechanism with floating-point arithmetic blatantly fail to ensure differential privacy, or any form of privacy at all. As a remedy, Mironov introduced the snapping mechanism, which serves as a safe replacement for the Laplace mechanism in the floating-point setting. The snapping mechanism performs rounding and truncation on top of the floating-point arithmetic. However, properly implementing and analyzing the snapping mechanism can be involved [Cov19], due to the idiosyncrasies of floating-point arithmetic. Furthermore, the snapping mechanism requires a compromise on privacy and accuracy, relative to what is theoretically achievable. Our methods avoid floating-point arithmetic entirely and do not compromise the privacy or accuracy guarantees.

maxnames [GMP16] gave an alternate and more general analysis of Mironov’s approach of rounding the output of an inexact sampling procedure.

maxnames [GRS12] proposed and analyzed a discrete version of the Laplace mechanism, which is also private on finite computers. However, this has heavier tails than the Gaussian and requires the addition of more noise (i.e., higher variance) than the Gaussian in settings with a high degree of composition (i.e., many queries). We provide a detailed comparison in Section 4. The US Census Bureau is planning to to protect the data collected in the 2020 US Census using discrete Laplace noise [KCKHM18, Abo18]. However, their prototype implementation [Cen18] does not use an exact sampling procedure.

A concurrent and independent work [Goo20a] analyzed what is, effectively, a truncated version of the discrete Gaussian. That work provides an almost identical sampling procedure, but a very different privacy analysis. In particular, it shows that the truncated discrete Gaussian is close to a Binomial distribution, which is, in turn, close to a rounded Gaussian. And the privacy analysis is based on this closeness. Our analysis is more direct.

Going beyond noise addition, it has been shown that private histograms [BV17] and selection (i.e., the exponential mechanism) [Ilv19] can be implemented on finite computers. (Both of these results are for pure $(\varepsilon,0)$-differential privacy.) We remark that our noise addition methods can also form the basis of an implementation of these methods. For example, instead of the exponential mechanism, we can implement the “Report Noisy Max” algorithm [DR14], which uses Laplace or Exponential [MS20] or Gumbel [Ada13] noise to perform the same task of selection.

To the best of our knowledge, there has been no work on implementing an analogue of Gaussian noise addition on finite computers. An obvious approach would be to round the output of some inexact sampling procedure. Properly analyzing this may be difficult, due to the fact that the underlying inexact Gaussian sampling procedure will be more complex than the equivalent for Laplace. Furthermore, in Section 3, we show empirically that our approach yields better utility than rounding.

In Proposition 12 and Corollary 13, we give a conversion from Rényi and concentrated differential privacy to approximate differential privacy.\AtNextCite\AtEachCitekey\@nocounterrmaxnames [ALCKS20] provide an optimal conversion from Rényi differential privacy to approximate differential privacy as well as some approximations that subsume ours. Their optimal result is, by definition, tighter than ours (but only slightly) at the expense of being more complicated and less numerically stable. See Section 2.3.

Another of our secondary contributions is a simple and efficient method for sampling from a discrete Gaussian or discrete Laplace; see Section 5.\AtNextCite\AtEachCitekey\@nocounterrmaxnames [Kar16, DFW20] also provide such algorithms. We consider our method to be simpler. In particular, our method keeps all arithmetic within the integers or the rational numbers, where exact arithmetic is possible. In contrast, Karney’s method still involves representing real numbers, but this can be carefully implemented on a finite computer using a flexible level of precision and lazy evaluation – that is, although a uniform sample from $$ requires an infinite number of bits to represent, only a finite (but a priori unbounded) number of these bits are actually needed and these can be sampled when needed. There are also methods for approximate sampling [ZSS19], but our interest is in exact sampling.

Finally, we remark that (a multivariate version of) the discrete Gaussian has been extensively studied in the context of lattice-based cryptography [GPV08, Reg09, Pei10, Ste17, etc.].

Privacy

For completeness, we state the definitions of differential privacy [DMNS06, DKMMN06] and concentrated differential privacy [DR16, BS16].

The special case of $(\varepsilon,0)$-differential privacy is referred to as pure or pointwise $\varepsilon$-differential privacy, whereas, for $\delta>0$, $(\varepsilon,\delta)$-differential privacy is referred to as approximate differential privacy.

Note that $(\varepsilon,0)$-differential privacy implies $\frac{1}{2}\varepsilon^{2}$-concentrated differential privacy and $\frac{1}{2}\varepsilon^{2}$-concentrated differential privacy implies $\left(\frac{1}{2}\varepsilon^{2}+\varepsilon\cdot\sqrt{2\log(1/\delta)},\delta\right)$-differential privacy for all $\delta>0$ [BS16].

In this section, we prove our main result on concentrated differential privacy (CDP), showing that the discrete Gaussian provides the same CDP guarantees as the continuous one.

Theorem 4 follows from Proposition 5 and Definition 3 .

Furthermore, this inequality is an equality whenever $\alpha\cdot(\mu-\nu)$ is an integer.

We remark that, like its continuous counterpart, the discrete Gaussian can also be analysed in the setting where the scale parameter $\sigma^{2}$ is data dependent [BDRS18]. This arises in the application of smooth sensitivity [NRS07, BS19].

2 Approximate Differential Privacy

In this section, we prove our main result on approximate differential privacy; namely, a tight bound on the privacy parameters achieved by the discrete Gaussian.

Furthermore, this is the smallest possible value of $\delta$ for which this is true.

This privacy guarantee matches that of the continuous Gaussian: If we replace all occurrences of the discrete Gaussian with the continuous Gaussian above, then the same result holds [BW18, Thm. 8]. Empirically, these guarantees are very close.

In Figure 1, we empirically compare the optimal $\delta$ (given by Theorem 7) to the bound attained by the corresponding continuous Gaussian, as well as this analytic upper bound (5), the standard upper bound entailed by concentrated differential privacy, and an improved upper bound via concentrated differential privacy (Corollary 13). We see that the upper bounds are reasonably tight. The discrete and continuous Gaussian attain almost identical guarantees for large $\sigma$, but the discretization creates a small difference that becomes apparent for small $\sigma$.

To prove Theorem 7, we introduce the privacy loss random variable [DRV10, DR16, BS16] and relate it to approximate differential privacy.In the information theory literature, the term “relative information spectrum” is sometimes used for the distribution of what we call the privacy loss random variable [SV16, Liu18].

Approximate differential privacy can also be characterized via the privacy loss as follows. This characterization is implicit in the work of Bun and Steinke [BS16, Lemma B.2] and is explicit in the work of Meiser and Mohammadi [MM18, Lemma 1] (see also [Goo20, Observation 2] and references therein).

Let $\varepsilon,\delta\geq 0$. Let $M\colon\mathcal{X}^{n}\to\mathcal{Y}$ be a randomized algorithm. Then $M$ satisfies $(\varepsilon,\delta)$-differential privacy if and only if

for all $x,x^{\prime}\in\mathcal{X}^{n}$ differing on a single element.

Observe that, by Markov’s inequality, for all $\alpha>1$, it suffices to set

This is the usual expression that is used to convert bounds on the privacy loss or Rényi divergence into approximate differential privacy. Lemma 9 and Proposition 12 represent an improvement on this.

Thus, for all $E\subset\mathcal{Y}$, we have

Now it is easy to identify the worst event as $E=\{y\in\mathcal{Y}:1-e^{\varepsilon-f(y)}>0\}$. Thus

Alternatively, since the worst event is equivalently $E=\{y\in\mathcal{Y}:f(y)>\varepsilon\}$, we have

We will use Lemma 9, Equation 8. Thus our main task is to determine the distribution of the privacy loss random variable.

We note that $Y$ and $-Y$ and $Y^{\prime}-q(x^{\prime})$ and $q(x^{\prime})-Y^{\prime}$ all have the same distribution. Hence

We now provide the proofs of the two aforementioned analytical bounds for $(\varepsilon,\delta)$-differential privacy our theorem readily implies.

In the setting of Theorem 7, for $\Delta=1$, we have $\delta\leq e^{-\lfloor\varepsilon\sigma^{2}\rceil^{2}/2\sigma^{2}}/\sqrt{2\pi\sigma^{2}}\,,$ where $\lfloor\cdot\rceil$ denotes rounding to the nearest integer. More generally, $\delta\leq\sum_{k=\lceil\varepsilon\sigma^{2}/\Delta-\Delta/2\rceil}^{\lceil\varepsilon\sigma^{2}/\Delta+\Delta/2\rceil}\frac{e^{-k^{2}/2\sigma^{2}}}{\sqrt{2\pi\sigma^{2}}}$

and the result now follows from the bound on the normalization constant from Fact 19. ∎

Conversely, by a comparison series-integral, we can easily show that, for any integer $m$,

which, combined with Fact 19 on the normalization constant of the discrete Gaussian, yields

The result then follows from Theorem 7. ∎

3 Converting Concentrated Differential Privacy to Approximate Differential Privacy

We have stated guarantees for both concentrated differential privacy (Theorem 4) and approximate differential privacy (Theorem 7). Now we show how to convert from the former to the latter (Corollary 13). This is particularly useful if the discrete Gaussian is being used repeatedly and we want to provide a privacy guarantee for the composition – concentrated differential privacy has cleaner composition guarantees than approximate differential privacy. We include this result for completeness; this result was recently proved independently [ALCKS20, Lem. 1, Eq. 20].

We start with a conversion from Rényi differential privacy to approximate differential privacy.

In contrast, the standard bound [DRV10, DR16, BS16, Mir17] is $\delta\leq e^{(\alpha-1)(\tau-\varepsilon)}$. Note that $\frac{e^{(\alpha-1)(\tau-\varepsilon)}}{\alpha-1}\cdot\left(1-\frac{1}{\alpha}\right)^{\alpha}=\frac{e^{(\alpha-1)(\tau-\varepsilon)}}{\alpha}\cdot\left(1-\frac{1}{\alpha}\right)^{\alpha-1}$. Thus Proposition 12 is strictly better than the standard bound for $\alpha>1$. Equation 12 can be rearranged to

Fix neighbouring $x,x^{\prime}\in\mathcal{X}^{n}$ and let $Z\leftarrow\mathsf{PrivLoss}\left(M(x)\middle\|M(x^{\prime})\right)$. We have

We identify the smallest possible value of $c$:

where $f(z)=e^{z-\alpha\cdot z}-e^{\varepsilon-\alpha\cdot z}$. We have

Clearly $f^{\prime}(z)=0\iff e^{z}=\frac{\alpha}{\alpha-1}e^{\varepsilon}\iff z=\varepsilon-\log(1-1/\alpha)$. Thus

maxnames [ALCKS20] provide an optimal conversion from Rényi differential privacy to approximate differential privacy – i.e., an optimal version of Proposition 12. Specifically, the optimal bound is

Clearly, the expression in Proposition 12 is simpler than this. Moreover, our expression is numerically stable, whereas the alternative is unstable for small values of $\delta$. We implemented both methods and found that they yield very similar results for the parameter regime of interest, but numerical stability was a significant practical issue.

By taking the infimum over all divergence parameters $\alpha$, Proposition 12 entails the following conversion from concentrated differential privacy to approximate differential privacy.

Let $M\colon\mathcal{X}^{n}\to\mathcal{Y}$ be a randomized algorithm satisfying $\rho$-concentrated differential privacy. Then $M$ is $(\varepsilon,\delta)$-differentially private for any $\varepsilon\geq 0$ and

Corollary 13 should be contrasted with the standard bound [DRV10, DR16, BS16, Mir17] of

which holds when $\varepsilon\geq\rho>0$. \AtNextCite\AtEachCitekey\@nocounterrmaxnames [BS16] prove an intermediate bound of

For the looser expression in Corollary 13, we can analytically find an optimal $\alpha$. However, we can efficiently compute a tighter numerical bound: The equality in Equation 15 is equivalent to

Since $g$ is a smooth convex function withHere we assume $\varepsilon>\rho$, which is the setting of interest.

it has a unique minimizer $\alpha_{*}\in\left(\frac{\varepsilon+\rho}{2\rho},\max\{\frac{\varepsilon+\rho+1}{2\rho},2\}\right)$. We can find the minimizer $\alpha_{*}$ by conducting a binary search over the interval $\left(\frac{\varepsilon+\rho}{2\rho},\max\{\frac{\varepsilon+\rho+1}{2\rho},2\}\right)$. That is, we want to find $\alpha_{*}$ such that $g^{\prime}(\alpha_{*})=0$; if $\alpha<\alpha_{*}$, we have $g^{\prime}(\alpha)<0$ and, if $\alpha>\alpha_{*}$, we have $g^{\prime}(\alpha)>0$.

4 Sharp Approximate Differential Privacy Bounds for Multivariate Noise

Next we consider adding independent discrete Gaussians to a multivariate function. We begin with a concentrated differential privacy bound:

Theorem 14 follows from Proposition 5, composition of concentrated differential privacy, and Definition 3. If $\sigma_{1}=\sigma_{2}=\cdots=\sigma_{d}$, then the concentrated differential privacy guarantee depends only on the sensitivity of $q$ in the Euclidean norm; if the $\sigma_{j}$s are different, then it is a weighted Euclidean norm. Note that we only consider multivariate Gaussians with independent coordinates.

It is possible to obtain an approximate differential privacy guarantee for the multivariate discrete Gaussian from Theorem 14 and Corollary 13. While this bound is reasonably tight, we will now give an exact bound:

Fix neighbouring $x,x^{\prime}\in\mathcal{X}^{n}$. Without loss of generality, we may assume $q(x)=0$. Following the proof of Theorem 7, we will apply Lemma 9, which requires understanding the privacy loss random variable.

Then the privacy loss $Z\leftarrow\mathsf{PrivLoss}\left(M(x)\middle\|M(x^{\prime})\right)$ is given by

Substituting this expression into Equations 7 and 9 yields Equations 23 and 25 respectively.

Next we look at $Z^{\prime}\leftarrow\mathsf{PrivLoss}\left(M(x^{\prime})\middle\|M(x)\right)$, which is given by

Noting that each $Y_{j}$ has a symmetric distribution, we see that $Z^{\prime}$ has the same distribution as $Z$. Substituting these expressions into Equation 8 yields Equation 24. ∎

Theorem 14 gives three equivalent expressions for the approximate differential privacy guarantee of the multivariate discrete Gaussian. All of these expressions are in terms of the privacy loss random variable $Z\leftarrow\mathsf{PrivLoss}\left(M(x)\middle\|M(x^{\prime})\right)$. We make some remarks about evaluating these expressions:

Direct evaluation of the expressions is often impractical. Computing the distribution of $Z$ entails evaluating an infinite sum. Fortunately the terms decay rapidly, so the sum can be truncated, but this still leaves a number of terms that grows exponentially in the dimensionality $d$. Thus we must find more effective ways to evaluate the expressions.

The equivalence of the second expression (28) follows from the Poisson summation formula. When $2\pi^{2}\sigma_{j}^{2}>1/2\sigma_{j}^{2}$, then the second expression converges more rapidly; otherwise the first expression converges faster. In either case, accurately evaluating the characteristic function of the discrete Gaussian is easy.

It is then possible to compute the characteristic function of the privacy loss:

It is possible to compute the probability mass function of the privacy loss from the characteristic function:

This can form the basis of an algorithm for computing the guarantee of Theorem 14: The characteristic function can be easily computed from Equations 27, 28, and 29 and then we numerically integrate it according to Equation 30 to compute the probability distribution of the privacy loss and finally we substitute this into Equation 23.

The downside of this approach is that (i) it requires numerical integration and (ii) it only gives us the probabilities one at a time. Both of these downsides could make the procedure quite slow.

We propose to use the discrete Fourier transform (a.k.a. fast Fourier transform) to avoid these downsides.

Effectively, we will compute the distribution of $Z$ modulo $m\gamma$ for some integer $m$. (For fast computation, $m$ should be a power of two.) Call this modular random variable $Z_{m}$, so that

Rather than taking $Z_{m}$ to be supported on $\{0,\gamma,\cdots,(m-1)\gamma\}$ as is usual, we will take $Z_{m}$ to be supported on $\{(1-m/2)\gamma,(2-m/2)\gamma,\cdots,(m/2-1)\gamma,(m/2)\gamma\}$.

The inverse discrete Fourier transform allows us to compute the probability mass of $Z_{m}$ from the characteristic function of $Z$ (which is identical to the characteristic function of $Z_{m}$ at the points of interest):

Now we can compute an upper bound on the approximate differential privacy guarantee (23) using the inequality

The value of $m$ should be chosen such that this error term is tolerable. For example, if the intent is to obtain an approximate $(\varepsilon,\delta)$-differential privacy bound with $\delta=10^{-6}$, then we should choose $m$ large enough such that Equation 34 is less than, say, $10^{-9}$.

We should set $m=\frac{1}{\gamma}\cdot\left(\sqrt{8\log(1/\delta^{\prime})\sum_{j}^{d}\mu_{j}^{2}/\sigma_{j}^{2}}+\sum_{j}^{d}\mu_{j}^{2}/\sigma_{j}^{2}\right)$, where $\delta^{\prime}>0$ is the error tolerance in our final estimate of $\delta$.

To obtain lower bounds on $\delta$, we would use

The algorithm we have sketched above should be relatively efficient and numerically stable. The fast Fourier transform requires $O(m\log m)$ operations. We must evaluate the characteristic function of $Z$ at $m$ points; each evaluation requires evaluating the characteristic function of $d$ discrete Gaussians and multiplying the results together. (Of course, we must only evaluate coordinates where $\mu_{j}\neq 0$.) The characteristic function of the discrete Gaussian has a very rapidly converging series representation, so this should be close to a constant number of operations.

The discrete Fourier transform is also numerically stable, since it is a unitary operation. (Indeed this is the advantage of the characteristic function/Fourier transform over the moment generating function/Laplace transform.)

The main problem for this algorithm would be if $\gamma$ is extremely small (as the space and time used grows linearly with $1/\gamma$) or if the assumption that $\gamma$ exists fails. This depends on the choice of the parameters $\sigma_{1},\cdots,\sigma_{d}$.

Utility

We now consider how much noise the discrete Gaussian adds. As a comparison point, we consider both the continuous Gaussian and, in the interest of a fair comparison, the rounded Gaussian – i.e., a sample from the continuous Gaussian rounded to the nearest integral value. In Figure 2, we show how these compare numerically. We see that the tail of the rounded Gaussian stochastically dominates that of the discrete Gaussian. In other words, the utility of the discrete Gaussian is strictly better than the rounded Gaussian (although not by much for reasonable values of $\sigma$, i.e., those which are not very small).

To obtain analytic bounds, we begin by bounding the moment generating function:

The bound on the moment generating function shows that the discrete Gaussian is subgaussian [Riv12]. Standard facts about subgaussian random variables yield bounds on the variance and tails:

Thus the variance of the discrete Gaussian is at most that of the corresponding continuous Gaussian and we also have subgaussian tail bounds. In fact, it is possible to obtain slightly tighter bounds, showing that the variance of the discrete Gaussian is strictly less than that of the continuous Gaussian. We elaborate in the following subsections, providing tighter variance and tail bounds. However, these improvements are most pronounced for small $\sigma$, which is not the typical regime of interest for differential privacy. Nonetheless, these facts may be of independent interest.

By the Poisson summation formula [Poi, Wei],

as the second sum in the numerator is zero. ∎

As for the upper bound, it follows from a standard comparison between series and integral:

The above bounds, albeit simple to obtain, are not quite as tight as they could be. We state below a refinement, which can be found, e.g., in [Ste17, Claim 2.8.1]:

The first set of bounds is better for $\sigma\geq\frac{1}{\sqrt{2\pi}}$, and the second for $\sigma<\frac{1}{\sqrt{2\pi}}$.

The bounds obtained in Fact 20 are depicted in the figure below.

2 Tighter Variance and Tail Bounds

We now analyze the variance of the discrete Gaussian, showing that it is stricty smaller than that of the corresponding continuous Gaussian (and asymptotically the same), with a much better variance for small $\sigma$.

To prove Proposition 21 we use the following lemma which relates upper bounds on the variance of a discrete Gaussian to lower bounds on it, and vice-versa.

By applying the Poisson summation formula to both numerator and denominator of the variance, we have

where we set $\tau\coloneqq\frac{1}{2\pi\sigma}$. ∎

Next we have a lower bound on the variance.

We emphasize that the lower bound of Proposition 23 is not specific to the discrete Gaussian. It applies to any distribution $X$ such that adding $X$ to a sensitivity-1 function provides $\frac{1}{2}\varepsilon^{2}$-concentrated differential privacy.

Combining Lemma 22 with Proposition 23 (specifically, Corollary 24) yields the first claim:

Now we establish the last part of the proposition. We have $(m+1)^{2}\geq 2m+1$ and, hence,

It only remains to show that $\sum_{m=0}^{\infty}(m+1)^{2}e^{-m/\sigma^{2}}\leq 3/2$ when $\sigma^{2}\leq 1/3$. For $x\in(-1,1)$, one can show that $\sum_{m=0}^{\infty}(m+1)^{2}x^{m}=\frac{1+x}{(1-x)^{3}}$. Set $x=e^{-1/\sigma^{2}}$ to conclude. ∎

Moreover, if $\sigma\geq 1/\sqrt{2\pi}$, we have

Note that the above proposition focuses on upper tail bounds, but by symmetry of the discrete Gaussian one immediately gets similar lower tail bounds. The upshot is that, up to a small shift or $(1+o(1))$ multiplicative factor, discrete and continuous Gaussians display the same tails.

One can actually slightly refine the above upper bound, by comparing the discrete Gaussian to the rounded Gaussian $\mathcal{N}_{\text{round}}(0,\sigma^{2})$, obtained by rounding a standard continuous Gaussian to the nearest integer:

On the one hand, by definition of a rounded Gaussian, we have

by Fact 19. Similarly as before, we can write

using monotonicity of $x\mapsto e^{-x^{2}/(2\sigma^{2})}$ on $[0,\infty)$. Combining the three equations above gives

We highlight the fact that comparing with the rounded Gaussian, as the above proposition does, is meaningful, since by postprocessing any differential privacy guarantee implied by adding rounded Gaussian noise to discrete data is at least as good as that implied by adding continuous Gaussian noise to the same discrete data.

3 Other Discretizations, and Convergence to the Continuous Gaussian

In applications where query values are not naturally discrete, it is necessary to round them before adding discrete noise. A finer discretization (i.e., smaller $\alpha$) entails less error being introduced by the rounding.

Discrete Laplace

We now compare the discrete Gaussian with the most obvious alternative – the discrete Laplace. But first we give a formal definition and state some relevant facts.

The discrete Laplace (also known as the two-sided geometric) was introduced into the differential privacy literature by \AtNextCite\AtEachCitekey\@nocounterrmaxnames [GRS12], who showed that it satisfies strong optimality properties.

We remark that the discrete Laplace can also be efficiently sampled. Indeed, it is a key subroutine of our algorithm for sampling a discrete Gaussian; see Section 5.

There are two immediate qualitative differences between the discrete Laplace and the discrete Gaussian.The entire discussion in this section applies equally well to the continuous analogues of these distributions. In terms of utility, the discrete Laplace has subexponential tails (i.e., decaying as $e^{-\varepsilon m}$), whereas the discrete Gaussian has subgaussian tails (i.e., decaying as $e^{-m^{2}/2\sigma^{2}}$). In terms of privacy, the discrete Gaussian satisfies concentrated differential privacy, whereas the discrete Laplace satisfies pure differential privacy; pure differential privacy is a qualitatively stronger privacy condition than concentrated differential privacy.

Thus neither distribution dominates the other. They offer different privacy-utility tradeoffs. If the tails are important (e.g., for computing confidence intervals), then the discrete Gaussian is to be favoured. If pure differential privacy is important, then the discrete Laplace is to be favoured.

We now consider a quantitative comparison. To quantify utility, we focus on the variance of the distribution. (An alternative would be to consider the width of a confidence interval.) For now, we will quantify privacy by concentrated differential privacy. Pure $(\varepsilon,0)$-differential privacy implies $\frac{1}{2}\varepsilon^{2}$-concentrated differential privacy; thus both distributions can be evaluated on this scale.

Thus, asymptotically (i.e., for small $\varepsilon$), the discrete Gaussian has half as much variance as the discrete Laplace for the same level of privacy. In this comparison, the Gaussian clearly is better.

However, the above quantitative comparison is potentially unfair. Quantifying differential privacy by concentrated differential privacy may favour the Gaussian. If instead we demand pure $(\varepsilon,0)$-differential privacy or approximate $(\varepsilon,\delta)$-differential privacy for a small $\delta>0$, then the comparison would yield the opposite conclusion. It is fundamentally difficult to compare algorithms satisfying different versions of differential privacy, as there is no level playing field.

There is another factor to consider: A practical differentially private system will answer many queries via independent noise addition. Thus the real object of interest is the privacy and utility of the composition of many applications of noise addition.

For the rest of this section, we consider the task of answering $k$ counting queries (or sensitivity-1 queries) by adding either discrete Gaussian or discrete Laplace noise. We will measure privacy by approximate $(\varepsilon,\delta)$-differential privacy over a range of parameters. The results are summarized in Figure 4.

Concentrated differential privacy has an especially clean composition theorem [BS16]:

Let $M_{1}\colon\mathcal{X}^{n}\to\mathcal{Y}_{1}$ satisfy $\frac{1}{2}\varepsilon_{1}^{2}$-concentrated differential privacy. Let $M_{2}\colon\mathcal{X}^{n}\times\mathcal{Y}_{1}\to\mathcal{Y}_{2}$ be such that, for all $y\in\mathcal{Y}_{1}$, the restriction $M_{2}(\cdot,y)\colon\mathcal{X}^{n}\to\mathcal{Y}_{2}$ satisfies $\frac{1}{2}\varepsilon_{2}^{2}$-concentrated differential privacy. Define $M_{*}\colon\mathcal{X}^{n}\to\mathcal{Y}_{2}$ by $M_{*}(x)=M_{2}(x,M_{1}(x))$. Then $M_{*}$ satisfies $\frac{1}{2}(\varepsilon_{1}^{2}+\varepsilon_{2}^{2})$-concentrated differential privacy.

In contrast, analysing the composition of multiple invocations of discrete Laplace noise addition is not as clean. We use an optimal composition result provided by Kairouz, Oh, and Viswanath [KOV17, MV16]: The $k$-fold composition of $(\varepsilon,\delta)$-differential privacy satisfies $(\varepsilon^{\prime},\delta^{\prime})$-differential privacy if and only if

In Figure 4, we compare the discrete Gaussian and the discrete Laplace in two ways. First (on the left), we fix the utility and compare the approximate differential privacy guarantees. Specifically, we fix the task of answering $k=100$ counting queries with the noise added to each value having variance $50^{2}$. Both distributions yield different curves of $(\varepsilon,\delta)$-differential privacy guarantees and there are many points to consider. We see that, for this task, the discrete Gaussian attains better $(\varepsilon,\delta)$-differential privacy guarantees except for extremely small $\delta$ – specifically, $\delta<10^{-45}$. For $\varepsilon=1$, the discrete Gaussian provides $(1,10^{-7})$-differential privacy for this task, whereas the discrete Laplace only provides $(1,206\times 10^{-7})$-differential privacy. If we demand pure differential privacy, then the discrete Laplace provides $(2.83,0)$-differential privacy, but the discrete Gaussian cannot provide pure differential privacy. The separation becomes more pronounced as the number of queries grows.

Second (on the right of Figure 4), we fix the privacy goal to approximate $(1,10^{-6})$-differential privacy. We vary the number of counting queries (from $k=1$ to $k=100$) and measure the variance of the noise that must be added to each query answer. For a small number of queries ($k\leq 10$), the discrete Laplace gives lower variance. However, as the number of queries increases, we see that the discrete Laplace requires higher variance; for $k=100$, the variance is $69\%$ more.

Overall, Figure 4 demonstrates that the discrete Gaussian provides a better privacy-utility tradeoff than the discrete Laplace, except in two narrow parameter regimes: Either a small number of queries or if we demand something very close to pure differential privacy. We only compared variances; if we compare confidence interval sizes instead, then this would further advantage the Gaussian, which has lighter tails.

Sampling

In this section, we show how to efficiently sample exactly from a discrete Gaussian on a finite computer given access only to uniformly random bits. Such algorithms are already known [Kar16, DFW20]. However, we include this both for completeness because we believe that our algorithms are simpler than the prior work. A sample Python implementation is available online [Dga].

For simplicity, we focus our discussion of runtime only on the expected number of arithmetic operations; each such operation will take time polylogarithmic in the bit complexity of the parameters (e.g., in the representation of $\sigma^{2}$ as a rational number). We elaborate on this at the end of the section.

In Algorithm 3, we present a simple and fast algorithm for discrete Gaussian sampling, with the following guarantees:

At a high level, the idea behind the algorithm is to first sample from a discrete Laplace distribution and then “convert” this into a discrete Gaussian by rejection sampling. In order to do so, we provide two subroutines, which we believe to be of independent interest: the first, to efficiently and exactly sample from a Bernoulli with parameter $e^{-\gamma}$, for any rational parameter $\gamma\geq 0$ (Proposition 33). The second, to efficiently and exactly sample from a discrete Laplace with scale parameter $t$, for any positive integer $t$ (Proposition 34).

On input (rational) $\gamma\geq 0$, the procedure described in Algorithm 1 outputs one sample from $\mathsf{Bernoulli}(\exp(-\gamma))$, and requires a constant number of operations in expectation.

First, consider the case where $\gamma\in$. For the analysis, we let $A_{k}$ denote the value of $A$ in the $k$-th iteration of the loop in the algorithm, and $K^{\ast}$ denote the final value of $K$ upon exiting the loop. Then, for all $k\in\{0,1,2,\cdots\}$, we have

2 Sampling from a Discrete Laplace

Now we show how to efficiently and exactly sample from a discrete Laplace distribution; see Section 4 for more about this distribution. Other methods for sampling from the discrete Laplace distribution are known [SWSZW19].

Fix $p\in(0,1]$. Let $G$ be a $\mathsf{Geometric}(1-p)$ random variable, and $n\geq 1$ be an integer. Then $\left\lfloor\frac{G}{n}\right\rfloor$ is a $\mathsf{Geometric}(1-q)$ random variable for $q=p^{n}$.

With this in hand, we analyze the distribution of $Z$ conditioned on $\top$ (success), i.e., conditioned on $D=1$ and $(B,Y)\neq(1,0)$.Note that this later condition is added to the algorithm to avoid double-counting the probability that $Z=0$. Let $\alpha\coloneqq s/t$ for convenience. Recalling $B$ is independent of $D$ and that $B$ and $Y$ (conditioned on $D=1$) are independent, we have

We then bound the probability that a fixed iteration of the loop succeeds:

3 Sampling from a Discrete Gaussian

Therefore the probability that the algorithms succeeds and outputs a value in any given iteration of the loop is lower bounded by a positive constant. Thus the number of iterations of the loop follows a geometric distribution and is constant in expectation. Since, for each iteration, the expected number of operations required to sample $Y$ and $C$ is constant (by Propositions 34 and 33) the overall number of operations is constant in expectation. ∎

4 Runtime Analysis

We have stated that our algorithms require a constant number of operations in expectation. We now elaborate on this.

The runtime of our algorithms is random. Beyond showing that the number of operations is constant in expectation, it is possible to show, for all of our algorithms, that it is a subexponential random variable. We give a precise definition of this term.

Our algorithms effectively consist of a constant number of nested loops and the number of times each of them runs is subexponential. For most of our loops, they have a constant probability of terminating in each run, which means the number of times they run follows a geometric distribution, which is a subexponential random variable.

It turns out that such nested loops also have a subexpoential runtime. Specifically, one can show that, if $X_{1},\dots,X_{n},\dots$ are independent subexponential random variables and $T$ is a stopping time that is subexponential, then $\sum_{n=1}^{T}X_{n}$ is still subexponential:

Let $\alpha,\beta>1$. Suppose $(X_{n})_{1\leq n\leq\infty}$ are independent non-negative $\alpha$-subexponential random variables and $T$ is a $\beta$-subexponential stopping time. Then $S\coloneqq\sum_{n=1}^{T}X_{n}$ is $\alpha\beta$-subexponential.

We will require the following simple result:

Thus $S$ is $\alpha\beta$-subexponential. ∎

In our case, $T$ corresponds to the number of times the loop runs and $X_{n}$ corresponds to the number of operations required inside the $n$-th run of the loop. Applying the above lemma to each nested loop shows that the overall runtimes of Algorithm 1, Algorithm 2, and Algorithm 3 are all subexponential random variables.

By terminating (and outputting 0) after a pre-specified time limit, our algorithms can be made to have a deterministic runtime. However, this comes at the expense of now only satisfying approximate $(\varepsilon,\delta+\delta^{\prime})$-differential privacy or $\delta^{\prime}$-approximate $\frac{1}{2}\varepsilon^{2}$-concentrated differential privacy [BS16], where $\delta^{\prime}$ is the probability of reaching the time limit. Since the running time is roughly subexponential, this failure probability $\delta^{\prime}$ can be made astronomically small with no cost in accuracy and very little cost in runtime (i.e., only milliseconds overall). Realistically, a far greater concern than this failure probability is that the source of random bits is not perfectly uniform [GL20].

We have implemented the algorithms from Algorithms 1, 2, and 3 in Python (using the fractions.Fraction class for exact rational arithmetic and using random.SystemRandom() to obtain high-quality randomnesss). Overall, on a standard personal computer, our basic (non-optimized) implementation is able to produce over 1000 samples per second even for $\sigma^{2}=10^{100}$. The source code is available online [Dga].

Acknowledgments

We thank Shahab Asoodeh, Damien Desfontaines, Peter Kairouz, and Ananda Theertha Suresh for making us aware of several related works.

References