Differentially Private Stochastic Gradient Descent with Fixed-Size Minibatches: Tighter RDP Guarantees with or without Replacement

Authors: jeremiah  birrell
reza  ebrahimi
rouzbeh  behnia
jason  pacheco
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-15 19:51:15

Abstract: Differentially private stochastic gradient descent (DP-SGD) has been instrumental in privately training deep learning models by providing a framework to control and track the privacy loss incurred during training.  At the core of this computation lies a subsampling method that uses a privacy amplification lemma to enhance the privacy guarantees provided by the additive noise. Fixed size subsampling is appealing for its constant memory usage, unlike the variable sized minibatches in Poisson subsampling. It is also of interest in addressing class imbalance and federated learning. Current computable guarantees for fixed-size subsampling are not tight and do not consider both add/remove and replace-one adjacency relationships. We present a new and holistic Rényi differential privacy (RDP)  accountant for DP-SGD with fixed-size subsampling without replacement (FSwoR) and with replacement (FSwR). For FSwoR we consider both add/remove and replace-one adjacency, where we improve on the best current computable bound by a factor of $4$. We also show for the first time that the widely-used Poisson subsampling and FSwoR with replace-one adjacency have the same privacy to leading order in the sampling probability. Our work suggests that FSwoR is often preferable to Poisson subsampling due to constant memory usage. Our FSwR accountant includes explicit non-asymptotic upper and lower bounds and, to the  authors' knowledge, is the first such RDP analysis of fixed-size  subsampling with replacement  for DP-SGD. We analytically and empirically compare fixed size and Poisson subsampling, and show that DP-SGD gradients in a fixed-size subsampling regime exhibit lower variance in practice in addition to memory usage benefits.

Main paper: 
 introduction: 
 Differentially private stochastic gradient descent (DP-SGD) (Abadi et al., 2016) (DP-SGD) has been one of the cornerstones of privacy preserving deep learning. DP-SGD allows a so-called moments accountant technique to sequentially track privacy leakage (Abadi et al., 2016) . This technique is subsumed by Rényi differential privacy (RDP) (Mironov, 2017) , a relaxation of standard differential privacy (DP) (Dwork and Rothblum, 2016) that is widely used in private deep learning and implemented in modern DP libraries such as Opacus (Meta Platforms, 2024) and autodp (Zhu and Wang, 2019; Zhu et al., 2022) . RDP facilitates rigorous analysis when the dataset is accessed by a sequence of randomized mechanisms as in DP-SGD. While other privacy accountant frameworks such as f -DP (Dong et al., 2022) , privacy loss distributions (PLD) (Koskela et al., 2020) , Privacy Random Variable (PRV) (Gopi et al., 2021) , Analytical Fourier Accountant (AFA) (Zhu et al., 2022) , and Saddle-Point Accountant (SPA) Alghamdi et al. (2023) have been proposed to improve conversion to privacy profiles, RDP is still one of the main privacy accountants used in popular deep learning tools and applications. Thus, we expect providing tighter bounds specific to RDP will have important practical implications for the community. Each iteration of DP-SGD can be viewed as a private release of information about a stochastic gradient. DP accounting methods bound the total privacy loss incurred by applying a sequence of DP-SGD mechanisms in training. Privacy of the overall mechanism is ensured through the application of two random processes: 1) a randomized mechanism that subsamples minibatches of the training data, 2) noise applied to each gradient calculation (e.g., Gaussian). Informally, if a mechanism M is (ϵ, δ)-DP then a mechanism that subsamples with probability q ∈ (0, 1) ensures privacy (O(qϵ), qδ)-DP; a result known as the "privacy amplification lemma" (Li et al., 2012; Zhu and Wang, 2019) . Computing the privacy parameters of a mechanism requires an accounting procedure that can be nontrivial to design. Abadi et al. (2016) compute DP parameters for the special case of Gaussian noise, which was later extended to track RDP parameters by Mironov (2017) . These analyses are limited to Poisson subsampling mechanisms which produce minibatches of variable size. This creates engineering challenges in modern machine learning pipelines and also does not allow for sampling with replacement. Fixed-size minibatches are preferable as they can be aligned with available GPU memory, leading to higher learning throughput. In practice, fixed size minibathces are of interest in practical applications such as private federated learning. To address the issue, Wang et al. (2019) and Zhu et al. (2022) propose the latest fixed-size accounting methods. Wang et al. (2019) provide computable RDP bounds for FSwoR that cover general mechanisms under replace-one adjacency, while Zhu et al. (2022) give computable bounds for FSwoR under add/remove adjacency. Specifically, the study in Wang et al. (2019) provides a general formula to convert the RDP parameters of a mechanism M to RDP parameters of a subsampled (without replacement) mechanism, but the bounds they provide are only tight up to a constant factor at leading order (i.e., comparing dominant terms in the asymptotic expansion in q). Moreover, while the FSwoR results under add/remove adjacency relationship in (Zhu et al., 2022) are in parallel to our results, their accountant is not readily applicable to FSwoR under replace-one adjacency relationship.In our work, we (1) present a new and holistic RDP accountant for DP-SGD with fixed-size subsampling without replacement (FSwoR) and with replacement (FSwR) that considers both add/remove and replace-one adjacency. We note that FSwoR is equivalent to shuffling the dataset at the start of every iteration and then taking the first |B| elements, which is conveniently similar to the default behavior of common deep learning libraries (e.g. Pytorch's DataLoader iterator object) while FSwR latter is implemented, e.g., in Pytorch's RandomSampler, and is widely used to address class imbalance in classifiers (Mease et al., 2007) and biased client selection in federated learning applications (Cho et al., 2022) . (2) Our FSwoR improves on the best current computable bound under replace-one adjacency (Wang et al., 2019) by a factor of 4. (3) Our FSwR accountant includes non-asymptotic upper and lower bounds and, to the authors' knowledge, is the first such analysis of fixed-size RDP with replacement for DP-SGD. (4) For the first time, we show that FSwoR and the widely-used Poisson subsampling have the same privacy under replace-one adjacency to leading order in the sampling probability. This has important practical benefits given the memory management advantages of fixed-size subsampling. Thus our results suggest that FSwoR is often preferable to Poisson subsampling. The implementation of our accountant is included in Supplementary Materials and is also publicly available to the community at https://github.com/star-ailab/FSRDP.2 Background and Related Work Mironov (2017) proposed a differential privacy definition based on Rényi divergence, RDP, that subsumes the celebrated moments accountant (Abadi et al., 2016) used to track privacy loss in DP-SGD. RDP and many other flavors of DP rely on subsampling for privacy amplification to ensure that privacy guarantees hold when DP-SGD is applied on samples (i.e., minibatches) of training data (Bun et al., 2018; Balle et al., 2018) . Later works such as (Mironov et al., 2019; Zhu and Wang, 2019) propose a generic refined bound for subsampled Gaussian Mechanism that improves the tightness of guarantees in RDP. However, these methods also produce variable-sized minibatches. Bun et al. (2018) , Balle et al. (2018) , and Wang et al. (2019) propose fixed-size subsampled mechanisms. Balle et al. (2018) provide a general fixed-size subsampling analysis with and without replacement but their bounds focus on tight (ϵ, δ)-DP bounds for a single step, and are not readily applicable to DP-SGD.The AFA developed by Zhu et al. (2022) is an enhancement of the moments accountant offering a lossless conversion to (ϵ, δ) guarantees and is designed for general mechanisms. However, in the context of DP-SGD, there are still technical and practical difficulties. For instance, in DP-SGD there is generally a large number of training steps, in which case the numerical integrations required by Algorithm 1 in Zhu et al. (2022) involve the sum of a large number of oscillatory terms. This can lead to numerical difficulties, especially in light of the lack of rigorous error bounds for the double (Gaussian) quadrature numerical integration approach as discussed in their Appendix E.1. In addition, the bounds in Zhu et al. (2022) are not readily applicable under replace-one adjacency. Due to these considerations, and the fact that RDP is still widely used, in this work we focus on improving the RDP bounds from Wang et al. (2019) . Hayes et al. (2024) uses the results in (Zhu et al., 2022) to provide a bound specific to DP-SGD when one is only concerned with training data reconstruction attacks rather than the membership inference attacks, which leads to a relaxation of DP to provide a bound that is only applicable to defending against data reconstruction and not membership inference. Our work can be viewed as a complimentary work to (Hayes et al., 2024) that is specific to DP-SGD and provides conservative bounds protecting against membership inference attacks (and thus any other type of model inversion attacks including training data reconstruction) in the RDP context.In summary, to our knowledge, Wang et al. (2019) provides the best computable bounds in the fixed-size regime for RDP that are practical for application to DP-SGD. While the computations in Wang et al. (2019) have the attractive property of applying to a general mechanism, in this work we show that there is room for obtaining tighter bounds specific to DP-SGD with Gaussian noise. Our FS-RDP bounds address this theoretical gap by employing a novel Taylor expansion expansion approach, which precisely captures the leading order behavior in the sampling probability, while employing a computable upper bounds on the integral remainder term to prevent privacy leakage. For convenience, we provide the definition of RDP below. Definition 2.1 ((α, ϵ)-RDP (Mironov et al., 2019) ). Let M be a randomized mechanism, i.e., M(D) is a probability distribution for any choice of allowed input (dataset) D. A randomized mechanism is said to have (α, ϵ)-RDP if for any two adjacent datasets D, D ′ it holds that D α (M(D)∥M(D ′ )) ≤ ϵ, where the Rényi divergence of order α > 1 is defined byEQUATIONThe definition of dataset adjacency varies among applications. In this work we consider two such relations: 1) The add/remove adjacency relation, where datasets D and D ′ to be adjacent if one can be obtained from the other by adding or removing a single element; we denote this by D ≃ a/r D ′ .2) The replace-one adjacency definition, denoted D ≃ r-o D ′ , wherein D is obtained from D ′ by replacing a single element. In the next section we derive new RDP bounds when M is DP-SGD using fixed-size minibatches, with or without replacement.3 Rényi-DP Bounds for SGD with Fixed-size SubsamplingIn this section we present our main theoretical results, leading up to Rényi-DP bounds for SGD for fixed-size subsampling done without replacement (Thm. 3.3 for add/remove adjacency and Thm. 3.4 for replace-one adjacency) and with replacement (Thm. 3.7). In Sec. 3.4, we also compare our results with the analysis in (Wang et al., 2019) and show that our analysis yields tighter Rényi bounds by a factor of approximately 4. 
fs-rdp: definition and initial bounds: 
 Given a loss function L, a training dataset D with |D| elements, and a fixed minibatch size, we consider the DP-SGD NN parameter updates with fixed-size minibatches,EQUATION)where t is the iteration number, the initial condition Θ D 0 is independent of D, η t are the learning rates, the noises Z t are Gaussians with mean 0 and covariance C 2 σ 2 t I, and the B D t are random minibatches with fixed size; we denote the fixed value of |B D t | by |B| for brevity. We will consider fixed-size subsampling both without and with replacement; in the former, B D t is a uniformly selected subset of {0, ..., |D| -1} of fixed size |B| for every t and in the latter, B D t is a uniformly selected element of {0, ..., |D| -1} |B| for every t. We assume the Gaussian noises and minibatches are all independent jointly for all steps. Here Clip denotes the clipping operation on vectors, with ℓ 2 -norm bound C > 0, i.e., Clip(x) := x/max{1, ∥x∥ 2 /C}. We derive a RDP bound (Mironov, 2017) for the mechanism consisting of T steps of DP-SGD (2), denoted:EQUATIONSpecifically, we consider the cases without replacement, M F S woR [0,T ] (D), and with replacement, M F S wR [0,T ] (D). To obtain these FS-RDP accountants we boundD α (M F S [0,T ] (D)∥M F S [0,T ] (D ′ )), where D ′ is any dataset that is adjacent to D; as previously stated, we will consider both add/remove and replace-one adjacency. By taking worst-case bounds over the input state at each step, as done in Theorem 2.1 in (Abadi et al., 2016) , we decompose the problem into a sum over stepsEQUATIONwhere the time-inhomogeneous transition probabilities are,EQUATIONµ t (θ t , b, D) := θ t -η t 1 |B| i∈b Clip(∇ θ L(θ t , D i )) .Here N µ, σ 2 (θ) is the Gaussian density with mean µ and covariance to distinguish the transition probabilities for these respective cases.σ 2 I, σ 2 t := C 2 η 2 t σ 2 t /|B|The computations thus far mimic those in (Abadi et al., 2016; Mironov et al., 2019) , with the choice of subsampling method and adjacency relation playing no essential role. That changes in this section, where we specialize to the case of subsampling without replacement. First, in Section 3.2.1 we consider the add/remove adjacency relation; the proof in this case contains many of the essential ideas of our method but is simpler from a computational perspective. The more difficult case of replace-one adjacency will then be studied in Section 3.2.2In Appendix A we derive the following Rényi divergence bound for FS woR -subsampled DP-SGD under add/remove adjacency. Theorem 3.1. Let D ≃ a/r D ′ be adjacent datasets. With transition probabilities defined as in (5) and letting q = |B|/|D| we havesup θt D α (p FSwoR t (θ t+1 |θ t , D)∥p FSwoR t (θ t+1 |θ t , D ′ )) ≤ D α (qN 1,σ 2 t /4 + (1 -q)N 0,σ 2 t /4 ∥N 0,σ 2 t /4 ) . (6)The key step in the proof consists of the decomposition of the mechanism given in Lemma A.1. We note that the r.h.s. of (6) differs by a factor of 1/4 in the variances from the corresponding result for Poisson subsampling in Mironov et al. (2019) . This is due to the inherent sensitivity difference between fixed-size and Poisson subsampling under add/remove adjacency; see ( 28) -(29). To bound the r.h.s. of (3.1) we will employ a Taylor expansion computation with explicit remainder bound. As the r.h.s. of (6) has the same mathematical form as the Poisson subsampling result of Mironov et al. (2019) , and therefore can be bounded by the same methods employed there, the Taylor expansion method is not strictly necessary in the add/remove case. However, we find it useful to illustrate the Taylor expansion method in this simpler case before proceeding to the significantly more complicated replace-one case in Section 3.2.2, where the method in Mironov et al. (2019) does not apply.The Rényi divergences on the r.h.s. of ( 6) can be writtenEQUATIONwhere for any α > 1, σ > 0 we defineEQUATIONEq. ( 8) does not generally have a closed form expression. We upper bound it via Taylor expansion with integral formula for the remainder at order m (c.f. Thm. 1.2, (Stewart, 2022) ):EQUATIONwhere the remainder term is given byEQUATION)Note that we do not take m → ∞ and so H α,σ is not required to be analytic in order to make use of ( 9). Also note that ( 9) is an equality, therefore if we can compute/upper-bound each of the terms then we will arrive at a computable non-asymptotic upper bound on H α,σ (q), without needing to employ non-rigorous stopping criteria for an infinite series, thus avoiding privacy leakage. The order m is a parameter that can be freely chosen by the user. To implement (9) one must compute the derivativesEQUATION0) and bound the remainder (10). We show how to do both of those steps for general m in Appendix B; in the following theorem we summarize the results for m = 3, which we find provides sufficient accuracy in our experiments. Theorem 3.2 (Taylor Expansion Upper Bound). For q < 1 we haveEQUATIONwhere the remainder has the boundEQUATIONand the M and B parameters are given by (41) and (44) respectively.The reason for the complexity of the formulas ( 11) -( 12) is the need to obtain a rigorous upper bound, and not simply an asymptotic result. Results for other choices of m are given in Appendix B.Combining Theorem 3.2 with equations ( 4) -( 8) we now arrive at a computable T -step FS woR -RDP guarantee.Theorem 3.3 (T -step FS woR -RDP Upper Bound: Add/Remove Adjacency). Assuming q < 1, the mechanism M FSwoR [0,T ] (D), defined in (3), has (α, ϵ FSwoR [0,T ] (α))-RDP under add/remove adjacency, whereEQUATIONwhere M and R are given by (41) and (12) respectively.More generally, using the calculations in Appendix B one obtains RDP bounds for any choice of the number of terms, m. The result (13) corresponds to m = 3, which we find to be sufficiently accurate in practice; see Figure 6 .Theorem 3.1 applies to add/remove adjacency, but our method of proof can also be used in the replace-one adjacency case, where it yields the following FS woR -RDP upper bounds. Theorem 3.4 (FS woR -RDP Upper Bounds for Replace-one Adjacency). Let D ≃ r-o D ′ be adjacent datasets. Assuming q := |B|/|D| < 1, for any integer m ≥ 3 we havesup θt D α (p FSwoR t (θ t+1 |θ t , D)∥p FSwoR t (θ t+1 |θ t , D ′ )) (14) ≤ 1 α -1 log 1 + q 2 α(α -1) e 4/σt 2 -e 2/σt 2 + m-1 k=3 q k k! F α,σt,k + E α,σt,m (q) ,where F α,σt,k is given by ( 80) and E α,σt,m (q) by (87).The derivation, found in Appendix C, follows many of the same steps as the add/remove-adjacency case from Theorem 3.1, though deriving bounds on the terms in the Taylor expansion is significantly more involved and the resulting formulas are more complicated. The decomposition from Lemma A.1 again constitutes a key step in the proof. For comparison purposes, in Appendix C.1 we derive the analogous result for Poisson-subsampling under replace-one adjacency; see Theorem C.9.A corresponding T -step RDP bound, analogous to Theorem 3.3, for replace-one adjacency can similarly be obtained by combining Theorem 3.4 with Eq. ( 4). Theorem 3.5 (T -step FS woR -RDP Upper Bound: Replace-one Adjacency). Assuming q < 1 and for any integer m ≥ 3, the mechanism M FSwoR [0,T ] (D), defined in (3), has (α, ϵ FSwoR [0,T ] (α))-RDP under replace-one adjacency, whereEQUATIONHere F α,σt,k is given by (80) and E α,σt,m (q) by (87).We emphasize that one strength of our approach is that it provides a unified method for deriving computable bounds in both the add/remove and replace-one adjacency cases. We also note that Theorem 11(b) of Zhu et al. (2022) , Theorem 6 of Balle et al. (2018) , and Theorem 5 in Mironov et al. (2019) combine to provide an alternative method for deriving the add/remove result of Theorem 3.1 but not the replace-one result of Theorem 3.4. Unlike in the corresponding add/remove theorems, using m > 3 in Theorems 3.4 and 3.5 is often required to obtain acceptable accuracy, though we generally find that m = 4 is sufficient; see Figure 2 .We similarly obtain RDP bounds for DP-SGD using fixed-size minibatches with replacement, i.e., (2) with minibatches B D t that are iid uniformly random samples from {0, ..., |D| -1} |B| . First we present upper bounds on the RDP of M FSwR [0,T ] (D). The derivation follows a similar pattern to that of our results for FS woR -subsampling, first using a probabilistic decomposition of the mechanism (see Lemma D.1) and then expanding (see Appendix D.1). Theorem 3.6 (Fixed-size RDP with Replacement Upper Bound). Assuming |D| > |B| and under add/remove adjacency, the mechanism M FSwR [0,T ] (D) has (α, ϵ FSwR [0,T ] (α))-RDP with (H α,σ given by (8))EQUATIONq := 1 -(1 -|D| -1 ) |B| , a n := q -1 |B| n |D| -n (1 -|D| -1 ) |B|-n .The functions H α,σ were previously bounded above in Theorem 3.2 (for order m = 3), and more generally in Appendix B (for general m). It is straightforward to show that q = O(q) and a n = O(q n-1 ). Therefore, recalling the asymptotics (11) of H α,σ , the behavior of ϵ FSwR [0,T ] (α) and ϵ FSwoR [0,T ] (α) is the same at leading order in q. We also present the following novel lower bounds on the Rényi divergence for one step of FS wR DP-SGD; see Appendix D.2 for the derivation.Theorem 3.7 (Fixed-size RDP with Replacement Lower Bound). Let the mechanism M FSwR t (D, θ t ) denote the t'th step of FS wR DP-SGD (2) when using the training dataset D and starting from the NN parameters θ t , i.e., the transition probabilities (5) where the average is over minibatches b that are elements of {0, ..., |D| -1} |B| . In this result we consider the choice of |B| to be fixed but will allow D and θ t to vary. Suppose there exists θ t , d, andd ′ such that Clip(∇ θ L(θ t , d)) = -Clip(∇ θ L(θ t , d ′ )) and ∥Clip(∇ θ L(θ t , d))∥ = ∥Clip(∇ θ L(θ t , d ′ ))∥ = C ,i.e., when using the NN parameters θ t , the clipped gradients at the samples d and d ′ are anti-parallel and saturate the clipping threshold. Then for any integer α ≥ 2 and any choice of dataset size N we have the following worst-case Rényi divergence lower bound sup The lower bound (17) provides a complement to one step of the upper bound from Theorem 3.6. We demonstrate the behavior of the lower bound (17) in Figure 1 , where we employed the computational method discussed in Appendix D.2.1. Specifically, we show the FS wR -RDP lower bound as a function of α for fixed q and several choices of |B|. For small α, the lower bounds for different |B|'s are approximately equal. This holds up until a critical |B|-dependent threshold where there is a "phase transition" to a regime with non-trivial |B| dependence (even for fixed q) wherein the RDP bound increases quickly. This is in contrast to FS woR -RDP which depends on |B| only through the ratio q. The critical α is smaller for larger |B|, which implies worse privacy guarantees for larger |B|, with all else being equal. Intuitively, a larger minibatch size provides a greater chance for the element that distinguished D from D ′ to be selected multiple times in a single minibatch. If it is selected too many times then it can overwhelm the noise and become "noticeable", thus substantially degrading privacy; this is the key property that distinguishes FS wR from FS woR . In practice, this behavior implies that FS wR -RDP requires much smaller minibatch sizes than FS woR -RDP in order to avoid this critical threshold and maintain privacy.(θt,D,D ′ ):|D|=N, D ′ ≃ a/r D D α M FSwR t (D, θ t )∥M FSwR t (D ′ , θ t ) ≥ 1 α -1 log   |B| n1,...,3.4 Comparison with (Wang et al., 2019) We compare our fixed-size RDP upper bounds with the upper and lower bounds obtained by applying the general-purpose method from Wang et al. (2019) to DP-SGD with Gaussian noise. First we compare the results asymptotically to leading order in q and 1/σ 2 t , which is the domain relevant to applications, and then we will compare them numerically using the full non-asymptotic bounds. Unless stated otherwise, from here on the phrase leading order implicitly refers to the parameters q and 1/σ 2 t . The method from Wang et al. (2019) applied to DP-SGD (App. E.1) gives the one-step Rényi boundsEQUATION)In contrast, to second order in q, one step of our FS-RDP results under add/remove adjacency from Theorem 3.3 or Theorem 3.6 gives while our FS woR -RDP bound under replace-one adjacency from Theorem 3.4 givesϵ ′ t (α) ≤ q 2 2 α(e 4/σ 2 t -1) + O(q 3 ) = 2q 2 α/σ 2 t + O(q 2 /σ 4 t ) + O(q 3 ) ,ϵ ′ t (α) ≤q 2 α e 4/σ 2 t -e 2/σ 2 t + O(q 3 ) (19) =2q 2 α/σ 2 t + O(q 2 /σ 4 t ) + O(q 3 ) .Therefore, to leading order, our methods all have the same behavior and give tighter RDP bounds by a factor of 4 than that of Wang et al. (2019) . We emphasize that the adjacency relation used by Wang et al. (2019) is replaceone and therefore, strictly speaking, should be compared only with our Theorem 3.4, which also uses replace-one adjacency. Figure 2 shows a non-asymptotic comparison between these RDP bounds, with ours being tighter significantly tighter than that of Wang et al. (2019) . In particular, our bound from Theorem 3.4 with m = 4 remains close to the theoretical lower bound from Wang et al. (2019) over the entire range of α's pictured, which is a sufficiently large range for practical application to DP, e.g., larger than the default range used in Meta Platforms (2024). Following the setting in (Abadi et al., 2016) , we provide corresponding (ϵ, δ)-DP guarantees in Figure 3 shows corresponding ϵ for δ ∈ {1e -4, 1e -5, 1e -6, 1e -7, 1e -8, 1e -9, 1e -10} after 250 training epochs in DP-SGD. To translate RDP bounds to (ϵ, δ)-DP, we used Theorem 21 in (Balle et al., 2020) , which is also used in the Opacus library. As shown in Figure 3 , our FS woR bounds with m ∈ {3, 4, 5} (solid lines) are close to the lower-bound provided in (Wang et al., 2019) (circles) over a range of δ's that depend on the value of m. These guarantees are significantly tighter than the upper bound from Wang et al. (2019) (dashed line). In our experiments, we observe that m = 4 and m = 5 yields essentially the same results and we find that m = 4 suffices to obtain good results in practice while outperforming (Wang et al., 2019) . 
comparing fixed-size and poisson subsampling: 
 In Sec. 3.4 we demonstrated that our Rényi-DP-SGD bounds with fixed-size subsampling in Thm. 3.3 are tighter than the general-purpose fixed-size subsampling method from Wang et al. (2019) . However, the most commonly used implementations of DP-SGD rely on Poisson subsampling, not fixed-size, and there the comparison is more nuanced. In this section we compare Poisson and fixed-size subsampling and point out advantages and disadvantages of each.In terms of privacy guarantees, all else being equal, when using replace-one adjacency we find that DP-SGD with fixed-size subsampling yields the same privacy guarantees as Poisson subsampling to leading order; specifically, when moving from add/remove to replace one-adjacency, our fixedsize subsampling DP bounds do not change at leading order but the Poisson subsampling DP bounds change by a factor of 2 to match the fixed-size results. In contrast, when using add/remove adjacency we find that Poisson subsampling has a natural advantage over fixed-size subsampling by approximately a factor of 2. Here we outline the reason for this intuitively.The derivation of Rényi-DP bounds when using Poisson subsampling, (see Abadi et al. (2016) ; Mironov et al. (2019) ) shares many steps with our derivation for fixed-size subsampling without replacement in Appendix A. In both cases, to leading order the bound on the Rényi divergence is proportional to ∥µ -µ ′ ∥ 2 , where µ and µ ′ are the sum of the gradients under the adjacent datasets D and D ′ . When using Poisson subsampling, under add/remove adjacency, µ ′ differs from µ by the addition or deletion of a single term in the sum with probability q, while under replace-one adjacency, µ ′ differs from µ by the replacement of one element with probability q. Therefore in the former case, ∥µ -µ ′ ∥ is the norm of a single clipped vector, hence it scales with the clipping threshold, C, while in the latter, ∥µ -µ ′ ∥ is the norm of the difference between two clipped vector, hence it scales with 2C, as in the worst case the two vectors are anti-parallel. In fixed-size subsampling and under both adjacency relations, µ ′ differs from µ by replacing one element with another (due to the fixed-size minibatch constraint) and therefore ∥µ -µ ′ ∥ scales like 2C. After squaring and at leading order, this leads to a factor of 4 difference in the Rényi divergences between Poisson and fixed-size cases under add/remove adjacency while under replace-one adjacency the Poisson and fixed-size results agree; see Appendix C.1 for a more precise analysis of the latter comparison.Translating this to (ϵ, δ)-DP then leads to an approximate factor of 2 difference in ϵ for the same δ under add/remove adjacency, as shown in Appendix G, while under replace-one adjacency the different subsampling methods lead to the same DP guarantees at leading order. When the higher order terms are taken into account, Poisson subsampling regains a slight privacy advantage under replace-one adjacency; see Figures 4 and 7 .Our numerical comparisons in Section 3.4 showed that FS woR -RDP yields significantly tighter guarantees than its counterpart from Wang et al. (2019) . It is also useful to compare the empirical privacy guarantees of fixed-size versus the Poisson subsampling RDP commonly used in DP implementations such as Opacus. To this end, we use CIFAR10 and a simple convolutional neural network (CNN) to compare the canonical Poisson subsampling RDP with FS-RDP. Our setup closely follows the procedure in (Abadi et al., 2016) . Note that the goal of this comparison is not to achieve the state-ofthe-art privacy guarantees on CIFAR10 with any specific neural network. Instead we aim to compare the performance and guarantees with the alternative fixed-size and Poisson subsampling methods.Following Abadi et al. (2016) , we used a simple CNN network that is pre-trained on CIFAR100 dataset non-privately. Our CNN has six convolutional layers 32, 32, 64, 128, 128, and 256 square filters of size 3 × 3, in each layer respectively. Since using Batch Normalization layers leads to privacy violation, we used group normalization instead as a common practice that is also used the Opacus library (Meta Platforms, 2024). Our network has 3 group normalization layers of size 64, 128, and 256 after the second, fourth, and sixth convolution layer. Similar to Abadi et al. (2016) , this model reaches the accuracy of 81.33% after 250 epochs of training in a non-private setting. For the Poisson subsampling RDP, we simply used the one in Opacus. Figure 4 indicates that for a fixed σ, FS woR -RDP surpasses the accuracy of Poisson subsampled RDP in Opacus, but with slightly higher epsilon guarantees, as shown in the top panel. This behavior matches the approximate calculation from Eq. ( 150) as well as the discussion in Section 4.1; specifically, Poisson and FS woR -RDP agree to leading order but the higher order contributions give Poisson subsampling a slight privacy advantage at the same level of added noise. We emphasize that these results all apply to the replace-one notion of adjacency, as used in (Wang et al., 2019) and also our Theorem 3.4; we did not find existing computable RDP bounds for Poisson-subsampling under replace-one adjacency and so we used the new result in Theorem C.9. We conjecture that the performance difference observed in the bottom panel of Figure 4 is due to the difference in variance between Poisson and fixed-size subsampling, which we derive in Appendix F. Specifically, there we show that fixed-size subsampling has reduced variance compared to Poisson when the parameters are away from a minimizer. Additional experiments, gauging the sensitivity of FS woR -RDP to key parameters including σ and |B|, are given in Appendix H.The RDP accountant used in common DP libraries assumes a Poisson subsampling mechanism that leads to variable-sized minibatches during the training. That is, the size of each mini-batch cannot be determined in advance. Allowing the minibatches to have variable sizes creates an engineering challenge: an uncharacteristically large minibatch can cause an out-of-memory error, and depending on how fine tuned the minibatch size is to the available GPU memory this could be a frequent occurrence. To tackle this issue, the Opacus implementation requires the additional complication of wrapping the Pytorch DataLoader object into a BatchMemoryManager object, to alleviate the memory intensive nature of Poisson subsampling. BatchMemoryManager requires privacy practitioners to define two minibatch sizes: one is a logical minibatch size for Poisson sampling and the other is a physical minibatch size that determines the actual space allocated in memory and is determined by the practitioner via a variable named max physical batch size (Meta Platforms, 2024) . This design follows the distinction between 'lots' and minibatches described in the moments accountant (Abadi et al., 2016) . Our proposed FS-RDP does not require this memory management and is much less memory intensive. Figure 5 depicts the GPU memory footprints during 100 epochs of privately training the above CNN on CIFAR10 with FS woR -RDP and three of the most common DP accountants: RDP (Mironov et al., 2019) , f -DP (Dong et al., 2022) , and PRV (Gopi et al., 2021) . As shown in Figure 5 , the memory usage of FS-RDP remains constant at 14,826 MB. However, the memory usage varies from 7,984 MB to 23,280 MB (almost twice larger than FS-RDP) for Poisson subsampled RDP in Opacus. Similar memory usage is observed for other Opacus accountants that use Poisson subsampling mechanism (i.e., f -DP and PRV). This fluctuating memory usage could result in undesirable outcomes in settings where resource-constrained devices are considered. For instance, in federated learning environments, erratic memory usage could result in a higher dropout rate and delay the global model's convergence (Liu et al., 2021) . It could also undermine model's accuracy and fairness, preventing lower-end devices from participating in the training process (Imteaj et al., 2021) . 
conclusion: 
 Differentially private stochastic gradient descent with fixed-size minibatches has attractive properties including reduced gradient estimator variance and simplified memory management. Our work presents a holistic RDP accountant for DP-SGD with fixed-size subsampling without replacement (FSwoR) and with replacement (FSwR) and, in the FSwoR case, consider both add/remove and replace-one adjacency. As we showed theoretically and empirically, since FSwoR under replace-one adjacency leads to the same leading-order privacy guarantees as the widely-used Poisson subsampling, we suggest using the former over the latter to benefit from the memory management and reduced variance properties. For subsampling without replacement under replace-one adjacency, we obtained significantly tighter RDP bounds (4 times improvement) over the most recent computable results (Wang et al., 2019) . For subsampling with replacement we obtained the first non-asymptotic upper and lower RDP bounds for DP-SGD. We also provided the first comparison of gradient estimator's variance and privacy guarantees between FS-RDP and Poisson subsampled RDP commonly used in DP libraries. Our analysis revealed that FS-RDP reduces the variance of the gradient estimator. We highlighted the memory usage advantages of FS-RDP over Poisson subsampled RDP, which makes it a practical choice for privately training large deep learning models. We made FS-RDP's implementation and its evaluations available online such that it can be added to common DP libraries.A Proof of Theorem 3.1: FS woR -RDP Bound under Add/Remove AdjacencyIn this appendix we prove a one-step Rényi bound for DP-SGD using subsampling without replacement under the add/remove adjacency relation, as stated in Theorem 3.1. To do this, we will show that the sum over minibatches b in (5) can be rewritten as an average of Gaussian mixtures with parameter q = |B|/|D|. This will then allow us to mimic several steps from the derivation in Mironov et al. (2019) , which applies to Poisson subsampling. In the Poisson subsampling case the Gaussian mixture structure is more or less apparent from the start, but here it requires some additional work to isolate. We first consider the case where D ′ is obtained from D by removing one element. Without loss of generality we can let D ′ = (D 0 , ..., D |D|-2 ). Later we will show that the bound we derive in this case also bounds the case where D ′ has one additional element.Noting that 1N Db is simply the expectation with respect to the distribution of B D t , we start by showing that the distribution of B D t can be obtained from the distribution of B D ′ t through the introduction of two auxiliary random variables J and B as follows: Let J be a Bernoulli(q) random variable where q = |B|/|D| and let (B ′ , B) be random variables independent of J such that B ′ is a uniformly selected subset of {0, ..., |D| -2} of size |B| (so EQUATIONP (B = b) =P (B ′ = b)(1 -q) = 1 -q |D|-1 |B| = (|D| -|B|)/|D| (|D| -1)!/(|B|!(|D| -1 -|B|)!) (22) = 1 |D|!/(|B|!(|D| -|B|)!) = 1 |D| |B| = P (B D t = b) as claimed. Now suppose |D| -1 ∈ b. Then b \ {|D| -1} is a subset of {0, ..., |D| -2} of size |B| -1, hence P (B = b) =qP ( B ∪ {|D| -1} = b) = |B| |D| P ( B = b \ {|D| -1}) (23) = |B| |D| E b ′ ∼B ′ P ( B = b \ {|D| -1}|B ′ = b ′ ) = |B| |D| E b ′ ∼B ′ 1 b\{|D|-1}⊂b ′ P ( B = b \ {|D| -1}|B ′ = b ′ ) = 1 |D| P (b \ {|D| -1} ⊂ B ′ ) = 1 |D| |D| -1 -(|B| -1) |D|-1 |B| = 1 |D| |B| = P (B D t = b) .This completes the proof.Using Lemma A.1, we can express the expectation with respect to B D t in terms of J, B ′ , and B. This will provide us with the desired decomposition of the sum in (5). The definition (20) is reminiscent of the constructions in Proposition 23 in Appendix C of Wang et al. (2019) , however our analysis will also result in tighter bounds for DP-SGD than one gets from the general-purpose method in Wang et al. (2019) .Using Lemma A.1 in the manner discussed above, the transition probabilities (5) for one step of DP-SGD with FS woR -subsampling can be rewritten as follows:p t (θ t+1 |θ t , D) =E b∼B D t N µt(θt,b,D), σ 2 t (θ t+1 ) (24) =E (b ′ , b)∼(B ′ , B) qN µt(θt, b∪{|D|-1},D), σ 2 t (θ t+1 ) + (1 -q)N µt(θt,b ′ ,D), σ 2 t (θ t+1 ) .Recalling that B ′ ∼ B D ′ t and noting that µ t (θ t , B ′ , D ′ ) = µ t (θ t , B ′ , D) (as D ′ agrees with D except for the last missing entry, which by definition cannot be contained in B ′ ), the transition probabilities p t (θ t+1 |θ t , D ′ ) can similarly be writtenp t (θ t+1 |θ t , D ′ ) =E b ′ ∼ B D ′ t N µt(θt,b ′ ,D ′ ) (θ t+1 ) (25) =E (b ′ , b)∼(B ′ , B) N µt(θt,b ′ ,D), σ 2 t (θ t+1) . The outer expectations in the expressions ( 24) and ( 25) are now both with respect to the same distribution, with the key difference being captured by the parameter q. The construction (20) was introduced precisely so that we might obtain such a decomposition. It allows the next step to proceed as in DP-SGD with Poisson subsampling Abadi et al. (2016) ; Mironov et al. (2019) . Namely, we use quasiconvexity of the Rényi divergences (i.e., the data processing inequality) Van Erven and Harremos (2014) to obtainEQUATION≤ max(b ′ , b) D α qN µt(θt, b∪{|D|-1},D), σ 2 t + (1 -q)N µt(θt,b ′ ,D), σ 2 t ∥N µt(θt,b ′ ,D), σ 2 t = max (b ′ , b) D α qN ∥µt(θt, b∪{|D|-1},D)-µt(θt,b ′ ,D)∥, σ 2 t + (1 -q)N 0, σ 2 t ∥N 0, σ 2 t ,where, as in Abadi et al. (2016) ; Mironov et al. (2019) , we used a change of variables to reduce the computation to 1-dimension in the last line.Remark A.2. We note that the above analysis can be applied to any divergence, D, that is quasiconvex in its arguments, not just the Rényi divergences. In the general case one similarly findsD(p(θ t+1 |θ t , D)∥p t (θ t+1 |θ t , D ′ ) (27) ≤ max (b ′ , b) D qN µt(θt, b∪{|D|-1},D), σ 2 t + (1 -q)N µt(θt,b ′ ,D), σ 2 t ∥N µt(θt,b ′ ,D), σ 2 t .For instance, this holds for all f divergences, including the important class of hockey-stick divergences that are important for other DP frameworks such as in Zhu et al. (2022) .Returning to (26), recalling that C is the clipping bound and that b ⊂ b ′ and |b ′ \ b| = 1, we can computeEQUATIONWe note that this result is 4 times as large as the corresponding calculation in the case of Poisson subsampling, i.e., the sensitivity in the case of fixed-sized subsampling is inherently twice that of Poisson subsampling when using add/remove adjacency. This is because in Poisson subsampling when the minibatch from D differs from that of D ′ it is due to the inclusion of a single additional element. However, in fixed-size subsampling, when the minibatches are not identical then they differ by a replacement; this contributes more to the difference in means by a factor of 2.Combining (29) with the fact that scaling the means in Gaussian mixtures by the same factor c ≥ 1 can only increase the Rényi divergence, see Section 2 in Mironov et al. (2019) , we obtain the following uniform Rényi bound.EQUATIONNow we show that the right-hand-side of (30) also bounds the case where D ′ has one additional element. Repeating the steps in the above derivation, in this case we obtainEQUATIONwhere q = |B|/(|D| + 1); the appearance of |D ′ | = |D| + 1 here is due to the interchange of the roles of D and D ′ in the decomposition (20) for this case. Next use Theorem 5 in Mironov et al. (2019) to obtainEQUATION)Quasiconvexity of the Rényi divergences implies the right-hand-side of ( 32) is non-decreasing in q. Therefore, as q ≤ q, we obtainEQUATION)The upper bound here is the same as the one we obtained in (30) when D ′ had one fewer element than D. Finally, by changing variables we can rescale the mean to 1 and therefore, noting thatEQUATIONwe obtainEQUATIONThis completes the proof. 
broader impacts: 
 Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?Answer: [Yes] Justification: The introduction and conclusion both highlight the positive societal impact of our work on privacy preserving machine learning by providing a new privacy accountant that advances the status quo in enhancing users privacy as machine learning and deep learning becomes more ubiquitous in individuals lives. As noted above, we do not envision any negative societal impact.Guidelines:• The answer NA means that there is no societal impact of the work performed.• If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact. • Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations. • The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is a direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster. • The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology. • If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML). 
safeguards: 
 Question: Does the paper describe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?Answer: [NA] Justification: The dataset is a publicly available common benchmark. Also the the proposed accountant helps privacy-preserving ML and is not envisioned to be used in a malicious manner.Guidelines:• The answer NA means that the paper poses no such risks.• Released models that have a high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters. • Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images. • We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make a best faith effort. • We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurIPS Code of Ethics and the guidelines for their institution. • For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review. 
 
 Appendix: In this appendix we provide details regarding the Taylor expansion (9) of H α,σ (q) (Eq. ( 8)), which can be writtenH α,σ (q) := q N 1,σ 2 /4 (θ) N 0,σ 2 /4 (θ) + (1 -q) α N 0,σ 2 /4 (θ)dθ .In Section 3.2.1 we obtained one-step RDP bounds in terms of the Rényi divergenceEQUATIONTo obtain the computable RDP bounds in Theorem 3.3 and also in Theorem 3.6 we must bound H α,σ (q). We proceed by computing the coefficients in its Taylor expansion, including an explicit upper bound on the remainder term (10). First note that using the dominated convergence theorem (see, e.g., Theorem 2.27 in Folland (1999)), it is straightforward to see that H α,σ (q) is smooth in q and can be differentiated under the integral any number of times:EQUATION)This justifies the use of Taylor's theorem with integral remainder to H α,σ for any order m ∈ Z + to arrive atEQUATIONwhereEQUATION(note that the M σ,1 = 0, which is why the summation in Eq. ( 38) starts at k = 2) and the remainder term is given byEQUATION)For integer k ≥ 0, the M σ,k 's can be computed using the binomial theorem together with the formula for the moment generating function (MGF) of a Gaussian:M σ,k = k ℓ=0 (-1) k-ℓ k ℓ N 1,σ 2 /4 (θ) N 0,σ 2 /4 (θ) ℓ N 0,σ 2 /4 (θ)dθ (41) = k ℓ=0 (-1) k-ℓ k ℓ e -ℓ/(2σ 2 /4) e ℓθ/(σ 2 /4) N 0,σ 2 /4 (θ)dθ = k ℓ=0 (-1) k-ℓ k ℓ e ℓ(ℓ-1)/(σ 2 /2) = k ℓ=2 (-1) k-ℓ k ℓ e 2ℓ(ℓ-1)/σ 2 + (-1) k-1 (k -1) .In particular, M σ,0 = 1 and M σ,1 = 0. If α is an integer then the expansion (38) truncates at a finite number of terms (terms with k > α are zero), but for non-integer α one must bound the remainder term. To enhance numerical stability of the computations, in practice we add terms starting from the highest order and proceeding to the lowest order.To bound the remainder term (40) we need to bound d m dq m H α,σ (sq) for t, q ∈ (0, 1). We will break the calculation into two cases, where different methods are appropriate. Here we assume q < 1.In these bounds, it will be useful to employ the following definition for integer j ≥ 0:EQUATION)For j even we have B σ,j = M σ,j and for j odd we can use the Cauchy-Schwarz inequality to boundEQUATIONas j ± 1 are both even. Therefore an explicitly computable upper bound for all integer j ≥ 1 is given byEQUATIONwhere M σ,k is given by (41).Case 1: α -m > 0 Using the bound x α-m ≤ 1 + x ⌈α⌉-m for all x > 0 we can computed m dq m H α,σ (sq) (45) ≤ m-1 j=0 |α -j| 1 + sq N 1,σ 2 /4 (θ) N 0,σ 2 /4 (θ) + (1 -sq) ⌈α⌉-m N 1,σ 2 /4 (θ) N 0,σ 2 /4 (θ) -1 m N 0,σ 2 /4 (θ)dθ = m-1 j=0 |α -j| sq N 1,σ 2 /4 (θ) N 0,σ 2 /4 (θ) -1 + 1 ⌈α⌉-m N 1,σ 2 /4 (θ) N 0,σ 2 /4 (θ) -1 m N 0,σ 2 /4 (θ)dθ + B σ,m = m-1 j=0 |α -j| ⌈α⌉-m ℓ=0 ⌈α⌉ -m ℓ (sq) ℓ N 1,σ 2 /4 (θ) N 0,σ 2 /4 (θ) -1 ℓ (46) × N 1,σ 2 /4 (θ) N 0,σ 2 /4 (θ) -1 m N 0,σ 2 /4 (θ)dθ + B σ,m ≤ m-1 j=0 |α -j|   ⌈α⌉-m ℓ=0 ⌈α⌉ -m ℓ (sq) ℓ B σ,ℓ+m + B σ,m   . Therefore |R α,σ,m (q)| (47) ≤q m m-1 j=0 |α -j|   ⌈α⌉-m ℓ=0 ⌈α⌉ -m ℓ q ℓ B σ,ℓ+m 1 0 (1 -s) m-1 (m -1)! s ℓ ds + 1 m! B σ,m   =q m m-1 j=0 |α -j|   ⌈α⌉-m ℓ=0 q ℓ (⌈α⌉ -m)! (⌈α⌉ -m -ℓ)!(m + ℓ)! B σ,ℓ+m + 1 m! B σ,m   .Case 2: α -m ≤ 0In this case we haveEQUATIONand thereforeEQUATIONThis leads to the following bound on the remainderEQUATIONwhere B σ,m is defined by (44).Explicit bounds for m = 3: We find that m = 3 generally yields sufficiently accurate bounds. Here we specialize the above results to this case. For q < 1 we haveEQUATION)This completes the proof of Theorem 3.2.To illustrate the use of these bounds, and the fact that m = 3 is often sufficient in practice, in Figure 6 we plot the FS woR -RDP bound from Theorem 3.3 (which uses the expansion and bound ( 52) -( 53)) for one step (i.e, T = 1) and with the parameter values σ t = 6, |B| = 120, and |D| = 50, 000. Note that the difference between the bound for m = 3 and m = 4 is negligible over a wide range of α (sufficient to cover the default α range needed by Meta Platforms ( 2024)). The value of q in this example is far from being unrealistically small; smaller values of q will only make the effect of m even less important. Therefore m = 3 is generally sufficient. C Proof of Theorem 3.4: FS woR -RDP Bound under Replace-One AdjacencyWe now prove the FS woR -RDP bound under replace-one adjacency as stated in Theorem 3.4 above. Without loss of generality, we can assume that D and D ′ share all elements except the last. The derivation shares several initial steps with the add/remove case from Appendix A. Specifically, we apply Lemma A.1 to the transition probabilities for one step of DP-SGD with FS woR -subsampling, both for D and for D ′ (as they now have the same number of elements), to obtainEQUATIONsimilarly to (24). Now we again use quasiconvexity, this time in both arguments simultaneously, to bound the effect of the mixture (expectation) over (b ′ , b) by the worst case, yielding the one-step Rényi boundEQUATIONQ b := N µt(θt, b∪{|D|-1},D), σ 2 t , Q ′ b := N µt(θt, b∪{|D ′ |-1},D ′ ), σ 2 t , P b ′ := N µt(θt,b ′ ,D), σ 2 t , P ′ b ′ := N µt(θt,b ′ ,D ′ ), σ 2 t .Using the fact D and D ′ agree on the index minibatch b ′ , and hence µt (θ t , b ′ , D) = µ t (θ t , b ′ , D ′ ), we see that P b ′ = P ′ b ′ .We can therefore change variables to center both of their means at zero, givingD α (p(θ t+1 |θ t , D)∥p t (θ t+1 |θ t , D ′ )) (56) ≤ max (b ′ , b) D α qN ∆µt(b ′ , b), σ 2 t + (1 -q)N 0, σ 2 t ∥qN ∆µ ′ t (b ′ , b), σ 2 t + (1 -q)N 0, σ 2 t , ∆µ t (b ′ , b) := µ t (θ t , b ∪ {|D| -1}, D) -µ t (θ t , b ′ , D) , ∆µ ′ t (b ′ , b) := µ t (θ t , b ∪ {|D| -1}, D ′ ) -µ t (θ t , b ′ , D) .The means of the Gaussians satisfy the boundsEQUATIONwhere r t was defined in (29); the first two inequalities follow directly from (28) while the proof of the third is almost identical due to D and D ′ differing only in the last element.Remark C.1. Note that, similarly to Remark A.2, Eq. ( 56) remains true when the Rényi divergences D α are replaced by any other divergence D that is quasiconvex in both of its arguments, e.g., the hockey-stick divergences. This fact is useful for other differential privacy paradigms, though using the quasiconvexity bound does not always lead the tightest possible bounding pair of distributions; see Lebeda et al. (2024) .The key difference between the bound (56) and the corresponding result (26) in the add/remove adajcency case is that both arguments of the Rényi divergence on the right-hand side of ( 56) are now Gaussian mixtures with mixing parameter q; in the add/remove adjacency case only one argument is a Gaussian mixture. This difference makes the following computations more involved, though the same Taylor expansion technique can be applied.Using the definition of Rényi divergences we can writeEQUATION)Next we Taylor expand the argument of the logarithm:EQUATIONwhere here the remainder term is given byEQUATION)Using the formulad k dx k (f g) = k j=0 k j d j dx j f d k-jdx k-j g along the dominated convergence theorem to justify differentiating under the integral, the derivatives can be computed as follows:d k dq k F α (q) (61) = k j=0 k j j-1 ℓ=0 (α -ℓ) k-j-1 ℓ=0 (1 -α -ℓ) qN ∆µt(b ′ , b), σ 2 t (θ) + (1 -q)N 0, σ 2 t (θ) α-j qN ∆µ ′ t (b ′ , b), σ 2 t (θ) + (1 -q)N 0, σ 2 t (θ) α+k-j-1 × N ∆µt(b ′ , b), σ 2 t (θ) -N 0, σ 2 t (θ) j N ∆µ ′ t (b ′ , b), σ 2 t (θ) -N 0, σ 2 t (θ) k-j dθ .Taylor expansion coefficients for k = 0, 1, 2: Evaluating (61) at q = 0 we haved k dq k F α (0) (62) = k j=0 (-1) k-j k j j-1 ℓ=0 (α -ℓ) k-j-1 ℓ=0 (α + ℓ -1) × N ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 j N ∆µ ′ t (b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 k-j N 0, σ 2 t (θ)dθ .Note that we have not yet been able to maximize over b ′ , b; doing so while incorporating the constraints ( 57) is nontrivial and constitutes one of the the primary difficulties in obtaining sufficiently tight bounds. We will make particular attention to the first few terms and then use looser (but still sufficiently tight in practice) approximations for the higher order terms. First note that F α (0) = 1 andd dq F α (0) = -(α -1) N ∆µ ′ t (b ′ , b), σ 2 t (θ) -N 0, σ 2 t (θ) dθ (63) + α N ∆µt(b ′ , b), σ 2 t (θ) -N 0, σ 2 t (θ) dθ = 0 .For k = 2 we can expand (62) and then evaluate the integrals by using the formula for the MGF of a Gaussian to obtaind 2 dq 2 F α (0) (64) =α(α -1) exp ∥∆µ t (b ′ , b)∥ 2 σ 2 t + exp ∥∆µ ′ t (b ′ , b)∥ 2 σ 2 t -2 exp ∆µ t (b ′ , b) • ∆µ ′ t (b ′ , b) σ 2 t .Using the third constraints in Eq. ( 57) we obtain the boundexp ∆µ t (b ′ , b) • ∆µ ′ t (b ′ , b) σ 2 t (65) = exp ∥∆µ t (b ′ , b)∥ 2 + ∥∆µ ′ t (b ′ , b)∥ 2 -∥∆µ t (b ′ , b) -∆µ ′ t (b ′ , b)∥ 2 2 σ 2 t ≥ exp ∥∆µ t (b ′ , b)∥ 2 + ∥∆µ ′ t (b ′ , b)∥ 2 -r 2 t 2 σ 2 t .Elementary calculus along with the bounds ∥∆µEQUATION)where we used (34). These terms are sufficient to precisely capture the leading-order behavior.Bounding the coefficients for k ≥ 3: For the higher-order terms we obtain an upper bound which relies on the following lemmas. It will be useful to phrase various bounds in terms of the χ β divergences. Definition C.2. For β > 0 define the χ β divergencesEQUATIONIt will also be convenient to adopt the notation D χ 0 := 1, despite it not defining a divergence.The following property is a imply consequence of Hölder's inequality.EQUATION)The next lemma will be used to obtain a uniform bound, i.e., not dependent on b ′ , b, or θ t . Lemma C.4. For k ≥ 0 an integer we haveEQUATIONmax b ′ , b D χ k N ∆µt(b ′ , b), σ 2 t ∥N ∆µ ′ t (b ′ , b), σ 2 t ≤ B σt,k ,where B was defined in (44). These bounds also hold with ∆µ t and ∆µ ′ t interchanged.Proof. D χ k is the f divergence corresponding to f β (y) = |y -1| k , therefore it is jointly convex in (Q, P ). Therefore the same argument as in Section 2 of Mironov et al. (2019) (which only relies on the convexity of the divergence) impliesD χ k N ∆µt(b ′ , b), σ 2 t ∥N 0, σ 2 t ≤ D χ k N rt, σ 2 t ∥N 0, σ 2 t (70) =D χ k N 1,σ 2 t /4 ∥N 0,σ 2 t /4 = B σt,k ≤ B σt,k, where we used ( 57) and ( 44). The second result in ( 69) is proven in a similar manner, where we must also change variables to center the second Gaussian at zero and then employ the third bound in (57).The second lemma we require is a slight improvement on Lemma 24 of (Wang et al., 2019) . Lemma C.5. Let p, q, r be strictly positive probability densities and k ≥ 1 be an integer. Thenp(θ) -q(θ) r(θ) k r(θ)dθ (71) ≤D χ k (q∥p) + D χ k (p∥q) + D χ k (p∥r) + D χ k (q∥r) if k is even 0 if k is odd .Remark C.6. Note that a key difference between our technique and that of Wang et al. (2019) is that we only apply Lemma C.5 to the higher order terms (i.e., when k ≥ 2) while we bounded the leading contribution (66) more precisely. In contrast, Wang et al. (2019) used Lemma C.5 (their Lemma 24) on all terms.Proof. We proceed by breaking the domain of integration into four regions, depending on the sizes of p, q, r and then use the appropriate bounds to replace one of the densities with another:p(θ) -q(θ) r(θ) k r(θ)dθ (72) = r≥p (p(θ) -q(θ)) k r(θ) k-1 dθ + r<p,r≥q (p(θ) -q(θ)) k r(θ) k-1 dθ + r<p,r<q,q≥p (p(θ) -q(θ)) k r(θ) k-1 dθ + r<p,r<q,q<p (p(θ) -q(θ)) k r(θ) k-1 dθ ≤ r≥p |p(θ) -q(θ)| k p(θ) k-1 dθ + r<p,r≥q |p(θ) -q(θ)| k q(θ) k-1 dθ + r<p,r<q,q≥p (q(θ)-r(θ)) k r(θ) k-1 dθ if k is even 0 if k is odd + r<p,r<q,q<p (p(θ) -r(θ)) k r(θ) k-1 dθ ≤ |p(θ) -q(θ)| k p(θ) k-1 dθ + |p(θ) -q(θ)| k q(θ) k-1 dθ + |q(θ)-r(θ)| k r(θ) k-1 dθ if k is even 0 if k is odd + (p(θ) -r(θ)) k r(θ) k-1 dθ .This completes the proof.Lemma C.7. Let β 1 , β 2 > 0 with β 1 + β 2 and integer. Thenmax b ′ , b N ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 β1 N ∆µ ′ t (b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 β2 N 0, σ 2 t (θ)dθ (73) ≤ B σt,β1+β2 ,where B σ,k was defined in (44).Proof. Let p, q > 1 be conjugate exponents, i.e., they satisfy 1/p + 1/q = 1. Using Hölder's inequality we can computeN ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 β1 N ∆µ ′ t (b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 β2 N 0, σ 2 t (θ)dθ (74) ≤D χ β 1 p N ∆µt(b ′ , b), σ 2 t ∥N 0, σ 2 t 1/p D χ β 2 q N ∆µ ′ t (b ′ , b), σ 2 t ∥N 0, σ 2 t 1/q.We now use a simple heuristic that allow for p, q to be chosen so as to minimize the right-hand side of this upper bound when β 1 , β 2 are large. Lemma C.4 together with ( 44) and ( 41) suggest that the two factors on the right-hand side approximately behave like exp(2pβ 2 1 /σ 2 t ) and exp(2qβ 2 2 /σ 2 t ) respectively. Therefore, to make the upper bound as tight as possible, we want to to minimize pβ 2 1 + qβ 2 2 subject to p > 1, q = p/(p -1). Simple calculus shows that the minimum occurs at p := 1 + β 2 /β 1 , q = 1 + β 1 /β 2 . Making this choice, using Lemma C.3 to convert to integer orders, and then employing Lemma 69 we findN ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 β1 N ∆µ ′ t (b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 β2 N 0, σ 2 t (θ)dθ (75) ≤D χ ⌈β 1 p⌉ N ∆µt(b ′ , b), σ 2 t ∥N 0, σ 2 t β1/⌈β1p⌉ D χ ⌈β 2 q⌉ N ∆µ ′ t (b ′ , b), σ 2 t ∥N 0, σ 2 t β2/⌈β2q⌉ ≤ B σt,⌈β1p⌉ β1/⌈β1p⌉ B σt,⌈β2q⌉ β2/⌈β2q⌉.In the case where β 1 + β 2 is an integer the above can be simplified viaB σt,⌈β1p⌉ β1/⌈β1p⌉ B σt,⌈β2q⌉ β2/⌈β2q⌉ = B σt,β1+β2 .We use these lemmas to bound the higher order derivatives as follows. First rewrite (62) asd k dq k F α (0) (76) =(α -1)α k-1 k j=0 (-1) k-j k j × N ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 j N ∆µ ′ t (b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 k-j N 0, σ 2 t (θ)dθ + (α -1)α k-1 k j=0 (-1) k-j k j α α -1 j-1 ℓ=0 (1 -ℓ/α) k-j-1 ℓ=0 (1 + (ℓ -1)/α) -1 × N ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 j N ∆µ ′ t (b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 k-j N 0, σ 2 t (θ)dθ =(α -1)α k-1 N ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -N ∆µ ′ t (b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) k N 0, σ 2 t (θ)dθ + (α -1)α k-1 k j=0 (-1) k-j k j α α -1 j-1 ℓ=0 (1 -ℓ/α) k-j-1 ℓ=0 (1 + (ℓ -1)/α) -1 × N ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 j N ∆µ ′ t (b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 k-j N 0, σ 2 t (θ)dθ .The first family of integrals can be bounded using Lemma C.5 and then Lemma C.4 to obtainN ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -N ∆µ ′ t (b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) k N 0, σ 2 t (θ)dθ (77) ≤    4M σt,k if k is even 3M 1/2 σt,k-1 M 1/2 σt,k+1if k is odd .The second family of integrals can be bounded by using Lemma C.7 (when j, k -j > 0) and Lemma C.4 (when j = 0 or j = k), which yieldsN ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 j N ∆µ ′ t (b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 k-j N 0, σ 2 t (θ)dθ (78) ≤ B σt,k ,Therefore we findEQUATIONwhere M and B are given by ( 41) and ( 44) respectively.Bounding the remainder term: We now proceed to bound the remainder term (60). Using (61), for q < 1 we can compute|E α,m (q)| ≤ q m 1 0 (1 -s) m-1 (m -1)! d m dq m F α (sq) ds (81) ≤q m m j=0 m j j-1 ℓ=0 |α -ℓ| m-j-1 ℓ=0 (α + ℓ -1) × 1 0 (1 -s) m-1 (m -1)! sqN ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) + (1 -sq) α-j sqN ∆µ ′ t (b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) + (1 -sq) α+m-j-1 × N ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 j N ∆µ ′ t (b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 m-j N 0, σ 2 t (θ)dθds ≤ q m m! m j=0 (1 -q) -(α+m-j-1) m j j-1 ℓ=0 |α -ℓ| m-j-1 ℓ=0 (α + ℓ -1) K α,m,j (q) , K α,m,j (q) := m 1 0 (1 -s) m-1 sqN ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) + (1 -sq) α-j × N ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 j N ∆µ ′ t (b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 m-jN 0, σ 2 t (θ)dθds , where to obtain the final inequality we used the fact that α + m -j -1 > 0 andEQUATION)to upper bound the reciprocal of the denominator.To bound the integrals K α,m,j we consider the following two cases.Case 1: α -j ≤ 0:In this case we usesqN ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) + (1 -sq) α-j ≤ (1 -q) α-j (83)to bound the first factor in the integrand. Then, using Lemma C.7 and Lemma C.4, we can compute K α,m,j (q) (84)≤(1 -q) α-j N ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 j N ∆µ ′ t (b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 m-j N 0, σ 2 t (θ)dθ ≤(1 -q) α-j B σt,m .Case 2: α -j > 0:In this case use the bound x α-j ≤ 1 + x ⌈α⌉-j for all x > 0, followed by the binomial formula and Lemma C.7 to computeK α,m,j (q) ≤ m 1 0 (1 -s) m-1 1 + sqN ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) + (1 -sq) ⌈α⌉-j (85) × N ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 j N ∆µ ′ t (b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 m-j N 0, σ 2 t (θ)dθds = N ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 j N ∆µ ′ t (b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 m-j N 0, σ 2 t (θ)dθ + m ⌈α⌉-j ℓ=0 q ℓ ⌈α⌉ -j ℓ 1 0 (1 -s) m-1 s ℓ N ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 ℓ × N ∆µt(b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 j N ∆µ ′ t (b ′ , b), σ 2 t (θ)/N 0, σ 2 t (θ) -1 m-j N 0, σ 2 t (θ)dθds ≤ B σt,m + ⌈α⌉-j ℓ=0 q ℓ (⌈α⌉ -j)!m! (⌈α⌉ -j -ℓ)!(m + ℓ)! B σt,ℓ+m .Putting the above two cases together we arrive at the remainder boundEQUATION)Combining the above results, we obtain the following FS woR -RDP bounds under replace-one adjacency:sup θt D α (p(θ t+1 |θ t , D)∥p t (θ t+1 |θ t ), D ′ ) (88) ≤ 1 α -1 log 1 + q 2 α(α -1) e 4/σt 2 -e 2/σt 2 + m-1 k=3 q k k! F α,σt,k + E α,σt,m (q) ,where m ≥ 3 is an integer, F α,σt,k is given by ( 80), and E α,σt,m (q) is given by ( 87); both of these quantities are expressed in terms of M and B, as given by ( 41) and ( 44) respectively. This completes the proof of Theorem 3.4.Remark C.8. The same M 's appear many times during in the formulas for F α,σt,k and E α,σt,m (q). Therefore, in practice, we compute all required values once and recall them as needed.In this appendix we provide a precise comparison between RDP for fixed-size and Poisson subsampling under replace-one adjacency, complimenting the intuitive discussion in Section 4.1.Here we let p Poisson t (θ t+1 |θ t , D) denote the transition probabilities for one step of DP-SGD with Poisson subsampling. By Mirroring the analysis in Section 2 of Mironov et al. (2019) , but now under the assumption that D and D ′ have the same size and differ in only their last elements, we arrive at the one-step Rényi divergence boundD α p Poisson t (θ t+1 |θ t , D)∥p Poisson t (θ t+1 |θ t , D ′ ) (89) ≤ max T D α qN ∆µt(T ), σ 2 t + (1 -q)N 0, σ 2 t ∥qN ∆µ ′ t (T ), σ 2 t + (1 -q)N 0, σ 2 t , ∆µ t (T ) := µ t (θ t , T ∪ {|D| -1}, D) -µ t (θ t , T, D) , ∆µ ′ t (T ) := µ t (θ t , T ∪ {|D| -1}, D ′ ) -µ t (θ t , T, D) ,where T denotes the collection of indexes from {0, ..., |D| -2} that were randomly selected during Poisson sampling.This is similar to the result ( 56) for FS woR -subsampling, except here the means satisfy the constraintsEQUATIONand these bounds are achieved when the differing elements from D and D ′ are anti-parallel and saturate the clipping threshold C. As we will show, this difference only contributes at higher order.We now mirror the analysis from earlier in this section, taking the new constraints (90) into account.We again have a Taylor expansion of the formEQUATION)The coefficient at k = 2 can be computed as in (62), which givesd 2 dq 2 F α (0) (93) =α(α -1) exp ∥∆µ t (T )∥ 2 σ 2 t + exp ∥∆µ ′ t (T )∥ 2 σ 2 t -2 exp ∆µ t (T ) • ∆µ ′ t (T ) σ 2 t ≤2α(α -1) e 1/σ 2 t -e -1/σ 2 t .This bound is optimal in that equality is achieved when the distinct elements in D and D ′ are anti-parallel and saturate the clipping bound. To leading order in 1/σ 2 t we therefore haveEQUATION)which agrees with the leading order behavior of the corresponding term in the FS woR -subsampling case (66). As the second order coefficient in the Taylor expansion controls the dominant behavior of the RDP bound, we can therefore conclude that Poisson and FS woR subsampling give the same DP bounds under replace-one adjacency to leading order.Our method for bounding the higher order terms, as in (80), and remainder, as in ( 86), only uses bounds on the norms of ∆µ ′ and µ ′ ; it does not incorporate the bound that couples these two quantities. Therefore the corresponding bounds on the higher order terms for Poisson subsampling are obtained by simply replacing r t with r t /2 in our prior calculations; this is equivalent to replacing σ t with 2σ t in ( 80) and ( 86). Therefore we arrive at the following non-asymptotic bound for Poisson subsampling under replace-one adjacency.Theorem C.9 (Poisson-RDP Upper Bounds for Replace-one Adjacency). Let D ≃ r-o D ′ be adjacent datasets. Assuming q := |B|/|D| < 1, for any integer m ≥ 3 we havesup θt D α (p Poisson t (θ t+1 |θ t , D)∥p Poisson t (θ t+1 |θ t , D ′ )) (95) ≤ 1 α -1 log 1 + q 2 α(α -1) e 1/σ 2 t -e -1/σ 2 t + m-1 k=3 q k k! F α,2σt,k + E α,2σt,m (q) ,where F α,2σt,k and E α,2σt,m (q) are computed via (80) and ( 87) respectively (with σ t replaced by 2σ t ).While our FS woR and Poisson-subsampled RDP bounds agree to leading order, the effect of the higher order terms leads to Poisson subsampling having a slight privacy advantage over FS woR -subsampling in practice, especially when larger α's are required. This proven in Figure 7 , where we see that the upper bound on Poisson-subsampled RDP (dashed line) is slightly below the theoretical lower bound for FS woR -subsampled RDP derived in In this Appendix we provide a detailed derivation of Rényi DP-SGD bounds when using fixed-size subsampling with replacement (FSwR). In addition to the lower bound stated in Theorem 3.7, we also derive an upper bound.As in case without replacement in Appendix A, a key step is a decomposition of the sampling distribution into a baseline (where D and D ′ agree) and a perturbation where they disagree. The derivation for subsampling with replacement is a non-trivial generalization of sampling without replacement, with Lemma D.1 below requiring a several key new ingredients, as compared to Lemma A.1.First let D ′ be obtained from D by removing one element. Without loss of generality, assume it is the last (i.e., at index |D| -1) and call that element d.Let B ′ = (B ′ 1 , ..., B ′ |B| )where the components are iid Uniform({0, ..., |D| ′ -1}) random variables, so that B ′ uniformly samples from {0, ..., |D| ′ -1} with replacement. Now let N ∼ Binomial(|B|, |D| -1 ) and Π be a uniformly random permutation of (1, ..., |B|). Assume that B ′ , N , and Π are independent. In the following lemma we show how sampling from D with replacement can be obtained from B ′ , N , and Π. In short, N will determine the number of entries of B ′ are replaced by |D| -1 and Π will determine the indices of the entries that are to be replaced. Lemma D.1. For i = 1, ..., |B| define EQUATIONP (B = b) =P (B = b, N = n, {Π(1), ..., Π(n)} = {i 1 , ..., i n }) (97) = τ ∈Sn P B i = b i for i ̸ ∈ {i 1 , ..., i n }, N = n, (Π(1), ..., Π(n)) = (i τ (1) , ..., i τ (n) ) = τ ∈Sn P B ′ i = b i for i ̸ ∈ {i 1 , ..., i n }, N = n, (Π(1), ..., Π(n)) = (i τ (1) , ..., i τ (n) ) .Using independence and that |D ′ | = |D| -1 and |S n | = n! we can then computeP (B = b) (98) = τ ∈Sn   i̸ ∈{i1,...,in} P (B ′ i = b i )   P (N = n)P (Π(1), ..., Π(n)) = (i τ (1) , ..., i τ (n) ) = τ ∈Sn |D ′ | -(|B|-n) |B| n |D| -n (1 -|D| -1 ) |B|-n (|B| -n)! |B|! =n!|D| -|B| |B| n (|B| -n)! |B|! =|D| -|B| .This holds for all b, therefore we can conclude that the distribution of B is the uniform distribution on {0, ..., |D| -1} |B| . All the remaining claims then immediately follow from this fact.Using Lemma D.1 we can decompose the transition probabilities for one step of DP-SGD with FS wR -subsampling as followsEQUATIONEQUATION= |B| n=0 |B| n |D| -n (1 -|D| -1 ) |B|-n N µt(θt,Φ(b ′ ,n,π),D), σ 2 t (dθ t+1 )P B ′ ,Π (db ′ dπ) . Noting that µ t (θ, b ′ , D ′ ) = µ t (θ, b ′ , D), as D ′ differs from D only by removal of the last element, we now use quasiconvexity to computeD α (p t (dθ t+1 |θ t , D)∥p t (dθ t+1 |θ t , D ′ )) (101) ≤ max b ′ ,π D α   |B| n=0 |B| n |D| -n (1 -|D| -1 ) |B|-n N µt(θt,Φ(b ′ ,n,π),D), σ 2 t N µt(θt,b ′ ,D), σ 2 t   .Remark D.2. As in Remarks A.2 and C.1, Eq. ( 101) remains true when the Rényi divergences D α are replaced by any other divergence D that is quasiconvex in both of its arguments, e.g., the hockey-stick divergences. This fact is useful for other differential privacy paradigms.The next stages in the calculation are more complex than in the case of sampling without replacement, as the first argument of the above Rényi divergence is a mixture of |B| + 1 Gaussians, as opposed to simply two as in the case without replacement. This is because the binomial random variable N in (D.1) is playing the role that the Bernoulli random variable J did in Lemma A.1.Start by rewriting Eq. ( 101) asEQUATIONand decompose the mixture of distributions as follows|B| n=0 a n N ∆µn, σ 2 t =a 0 N 0, σ 2 t + K n=1 a n N ∆µn, σ 2 t + |B| n=K+1 a n N ∆µn, σ 2 t (104) = a 0 N 0, σ 2 t + K n=1 a n qN ∆µn, σ 2 t + (1 -q)N 0, σ 2 t + |B| n=K+1 a n N ∆µn, σ 2 t , a 0 := a 0 -(1/ q -1) K n=1a n , a n := a n / q for n ≥ 1 .The assumption (103) implies a n ∈ [0, 1] and a 0 + K n=1 a n + |B| n=K+1 a n = 1 and therefore we can use convexity to bound (102) as follows:D α (p t (dθ t+1 |θ t , D)∥p t (dθ t+1 |θ t , D ′ )) (105) ≤ max b ′ ,π 1 α -1 log a 0 + K n=1 a n q N ∆µn, σ 2 t (θ) N 0, σ 2 t (θ) + (1 -q) α N 0, σ 2 t (θ)dθ + |B| n=K+1 a n N ∆µn, σ 2 t (θ) N 0, σ 2 t (θ) α N 0, σ 2 t (θ)dθ   , ≤ 1 α -1 log a 0 + K n=1 a n q N 1,σ 2 t /(4n 2 ) (θ) N 0,σ 2 t /(4n 2 ) (θ) + (1 -q) α N 0,σ 2 t /(4n 2 ) (θ)dθ + |B| n=K+1 a n e 2α(α-1)n 2 /σ 2 t   = 1 α -1 log   a 0 + K n=1 a n H α,σt/n ( q) + |B| n=K+1 a n e 2α(α-1)n 2 /σ 2 t   ,where H was defined in Eq. ( 8) and can be bounded using the techniques from Appendix B. D α (p t (dθ t+1 |θ t , D)∥p t (dθ t+1 |θ t , D ′ )) (106) ≤ 1 α -1 log   a ′ 0 + K n=1 a ′ n H α,σt/n ( q) + |B| n=K+1 a ′ n e 2α(α-1)n 2 /σ 2 t   , a ′ n := |B| n |D ′ | -n (1 -|D ′ | -1 ) |B|-n , a ′ 0 := a ′ 0 -(1/ q -1) K n=1 a ′ n , a ′ n := a ′ n / q for n ≥ 1 ,whenever q satisfies (103). Note that (103) also impliesq ≥ 1 + a ′ 0 K n=1 a ′ n -1 (107)due to the fact that a 0 is increasing in |D| and a n is decreasing in |D| for n ≥ 1. These properties also implya ′ 0 + K n=1 a ′ n H α,σt/n ( q) + |B| n=K+1 a ′ n e 2α(α-1)n 2 /σ 2 t (108) ≤ a 0 + K n=1 a n H α,σt/n ( q) + |B| n=K+1a n e 2α(α-1)n 2 /σ 2 t , and therefore we can conclude thatEQUATION≤ 1 α -1 log   a 0 + K n=1 a n H α,σt/n ( q) + |B| n=K+1a n e 2α(α-1)n 2 /σ 2 t   whenever D ′ ≃ a/r D and any choice of q that satisfiesEQUATION)Note that if one chooses q = 1 +a0 K n=1 an -1 then a 0 = 0. If one also takes K = |B| then q = 1 -(1 -|D| -1) |B| . This is the version of the result we present above in Theorem 3.6.Here we derive a lower bound on the one-step RDP of SGD with fixed-size subampling with replacement. To do this, suppose that all entries of D ′ are the same, equal to d ′ , and that the extra entry d in D satisfies d ̸ = d ′ . Moreover, suppose there exists θ t such that the clipped gradients at d and d ′ are anti-parallel with both having norm equal to the clipping bound, C. In this case, µ t (θ t , b ′ , D) and µ t (θ t , Φ(b ′ , n, π), D) are independent of both π and b ′ andEQUATIONThis also implies that ( 101) is an equality (and the argument of the divergence is independent of b ′ and π). By changing variables and using (111) we can therefore writeD α (p t (dθ t+1 |θ t , D)∥p t (dθ t+1 |θ t , D ′ )) (112) = max b ′ ,π D α   |B| n=0 |B| n |D| -n (1 -|D| -1 ) |B|-n N ∥µt(θt,Φ(b ′ ,n,π),D)-µt(θt,b ′ ,D)∥, σ 2 t N 0, σ 2 t   =D α   |B| n=0 |B| n |D| -n (1 -|D| -1 ) |B|-n N nrt, σ 2 t N 0, σ 2 t   =D α   |B| n=0 |B| n |D| -n (1 -|D| -1 ) |B|-n N n,σ 2 t /4 N 0,σ 2 t /4   .This gives the one-step RDP lower boundEQUATIONThe integral over θ can be computed as follows. First note summation in ( 113) is the average of N n,σ 2 t /4 (θ)/N 0,σ 2 t /4 (θ) under the Binomial distribution Binomial(|B|, |D| -1 ). As a convenient shorthand we write this average as ...dn. With this notation, for integer α we can rewrite the argument of the logarithm in the lower bound (113) an iterated integral, then exchange the order of the integrals and evaluate the integral over θ using the formula for the MGF of a Gaussian as follows EQUATIONTherefore we arrive at the lower boundEQUATIONRecalling the definition of the dn i 's, this completes the proof of Theorem 3.7.The lower bound (115) behaves somewhat differently from the one derived for subsampling without replacement in Wang et al. (2019) and restated in Appendix E.2. In particular, note that the explicit |B| and |D| dependence in the lower bound ( 113) cannot be removed simply by reparameterizing the expression in terms of q = |B|/|D|. In contrast, the FS woR -RDP lower bound in Appendix E.2 depends on |B| and |D| only through q. This fact can be further illustrated through a simple loosening of the bound (113), where we bound the sum below by the term at n = |B| to obtainϵ ′ t (α) ≥ 1 α -1 log |D| -α|B| N |B|,σ 2 t /4 (θ) N 0,σ 2 t /4 (θ) α N 0,σ 2 t /4 (θ)dθ (117) = α|B| α -1 2|B|(α -1)/σ 2 t -log |B| -log(1/q) .For fixed ratio q = |B|/|D|, the right-hand-side of (117) approaches ∞ as |B| → ∞. Therefore we can conclude that ϵ ′ t (α) has nontrivial |B| dependence, even for fixed q as claimed. This accounts for some of the difference in behavior between fixed-size subsampling with and without replacement. This behavior is also illustrated numerically in Figure 8 , where we show the exact value of the bound (113) as a function of |B| for α = 2 (this case is straightforward to compute via the binomial theorem and MGF of a Gaussian) and for fixed q = 0.001; we show results for several values of σ.Note that initially each curve starts out horizontal, mirroring the case of FS woR -RDP, which has no |B| dependence outside of q. However, for each σ there is a threshold value of |B| as which there is a "phase transition", i.e., after which the explicit |B| dependence becomes important and the RDP value increases rapidly.When α is an integer the lower bound (116) can, in principle, be computed exactly by recalling that the integrals with respect to dn i are shorthand for the expectation with respect to the distribution Binomial(B, |D| -1 ). Therefore the lower bound (116) consists of α nested sums of |B| + 1 terms each. Clearly, when α is larger than 2 and |B| is not small iterated summation this quickly becomes computationally impractical. The looser lower bound ( 117) is easy to compute, but is much too inaccurate for most purposes. Here, for integer α, we will show how obtain lower bounds that are both practical and significantly more accurate.Recall that the integrals with respect to the dn's are actually finite sums, each over |B| + 1 elements, and so in principle the above expression can be computed exactly but it is still not practical. We next show how to obtain recursively defined lower bounds that are significantly less computationally expensive than the O(|B| α ) computations that the more naive formula (116) requires. To begin, for any integer k ≥ 2 defineEQUATION)For k = 2 we can reduce the computation to a single summation by using the formula for the MGF of the Binomial distribution:EQUATION)For k > 2 we explicitly write out the integral over n k as a summation and then simplify to obtainF k (c, d) (120) = ... e c k-1 j=2 j-1 i=1 ninj +d k-1 i=1 ni e (c k-1 i=1 ni+d)n k dn k dn k-1 ...dn 1 = |B| n=0 |B| n |D| -n (1 -|D| -1 ) |B|-n ... e c k-1 j=2 j-1 i=1 +d k-1 i=1 ni e (c k-1 i=1 ni+d)n dn k-1 ...dn 1 = |B| n=0 |B| n |D| -n (1 -|D| -1 ) |B|-n e dn ... e c k-1 j=2 j-1 i=1 ninj +(d+cn) k-1 i=1 ni dn k-1 ...dn 1 = |B| n=0 |B| n |D| -n (1 -|D| -1 ) |B|-n e dn F k-1 (c, d + cn) .As stated, the recursive formula provided by ( 119) and ( 120) is still too computationally expensive due to each recursive call requiring the summations over |B| + 1 elements. However, noting that terms in the summations are all positive, we can obtain a practical lower bound by retaining only a small number of terms. More specifically, given choices of T k ⊂ {0, ..., |B|} for each k, if we recursively defineEQUATIONwithEQUATIONThen, by induction, we have F k (c, d) ≥ F T,k (c, d) for all k, T, c, d and thereforeEQUATIONfor all α and all choices of T k , k = 2, ..., α. By appropriately choosing the number of terms T k for each k one trades between accuracy and computational cost. In particular, if each |T k | = 1 then the computation of (123) only requires O(α) computations. If |T k | = ℓ > 2 for all k then the computation of (123) required O(ℓ α-1 ) computations. In contrast, using the exact recursive formula (120) requires O(|B| α-1 ) computations, which is intractable unless α or |B| is sufficiently small. However, despite this improvement, we should note that computational difficulties persist for cases where the number of nontrivial terms is too large.In our experiments in Figure 9 we find that that using T k = {|B|} for k > 2 and T 2 = {0, ..., |B|} gives an accurate and computationally tractable lower bound. In the example in Figure 1 more terms were required; we show results obtained by using T k = {0, 1, 2, |B|} for k > 2 and T 2 = {0, ..., |B|}.Now we will compare our lower bound for FS wR -RDP from Theorem 3.7 with the upper bound derived in Appendix D.1. In Figure 9 we show the FS wR -RDP upper bound (solid line) and lower bound (circles) and, for comparison, we also show our bound on FS woR -RDP (dashed line). Note that all three are very similar when α is close to 1; this is to be expected, due to them having the same leading order behavior in q, given by ( 19). However, as α increases the non-leading order behavior becomes more important and we observe that, after some threshold α, both the FS wR -RDP upper and lower bounds increase rapidly. A similar "phase transition" was observed earlier in Wang et al. (2019) for subsampling without replacement, though for these parameters it occurs at much higher values of α. A key difference between these two cases is that, for FS wR -RDP, the phase transition can be brought on by increasing |B| while leaving q = |B|/|D| fixed, as demonstrated in Appendix D.2; see (117) and Figure 8 . In contrast, the lower bound for FS woR -RDP from Wang et al. (2019) depends on |B| and |D| only through their ratio q. We anticipate that the gap between the upper and lower bounds on FS wR -RDP in Figure 9 could be reduced by using a higher order expansion in q, but we leave that to future work. We also emphasize that the dashed FS woR -RDP upper bounds is only for comparison purposes as it does not apply to the FS wR case considered in this section. E Comparison with Wang et al. (2019) In Wang et al. (2019) , Rényi-DP enhancement of fixed-size subsampling was analyzed in a framework that applies to arbitrary mechanisms. In this appendix we specialize their results to fixed-size subsampling DP-SGD with Gaussian noise. This will facilitate comparison with our results. Specifically, we will show that our specialized method results in a factor tighter bounds by approximately a factor of 4.E.1 Upper Bound from Wang et al. (2019) , Theorem 9For integer α ≥ 2, the following the general fixed-size subsampling Rényi-DP enhancement bound (for one step) was derived in Wang et al. (2019) :EQUATIONK t (α) := log   1 + 2q 2 α(α -1)(e ϵt(2) -1) + α j=3 2q j j-1 ℓ=0 (α -ℓ) j! e (j-1)ϵt(j)   (if α = 2 then the empty summation is interpreted to be zero), where ϵ t (j) denotes a Rényi-DP bound for the mechanism without random subsampling. For DP-SGD we haveEQUATION)and therefore we can useEQUATIONNote that (125) uses the same bound on the means that we use in our method, and therefore (126) will yield a fair comparison. Comparing the leading order behavior of (124) with the leading order behavior of a single step of our result in Theorem 3.3, we see that our bound is smaller by a factor of 4. In practice we find similar behavior even when including the higher order terms; see Figure 2 .For non-integer α ≥ 2 they use the following bound, which follows from convexity of K t (α) in α:EQUATIONwhere ⌊x⌋ denotes the floor of x. For α ∈ (1, 2) the convexity bound isEQUATIONwhich simplifies toEQUATIONE.2 Lower Bound from (Wang et al., 2019) , Proposition 11A general lower-bound on fixed-size subsampling enhancement was also derived in (Wang et al., 2019) ; we repeat it here for convenience. For integer α ≥ 2 they showEQUATIONwhere ϵ t (j) is again given by (126). In Figure 2 we demonstrate that our result is very close to this theoretical RDP lower bound.The analysis in Wang et al. (2019) is done under the replace-one adjacency definition and is therefore directly compatible to our Theorem 3.4. In Section 3.4 we showed that our replace-one and our add/remove adjacency results, Theorem 3.4 and 3.1 respectively, lead to the same RDP bounds to leading order in q and 1/σ 2 t , including the factor of 4 improvement over Wang et al. (2019) ; this confirms that the improvement made by our approach is not reliant on the adjacency definition.In this appendix we provide a conjecture as to the origin of this difference. The leading order behavior comes from the calculations ( 64)-( 66), where we used the constraints on the means from (57):EQUATIONHowever, if one doesn't use the constraint on the difference and only uses the first two upper bounds on the means then the best bound on the dot product becomesEQUATION)and so the bound ( 66) is weakened toEQUATIONwhich at leading order is a factor of 4 larger than our result (66) which was derived using all three constraints in (57). The strategies we employ differ from Wang et al. (2019) in a way that makes it difficult to pinpoint exactly what the essential differences are, however the above calculation leads us to conjecture that the difference originates in the failure to account for the constraint that couples ∆µ t and ∆µ ′ t . See also Remark C.6.In this appendix we provide detailed calculations of the variance comparison that was discussed in Section 4.2. Let a i ∈ R, i = 1, ..., |D| (for application to DP-SGD, one can take a i = ∇ θ L(d i ) • v,Z P := 1 |B| |D| i=1 a i J i (134)with that of fixed-size subsampling with replacement,EQUATIONand with fixed-size subsampling without replacement,EQUATIONThey all have the same expected value, equaling the average of the samples:EQUATIONEQUATIONEQUATIONThe variances areVar[Z P ] = 1 |B| 2 |D| i,j=1 a i a j E[J i J j ] -a 2 (140) = 1 |B| 2 i a 2 i E[J i ] + 1 |B| 2 i̸ =j a i a j E[J i ]E[J j ] -a 2 = 1 |B||D| i a 2 i + 1 |D| 2 i̸ =j a i a j - 1 |D| 2 i,j a i a j = 1 |B| - 1 |D| 1 |D| |D| i=1 a 2 i , Var[Z R ] = 1 |B| 2 |B| ℓ,j=1 E[a Rj a R ℓ ] -a 2 = 1 |B| 2 |B| j=1 E[a 2 Rj ] + 1 |B| 2 ℓ̸ =j E[a Rj ]E[a R ℓ ] -a 2 = 1 |B| 2 |B| j=1 E[a 2 Rj ] + 1 |B| 2 ℓ̸ =j a 2 -a 2 = 1 |B| 2 |B| j=1 E[ i 1 Rj =i a 2 i ] + 1 |B| 2 (|B| 2 -|B|)a 2 -a 2 = |D| i=1 a 2 i 1 |B| 2 |B| j=1 E[1 Rj =i ] + (1 -1/|B|)a 2 -a 2 = |D| i=1 a 2 i 1 |B| 2 |B| j=1 P (R j = i) - 1 |B| a 2 = 1 |B||D| |D| i=1 a 2 i - 1 |B| a 2 =(1 -|B|/|D|) -1 Var[Z P ] - 1 |B| a 2 ,andVar[Z B ] = 1 |B| 2 |D| i,j=1 a i a j P (i ∈ B, j ∈ B) -a 2 (141) = 1 |B| 2 i a 2 i P (i ∈ B) + 1 |B| 2 i̸ =j a i a j P (i ∈ B, j ∈ B) -a 2 = q |B| 2 i a 2 i + 1 2 i̸ =j a i a j |D|-2 |B|-2 |D| |B| -a 2 = 1 |B| - 1 |D|   1 |D| i a 2 i - 1 |D|(|D| -1) i̸ =j a i a j   = |D| |D| -1 Var[Z P ] - 1 |B| - 1 |D| a 2 .Therefore, defining a 2 = 1 |D| |D| i=1 a 2 i , when |B| < |D| and a 2 ̸ = 0 we obtain the following variance relationsEQUATIONEQUATIONandEQUATIONIn particular, we see that fixed-size subsampling with replacement always has larger variance than fixed-size subsampling without replacement by the constant multiple 1-1/|D| 1-q > 1; however, this factor is close to 1 when q is small. Therefore, for practical purposes we can consider them equivalent from the perspective of the variance. When |B| < |D| the ratio of fixed-size and Poisson variances isEQUATION)The Cauchy-Schwarz inequality implies 0 ≤ a 2 /a 2 ≤ 1 and so we see that Var[Z B ] < Var[Z P ] when a ̸ = 0 and |D| is sufficiently large; in DP-SGD the condition a ̸ = 0 corresponds to being away from a (local) minimizer or critical point of the loss. This calculation suggests that fixed-sized minibatches are better to use when one is away from the optimizer, as they lead to a (potentially substantial) reduction in variance, but they can become (slightly) worse than Poisson-sampled minibatches when one is close to a local minimizer (i.e., a ≈ 0). However, for a typical training-set size |D| ≫ 1 we have (1 -1/|D|) -1 ≈ 1 and so the variance advantage of Poisson subsampling when a ≈ 0 is very minor. In addition, due to the effect of noise, DP-SGD cannot reach a = 0. Therefore, from the perspective of the variance, we conclude that fixed-size subsampling (with or without replacement) is generally preferable over Poisson subsampling.Using the conversion from Proposition 3 in Mironov (2017) , if a mechanism satisfies ϵ ′ (α)-RDP then it satisfies (ϵ, δ)-DP withEQUATION)Under add/remove adjacency, to leading order, both Poisson and fixed-size subsampling satisfy RDP bounds with ϵ ′ (α) = αR for some R > 0 (see Abadi et al. (2016) ; Mironov et al. (2019) ) and ( 19).Assuming that the optimal α is not close to one, and hence we can approximate α with λ := α -1, minimizing over α givesEQUATION)For t steps of fixed-size subsampling (with or replacement), and with constant variance σ, Eq. ( 19) impliesEQUATION)In contrast, for Poisson subsampling the results in Abadi et al. (2016) ; Mironov et al. (2019) implyEQUATION)Therefore R F S ≈ 4R P at leading order and so (147) impliesEQUATIONFrom a privacy perspective, this suggests that ϵ for Poisson subsampling has an intrinsic advantage over fixed-size subsampling by approximately a factor of 2 when using the add/remove adjacency relation. Noting the scaling of ( 150) with σ, we see that fixed-size subsampling requires noise 2σ to obtain the same privacy guarantee (at leading order) as Poisson sampling with noise σ. Tighter translation between RDP and (ϵ, δ)-DP, such as Theorem 21 in Balle et al. (2020) which is used by Meta Platforms (2024) , alters this story somewhat but the general picture we have outlined here is corroborated by the empirical comparison provided in Section 4.2. We emphasize that this factor of two gap disappears when using replace-one adjacency; in that case Poisson subsampling and FS woR subsampling lead to the same DP guarantees, as discussed in Section 4.1.To assess the sensitivity of our proposed method to important parameters such as σ and batch size |B|, we conducted a series of additional experiments to compare our proposed fixed-size method with the non-fixed size setting while varying values of B ∈ {50, 100, 150, 200} and σ ∈ {3.0, 4.5, 6.0, 12.0}. Justification: We noted in the main body that Figure 4 was plotted with 3-sigma standard deviation around mean based on five runs each with different random seeds for sampling: 0, 1, 364, 2, 560, 3, 000, and 4, 111. Also, the code that generates the plots were provided in our provided source code. Guidelines:• The answer NA means that the paper does not include experiments.• The authors should answer "Yes" if the results are accompanied by error bars, confidence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper. • The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions). • The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.) • The assumptions made should be given (e.g., Normally distributed errors).• It should be clear whether the error bar is the standard deviation or the standard error of the mean. Justification: We envision no negative impacts of our work. Not human subjects was used. Datasets are publicly available. Also, the code associated with this work is not envisoned to be used for malicious purposes. Guidelines:• The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.• If the authors answer No, they should explain the special circumstances that require a deviation from the Code of Ethics.• The authors should make sure to preserve anonymity (e.g., if there is a special consideration due to laws or regulations in their jurisdiction). 


Decision: The paper studies the privacy accounting of DP-SGD under the sampling-without-replacement model but unlike in previous work that worked with a "replace" neighbor relationship, this paper uses an "add/remove" neighbor relationship. 

The main claim is that for DP-SGD, they can improve the best known Renyi-DP bound by a factor of 4 when the order $\alpha$ is smaller than a threshold and that (at least numerically) the restricted range $\alpha$ still supports stronger composition in DP-SGD.

Reviewers acknowledged the importance of the problem and the claimed improvement but found the following aspects unsatisfying.

1. the use of "add/remove" relationship in amplification by WOR sampling was studied in a prior work [ZDW'22], with a direct analysis of $(\epsilon,\delta)$-DP by analyzing the hockey-stick divergence of a dominating pair.

2. One reviewer is skeptical of the claim that the particular Gaussian mixture being a dominating pair for the WOR sampling case in the Renyi ordering, as it was shown in the recent work of Lebeda et al that it isn't in the "Blackwell ordering".  The same reviewer also comments on other issues in Appendix D.4. such as:

> "I don't see where the bound of Eq. (108) comes from, why is there no additional factor 2 (or 4 after squaring)."
 
3.  The techniques used are similar to those in Mironov et al. 2019, and Wang et al. 2019.

I only have a limited amount to time to look into these issues. For Point 1,  I do believe that a tighter RDP analysis has its own value even if it does not give a tight epsilon-delta DP bound as the Hockey-Stick Divergence. 

- If the authors add a comparison of the RDP-converted approx-DP privacy profile to the latest FFT-based numerical account as in Lebeda et al., it should be sufficient in addressing this concern.

Point 3 was acknowledged in the paper. The particular coupling used is similar to Wang et al. 2019 and the Taylor expansion of the mixture of Gaussian into quantities one can evaluate with Gaussian CDF calls is due to Mironov et al. 2019.  The authors however clearly discussed the differences in the proof.

In the discussion period, the reviewers and the AC tried to resolve Point 2, but it requires more time than we had. However, in general, dominance in Blackwell ordering and dominance in Renyi ordering are different. It is perfectly possible that the mixture of gaussian pair does not dominate in Blackwell ordering (as was shown in prior work) but it dominates in Renyi ordering.  We did not catch any obvious error in the proof of Theorem 3.1, so it should check out! 

I wanted to add that the authors also proved new results about the "replace" neighbor in appendix D.4, which the authors wrote in the rebuttal that they should have focused a bit more on it.  I agree that "replace" neighbor is more natural for random subset sampling.   The paper ran into the same issue of having mixture of two Gaussians in both the numerator and denominator as in Wang et al. 2019 when applying joint convexity directly.   Wang et al., 2019 went through a ternary-type divergence and then the binomial expansion, which resulted in a general result with more precise control of the higher order terms, but incurred a factor of 4 in the leading term.  The current paper directly applied Taylor expansion twice (first over q near q = 0, then for e^x near x = 0) in a way similar to the approach Abadi et al (2016) took. This argument ended up resulting in an expression that allows one to complete the square in Eqn (107), hence reducing the factor of 2 ( factor of 4 after squaring) in the leading term.   The reviewer's confusion of how (108) avoids the factor of 2 stems from that the sensitivity is now defined to be 2x the clipping threshold. This got clarified in the AC-reviewer discussion.

I think this is a nice step forward on the problem introduced 5 years ago by Wang et al 2019.  For this reason, I think the paper deserves to be published.

------------
Other comments the authors could consider:

1.  For the Add/Remove case, q is different for the two neighbor datasets.   |B| / n vs |B| / (n-1).  For the definition to make sense, one needs to either make n known or define the space of the dataset to be n > certain public threshold, then state the results according to that pubic threshold.  This was discussed in [ZDW'22].

2. Because of the above, it is more natural to use the Replace neighbors for WOR, and the result hidden in Appendix D.4. could be emphasized and moved to the main paper?

3. Results in Appendix D.4 are worth developing a bit further so that the higher-order terms in the Taylor expansions are replaced with valid upper bounds.  I also suspect the bound to be weaker than that of Wang et al 2019 when we are outside the regime of \sigma being large and q being very small.

4. As I mentioned before, numerically adding comparison to the FFT-based numerical accountants that skips RDP would make it clear the practical impact of the result here.

Reviews: review: paper_summary: paper_summary: This paper studies the Renyi DP guarantees of a with- or without replacement subsampled Gaussian mechanism. Authors present a privacy analysis tailored for the subsampled Gaussian mechanism, which improves the earlier general bound for $\epsilon(\alpha)$ by Wang et al. 2019 by a factor of four. Authors show analytically, that the subsampling induced variance (i.e. the variance of the noise arising from using minibatches instead of full data) is smaller for without replacement subsampling than for the more commonly used Poisson subsampling. Authors also give a theoretical analysis of the differences between WOR and Poisson subsampling, showing that the Poisson subsampling leads to approximately half the $\epsilon$ of the WOR sampling. Authors demonstrate empirically that for fixed noise level, the WOR leads to better accuracy than Poisson subsampling in a CIFAR10 based deep learning task, suggesting that the difference is due to the lower subsampling noise. Finally, authors show that the WOR sampling leads to more stable memory usage.
strength_weakness: strength_weakness: While the presented analysis provides important insights and improvements over the previous RDP analysis for WOR sampled Gaussian mechanism, I wonder if the problem is already solved with the modern privacy accounting tools. For example, Zhu et al. 2022 solve the WOR sampled Gaussian mechanism in their characteristic function formalism. Given that the conversion from RDP to approximate DP is lossy, I would imagine the bounds presented in this work are more loose than the bounds of Zhu et al when converted to approx DP domain. Since the approx DP bounds are more commonly used than RDP bounds, I'm not sure if the tight RDP analysis is really needed. Or is there some practical reason, e.g. an implementation difficulty, that would not allow using the characteristic function accounting for this problem?

Zhu, Yuqing, Jinshuo Dong, and Yu-Xiang Wang. "Optimal accounting of differential privacy via characteristic function." International Conference on Artificial Intelligence and Statistics. PMLR, 2022.
questions: questions: - As you acknowledge towards the end of Section 3.5, Wang et al. 2019 use the substitute neighbourhood relation which differs from the add/remove used in this paper. Can you clarify, is this difference taken into account in Figure 3?
- Fig. 2: I'm a bit confused on the upper bounds show in this Figure. It seems that after some $\alpha$, your proposed upper bound exceeds the one from Wang et al. Does this suggest that the Wang et al. bound is tighter in some regime of $\alpha$?
-  In the Appendix D, you derive the RDP bounds for the substitute relationship using mixture of Gaussians with weight $1-q$ on Gaussian centered at $0$ and weight $q$ for the "adjacent" point. A recent work by Lebeda et al. 2024, suggest that this might not be the worst-case distribution for woR-sampling (see Section 7 in Lebeda et al.). While their analysis is tailored for approximate-DP and not RDP, I wonder if the same holds RDP.

**Typos and other minor things**
- Eq. 57, extra parenthesis in the expression for $M_{\sigma, 4}$
- I guess the NN abbreaviation is never explained. However, I don't think using it is necessary for the paper to begin with as you analysis applies to any learning task using DP-SGD.
- Fig. 2: the caption is overlapping with the axis label. Also, I think this figure is never referred to in text.
- "... even without using the convexity technique; ...": which convexity technique are you talking about?
- "$a_i = \nabla_\theta L(d_i) \cdot v$", what is the $v$ here? Also, since the gradient is multidimensional, are you talking about dimension-wise variance?
- "addtional"

Lebeda, Christian Janos, et al. "Avoiding Pitfalls for Privacy Accounting of Subsampled Mechanisms under Composition." arXiv preprint arXiv:2405.20769 (2024).
limitations: limitations: I think authors should still address whether other accounting tools can solve this problem more efficiently that their method in the approximate DP domain. Other than that I believe the limitations are well addressed.
paper_summary: paper_summary: The paper first proves new privacy bound for the subsampled Gaussian mechanism under fixed-size sampling with and without replacement, which improves over the tightest known prior results in (Wang et al. 2019) by a constant of four. The proofs rely on a careful coupling of the sampling processes on neighboring datasets (similar to the design in Wang et al. 2019) to reduce the problem of analyzing divergence between Gaussian mixtures, for which the analyses in (Mironov et al. 2019) are applicable. 

- For fixed-size sampling without replacement, the paper proves integral upper bounds for the remainder term of the Taylor expansion of the privacy bound, to obtain tighter constants compared to the privacy bound in Wang et al. 2019. 
- The method extends to fixed-size sampling with replacement. It allows an upper bound that is similar to the upper bound for sampling without replacement in the leading term to sampling probability q. The authors also prove a lower bound under such settings and numerically investigate its dependence on the batch size and the Renyi divergence order.
- The authors analytically showed that the empirical gradient variance is larger under Poisson sampling compared to fixed-size sampling, but the privacy bound is smaller under Poisson sampling compared to fixed-size sampling, indicating an interesting privacy-utility trade-off depending on the choice of sampling schemes, when fixing the sampling probability the same.
strength_weakness: strength_weakness: - The reason for the improved constant factor compared to (Wang et al.) is not crystal clear. Is the integral upper bound for the remainder term of Taylor expansion of the DP bound contributing to the tighter constant? Unfortunately, the Taylor expansion is not explained in much detail in the main paper (e.g., in Theorem 3.2, the crucial terms related to $A$ and $M$ are not presented nor explained).
- Although it is interesting that the paper shows tighter constant in leading order term of the DP bound (via analytical approach), the value of this contribution in practice needs more explanation. As the divergence between Gaussian mixtures could be tightly computed numerically following (Mironov 2019.), it is unclear why we need an tighter analytical bound given by integral upper bound for the remainder term of Taylor expansion of the DP bound.
questions: questions: 1. The reason for improved tightness of DP bound by a constant factor (see weakness for more details)
2. It is interesting that the paper shows that Poisson sampling enables smaller RDP bound than fixed-size sampling when keeping the sampling ratio the same. However, I wonder if it is a phenomenon that is unique to the add-or-remove-one notion of differential privacy. Could the authors comment on whether there will be a similar gap between Poisson sampling and fixed-size sampling when considering the replace-one notion of differential privacy?
limitations: limitations: N/A
paper_summary: paper_summary: This paper analyzes Differentially Private-SGD with fixed batch size (with and without replacement), through the lens of Rényi-DP. The bounds without replacement have a much better constant than previous ones; the ones with replacement are brand new.
strength_weakness: strength_weakness: It’s a bit of a nitpick, but the paper keeps claiming that this analysis is DP-SGD specific, and that’s what enables the tight bounds. It seems like it is specific to:
- the Gaussian mechanism (with specific subsampling approaches)
- and I guess additive contributions of datapoints to the function?
Anyway, it might be good to give the specifics in a slightly more general/abstract fashion, and then discuss why it applies to well to SGD? Just so that the reader doesn’t expect tight couplings with the optimization or something, and also because it may be slightly more broadly useful.

Minor:
- several citations seem off: the “standard differential privacy (DP)” citation is CDP, poisson sampling from RDP cites the RDP paper instead of this one I assume https://arxiv.org/pdf/1908.10530.
- "even without using the convexity technique” p. 6 -> this has never been mentioned before. You need more details here for context.
questions: questions: Since you reduce to the same intermediary quantity as https://arxiv.org/pdf/1908.10530 (in Eq. 6), does it mean that:
- the bound of that paper applies as-is to fixed batches using your proof (i.e., if one were using the Opacus accountant with fixed batches as a heuristic, it’s actually not a heuristic)?
- if not, why is that?
- if so, how does your new bound compares? (other than enabling the analysis with replacement, which is cool!)
I think that the paper would benefit from discussing those to give a bit more context. S4.2 is related, as it gives an interesting discussion of where the different factor comes from (so the quantity is not identical I guess), but not all of this, like whether you could use the previous approximation as is, and how the two numerical procedures compare.

After reading S4.2, I am wondering if the difference is a bit artificial, or at least tied to the notion of sensitivity. For instance, Poisson sampling would also pay a factor two under the change one definition (I think?) whereas your approach would mostly not change I believe (I only skimmed the appendix on sensitivity to I’m not completely sure, but intuitively it’d make sense—the section doesn’t seem to discuss the implications on Poisson vs. fixed batch). Is it the case that your approach pays a 2x factor, but also gets the stronger change-one for free, whereas Poisson gets the better guarantee under add/remove, but its equivalent under change-one? It would be an interesting thing if that’s the case I think (change-one is pretty popular in other applications, I don’t think it’s really less conventional as the paper claims).
limitations: limitations: N/A
paper_summary: paper_summary: The paper studies the Rényi differential privacy (RDP) guarantees of the subsamplied Gaussian mechanism with fixed-sized random minibatches, when the so called add/remove neighborhood relation of datasets is considered. It uses similar techniqes as were used in the paper 

Mironov, Ilya, Kunal Talwar, and Li Zhang. "R\'enyi differential privacy of the sampled gaussian mechanism." arXiv preprint arXiv:1908.10530 (2019),

which augmented the seminar paper Deep Learning with DP (Abadi et al., 2016) and gave rigorous RDP guarantees for the Poisson subsampled DP-SGD.

The paper also studies both analytically and numerically the benefits of carrying out subsampling with fixed-sized minibatches instead of Poisson subsampling.
strength_weakness: strength_weakness: - Clearly the biggest deficit of the paper is that it does not take into account some of the recent research in this area. Most importantly, it overlooks the work

Zhu, Yuqing, Jinshuo Dong, and Yu-Xiang Wang. "Optimal accounting of differential privacy via characteristic function." International Conference on Artificial Intelligence and Statistics. PMLR, 2022. https://proceedings.mlr.press/v151/zhu22c/zhu22c.pdf

The main result of this paper is a special case of the Thm. 11 by Zhu et al. (2022) : If $(P,Q)$ are a dominating pair of distributions for the base mechanism in substitute relation, then the fixed-size random subsampling without replacement give a dominating pair $(q P + (1-q) Q, Q)$ for removal neighbours and $(P, (1-q) P + q Q)$ gives a dominating pair for add neighbors. Then, we know the dominating pair $(P,Q)$ of the Gaussian mechanism in case of substiture relation of datasets (pair of one-dimensional Gaussians) and the result of this paper's mmain result (Thm. 3.1) follows since if a pair if dominating pair for the hockey-stick divergence, by the Blackwell theorem it dominates also for other $f$-divergences. 

Since this "without replacement" result consitutes such a big part of the paper, I think this is a major deficit and there should be a major revision before accepting this paper.

This with replacement upper bound seems interesting, however my impression is that it is quite conservative (see e.g. Fig. 7 of the appendix). Only the experimental results for the lower bound in case of "with replacement" are given in the main text.
questions: questions: Can you comment on the accuracy of the upper bound in case of "with replacement"? Is it still open/unclear, how tight those bounds are?

Could you provide "with replacement" upper bounds also for the hockey-stick divergence, i.e., do you obtain a dominating pair of distributions for the hockey-stick divergence?
limitations: limitations: Yes.

Comments: comment: Thanks for the clarifications. I will keep the score -- I now understand better about the efficiency value of the work (compared to numerical computation for Gaussian mixtures). However, the reason for constant improvement is still not clear to me analytically, and I encourage the authors to add more discussions of the Taylor expansion terms in the main paper.
Thanks for the question. We want to clarify the notion of "worst-case" in this context that might have created some confusion: Our method gives "an" upper bound vs. the "tightest" bound.  Based on your comments and also those of uSpK, we have added additional details to Appendix D to make our derivation clearer.  We hope that the following outline clarifies it:

To show that the Gaussian mixtures provide an upper bound on the worst case, we use the probabilistic Lemma A1 that decomposes the transition probabilities of mechanism.  In the replace-one case the role of $D$ and $D'$ are symmetric and so they both satisfy a decomposition of the form (31).  Then using convexity, the Renyi divergence between Eq.31(with $D$) and Eq.31(with $D'$) is bounded above by maximum over $b^\prime,\widetilde{b}$ of the Renyi divergence between Gaussian mixtures, as in Eq.102.  The only property of the divergence used in these steps is quasiconvexity, and so they apply to the hockey-stick divergence as well, though our subsequent Taylor expansion calculations that lead to explicit computable bounds are specific to Renyi. We emphasize that we only compute an upper bound and not the tightest upper bound, so there is no conflict with Lebeda et al, though we still obtain tighter bounds than the previous state-of-the-art explicit computable bound in Wang et al. The examples in Lebeda et al. seem to highlight the difficulty of obtaining the tightest possible bound on the hockey-stick divergence, which is an interesting but different theoretical problem than what we achieved in this work. 

As to your other question, the factor of 4 is contained in our definition of $r_t$; see eq. (34).
Thanks for the question. For FS-woR under add/remove adjacency the worst case bound is stated in eq. (6) in Theorem 3.1.
For FS-woR under replace-one the worst case bound is stated in eq (102) in appendix D4.  Both are derived by the method in Appendix A, wherein the mechanisms are decomposed into mixtures via Lemma A1: see eq (31)-(32).  Then  quasiconvexity is used to bound the mixture over $b',\widetilde{b}$, leaving only the mixture over q, which is the probability that the last sample (the one that can differ) is included.  The essential difference that the adjacency relation makes is as follows:

- In the add/remove case the first argument of Renyi divergence is decomposed via the mixture (31) while the second via (32) (as D' has one less element).  After using quasiconvexity to bound the expectation over $b',\widetilde{b}$ by the maximum one obtains (6). where the Gaussian mixture is only in the first argument of the Renyi divergence.

- In the replace-one case, D and D' have the same number of elements and so both arguments of the Renyi divergence are decomposed as in (31) but with different means (coming from D and D'). After using quasiconvexity, this leaves a Gaussian mixture in both arguments of the Renyi divergence.
In both cases we then go on to obtain computable worst-case bounds using Taylor expansion with remainder.

We have added additional details to Appendix D to make this derivation clearer. Also, we agree with you that investigating the practical boundaries of numerical methods in this context is a fruitful but separate research endeavor.
Thanks for your response, and apologies for my very late response! I have now read the rebuttal, and authors have addressed many of my concerns. 

In your common response you state that "AFA is not specific to DP-SGD". This is of course true, since the particular accounting technique is not tailored for DP-SGD. However, my point was that the accounting of DP-SGD can be carried out with the results of Zhu et al. using the AFA (or any other accountant like the PRV accounting). However, I do acknowledge the possible numerical challenges these methods pose, whereas your method can give at least an upper bound for the privacy cost (when converted to approximate DP). Since the proposed bound is more tight than the one proposed by Wang et al. I believe the proposed method can provide useful insights on the privacy accounting. However, it would be very interesting to see, if the numerical accounting indeed fails for this privacy analysis, or how far it is from your bounds.

One more question related to your response **Worst-case bounds for woR-sampling**. I'm not sure where in Appendix D you actually show that the worst-case is what is presented. It seems to me that you just plug in the Gaussian mixtures as the upper bound after eq. 102.
Thank you for the reply. I confused your Section C.2 with your Section D.4. There, in Section D.4, it is diffiult for me to see why that RDP bounds should hold. Namely, you have a bound with a Rényi divergence between the mixtures $q \cdot \mathcal{N}( \Delta_1, \sigma^2) + (1-q) \cdot \mathcal{N}( 0, \sigma^2)$ and $q \cdot \mathcal{N}( \Delta_2, \sigma^2) + (1-q) \cdot \mathcal{N}( 0, \sigma^2)$. The derivation of this bound is omitted. Recently, [Lebeda et al.](https://arxiv.org/pdf/2405.20769) showed that the bound for the subsampling without replacement under the substiture relation does not at least have a hockey-stick upper bound with mixtures $q \cdot \mathcal{N}( -1, \sigma^2) + (1-q) \cdot \mathcal{N}( 0, \sigma^2)$ and $q \cdot \mathcal{N}( 1, \sigma^2) + (1-q) \cdot \mathcal{N}( 0, \sigma^2)$, one has to use bounds such as Prop. 30 by [Zhu et al.](https://proceedings.mlr.press/v151/zhu22c/zhu22c.pdf). Thus, I am a bit suprised to see this result given in Section D.4, and I think a detailed derivation would be needed. Also, I don't see where the bound of Eq. (108) comes from, why is there no additional factor 2 (or 4 after squaring). 

I do acknowledge the contribution in making the RDP bound calculations faster, but I think the paper would benefit some polishing and by piutting the results in a context.
We thank reviewer kPGW again for their time and useful insights.  As the discussion period comes to a close we ask that you consider our responses and whether your concerns have been adequately addressed.  We tried to pay particular attention to clarifying remarks about improvement of the constant factor and the value of the contribution.  Thanks again!
We thank reviewer uSpK again for their time and helpful feedback.  As the discussion period comes to a close we ask the reviewer to consider our responses and whether their concerns are adequately addressed.  In particular we have provided extensive discussion relating our work to the work of Zhu et al. in the general comments, as well as in response to reviewer Rv66.  Thanks again!
Thanks for the answer. It is not immediately clear to us that the calculation you propose can be carried out explicitly in that case, but either way it requires substantial analysis/computation which is not done in Zhu et al. As we noted, the result in Zhu et al, as presented, is not amenable to implementation in practice, while our results are. Specifically, a significant aspect of our proposed method is the use of Taylor series expansions to capture the leading order behavior (when $q$ is small) while providing explicit computable bounds on the error terms to prevent privacy leakage. This aspect of our work is not considered in Zhu et al or Mironov et al.
Thank you for the replies. In fact, I think the replace one result could also be derived directly from the results by [Zhu et al.](https://proceedings.mlr.press/v151/zhu22c/zhu22c.pdf). If you use their proposition 30, and combine it with the results of Section 3 of [Mironov et al.](https://arxiv.org/pdf/1908.10530) (i.e., onlly on of the cases needs to be considered), you get the RDP bound for that case as well.
Thanks for asking for clarification on this. We get to a similar, but not identical expression, to that of Mironov et al. The difference is in the right-hand side of our eq.6, which differs from Mironov by the factor of 1/4 in the variances. This is due to the effects of fixed-size vs. Poisson subsampling. The key difference in the derivation in our Appendix A is contained in our Lemma A.1, as compered to Theorem 4 in Mironov et al. There we construct the random variables that lead to the appropriate decomposition of the mechanism. In particular, compare the construction of our $J$, $B'$, and $\widetilde{B}$ in the paragraph about Lemma A.1 with the random variable T defined at the start of the proof of Theorem 4 in Mironov, along with (implicitly) the random variable that selects whether the additional element x is included, leading to the mixture at the bottom of page 3 of Mironov. The construction of our $J$, $B'$, and $\widetilde{B}$ are the analogues that lead to the appropriate decomposition of the mechanism in the case of fixed-size subsampling, as proven in our Lemma A.1. This difference then propagates through the derivation, eventually leading to the difference by a factor of 1/4 in the variance as noted above. The key reason for this is that our eq 34-35 shows that the mechanism in the case of fixed-size subsampling has a sensitivity of twice that in the Poisson case. This is because in Poisson subsampling when the minibatch from $D$ differs from that of $D^\prime$ it is due to the inclusion of a single additional element. However, in fixed-size subsampling, when the minibatches  are not identical then they differ by a replacement; this  contributes more to the difference in means  by a factor of $2$. Furthermore, in our analysis of fixed-size subsampling with replace-one adjacency in Appendix D.4 there are even greater differences between our analysis and that of Mironov et al.

We added the text to highlight this difference with Mironov et al. and directed readers to our Appendix A and Lemma A.1 in this regard. Also, in our Appendix A, we added the intuition described above.
Thanks for asking for clarification on this. We agree with you that the steps you outline could provide an alternative path to the add/remove result in our Theorem 3.1. What we were highlighting in our last response is that our method of proof leads to a unified method for treating both add/remove and replace-one. The replace-one result is also of significant application value as was pointed out by JRin. Originally, we had the add/remove result in Theorem 3.1 in the main text while discussing replace-one in Appendix D.4. Based on your and JRin's comment, we have changed the main text to highlight the replace-one as well. We have also pointed out that the add/remove result could alternatively be obtained via the method you outline.
Thank you. My impression looking at how the theory works in your paper was that it did reduce to the same accounting at Mironov et al. 2019 (Since you reduce to the same intermediary quantity as https://arxiv.org/pdf/1908.10530 (in Eq. 6)). Could you give an intuition / place of divergence in your proof for why that is not the case? Thanks!
Thanks for the answers. But why would you need a result for replace-one adjacency? You define the neighborhood relation to be the add/remove neighborhood relation. On lines 80-81:

> We define datasets $D$ and $D′$ to be adjacent if one can be obtained from the other by adding or removing a single element.

Also, the proof of Thm 3.1 clearly shows that $D$ and $D'$ are of different sizes.

Then, I think [Zhu et al., Thm. 11b](https://proceedings.mlr.press/v151/zhu22c/zhu22c.pdf) will give your result. On the right hand side of the bound of your Thm 3.1, you exactly have the dominating pair for the Gaussian mechanism under the replacement relation, just as required by [Zhu et al., Thm. 11b](https://proceedings.mlr.press/v151/zhu22c/zhu22c.pdf). Then you can use the results of Section 3 of [Mironov et al.](https://arxiv.org/pdf/1908.10530) to show that one of the Rényi divergences dominates the other which exactly gives your Thm. 3.1 (you also seem to use those results by Mironov et al.).
Thanks for further clarification. We investigated the extension of Zhu et al.'s Thm 11 using Blackwell Thm as you suggested to validate if one can obtain our bounds in that manner.

In the case of replace-one adjacency relations, we do not believe that the general Blackwell theorem ideas along with the results in Zhu et al are sufficient to reproduce our results. Theorem 11(b) in Zhu et al applies to add or remove adjacency for the subsampled mechanism. Proposition 30 of Zhu et al does address replace-one adjacency, but it does not derive a full dominating pair, as it must break the analysis into two cases based on alpha>=1 and alpha<1. Corollary 32 in Zhu et al does give general theory for constructing a dominating pair, but it is not explicit or readily implementable in practice due to the need to compute a Legendre of transform (Fenchel conjugate) as a function of $x$ for `every' $x\in[0,1]$.  Therefore to the best of our knowledge, the result in Zhu et al. cannot be directly used to reproduce our new RDP bound for replace-one adjacency as found in our equation 110.
Thank you for the replies! Thm. 11 by Zhu et al. 11 is independent of any accountant used to compute the privacy profiles. As I have pointed out in my review, from that theorem follows the Blackwell dominance as well and an existence of a post-processing function and furthermore the RDP bound you state in your paper. This unfortunately reduces the relevance of that particular result I think, and therefore I plan to keep my score.