Banded Square Root Matrix Factorization for Differentially Private Model Training
Authors: Nikita Kalinin
Christoph H. Lampert
Venue: NeurIPS 2024 poster (2024)
Keywords: Differential Privacy
Machine Learning
Federated Learning
Primary area: privacy
TL;DR: A more efficient technique for matrix-factorization-based differentially private model training.

Abstract: Current state-of-the-art methods for differentially private model training are based on matrix factorization techniques. However, these methods suffer from high computational overhead because they require numerically solving a demanding optimization problem to determine an approximately optimal factorization prior to the actual model training. In this work, we present a new matrix factorization approach, BSR, which overcomes this computational bottleneck. By exploiting properties of the standard matrix square root, BSR allows to efficiently handle also large-scale problems. For the key scenario of stochastic gradient descent with momentum and weight decay, we even derive analytical expressions for BSR that render the computational overhead negligible. We prove bounds on the approximation quality that hold both in the centralized and in the federated learning setting. Our numerical experiments demonstrate that models trained using BSR perform on par with the best existing methods, while completely avoiding their computational overhead.

Banded Square Root Matrix Factorization for Differentially Private Model Training
Nikita Kalinin
Christoph Lampert
Institute of Science and Technology (ISTA)
Institute of Science and Technology (ISTA) Klosterneuburg, Austria ch1 mil@ist... ac.at
Klosterneuburg, Austria nikita. kaltta.kaliin inin@ist.....mm ac.at
Abstract
Current state-of-the-ar methods for differentially private model training are based matrix factorization techniques. However, these methods suffer from high computational overhead because they require numerically solving demanding optimization problem determine approximately optimal factorization prior to the actual model training. In this work, we present a new matrix factorization approach, BSR, which overcomes this computational bottleneck. By exploiting properties of the standard matrix square root, BSR allows to efficiently handle also large-scale problems. For the key scenario of stochastic gradient descent with momentum and weight decay, we even derive analytical expressions for BSR that render the computational overhead negligible. We prove bounds on the approximation quality that hold both in the centralized and in the federated learning setting. Our numerical experiments demonstrate that models trained using BSR perform on par with the best existing methods, while completely avoiding their computational overhead.
1 Introduction
We study the problem of differentially private (DP) model training with stochastic gradient descent (SGD) in the setting of either federated or centralized learning. This task has recently emerged as one of the most promising ways to train powerful machine learning models but nevertheless guarantee the privacy of the used data, which led to number of studies, both theoretical as well as application-driven [Abadi et al. 2016 Yu et 2020 Zhang al., 2021 Kairouz al. 2021
Denisov et al. 2022 The state of the art in the field are approaches based on the mairix factorization
(MF) mechanism [i al. 2015 Henzinger et al. 2024
which combines theoretical guarantees 2024 It is based on the observation
with practical applicability (Choquette- Choo et al. 2023a that all iterates of SGD are simply linear combinations of model gradients, which are computed at intermediate time steps. Consequently. the iterates can be written formally as the result of multiplying the matrix coefficients, called workload matrix, with the row-stacked gradient vectors. To preserve the privacy of the training data in this process one adds suitably scaled Gaussian noise at intermediate steps of the computation. The MF mechanism provides a way to select the noise covariance structure based on factorization of the workload matrix into two matrices.
Identifying the minimal amount of noise necessary achieve desired privacy level requires solving an optimization problem over all possible factorizations, subject to data participation constraints. For some specific settings, the optimal solutions have been characterized: for streaming learning, when each data batch contributes at most once to the gradients, [i et al. [2015] presented formulation of
this problem as a semi-definite program. Henzinger al. 20044 proved that square root factorization
Note that this specific type MF should be confused with other occurrences matrix factorization in,
potentially private, machine learning, such as recommender systems [Shin et al. 2018 ct al. 2021
38th Conference on Neural Information Processing Systems (NcurIPS 2024).

---

of the workload matrix is asymptotically optimal for different linear workloads, including continual summation and decaying sums.
In this work, our focus lies on the settings that are most relevant to machine learning tasks: workload matrices that reflect SGD-like optimization, and participation schemes in which each data batch can potentially contribute to more than one gradient vector, as it is the case for standard multi-epoch training. Assuming that this happens at most once every b steps, for some value > leads the problem of optimal matrix factorization in the context of b-min-separated participation sensitivity. Unfortunately, as shown in Choquette- Chhoo et al. [023a] finding the optimal matrix factorization in this setting is computationally intractable. Instead, the authors proposed an approximate solution by posing additional constraints on the solution set. The result semi-definie program that is tractable, but still has high computational cost, making it practical only for small to medium-sizen problem settings.
Subsequent work concentrated on improving or better understanding the factorizations for specific algorithms, such as plain SGD without gradient clipping, momentum, or weight decay [Koloskoval et all [203]] or specific, e.g. convex, objective functions [Choquette-Choo et al. 2024]. Often, streaming data was asssumed, i.e. each data item can contribute at most to one model update, which is easier to analyze theoretical, but further removed from real-world applications (Dvijotham et al. 2024]. A concurrent line of works has also focused on scaling matrix factorizations for large-scale training McMahan et al. [2024) extended the Buffered Linear Toeplitz (BLT) mechanism to support participaatiiiiioonnnn in federated learning, improving privacy-utility tradeoffs, while ensuring memory efficiency. [McKenna (2024) improved DP Banded Matrix Factorization's 's scalability, enabling it to handle millions iterations and large models efficiently. These advancements enhance the practicality of differentially private training in real-world applications.
Our work aims general-purposs. solution that covers as many realistic scenarios as possible. Our ultimate goal is to make general-purposs differentially private model training as simple and efficient use as currently dominating non-private technique. Our main contributions are:
1. We introduce a new factorization, the banded squared root (BSR), which is efficiently computable even for large workload matrices and agnostic to the underlying training objective. For matrices stemming from SGD optimization potentially with momentum and/or weight decay, even provide closed form expressions.
2. We provide general lower bounds on the approximate error for any factorization, along with specific upper and lower bounds for the BSR, Square Root, and baseline factorizations, in the contexts of both single participation (streaming) and repeated participation (e.g..
multi-epoch) training.
3. We demonstrate experimentally that BSR's 's approximation error is comparable to the state-of-the-art method, and that both methods also perform comparably in real-world training tasks.
Overall, the proposed BSR factorization achieves training high-accuracy models with provable privacy guarantees while staying computationally efficient even for large-seale training tasks.
Background
Our work falls into the areas of differentially private (stochastic) optimization, of which we remind the reader here, following mostly the description of Denisov et al. [2022] and CnohCeeettte-hhee et al. [2023a]. The goal is to estimate sequence of parameter vectors, θ (61...., Hn) E R", where each 0g is a linear combination of update vectors, Ii,..., dt E Rdd that were computed in previous steps, typically as gradients of model with respect to training data that is meant to stay private. We assume that all x, have bounded norm, || S Ç. Compactly, we write θ AX, where the lower triangular workload matrix A contains the coefficients and X E Raxd is formed by stacking the update vectors as rows. With different choices of A, the setting then reflects many popular first-order optimization algorithms, in particular stochastic gradient descent (SGD), potentially with momentum and/or weight decay. Depending on how exactly Zi...., In are obtained, the setting can express different centralized as well as federated training paradigms. To formalize this aspect, we adopt the concept of amin-separated participation [Choquette- Chho et al. 2023a]. For some integer 2 1, it states that if data item (e.g. a single training example in central training, client batch in
2

---

Algorithm Differentially Private SGD with Matrix Factorization
Input: Initial model fg E Rd, dataset D, batchsize b, matrix C' E R R"x", model loss ((0,, d), clipnorm (, noise matrix IRaxd with i.i.d. entries NO0 s), where = a senst for i = 1,2,...., do
S - fdi,...,de dm}
C
D select a data batch, respecting the data participation constraints
9i V,l(0,-1,, dg))
j= = I...., where clip,(d) = min(li<//ld[)d
I; + -j=1 clip,(g)) (g;)
á + x; + S(C- (6.l
0: + update 10 R;),
// SGD model updates
Output: θ (01.... ........ ,p. 0n)
federated learning) contributed to update T;, the carliest it can contribute again is the update 24+b- Additionally, let & I g be the maximal number any data point can contribute. In particular, this notion also allows to treat in unified way streaming data (b = n or k = 1), well as unrestricted access patterns (k- = n with = 1). but also intermediate settings, such as multi-epoch training on a fixed-size dataset.
The matrix factorization approach [L et al. 2015) adopts factorization = BC of the workload matrix and computes O BCC + Z), where Z is Gaussian noise that is chosen appropriately to make the intermediate result CX + Z private to the desired level. Algorithm [ shows the resulting algorithm in pseudocode. It exploits the fact that instead of explicit multiplication by C' and B, standard optimization toolboxes can be employed with suitably modified update vectors, because also MF A(X + C-' Z), and multiplication by A corresponds to performing the optimization.
Different factorizations recover different algorithms from the literature. For example, B A, Id recovers DP-SGD [Abadi et al. [016], where noise is added directly to the gradients. Conversely, = Id, C A simply adds noise to each iterate of the optimization [Dwork et al.] [0066]. However, better choices than these baselines are possible, in the sense that they can guarantee the same levels of privacy with less added noise, and therefore potentially with higher retained accuracy. The reason lies in the fact that B and C' play different roles: B acts as a post-processing operation of already private data. Hence, has no further effect on privacy, but it influences to what amount the added noise affects the expected error in the approximation of θ. Specifically, for Z (O; Id),
Ez - Ez BZ s² |B]
(1)
In contrast, CX is the quantity that is meant to be made private. Doing so requires noise of a strength proportional to C's sensitivity, sCns(C) == supx-x' |CX CX'I|F,, where the neighborhood relation, X ~ X', indicates that the two sequences of update vectors differ only in those entries that correspond to single data itemr] As shown Choquette- et al. [2023a], in the setting of b-min-separated repeated participation, it holds that
senss,,(())
(C) S max #Elle,6
Σ E(CTraaa J.j
where 11 (T C (1,...., n)
C ii
b) J. is the set of
possible b-min-separated index sets with at most participation. Furthermore, 2 holds even with equality if all entries of CC are non-negative.
Combining I with = sensk, (C)) yields a quantitative measure for the quality of a factorization.
Definition 1. For any factorization A BC, its expected approximation error
1 Vn.
E(B == V Ez e - 1MF
where the 1//n factor is meant make the quantity comparable across different problem sizes.
As proved in Denisov et al. 2022 Theorem 2.1]. establishing privacy in this non-adaptive setting suffices to guarantee also privacy in the adaptive setting, where the update vectors depend not only on the data items but also the intermediate estimates of the model parameters, is the case for private model training.

---

The optimal factorization by this reasoning would be the one of smallest expected approximation error. Unfortunately, minimizing (3) across all factorizations it is generally computationally intractable.
Instead, Choquette-Chov et al. 2023a propose an approximately optimal factorization.
Definition 2. For a workload matrix A, let S be the solution to the optimization problem
arg min trace[AAASS SE84
subject to diag(S) = and Ska)== =0 ffor l - jl 2 b,
(4)
where S4 is the cone of positive definite n x matrices. Then, A = BC is called the approximately optimal factorization (AOF), C is lower triangular and fulfills C-C=5
The optimization problem is semi-definite program (SDP), and can therefore be solved numeri- cally using standard packages. However, this is computationally costly, and for large problems (e.g. 5000) computing the AOF solution impractical. This poses a problem for real-world training tasks, where the number of update steps are commonly thousands or tens of thousands.
mitigate this problem to some extent, but as we will discuss in Section this can lead to robustness problems, especially because the recovery of C from S in Definition e.g. by Cholesky decomposition, tends be sensitive numerical errors.
Solving (4 itself only approximately
3 Banded Square Root Factorization
In the following section, we introduce our main contribution: a task of differentially private stochastic optimization that can be computed efficiently even for large problem sizes.
factorization for the
Definition 3 (Banded Square Root Factorization). Let A E Ruxn be a lower triangular workload matrix with strictly positive diagonal entries. Then, we call = CC the square root factorization (SR), when C denotes the unique matrix square root that also has strictly positive diagonal entries. Furthermore, for any bandwidth E [I,...,n)) we define the banded square root factorizaaiinn of bandwidth p (p-BSR) of A. as
A = Blll clrl
(5)
where Clp is created from C by setting all entries below the p-th diagonal O,
Cjal
i
C.Pl
and Birl_ A(Cirl)-
(6)
0
otherwise.
Note that determining the SR, and therefore any p-BSR, is generally efficient even for large workload matrices, because explicit recursive expressions exist for computing the square root of a lower triangular matrix [Böörck and Hammarling 1983 Deadman et al. 2012).
In the rest of this work, we focus on the case where the workload matrix stems from SGD with momentum and/or weight decay, and we show that then even closed form expressions for the entries of Clpl exist that renders the computational cost negligible.
3.1 Banded Square Root Factorization for SGD with Momentum and Weight Decay
We recall the update steps of SGD with momentum and weight decay:
A: at;-1 171; for m, = 1m,mmm1- 24
(7)
where 3i,..., Fn are the update vectors, v > is the learning and 0 < B 1 is the momentum strength and 0 < a < 1 is the weight decay parameter. Note that our results also hold for 3 = 0, i.e. without momentum, and for a 1, i.c., without weight decay. In line with real algorithms and to avoid degenerate cases, we assume 3 < a throughout this work. The update vectors are typically gradients of the model with respect to its parameters, but additional operations such normalization or clipping might have been applied.
Unrolling the recursion, we obtain an expression for 0 as linear combination of update vectors
0: =
33
(8)
j=1
k=j

---

Consequently, the workload matrix has the explicit form A = 7A0, for
ag
0
a1
4g
0
a Bd+1
Aa,3 =
02
01
dg
0
with aj = =aa
(9)
a - B
:
1=0
On-1
0n
00
As one can see, Aa.3 is lower triangular Toeplitz-matrix, so it is completely determined by the entries of its first column. In the following, we use the notation I.DToep(m1.)) mmeem....... m.)) to denote lower triangular Toeplitz matrix with first column mi,...., Men, i.e. Aa, = LDToep(eo,.....
Our first result is an explicit expression for the positive square root of A..3 (and thereby its p-BSR). Theorem (Square-Rood Root SGD Workload Matrix), Let An.3 be the workloud matrix 9 Then Aa.3 C. for Ca.3 = LDTooep(c.....,CC Cn-1), with = and Cj = Ei=o a -- for n,], the p-banded BSR matrix clel is obtained from this by setting all coefficients cy 0. 2
j =1, 1,...,n- with coefficients r; = [(- For any [1,...,m]]
Proof sketch. The proof be found in Appendix F.1 Its main idea is to factorize A.,s into a product of two simpler lower triangular matrices, each of which has closed-form square root. We show that the roots commute and that the matrix Ca. is their product, which implies the theorem. O
3.2 Efficiency
We first establish that the p p-BSR for SGD can be computed efficiently even for large problem sizes.
Lemma (Efficiency of BSR). The entries of cli), can be determined in runtime O(p log p), i.e., in a.
particular independent of
Proof sketch. As a lower triangular Toeplitz matrix, clpl is fully determined by the values of its first column. By construction cptl....,ee G2 0, so only the complexity of computing C1, Cp-1 mat- ters. These can be computed efficiently by writing them as the convolution of vectors (a r s=........ and (3"r;) (B0..)) and, e.g.. employing the fast Fourier transform.
Note that for running Algorithm the matrix B of the factorization is not actually required. However, needs to know the sensitivity of c.Pl as this determines the necessary amount of noise. The following theorem establishes that for large class of matrices, including the BSR in the SGD setting, this is possible exactly and in closed form.
Theorem 2 (Sensitivity for decreasing non-negative Toeplitz matrices). Let M I.Doopp(mo..... mm-1) ) be a lower triangular Toeplitz. matrix with decreasing non-negative entries, ie.
mo 2 mi m2 2 ... ...m-11 0. Then its sensitivity ② in the setting of b-min-seppration is
1/2
2
k-1
n-1
miink-1.i/45
-1,1bb Σ
sensk,((M)
Σ
Mf. Mr,1
mi-jb
(10)
j=0
j=
where M denotes the (1 jb)-th column of M.
Proof sketch. The proof can be found in Appendix F2 It builds the identity 2 which holds with equality because of the non-negative entries of M. Using the fact that the entries of M are non-increasing one establishes that optimal b-separated index set is {1, 1 b.... (k 1)6). From this, the identity (T0) follows.
Corollary 1. The sensitivity of the p-BSR for SGD can be computed using formula 1O
Proof sketch. It suffices show that the coefficients Cp- of Theorem are monotonically decreasing, We do by explicit computation, see Appendix R3
5

---

3.3 Approximation Quality Single Participation
Having established the efficiency of BSR, we now demonstrate its suitability for nigh-quality model training. To avoid corner cases, we assume that is an integer, which does not affect the asymptotic behavior. We also discuss only the case in which the update vectors have bounded norm = 1. Results for general can readily be derived using the linearity of the sensitivity with respect to Ç.
We first discuss the case of model training with single participation (k 1); where more precise results are possible than the general case. Our main result are bounds on the expected approximation error of the square root factorization that, in particular, prove its asymptotic optimality.
Theorem 3 (Expected approximation error with single participation). Let Aa.ß E [Rax be the workload matrix (9) of SGD with momentum 0 S B and weight decay parameter < a S 1, where a > B. Assume that each data item can contribute at most once to update vector (e.g- single participation, k = 1). Then, the expected approximation error of the square root factorization, Aa,3 C o. fulfills
1 log
S E(Ca,B Ca.p) S (a- B)2
(11)
1- a²
for and
log(n + +1 +1)-1} < &E(C,,, C1,0)
1 log(n) (1 8)
max
(12)
4
Proof sketch. For the proof, we establish relations between sens (C) and || flle, and then we bound the resulting expressions by an explicit analysis of the norm. For details, see Appendix F.5 n
The following two results provide context for the interpretation of Theorem Theorem 4. Assume the setting of Theorem 3
Then, for any factorization Aa.3 BC with CT C 2 O, the expected approximation error fulfills
(111)
for < 1, for ===a 1.
E(B,C) = [2(logn)
(13)
Proof sketch. The theorem the special case k of Theorem which we state in the next section and prove in Section F9
Theorem 5. Assume the setting of Theorem Then, the baseline factorizations A.,3 = Aa. Id and Aa,ß ld Aa, fulfill, for < 1,
v1 + aB
VT
E(Aa,g, Td)
o(1) and E(A1,a; V2(1 - B) +o(vn)
v (1 - aB) a²) 32
(14)
a
Vn
E(Id, (a. = v(1 a11-a a²)
+ o(1 and E(dd, < 1- + +olvn)
82)
(15)
Proof sketch. The result follows from an explicit analysis of the coefficients, see Appendix
F.6
Discussion. Theorems to provide a full characterization of the approximation quality of the square root factorization as well as its alternatives: 1) the square root factorization has asymptotically optimal approximation quality, because the upper bounds in Equation (12) match the lower bounds in Equation (3): 2) the AOF from Definition 2 also fulfills the conditions of Theorem 4 Therefore, it
must also adhere the lower bound 333 and cannot be asymptotically better than the square root
factorization: 3) the approximation qualities of the baseline factorizations in Equation (14) and (15)
are asymptotically worse than optimal in the a = setting, and worse by a constant factor for a The BSR factorization can be applied even in more general scenarios, such as with varying learning rates. However, in this case, the workload matrix will no longer be Toeplitz. This makes difficult to provide analytical guarantees for the matrix, but it can still be applied numerically.

---

3,4 Approximation Quality Repeated Participation.
We now provide mostly asymptotic statements about the approximation quality BSR and baselines in the setting where data items can contribute more than once to the update vectors.
Theorem 6 (Approximation error of BSR). Let Aa,3 Rux= be the workload matrix 9 of SGD
with momentum 0 < 3 < and weight decay 0 < S 1, with a > B. Let A,,3 Ba, P
0.3 be
its banded square root factorization as in Definition Then, for any b E nI. S b and E (1....: #) holds:
Og (v nk logp
+ +0ap(Vh)
for a 1,
E(Blr) Clel a,s a.3)
P
(16)
O OR.ni (Vh)
for a <1
Proof sketch. For the proof, separately bound the sensitivity of C.P. and the Frobenius of BlP
a- The former straightforward because of the matrix's band structure. The latter requires an in-depth analysis of the inverse matrix' coefficient. Both steps are detailed in Appendix F.
The following results provide context for the interpretation of Theorem
Theorem 7 (Approximation error of Square Root Factorization), Let Aa.3 Rxxx be the workload matrix (9) of SGD with momentum 0 < B < 1 and weight decay 0 < a S 1, with > B. Let be its square root factorization. Then, for any E {I...., ,n] and k = 3 holds:
An.3 C2 o,p
Og (lvlogn + vlllog))
for a = 1,
E(Ca,gs Ca.B)
(17)
+ O (VA)
for 1.
Proof sketch. We bound sensæ & (C and Co.,3 || using the explicit entries for Ca.g from Theo- rem Details are provided in Appendix F.8
Theorem 8. Assume the setting of Theorem CTC C 0, the approximation error fulfills
Then, for any factorization Aa,3 BC with
log
for =
8(B,))
(18)
Vk
for a < 1,
Proof sketch. The proof is based on the observation that || EX|rIYI 2 ||X for any matrices X,Y, where || ||- denotes the nuclear norm. To derive (18). we show that sensk,u(C) is lower bounded by //k ||C'|] and derive explicit bounds the singular values of Aa,3.
Theorem 9. Assume the setting of Theorem Then, the baseline factorizations Aa, Aa, Id and Aa, Id Aa,8 fulfill
kVn
nk
V3 for a =
for a =
E( E(An,,,
>
2
E(Id, An,a) 2
(19)
V/
for 1
Vk for a 1.
Proof sketch. The proof relies on the fact that the workload matrices can be lower bounded compo- nentwise by simpler matrices: Aa,9 A..o and A.,0 2 Id. For the simpler matrices, the bounds 19 can then be derived analytically, and the general case follows by monotonicity.
7

---

Discussion. Analogously to the case of single participation, Theorems 6 to 9 again establish that the proposed BSR is asymptotically superior to the baseline factorizations if 1. comparison of Theorems 6 and suggests that, at least for maximal participation, h ? and p = b, the bandedness of the p-BSR improves the approximation quality, specifically in the practically relevant regime
where b < n. While none of the methods match the lower bound of Theorem [6] we conjecture that this is not because any asymptotically better methods would exist, but rather a sign of Equation (18) is not tight. Both theoretical consideration and experiments suggest that a term linear in & should appear there. For a < 1, all studied methods are asymptotically identical and, in fact, optimal.
4 Experiments
To demonstrate that BSR can achieve high accuracy not only in theory but also in practice, we compare it to AOF and baselines in numerical experiments. Our results show that BSR achieves quality comparable the AOF, but without the computational overhead, and it clearly outperforms the baseline factorizations. The privacy guarantees are identical for all methods, so we do not discuss them explicitly.
Implementation and computational cost. of Theorem ( For single data participation, repeated data participation we use p-BSR with p= b. Using standard python/numpy code, computing the BSR as dense matrices are memory-bound rather than compute-bound Even sizes of = 10, 000 or more take at most a few seconds. Computing only the Toeplitz coefficients is even faster, of course.
We implement BSR by the closed-form expressions
the square root decomposition directly. For
To compute AOF, solve the optimization problem using the cvxpy package with SCS backend, see Algorithm B for the source code") With the default numerical tolerance, 10 -4- each factorization took a few minutes (n I 100) to hours (n 500) to several days (n 2 700) of CPU time. Note that this overhead reappears for any change in the number of update steps, n, weight decay, o, or momentum, 33 as these induce different workload matrices. In our experiments, when the optimization for AOF did not terminate within 10 days, we reran the optimization problem with the tolerance increased by factor of 10. The runtime depends not only on the matrix size but also on the entries. In particular, we observe matrices with momentum to be harder to factorize than without. For large matrix sizes we frequently encountered numerical problems: the intermediate matrices, S. in 4 often did not fulfill the positive definiteness condition required to solve the subsequent Cholesky decomposition for C. Unfortunately, simply projecting the intermediates back the cone of positive semi-definins matrices is not enough, because the resulting C' matrices also have to be invertible and not too badly conditioned. Ultimately, we adopted a postprocessing step for S that ensures that all its eigenvalues were at least of value /I/n,, which we find to be a reasonable modification to ensure the stability of the convergence. Enforcing this empirically found value leads to generally good results, as our experiments below show, but it does add an undesirable extra level of complexity to the process. In contrast, due to its analytic expressions, BSR does not suffer from numerical problems. also does not possess additional hyperparameters, such as numeric tolerance or the number of optimization steps.
Apart from the factorization itself, the computational cost of BSR and AOF are nearly identical. Both methods produce (banded) lower triangular matrices, so computing the inverse matrices or solving linear systems can be done within milliseconds to seconds using forward substitution. Note that, in principle, one could even exploit the Toeplitz structure of p-BSR, but we found this not to yield any practical benefit in our experiments. Computing the sensitivity trivial for p- p-BSR using Corollary and it is still efficient for AOF by the dynamic program proposed in Choquette-Choo al. [203aa]
Expected Approximation Error. As a first numeric experiment, approximation error for workload matrices that reflect different SGD settings. Specifically, we use workload matrices (9) for n e {100, 200,... .... 1000, 20000, with = {0.99, 0.999, 0.9999, 1}. and / E {0, 0.9}. either with single participation, k = 1, or repeated participation, b = 100, k n/ 100. Figure [ shows the expected approximate error, EBB C), of the proposed BSR, AOF, as well as the baseline factorizations, A = A Id and A = Id .A in two exemplary cases. Additional results for other privacy levels can be found in Appendix
evaluate the expected
Additional experiments with gradient-basse optimizers can be found in Appendix E
8

---

a=
0.999
[ 0.0
b 100
1.0
8 009.9
b = 100
approximation quality. 3
BSR
Sanooiimaaa quality.
BSR
80
ADF
AOF
60
1000
*6 id
40
500
20
600 800 1000
600 800 1000 1500
workload matrix size,
workload matrix size,
Figure 1: Expected approximation error of BSR, AOF and baseline factorizations for two different hyperparameter settings (left: 0.999, B = O, right: a = 8 = 0.9) with repeated participation 100, k n/100). See Section for details.
CIFAR-10 (single participation)
CIFAR- 10 (repeated participation)
accuracy
30
accuracy
40
AOF
30
20
25 200
500
1000
2000
10
15
20
number update steps
maximum number epochs
Classification accuracy (mean and standard deviation over runs with different random seeds) on CIFAR-10 10 for BSR, AOF, and baselines for (c, 6) = (4, 10-5) for independent training runs. Left: one epoch, different batch sizes. Right: different number of epochs, constant batch size.
Figure 2:
The results confirm expectations from the theoretical analysis: in particular, BSR's 's expected approximation error is quite close to AOF's, typically within a few percent (left plot). Both methods clearly superior to the naive factorizations. For large matrix sizes, BSR sometimes even yields slightly better values than AOF (right plot). However, we believe this to be a numeric artifact of us having to solve AOF with less-than-perfee precision.
Private Model Training on CTFAR R-10. To demonstrate the usefulness of BSR in practical settings, we follow the setup of Kairouz et al. [022]] and report results for training simple ConvNet on the CIFAR-10 dataset (see Table in Appendix Cfor the architecture). Specifically, we adapt Google's reference implementation of DP-SGD in jax Bradbury et al. [2018) to work with the different matrix factorizations: BSR, AOF, and the two baselines. To reflect the setting of single- paagle-paticippaioon training, split the 50,000 training examples into batches of size E {1000, 500, 250, 200, 100, 50, 25), resulting in n e { (100, 200, 400, 500, 1000, 2000} update steps. For repeated participation, we fix the batch size to 500 and run k: (1,2,.... 10, 15, 20) epoch of training, i.e. 100% and b = 100. In both cases, 20% of the training examples are used as validation sets to determine the learning rate n E (0.01, 0.05. 0.1, 0.5. 1). weight decay parameters a E (0.99, 0.999, 0.9999, 1). and momentum E (0. 0.9). Figure shows the test set accuracy of the model trained with hyperparameters that achieved the highest validation accuracy One can see the expected effect that in DP model training, more update steps/epochs do not necessary lead to higher accuracy due to the need to add more noise. The quality of models trained with BSR mostly identical to AOF, When training for large number of epochs it achieves even better slightly results, but this could also be an artifact of us having to solve AOF with reduced precision in this regime. Both methods are clearly superior the basclines.
*Such a setting would not optimal for world private training. because the many repeated experiments reduce the privacy guarantees [Papernot and Steinke 202 Kurakin et al. 2022 Ponomareva et al, 2023). We nevertheless adopt here allow for simpler and fair comparison between methods.

---

Conclusion and Discussion
We introduce an efficient and effective approach to the matrix factorization mechanism for SGD- based model training with differential privacy. The proposed banded square root factorization (BSR) factorization achieves results on par with the previous state-of-the-art, and clearly superior to baseline methods. At the same time, it does not suffer from the previous method's 's computational overhead, thereby making differentially private model training practical even for large scale problems.
Despite the promising results, some open questions remain. On the theoretical side, the asymptotic optimality of BSR without weight decay still unresolved because the current upper bounds on the expected approximation error do not match the provided lower bounds. Based on the experimental results, we believe this discrepancy lies with the lower bounds, which we suspect should be linear in the number of participations. We observe that BSR achieves results comparable to AOF, although we cannot currently prove this due to the insufficient understanding of AOF's theoretical properties; nonetheless, we consider it promising research direction. On the practical side, it would be interesting to extend the guarantees to even more learning scenarios, such as variable learning rates.
Acknowledgments and Disclosure of Funding
This research was supported by the Scientific Service Units (SSU) of ISTA through resources provided by Scientific Computing (SciComp). We thank Monika Henzinger and Jalaj Upadhyay for their valuable comments on the earlier versions of this manuscript.
References
M. Abadi, A. Chu, I. Goodfellow, H. B. McMahan, I. Mironov, K. Talwar, and L. Zhang. Deep learning with differential privacy. In ACM Special Interest Group on Security, Audit and Control (SIGSAC), 2016.
N. Batir, H. Kücük, and S. Sorgun. Convolution identities involving the central binomial coefficients and Catalan numbers. Transactions Combinatorics, 2021.
A. Björck and S. Hammarling. A Schur method for the square root a matrix. Linear Algebra and its Applications, 1983.
J. Bradbury, R. Frostig. P. Hawkins, M. J. Johnson, C. Leary, D. Maclaurin, G. Necula, A. Paszke, J. VanderPlas, S. Wanderman-Milne, and Q Zhang. JAX: composable transformations of Python+NumPy programs, 2018. URL http: http::///ttbbb www.google////oooooooooooooooolo
C. A. Choquette-Choo, A. Ganesh, R. McKenna, H. B. McMahan, J. K. Rush, A. G. Thakurta, and X. Zheng. (Amplified) banded matrix factorization: unified approach to private training. In Conference on Neural Information Processing Systems (NeurIPS), 2023a.
C. A. Choquette-Choo, A. Ganesh, T. Steinke, and A. Thakurta. Privacy amplification for matrix mechanisms. In International Conference on Learning Representations (ICLR), 2023b.
C. A. Choquette-Choo, H. B. McMahan, K. Rush, and A. Thakurta. Multi-epoch matrix factorization mechanisms for private machine learning. International Conference on Machine Learing (ICML) 2023c.
C. A. Choquette- Choo, K. Dvijotham, K. Pillutla, A. Ganesh, T. Steinke, and A. Thakurta. Corre- lated noise provably beats independent noise for differentially private learning. In International Conference on Learning Representations (ICLR), 2024.
E. Deadman, N. J. Higham, and R. Ralha. Blocked schur algorithms for computing the matrix square root. In International Workshop on Applied Parallel Computing (PARA), 2012.
S. Denisov, H. B. McMahan, J. Rush, A. Smith, and G. A. Thakurta. Improved Differential Privacy for SGD via optimal private linear operators on adaptive streams. In Conference on Neural Information Processing Systems (NeurIPS), 2022.
10

---

K. Dvijotham, H. B. McMahan, K. Pillutla, T. Steinke, and A. Thakurta. Efficient and near- near-optimal noise generation for streaming differential privacy, 2024. 1rX1v:2404.16706 [es.DS].
C. Dwork. Differential privacy. In International colloquium on automata, languages, and program- ming (ICALP). 2006.
C. Dwork and A. Roth. The algorithmic foundations of differential privacy. Foundations and Trends in Theoretical Computer Science, 2014.
C. Dwork, F. McSherry, K. Nissim, and A. Smith. Calibrating noise to sensitivity in private data analysis. In Theory of Cryptography Conference (TCC), 2006.
F. Elliott. The characteristic roots of certain real symmetric matrices, 1953. M. thesis, University of Tennessee,
S. Gershgorin. Über die Abgrenzung der Eigenwerte einer Matrix. Proceedings of the USSR Academy of Sciences. Mathematics Series, (6):749-755) -544, 1931.
Granqvist, C. Song, Åine Cahill, R. van Dalen, M. Pelikan, S. Chan, X. Feng, N. Krishnaswami, V. Jina, and M. Chitnis. pfl-research: simulation framework for accelerating research in private federated learning, 2024. Uttthtts:://aat http///rixv..o/aabaaa/a/a 2404.06430
D. H. Greene and D. E. Knuth. Mathematics for the Analysis of Algorithms. Springer, 1990.
M. Henzinger, Upadhyay, and S. Upadhyay. A unifying framework for differentially private sums under continual observation. In Symposium on Discrete Algorithms (SODA), 2024.
Kairouz, B. McMahan, S. Song, O. Thakkar, A. Thakurta, and Z. Xu. Practical and private (deep) learning without sampling shuffling. In International Conference on Machine Learing (ICML), 2021.
A. Koloskova, R. McKenna, Z. Charles, J. Rush, and H. B. McMahan. Gradient descent with linearly correlated noise: Theory and applications to differential privacy. In Conference on Neural Information Processing Systems (NeurIPS). 2023.
A. Kurakin, S. Song, S. Chien, R. Geambasu, A. Terzis, and A. Thakurta. Toward training at imagenet scale with differential privacy, 2022. arXiv:2201.12322
C. Li, G. Miklau, M. Hay, A. McGregor, and V. Rastogi. The matrix mechanism: Optimizing linear counting queries under Differential Privacy. International Conference on Very Large Data Bases (VLDB), 2015.
Z. Li, B. Ding, C. Zhang, N. Li, and Zhou. Federated matrix factorization with privacy guarantee. International Conference Very Large Data Bases (VLDB), 2021.
MacWilliams and J. Sloane. The theory of error-correcting codes. Elsevier, 1977.
R. McKenna. Scaling up the banded matrix factorization mechanism for differentially private ml, 2024. urXiv:2405.1591. DS).
H. B. McMahan, Z. Xu, and Y. Zhang. hassle-free algorithm for private learning in practice: Don't use tree aggregation, use blts, 2024. arXiv:2408.08888 [cs.DS].
T. V. Nguyen, Y. Mori, and T. Mori. Relaxed monotonic conditions for Schur stability of real polynomials. IEICE Transactions Fundamentals of Electronics, Communications and Computer Sciences, (0(10))3326233) 2328, 2007.
N. Papernot and T. Steinke. Hyperparameter tuming with renyi differential privacy. In International Conference on Learning Representations (ICLR), 2021.
N. Ponomareva, H. Hazimeh, A. Kurakin, Z. Xu, C. Denison, H. B. McMahan, S. Vassilvitskii, S. Chien, and A. G. Thakurta. How to DP-fy MI: practical guide to machine learning with differential privacy. Journal of Artificial Intelligence Research (JAIR), 11:111- 2023.

---

SebastienB.
Trace
norm of a triangular matrix with only ones above the diag-
onal, 2017.
URL otw.
H. Shin, S. Kim, J. Shin, and X. Xiao. Privacy enhanced matrix factorization for recommendation with local differential privacy. IEEE Transactions on Knowledge and Data Engineering. 30(9), 2018.
S. Vadhan. The complexity of differential privacy. Tutorials on the Foundations of Cryptography, 2017.
D. Yu, H. Zhang, W. Chen, J. Yin, and T.-Y. Liu, Gradient perturbation is underrated for differentially private convex optimization. In International Joint Conference on Artificial Intelligence (IJCAI), 2020.
X. Zhang, J. Ding, M. Wu, S. T. Wong, H. Van Nguyen, and M. Pan. Adaptive privacy preserving deep learning algorithms for medical data, In Winter Conference on Applications of Computer Vision (WACV), 2021.
12

---

Appendix
A
General introduction to differential privacy.
Differential privacy [Dwork, 2006] is robust framework designed to provide strong privacy guaran- lees for statistical analyses and data sharing. It aims protect individual data points dataset while still allowing meaningful aggregate information to be extracted. Unlike traditional data anonymization techniques, which might involve removing identifiers or aggregating data, differential privacy offers mathematical definition of privacy that quantifies the amount of privacy loss and ensures that the risk of identifying any individual's data remains low, even when combined with other data sources. To formalize this concept, randomized mechanism M is said to provide (e, -)-differential privacy if, for all data sets and D' that differ in one element, and for all subsets of the mechanism's output space S:
PD(M(D) S E Pr(M(D') E S] +
For detailed introduction to differential privacy, recommend the books "The Algorithmic Foundations of Differential Privacy" by Dwork and Roth 2014] and "The Complexity of Differential Privacy" by Vadhan 2017
B Source code for computing AOF
Algorithm 2 Source code for computing AOF using cvxpy.
import cvxpy as import numpy as
def banded_factorizationd b): n len X cp. Variabla ((n, n), PSD=True)
сp pp..maaai matrix_frac(4 X) tr e X--1 A)
objective cp. constraints Minimize (ep.iiizeeppmmarii.. frac (A. mp. (n b))
10
(X) == 1]
11
for i in range n): constraints [ep. diag i) 0, cp diag 0]
13
14
15
prob = cp. Problem (objective, constraints)
return prob. solve (solver=/ SCS')) ; value
C Network architecture for CIFAR-10 experiments
Table ConvNet architecture for CIFAR-10 10 experiments
Conw2DDchammels=322.kk kernel=(3, 3), strides= 1), SAME activation='relu))
Co kernel=(3, 3), strides- 1), PSdding=SSAME: activation=
MaxPool/kernel=:22 2), strides=(2, 2))
Conv2D[channels=64,, kernel=(3, 3), strides= 1), padding='SAME') activation='relu")
Conv2DDchanmels-644 kernel=(3, 3), strides- 1), padding= activation='relu)
MaxPool/kernell122 2), strides=(2, 2))
Conv2Dichannels=122 128, kernel= 3), strides- 1), SAME), activation=Irelu)
Conv2D[channels=1222 128 kernel=(3, 3), strides= padding=SSAEE'' activation=Trelu)
Max 2), strides=(2, 2))
Flatten()
(onse(Outputs=10)
13

---

D Additional Experimental Results
In this section we provide additional experiments comparing BSR, AOF and baselines: Figures and following tables show their expected approximation error (lower is better) for workload matrices stemming from SGD with different hyperparameter settings. Figure(5) and following tables show the accuracy of resulting classifiers CIFAR-10 (higher is better) for different privacy levels.
The results show the same trends as the one in Section BSR achieves almost identical expected approximation error AOF, and results in equally good classifiers. In some cases, results for BSR even improve over sOF's. Presumably this is because of numerical issues in solving the optimization problem for AOF.
a=
1.0
B B=
100
0= 1.0
0= b
150
approximation quality.
approximation quallty:
AOF
ADF
100
1000
id
A
200
400 600 800 1000
1500
2000
200 400 600 800 1000
1500
2000
workload matrix size,
workload matrix size,
0.9999
3= b = 100
0.9999
B= 0.9 b = 100
approximation qualtty. E
approximation quality, 8
BSR
100
1000
50
500
200
400 600 800 1000
1500
2000
200
400 600 800 1000
1500
2000
workload matrix size,
workload matrix size,
0.999
B
b= 100
0.999
B
100
approximation quality, e
spproximation quality, e
80
AOF X
BSR
800
ADF
BSR
60
600
id
40
20
200
0
200
400 600
800 1000
1500
2000
200
400
600 800 1000
1500
2000
workload matrix
workload matrix
a= 0.99
0=0. b 100
0.99 0= b 100
approximatoo quality, 2
30
BSR
wwwroximaaaoo.. quality
300
+
BSR
AOF
ADF
20
200
id
A
10
100
200 400
600 800 1000
1500
2000
200
400
600
800 1000
1500
2000
worklead
Figure 3: Expected approximation error of p-SSR, AOF and baseline factorizations with Imultiple participations and p = b 100.
14

---

1.0 ==
0.0
k=
1.0
B 0.9 k=
approxinatoon quality.
30
approximation quality.
300
BSR
ADF
ADF
id
id
10
100
04-24
500
workload ma trix size,
workload matrix size, n
0.9999
B 0.0.0 k= 1
0.9999 0.9 k= 1
approximation quality
30
approximation quality.
300
ADF
20
200
id
10
100
0
200
500
1000
2000
200
1000
2000
workload matrix size,
workload matrix size,
0= 0.999 8=0. k=1
d= 0= k=1
approximation quality, 2
BSR
approximation quality.
BSR
AOF
15
150
10
100
0
200
500
1000
2000
200
1000
2000
workload matrix size,
workload matrix size,
0.99 8=0.0 k=1
0.99 0= k=1
approkimatoo qualtyy e
approximation quality, E
BSR
BSR
6
AOF
60
AOF
Id
4-
40
20
200
500
1000
2000
200
500
1000 workkooad matrix size,
2000
workload matrix
Figure 4: Expected approximation error of BSR, AOF and baseline factorizations with single participation (k: 1, p n)).
15

---

expected factorization error
π
BSR
sqrt
AOF
ld
A
100
2.4
2.4
2.2
7.1
10.0
200
3.6
4.0
3.5
14.2
22.4
300
4.8
5.5
4.6
21.2
37.4
400
5.9
6.9
5.7
28.3
54.8
500
6.9
8.4
6.7
35.4
74.2
600
8.0
9.9
7.7
42.5
95.4
700
9.0
11.3
8.6
49.5
118.3
800
10.1
12.8
9.5
56.6
142.8
900
11.1
14.2
10.4
63.7
168.8
1000
12.1
15.7
11.3
70.7
196.2
1500
17.2
23.1
15.7
106.1
352.1
2000
22.2
30.6
19.9
141.5
535.7
Table Numeric results for Figure well
plain square root decomposition: a 1, B
k 1, = k/n
expected factorization error
BSR
sqrt
AOF
Id
100
13.9
13.9
13.3
61.8
92.9
200
22.9
25.7
22.4
132.3
213.2
300
31.4
37.4
31.0
202.9
361.2
400
39.7
49.1
39.3
273.6
532.6
500
48.0
61.0
47.4
344.3
724.7
600
56.2
73.0
55.4
415.0
935.3
700
64.3
85.1
63.2
485.7
1163.0
800
72.5
97.2
71.1
556.4
1406.6
900
80.6
109.5
78.8
627.1
1665.2
1000
88.7
121.8
90.6
697.8
1937.9
1500
129.2
184.2
137.2
1051.3
3491.5
2000
169.7
247.8
196.9
1404.9
5322.6
Table 3: Numeric results for Figure as well as plain square root decomposition: a: 1, 009.
k 1, k/n
expected factorization error
BSR
sqrt
AOF
Id
A
100
2.4
2.4
2.2
7.1
10.0
200
3.6
3.9
3.5
14.1
22.2
300
4.8
5.4
4.6
21.0
36.9
400
5.8
6.9
5.7
27.9
53.9
500
6.9
8.3
6.6
34.8
72.7
600
7.9
9.7
7.6
41.6
93.1
700
8.9
11.1
8.5
48.4
115.1
800
9.9
12.5
9.4
55.1
138.4
900
10.9
13.9
10.3
61.8
163.0
1000
11.8
15.3
11.1
68.5
188.7
1500
16.5
22.3
15.2
101.1
332.6
2000
21.0
29.1
19.7
132.6
496.8
Table Numeric results for Figure]
as well as plain square root decomposition: a 0.9999.
8 = 0, 1: 1, k/n
16

---

expected factorization error
BSR
sqrt
AOF
ld
A
100
13.9
13.9
13.2
61.6
92.4
200
22.8
25.5
22.3
131.4
211.4
300
31.2
37.0
30.8
200.9
356.7
400
39.3
48.6
38.9
270.0
524.0
500
47.3
60.1
46.8
338.6
710.3
600
55.3
71.7
54.5
406.8
913.3
700
63.1
83.3
62.1
474.6
1131.5
800
70.9
95.0
69.6
541.9
1363.5
900
78.6
106.6
77.0
608.8
1608.1
1000
86.2
118.3
88.1
675.3
1864.6
1500
123.7
176.5
131.0
1001.3
3298.4
2000
159.9
234.4
187.0
1317.1
4937.9
Table 5: Numeric results for Figure
well
plain square root decomposition: 0.9999.
= 0.9, k 1, b,b= k/n
expected factorization error
BSR
sqrt
AOF
Id
A
100
2.3
2.3
2.2
6.9
9.5
200
3.5
3.8
3.4
13.3
20.5
300
4.6
5.1
4.4
19.3
33.1
400
5.5
6.4
5.4
25.0
46.7
500
6.4
7.6
6.2
30.4
61.1
600
7.2
8.7
7.0
35.4
76.0
700
8.0
9.8
7.7
40.2
91.2
800
8.7
10.8
8.4
44.8
106.5
900
9.4
11.8
9.1
49.2
122.0
1000
10.0
12.8
53.3
137.4
1500
13.0
17.2
12.5
71.6
212.7
2000
15.4
14.8
86.9
282.9
Table 6: Numeric results for Figure 3 as well as plain square root decomposition: =0.999,
0, k 1, k/n
expected factorization error
BSR
sqrt
AOF
ld
A
100
13.5
13.5
12.9
59.8
88.6
200
21.8
24.2
21.4
124.0
195.8
300
29.2
34.3
28.9
184.5
319.9
400
36.1
44.0
35.8
241.4
455.3
500
42.5
53.3
42.2
295.2
598.5
600
48.6
62.3
48.1
346.0
746.7
700
54.3
70.9
53.7
394.2
808.2
800
59.8
79.2
59.0
440.0
1051.7
900
65.0
87.2
64.1
483.6
1205.9
1000
70.0
95.0
71.8
525.1
1360.2
1500
91.9
130.3
95.1
708.4
2113.9
2000
110.3
161.0
115.0
861.2
2816.2
Table Numeric results for Figure as well plain square root decomposition: 0.999,
8 = 0.9, k = 1, = k/n
17

---

expected factorization error
T2
BSR
sqrt
AOF
Id
A
100
2.1
2.1
2.0
5.4
6.6
200
3.0
3.1
2.9
8.7
11.2
300
3.7
3.8
3.6
11.2
14.9
400
4.3
4.5
4.2
13.3
18.1
500
4.8
5.0
4.7
15.1
20.8
600
5.2
5.5
5.2
16.6
23.3
700
5.7
6.0
5.6
18.1
25.5
800
6.1
6.4
6.0
19.4
27.5
900
6.4
6.8
6.4
20.7
29.4
1000
6.8
7.2
6.8
21.9
31.2
1500
8.3
8.9
8.3
27.0
38.9
2000
9.6
10.3
9.6
31.3
45.4
Table 8: Numeric results for Figure
as well
plain square root decomposition: 0.99, 3
k 1, b- = k/n
expected factorization error
m
BSR
sqrt
AOF
Id
100
10.9
10.9
10.5
46.2
61.4
200
16.2
17.0
15.9
79.9
106.7
300
20.3
21.7
20.0
104.5
144.0
400
23.6
25.7
23.4
124.5
175.4
500
26.6
29.1
26.4
141.7
202.7
600
29.2
32.1
29.1
157.1
226.9
700
31.7
34.9
31.5
171.1
248.8
800
33.9
37.5
33.8
184.0
268.9
900
36.0
40.0
35.9
196.1
287.7
1000
38.0
42.3
38.1
207.4
305.3
1500
46.7
52.2
46.8
256.9
381.4
2000
54.1
60.6
56.0
298.2
444.6
Table 9: Numeric results for Figure
as well as plain square root decomposition:
0.99,
0.9, 1, b = k/n
expected factorization error
n
BSR
sqrt
AOF
Td
50
2.2
2.2
2.0
5.0
7.1
100
2.4
2.4
2.2
7.1
10.0
200
2.6
2.6
2.4
10.0
14.1
400
2.8
2.8
2.7
14.2
20.0
500
2.9
2.9
2.7
15.8
22.4
1000
3.1
3.1
2.9
22.4
31.6
2000
3.3
3.3
3.2
31.6
44.7
Table 10: Numeric results for Figure well plain square root decomposition:
1,
b = 0
18

---

expected factorization error
π
BSR
sqrt
AOF
Id
A
50
11.4
11.4
10.6
38.2
60.3
100
13.9
13.9
13.3
61.8
92.9
200
16.3
16.3
15.7
93.5
136.5
400
18.6
18.6
17.8
136.8
196.5
500
19.3
19.3
18.6
154.0
220.5
1000
21.6
21.6
23.9
220.7
314.0
2000
23.8
23.8
27.1
314.1
445.7
Table 11: Numeric results for Figure as well as a plain square root decomposition: a = 1, 0.9,
b = 0
expected factorization error
BSR
sqrt
AOF
Id
50
2.1
2.1
2.0
5.0
7.1
100
2.4
2.4
2.2
7.1
10.0
200
2.6
2.6
2.4
10.0
14.0
400
2.8
2.8
2.7
14.0
19.6
500
2.9
2.9
2.7
15.6
21.8
1000
3.1
3.1
2.9
21.7
30.1
2000
3.2
3.2
3.1
29.7
40.6
Table 12: Numeric results for Figure well plain square root decomposition: 0.9999.
0,
expected factorization error
7
BSR
sqrt
AOF
Id
50
11.4
11.4
10.6
38.2
60.2
100
13.9
13.9
13.2
61.6
92.4
200
16.2
16.2
15.6
92.9
135.2
400
18.4
18.4
17.7
135.0
192.7
500
19.1
19.1
18.4
151.4
215.2
1000
21.1
21.1
23.4
213.5
299.1
2000
23.0
23.0
26.0
294.5
404.7
Table 13: Numeric results for Figure as well plain square root decomposition: 0.9999,
3 0.9, b = 0
expected factorization error
BSR
sqrt
AOF
Id
A
50
2.1
2.1
2.0
5.0
6.9
100
2.3
2.3
2.2
6.9
9.5
200
2.5
2.5
2.4
9.4
12.8
400
2.7
2.7
2.6
12.5
16.6
500
2.7
2.7
2.6
13.6
17.8
1000
2.8
2.8
2.7
16.9
20.8
2000
2.8
2.8
2.8
19.4
22.2
Table 14: Numeric results for Figure 0,
as well as a plain square root decomposition: a 0.999,
19

---

expected factorization error
π
BSR
sqrt
AOF
Id
A
50
11.2
11.2
10.4
37.6
59.0
100
13.5
13.5
12.9
59.8
88.6
200
15.5
15.5
14.9
87.7
124.2
400
17.0
17.0
16.5
120.7
163.3
500
17.4
17.4
17.0
132.0
175.6
1000
18.4
18.4
20.0
166.1
206.6
2000
18.8
18.8
20.8
192.6
220.5
Table 13: Numeric results for Figure 4 as well as a plain square root decomposition: a = 0.999,
0.9, b = 0
expected facttorization error
7
BSR
sqrt
AOF
Id
50
2.0
2.0
1.9
4.3
5.6
100
2.1
2.1
2.0
5.4
6.6
200
2.1
2.1
2.1
6.2
7.0
400
2.1
2.1
2.1
6.6
7.1
500
2.1
2.1
2.1
6.7
7.1
1000
2.1
2.1
2.1
6.9
7.1
2000
2.1
2.1
2.1
7.0
7.1
Table 16: Numeric results for Figure as well as plain square root decomposition: 0.99,
0,
expected factorization error
n
BSR
sqrt
AOF
Id
A
50
9.8
9.8
9.2
32.8
48.7
100
10.9
10.9
10.5
46.2
61.4
200
11.5
11.5
11.2
56.5
66.9
400
11.7
11.7
11.6
62.3
67.7
500
11.8
11.8
11.7
63.4
67.7
1000
11.9
11.9
12.2
65.6
67.7
2000
11.9
11.9
12.3
66.7
67.7
Table 17: Numeric results for Figure as well as a plain square root decomposition: a 0.99, 0.9, b = 0
20

---

CIFAR-10 (single participation)
CIFAR-10 (repeated participation)
40
accuracy
BSR
BSR
20
AOF
ADF
20:
id
A
10
25
200
500
1000
2000
10
20
number update steps
maximum number of epochs
CIFAR-10 (single participation)
CIFAR-10 (repeated participation)
50-
accurapy
30
40
AOF
id
20
25 200
500
1000
2000
5
10
20
number update steps
maximum number efeesochs
Figure 5: Classification accuracy (mean and standard deviation over runs with different random seeds) on CIFAR-10 for BSR, AOF, and baselines for independent training runs. Top row: classi- fication accuracy on CIFAR- with (e, d) = (2, 10-5). Bottom row: classification accuracy on CIFAR-10 with (e, 8) = (8, 10 10-5). Left plots: one epoch, different batch sizes. Right plots: different number epochs, constant batch size.
accuracy
number of updates
BSR
AOF
Id
50
29.5
2.1
31.5+
2.3
21.8 + 2.1
100
34.8
+ 2.2
3.0
30.8 0.8
22.7
200
34.5
14.5+0. 0.8
3.2
18.7 1.6
400
3.78
3..4+1.1
27.5
17.4
500
36.3+200
3.7
2.7
17.2 1.8
1000
36.2
3..7+1.11
2..2
13.6
2000
3.7
3.14
2.19
12.7 1.5
Table 18:
Numeric values for results Table
left plot (CIFAR-10, single participation, (c,
(4, 10-9).
21

---

accuracy
number of epochs
BSR
AOF
1d
A
34.8
36.2 1.0
0.8 3.0
2.
0.6
1.9
21.3
4.0
4.7
3.5
191
4.7
0.6
36.6
2.0
4.90
4.9
0.6
21.3=
4..0
4.9
1.4
20.1
4..6++077
4.9
3.1
19.6
4.79
4.8
3.5
18.8
4..7
0.8
1..6 1..
10
0.5
0.8
14
2.7
15
5.7
47.9
3.8
17.8 = 1.0
20
53.1
47.0 + 0.2
0.7
17.7 = 1.2
Table 19:
Numeric values for results Table
right plot (CIFAR-10, 10, repeated participation, (c,
(4, 10 10-
accuracy
number of updates
BSR
AOF
Id
50
29.8 = 2.7
1.4
30.2 1.5
20.8 0.8
100
31.2 0.9
3..0
1.2
1.9
200
30.9
1.9
1.. 1.60
400
32.1
0.9
3.8
24.8
14.7
500
31.6 11.6222 2.3
3.8
2
14.9 + 2.0
1000
30.2 = 1.6
3.11
17.5 0.7
12.5
2000
29.9
2.0
30.5
1.01
1.13
Table 20: Numeric values for results in Table top left plot (CIFAR-10, 10, single participation,
(c, 6) = (2, 10-
accuracy
number of epochs
BSR
AOF
ld
A
3..2
34.0
29.8 1.2
16.5
2
34.5
0.6
3.9
18.
2.3
3.2
1..
17.6
2.1
40.1
2.9
17.4
4..4
4.9
2.04
16.1=211
4.7
1.3
2.9
1.3
4..2+1.6
4.9
15.4+ 28.4+1.0
8
4.2
4.14+
16.1 1..9
4..1
4.5
15.6 28.6
10
4.8
0.3
2.5
15.8 - 0.9
15
4..911.1
0.8
2.1
14.9 1.0
20
4.2
43.0
14.2 2.5
Table 21: Numerie values for results in Table
top right plot (CIFAR-10, repeated participation,
(c. 5) = (2, 10-5..
22

---

accuracy
number of updates
BSR
AOF
1d
A
50
31.4 1.8
3.4 32.6 2.0
25.3 + 1.4
100
37.7
1.9
1.3
26.3 2.7
200
1.8
41.
0.5
2.7
400
1.1
1.0
22.1 34.4144
500
4.14
0.6
3.5
2.5
1000
40.7 = 1.0
0.6
28.2
16.1 + 4.2
2000
39.2
1.5
37.6 1.5
23.5
16.3 + 1.5
Table 22: Numeric values for results in Table
bottom left plot (CIFAR-10, 10, repeated participation,
(c, 6) (8, 10-
accuracy
number of epochs
BSR
AOF
Id
A
37.7 2.6
39.7 + 1.9
38.3 + 1.3
26.3 + 2.7
1.4
4.2
1.0
2.0
1..
4.1 4.8
24.6=
1.0
1.6
4.7
24.1
1..
5.2
4.3
24.0=200
5.6
0.2
4..
24.0
5.8
0.9
1.3
24.5
55.7
0.5
4.
25.1
0.5
1.1
4.6
23.7
10
56.2+
5..1
4.1
+
15
16.0++00
52.3 + 0.5
3.0
2.3
20
55.2
- 1.2
53.8 + 0.5
20.6 1.
Table 23: Numeric values for results in Table
bottom right plot (CIFAR-I0, repeated participation,
(e, 6) = (8, 10-9).
E Experimental Results for Different Optimizers
In this section we report on experimental results when different optimizers are used to (approxi- mately) solve the AOF optimization problem (4). Besides cvxpy (CVX) these are standard gradient descent (GD) and the Limited-Memory The latter two implement in jax using the optax toolbox. Similar [Granqvist et al., 2024 ftrl_me sarr__mechaiise sm. s..ll, we use an adaptive line-search for the step size of the gradient-based meth- ods, which at the same time ensures the positive definiteness constraints of the optimization problem. Our implementation differs from theirs, however, in that our learning rate is not restricted to shrink monotonically, thereby avoiding premature termination.
Fepulen--letchhe--oollfrbbbbbbb Sahifaorb-Shnoo algorithm (LBFGS).
F.1 Runtime
We report the runtimes for the different methods in Tables 24 to 39 For comparison, we also include results for BSR and the CVX optimizers with three tolerance levels in the same settings, where practically feasible. Note that while the experiments for BSR and CVX used single-core CPU-only environment, the experiments for GD and LBFGS were run on an NVIDIA H100 GPU with 16 available CPU cores. As consequence, the absolute runtimes not directly comparable between the methods, but they should rather be seen as illustrations the scaling behavior the method for different workload types and problem sizes.
Indeed, the results show a clear trend: BSR is the fastest, with almost no overhead. Even for largest problem sizes of 10 BSR never took more than 2.58 to despite running in the single-core CPU-only setup. GD and LBFGS benefit strongly from the GPU hardware. In the multiple participation setting (p = 100, k n/p), they solve most workload sizes within a few
23

---

Table 24: 1.0, 0.9, p= 100, k n/pp
BSR
GD
LBFGS
CVX((ol=0.01)
CVX(tol=0.001)
CVX(tol=0.0001)
100
< 1s
28.5s
Im39s
4.58
7m18s
1h27m30s
200
1s
Iml0s
2m31s
37.1s
21m00s
10h35m00s
300
1s
1m44s
3m
2m16s
1h12m40s
22h36m40s
400
< 1s
2m35s
3m46s
6m06s
2h14m50s
53h03m20s
500
< 1s
3m47s
4m47s
1m445
5h31m40s
90h50m00s
600
1s
4m27s
Smlls
26m50s
17h36m40s
40h16m40s
700
1s
5ml3s
6ml2s
47m50s
22h10m00s
66h23m20s
800
1s
5m52s
7m30s
1h29m40s
38h03m20s
164h26m40s
900
< 1s
6m20s
7m29s
1h45m30s
62h30m00s
253h53m20ss
1000
< 1s
6m55s
8m01s
1h59m40s
83h36m40s
245h00m00s
1500
1s
10mlls
1l449s
6h08m20s
121h23m20s
timcout
2000
1s
13m39s
13m21s
15h10m00s
297h13m20s
timcout
5000
1.1s
1h09m55s
33m00s
10000
1.6s
6h07m10s
1h47m39s
Table 25: a = 1.0, 3 = 0.0, 100, n/p
n
BSR
GD
LBFGS
CVX(tol=0.01)
CVX(tol=0.001)
CVX(tol-0.000)
100
<
1.4s
12.3s
5.6s
48.6s
33.4s
200
< 1s
2.4s
20.0s
1m09s
2m38s
3m39s
300
< 1s
3.4s
25.0s
2m25s
1 1123s
19m20s
400
< 1s
6.4s
30.2s
12m15s
26m10s
1h14m40s
500
<
6.8s
36.6s
33m10s
1h07m00s
2h53m20s
600
< 1s
9.8s
44.1s
52m00s
7h45m00s
4h53m20s
700
< 1s
12.8s
51.5s
1h34m50s
14h30m00s
17h48m20s
800
< 1s
13.2s
54.3s
3h33m20s
17h26m40s
26h55m00s
900
17.9s
59.1s
4h43m20s
27h 8m20s
41h06m40s
1000
<
22.9s
Im09s
7h31m40s
45h33m20s
57h30m00s
1500
< 1s
47.3s
Im36s
29h43m20s
84h43m20s
200h50m00s
2000
< 1s
Im24s
Im45s
68h53m20s
258h36m40s
timeout
5000
2.4s
16m35s
4m18s
10000
1.6s
2h29m46s
14m23s
minutes, except the largest ones, which for GD can take few hours. In the single participation setting (k 1), LBFGS also occasionally need several hours converge. In general, stronger weight decay (smaller a) tends to lead to lower runtimes, while the use of momentum (3 0.9) to higher times until convergence. CVX (on weak hardware) orders of magnitude slower than the other methods. Furthermore, its runtime grow approximately cubic with the problem size, whereas for GD and LBFGS the relation is not too far from linear. Note that despite the stable patterns described above, all runtime results should be taken with caution, because internal parameters of the optimization, such as the convergence criterion and the specific implementation of the line search can substantially influence the overall runtime well.
E.2 Expected Approximation Error
Figures 6 and report the expected approximation errors achieved by the different optimizers of AOF 4 and by BSR. For CVX, we report the smallest error value across all tolerance levels for which the optimization converged.
The curves show several clear trends. GD and LBFGS generally perform similarly, and achieve expected approximation errors slightly (at most few percent) lower than BSR. An exception the problems with momentum (3 = 0.9) in the single participation setting, where appears that SGD) occasionally fails to find the optimum for large problem sizes (n > 1500). CVX performs comparably to the other methods for problems without momentum (3 0). With momentum, however, the solutions it find are often worse than the other methods, especially in the single-participation setting and for medium to large problem sizes (n 2 500). Presumably, even smaller tolerance values would be require here, which, however, would result in even longer runtimes.

---

Table 26: 0.9999, 3 0.9, p 100, k: n/p
BSR
GD
LBFGS
CVX(tol=0.01)
CVX(t0)=0.001)
CVX(tol=0.0001)
100
1s
29.9s
Im35s
6.1s
3m42s
47m00s
200
1s
Im10s
2m33s
36.3s
18m00s
6h56m40s
300
< 1s
1m45s
3m04s
1m48s
40m30s
25h46mm4ss
400
< 1s
2m32s
3m46s
6m21s
2h35m20s
57h13m20s
500
< 1s
3m32s
4m29s
11m05s
5h41m40s
24h08m20s
600
1s
4m21s
4m56s
10m41s
13h36m40s
41h56m40s
700
1s
5m
5m37s
54m20s
20h33m20s
76h56m40s
< 1s
5m39s
6m20s
m40s
38h20m00s
158h53m20s
900
< 1s
6m25s
7m05s
2h1 1m00s
51h06m40s
20m0s
1000
1s
6m166
7m46s
2h1 1m00s
68h03m20s
223h20m00s
1500
1s
10m09s
10m15s
8h08m20s
115h50m00s
timeout
2000
1s
13m39s
12m09s
12h46m40s
247h1 3m22s
timeout
5000
2.4s
1h10ml 18s
26m09s
10000
2.6s
6h06m46s
1h10m I0s
Table 27: 0.9999, 3 0.0, p 100, h: = nnpp
BSR
GD
LBFGS
CVX(tol=0.01) CVX(tol=0.001)
CVX(tol=0.0001)
100
< 1s
1.2s
13.7s
6.0s
14.1s
20.3s
200
< 1s
2.3s
18.5s
1m04s
4m15s
4m00s
300
< 1s
3.5s
23.8s
4m03s
15m10s
10m51s
400
<
5.3s
29.7s
7m45s
14m32s
19m50s
500
<
6.1s
38.5s
22m10s
1h23m50s
2h44m20s
600
< 1s
8.9s
43.3s
29m40s
3h46m40s
8h23m20s
700
< 1s
12.0s
44.5s
1h56m50s
12h03m20s
12h1 1m40s
800
< 1s
12.7s
51.6s
3h43m20s
16h40m00s
23h56m40s
<
16.1s
1m02s
1m40s
26h30m00ss
30h33m20s
1000
< 1s
20.3s
lm05s
8h23m20s
40h33m20s
41h23m20s
1500
< 1s
39.8s
1m19s
35h50m00s
78h53m20s
170h00m00s
2000
< 1s
1m10s
Im27s
72h30m00s
239h26m40s
timcout
5000
1.1s
9m41s
3m15s
10000
2.5s
1h00m37s
10m51s
Table 28: 0.999, 8 0.9, 100, h: n/p
BSR
GD
LBFGS
CVX(tol=0001)
CVXitol=0.00L)
CVX(tol=000001)
100
<
32.6s
Im52s
3.7s
3m36s
58m10s
200
< 1s
59.0s
2m08s
36.1s
20m50s
8h10m00s
300
< 1s
Im25s
2m36s
3m38s
37m30s
25h48m20s
400
< 1s
Im54s
3m02s
4m48s
1h30m00s
56h56m40s
500
< 1s
2m13s
3m33s
13m45s
2h58m20s
84h43m20s
600
<
3m07s
3m57s
26m20s
9h51m40s
41h56m40s
700
< 1s
3m04s
4m09s
52m30s
13h10m00s
85h33m20s
800
< 1s
3m28s
4m18s
59m00s
28h03m20s
164h43m20s
900
< 1s
3m59s
4m33s
57m00s
39h1 0m00s
280h33m20s
1000
<
4m34s
5m21s
1h27m00s
59h43m20s
258h03m20s
1500
<
5m48s
5m21s
4h41m40s
81h40m00s
timeout
2000
< 1s
7m26s
5m28s
14h16m40s
119hl0m00s
timeout
5000
< 1s
34m57s
9m16s
10000
2.5s
2h42m14s
28m36s
25

---

Table 29: 0.999, 0.0, 100, h n/p
BSR
GD
LBFGS
CVX(tol=0.01)
CVX(tol=0.001) CVX(tol=0.0001)
100
1s
1.1s
11.5s
7.3s
28.0s
27.6s
200
2.0s
15.6s
Im10s
46.8s
3m06s
300
3.1s
21.3s
3m39s
6m32s
8m10s
400
< 1s
3.8s
26.0s
8m19s
30m30s
29m50s
500
< 1s
4,6s
30.3s
9m47s
44m50s
2h17m50s
600
1s
5.3s
29.3s
48m30s
2h01m40s
1h59m50s
700
6.1s
39.58
1h10m20s
6h50m00s
7h01m40s
800
< 18
6.78
38.6s
3h16m40s
13h 13h10m00s
13050m00s
900
< 18
7.9s
37.4s
5h26m40s
21h51 1
28h03m20s
1000
< 1s
9.8s
45.1s
6h26m40s
30h50m00s
32h46mm40s
1500
12.8s
52.1s
29h26m40s
77h13m20s
49h43m20s
2000
50.3s
70h00m00s
91h23m20s
174h26m40s
5000
1.3s
Im19s
1m29s
10000
2.4s
6m40s
3m23s
Table 30: 0.99, 3 0.9, p 100, k= n/p
BSR
GD
LBFGS
CVX(tol=0.01) CVX(tol=0.001)
CVX(tol=0.0001)
100
< 1s
15.8s
Iml2s
3.5s
1m28s
1h19m10s
200
< 1s
19.0s
Im09s
14.4s
22m00s
1 11h23m20s
300
1s
19.6s
Iml5s
Im10s
2h29m50s
30446m40s
400
24.5s
Im23s
3m07s
4h55m00s
59h43m20s
500
1s
23.8s
Im21s
6m55s
13h06m440s
67h30m00s
600
< 1s
28.7s
Im20s
24m10s
3 31h40m00s
44h43m20s
700
1s
33.4s
Im24s
36m10s
61h23m20s
78h53m20s
800
1s
33.7s
Im29s
48m20s
73h36m40s
146h06m40s
39.2s
Im37s
1h06m50s
50h00m00s
280h33m20s
1000
< 1s
37.3s
Im32s
1h53m30s
66h06m40s
274h43m20s
1500
< 1s
53.9s
Im44s
7h36m40s
302h46m40s
timcout
2000
1s
1m07s
Im50s
22h45m00s
timcout
timcout
5000
1s
5m56s
2m51s
10000
2.4s
29m34s
7m45s
Table 31: a 0.99, 3 100, k: n/p
BSR
LBFGS
CVX(tol=0.01)
CVX(tol=0.001)
CVX(tol=0.0001)
100
1.0s
9.9s
3.9s
9.0s
15.9s
200
<
< 1s
9.8s
30.3s
1m40s
m338s
300
< 1s
< 1s
10.2s
2ml Is
2ml Is
Im50s
400
< 1s
1.2s
11.4s
1m166s
10m36s
4m16s
500
1.3s
14.6s
24m10s
17m40s
22m20s
600
1.3s
10.6s
32m00s
39m00s
1h44m00s
700
< 1s
1.4s
10.7s
lh06m106
1h33m20s
1h46m30s
800
< 1s
1.5s
11.2s
1h43m30s
3h36m40s
2m21m20s
900
< 1s
1.68
11.3s
3h15m00s
7h23m20s
8h21m40s
1000
1.8s
13.0s
5h11m400
8h31 m40s
9h26m40s
1500
2.7s
12.7s
20h10m00s
29h 10000s
30h00m00s
2000
< 1s
3.2s
12.5s
46h40m00s
11566m40s
33h03m20s
5000
< 1s
17.0s
18.3s
10000
1.4s
1m26s
42.7s
26

---

Table 32: 1.0, 0.9, =1
BSR
GD
LBFGS
CVX((ol=0.01)
CVX(1ol-0.001)
CVX(tol=0.0001)
100
< 1s
28.5s
Im39s
2.38
Im46s
50m20s
200
1s
1m02s
2m08s
25.5s
31m00s
7h15m00s
300
1s
1m37s
3m07s
Im45s
1h53m40s
21h03m20s
400
< 1s
2m25s
3m34s
5m48s
5h00m00s
50h00m00s
500
< 1s
2m39s
4m05s
7m38s
12h28m20s
71h40m00s
600
1s
5m10s
5m32s
26m30s
25h48m20s
43h20m00s
700
1s
4m22s
5m30s
37m10s
49h26m40s
94h26m40s
800
1s
4m22s
5m42s
1h09m40s
57h30m00s
175h 6mm0s
900
< 1s
5m37s
6m25s
2h04m50s
36h23m20s
305h33m20s
1000
< 1s
6m53s
6m32s
2h19m20s
53h03m20s
260h00m00s
1500
1s
10m24s
10m08s
13h13m20s
242h30m00s
timcout
2000
1s
13m39s
12m45s
33h03m20s
64hl 10000s
timeout
5000
1.1s
1h10m17s
44m58s
10000
< 1s
6h06m27s
2h30m09s
Table 33: a 1.0, 8 = 0.0, k =
n
BSR
GD
LBFGS
CVX(tol=0.01)
CVX(tol=0.001)
CVX(tol=0.0001)
100
<
12.3s
6.8s
6.4s
19.6s
200
< 1s
3.7s
21.8s
45.3s
2m53s
4m32s
300
< 1s
3.3s
25.1s
3m00s
11m09s
7m37s
400
< 1s
9.5s
33.8s
6m34s
19m30s
52m00s
500
19.3s
50.4s
16m10s
41m40s
22m50s
600
< 1s
9.68
42.9s
51m20s
2h26m440s
2h16m20s
700
< 1s
20.3s
54.3s
1h47m10s
4h40mm0ss
5h45m00s
800
< 1s
33.0s
1m00s
3h06m40s
9h40m00s
12h 8m20s
900
<
57.2s
1m22s
5h06m40s
14h05m00s
16h28m20s
1000
<
1ml2s
Im39s
7h03m20s
21h13m20s
28h20m00s
1500
< 1s
1m54s
Im39s
24h5 Im40s
80h 6mm0s
53h20m00s
2000
< 1s
5m39s
2m56s
61h56m40s
10h33m200
146h23m20s
5000
2.5s
35m20s
4m19s
10000
1.6s
5m15m15s
13m51s
Table 34: = 0.9999, 3 0.9,
BSR
GD
LBFGS
CVX(tol=0.01)
CVX(tol=0.001)
CVX(to==00001))
100
1s
29.9s
Im35s
6.5s
5ml2s
1h14m00s
200
1s
Im05s
2m13s
34.4s
33m20s
10h03m20s
300
< 1s
Im38s
2m49s
2m52s
1m41m10s
29h43m20s
400
< 1s
2m30s
3m29s
2m35s
6h30m00s
30h00m00s
500
1s
2m29s
4m03s
7m 7s
11h55m00s
75h50m00s
600
1s
4m57s
4m58s
19m00s
26h40m00s
41h40m00s
< 1s
5m13s
7m24s
42m20s
43h36m40s
99h 10m00s
800
< 1s
3m59s
5m13s
53m40s
56h23m20s
167h30m00s
900
< 1s
5m04s
5m47s
1h27m30s
34h10m00s
302h46m40s
1000
1s
6m53s
6m
2h15m10s
54h10m00s
1mm0ss
1500
1s
10m05s
13m04s
12h23m20s
308h20m00s
timcout
2000
1s
13m37s
9m39s
35h00m00s
timcout
timcout
5000
2.5s
1h10m23s
28m03s
10000
2.5s
6h06m15s
1h54m24s
27

---

Table 35: 0.9999, 0.0, k=
BSR
GD
LBFGS
CVX(10E=0.0))
CVX(101=0.001)
CVX(tol=0.0001)
100
< 1s
1.2s
13.7s
5.1s
21.5s
17.3s
200
< 18
3.1s
20.5s
55.7s
3m45s
7m08s
300
<
2.9s
22.7s
2m58s
4m 10s
8m32s
<
9.38
34.4s
6m59s
31m30s
26m20s
500
< 1s
19.3s
48.9s
30m40s
1h03m00s
41m40s
600
< 1s
8.2s
41.8s
1h24m40s
2h06m40s
1h22m10s
700
< 1s
17.7s
47.9s
1h54m40s
5h25m00s
4h38m20s
800
<
30.3s
57.8s
2h26m50s
8h36m40s
7h43m20s
900
< 1s
45.0s
Iml5s
4h36m40s
13h58m20s
15h35m00s
1000
< 1s
1m01s
Im23s
8h23m20s
18h46m40s
22h43m20s
1500
< 1s <
lm28s
Im33s
28h03m20s
80h33m20s
91h23m20s
2000
3m56s
2mlls
64h26m40s
98h36m40s
126h06m40s
5000
2.5s
58m03s
4m57s
10000
2.5s
6h00m1.3
11m30s
Table 36: 0.999, 8 0.9, k
n
BSR
GD
LBFGS
CVX(tol=0.01)
CVX(tol=0.001)
CVX(tol=000001)
100
<
Im52s
3.5s
6m20s
1h28m10s
200
< 1s
mm00s
2ms33
32.2s
26m30s
8h36m40s
300
< 1s
59.6s
2m20s
2m2ls
2h50m00s
28h20m00s
400
< 1s
2m08s
2m56s
5m53s
3h36m40s
58h03m20s
500
<
4m05s
4m19s
9m40s
13h16m 6mh0s
72h30m00s
600
< 1s
2m03s
3m35s
12m33s
22h50m00s
44h10m00s
700
< 1s
2m30s
3m51s
47m50s
46h06m40s
101h40m00s
800
< 1s
3m33s
4m01s
1h06m20s
78h03m20s
1 145h00m00s
900
<
4m08s
4m39s
1h40m20s
88h36m40s
305h33m20s
1000
<
5m35s
5m31s
2h42m00s
50h00m00s
timeout
1500
< 1s
10m10s
7m29s
10h10m00s
280h33m20s
timeout
2000
< 1s
13m47s
9m12s
26h08m20s
timeout
timeout
5000
< 1s
1h10ml9
19s
21m05s
10000
2.5s
6h08m30s
45m33s
Table 37: 0.999, 3 0.0,
BSR
GD
LBFGS
CVX(tol=0.01)
CVX(tol=0.001)
CVX(tol=0.0001)
100
< 1s
1.1s
11.5s
5.9s
6.3s
32.4s
200
2.4s
16.8s
47.4s
1m30s
1m54s
300
< 1s
6.3s
31.3s
2m48s
7m17s
17m40s
400
< 1s
5.0s
27.8s
8m10s 10s
17m30s
42m50s
500
< 1s
9.4s
33.0s
13m49s
45m50s
1h06m00s
600
<
13.6s
40.9s
50m20s
1h04m20s
3h21m40s
700
<
19.3s
1h53m00s
3h06m40s
3h21m40s
800
< 1s
25.58
52.3s
2h45m50s
7h03m20s
8h45m00s
900
< 1s
34.5s
Im07s
3h55m00s
9h23m20s
1m015m00s
1000
< 1s
42.4s
Iml0s
5h46m40s
17h08m20s
18h21m44
1500
<
Im35s
Im40s
26h55m00s
1m000
56h40m00s
2000
< 1s
2m45s
2m14s
60h33m20s
52h1 3m20s
82m3000ss
5000
< 1s
21m
4m32s
10000
2.4s
1h59m54s
12m44s
28

---

Table 38: a 0.99, 8 = =1
BSR
GD
LBFGS CVX(tol=0.01)
CVX(to==0.001)
CVX(tol=0.0001
100
<
15.8s
lm12s
3.6s
1m53s
35m00s
200
< 1s
15.6s
1m14s
31.8s
14m57s
9h00m00s
300
< 1s
22.5s
Iml5s
Im45s
51m30s
22h15m00s
400
<
33.6s
Im33s
3m166
2h19m00s
61h40m00s
500
41.1s
Im33s
17m20s
3m40s
73h53m20s
600
< 1s
47.8s
1m40s
31m50s
8h40m00s
41h40m00s
700
< 1s
54.5s
lm50s
51
20h26m40s
89h26m40s
800
< 1s
1m03s
Im53s
1h12m40s
39h43m20s
140h00m00s
900
1m12s
Im55s
2h09m30s
56h06m40s
270h00m00s
1000
<
Im05s
2m02s
3h20m00s
61h06m40s
265h 66m40s
1500
< 1s
1m53s
2m05s
13h46m40s
124h26m40s
timeout
2000
< 1s
2m16s
2m14s
32h30m00s
319h26m40s
timeout
5000
<
11m48s
3m20s
10000
2.5s
1h02m59s
7m50s
Table 39: a 0.99, 8 = k=
BSR
GD
LBFGS
CVX(tol=0.01)
CVX(tol=0.001)
CVX(to=00000)
100
1.0s
9.9s
3.8s
8.1s
22.1s 1s
200
< 1s
1.7s
14.0s
30.8s
Im15s
2m02s
300
< 1s
2.5s
18.8s
Iml8s 8s
5m23s
8m00s
400
3,0s
19.6s
3m45s
7m41s
7m43s
500
3.5s
23.9s
9m51s
14m17s
38m30s
600
< 1s
4.4s
26.4s
33m20s
35m00s
45m20s
700
< 1s
4,3s
24.7s
1h02m20s
1h22m40s
1h42m Os
800
< 1s
4,8s
27.8s
1h3 1m40s
2h27m50s
2h56m40s
900
5,4s
28.7s
2h1 1m10s
4h01m40s
4h26m40s
1000
5,6s
30.3s
2h5 1m40s
5h45m00s
6h45m00s
1500
< 1s
8.5s
39.38
11h21m40s
17h55m00s
17h26m40s
2000
< 1s
11.5s
42.9s
26h46m40s
45h00m00s
42h46m40s
5000
57.1s
Im07s
10000
1.38
4m58s
2m42s
29

---

a= 1.0,
B=0.9,
p=b=100 100, k= n/b
1.0,
8=
p=b= 100, =//b
10a
approximaton quality,
BSR
watimmmmbbooo quality,
BSR
AOF CVX
AOF-CVX
AOF-GD AOF-LBFGS
AOF-GD
20a
M
AOF-LBFGS
10
100
1000
10000
100
1000
10000
matrix size,
workload size,
a 0.9999, 0=9.9 p=b 100, k n/b
0.9999,
8=0.0, p=b 100, n/b
approximation quality.
BSR
approximatte quality,
BSR
AOF-CVX
AOF-CVX
AOF-GD
AOF-GD
102
*.
AOF-LBFGS
AOF-LBFCS
10
100
1000
10000
100
1000 workload matrix size, n
10000
workload matrix size,
a = 0.999,
0.9, p=b= 100,
k= n/b
a= 0.999,
8=0.0, p=b= 100,
n/b
approximation quality,
BSR
approxinn Eporoximation quality.
BSR
AOF-CVX
AOF-CVX
10?
AOF-GI
AQF-GD X
36
AOF-LBFGS
101
AOF-LBFGS
1000
10000
100
1000
10000
warkload matrix size,
workload matrix size, n
0.99,
0=9.9 p=b= 100, k n/b
0.99,
0=0.0, p=b= 100, k= n/b
approximation quality.. lit. E
10a
BSR
approximation quallty.
BSR
AOF-CVX
AQF-CVX
AOF-GD
100
AOF-GD 96
AOF-LBFGS
AOF-LBFGS
101
10000
1000
10000
workload matrix size, n
workload matrix size,
Figure 6: Expected approximation error for AOF with GD, LBFGS and CVX optimizers as well as for BSR in the multiple participations setting.
30

---

a=1,0, 0=0.9,
a= 1.0,
B 0.0.0
k=
approximaion quality
10
BSR
Saprokina quality
4 10
BSR
AOF-CVX
AQF-CVI
3 10
AOF-GD
AOF-GD AOF-LBEGS
AOF-LBFGS
3 10
2 10
100
1000
10000
100
1000
10000
matrix size,
n
0.9999,
0=9,
0.9999,
8 0.0,
A approximation quality q
10
BSR
3.5 109
BSR
AOF-CVX
*
AOF-CVX
AOF-GD
3x
*
ADF-GD
2
LBFGS
x
*
ADF-LBFGS
approxin
2.5x
10
100
1000
10000
100
1000
10000
workload matrix size,
workload matrix size,
a= 0.999,
== k=1 1
a 0.999,
==0.0, k=1
approximation quality
2.4x10
x.2
BSR
cenoooiimattee quality.
3 109
BSR
2x
-
AOF-CVX
AOF-CVX
1.8
X
AOF GD
ADF ADF-LBFGS
X
AOF-LBFGS
104
apgroxini
100
10
10
100
1000 workload matrix size, A
10000
1000 workload matrix size, n
10000
0.99, B=0.9, k=1
a= 0.99,
8=0.0, k= 1
approximation quality lky....
1.5 x 10
BSR
toproximation Ei
BSR
1.4 10
AQF-CVX
12x
AOF-CVX
1.3 10
AOF-GD
1.0
AOF-GD
96
AOF-LBFGS
AOF-LBFGS
1
109
10a
1000
10000
workload matrix size,
workload matrix size,
Figure 7: Expected approximation error for AOF with GD, LBFGS and CVX optimizers well as for BSR in the single participation setting.
31

---

Complete Proofs
In this section, we provide the complete proofs for our results from the main manuscript. For the convenience of the reader, we also restate the statements themselves.
F.1 Proof of Theorem
Theorem 1 (Square-Root of SGD Workload Matrix), Let Aug3 be the workload matrix ② Then A Ca for Co.3 = LDToee((o,..... ..,, with co and Cj = Ei=o a² -- for j = ,,..., with coefficients r; [(- For any p n.). the p-banded BSR matrix Carl is obtained from this by setting all coefficients Cj 0 for 2 p.
Proof. We observe that A can be written
1
0
0
a
A =
X
=: Ea
X
Eg
(20)
a-1
an
Bn-1
Bm-2
1
0
0
r1
Relying on the result from Henzinger et [2024], that E2/2
with
Tn-1
Tn-2
Tk [-1/9/
one can check that the square roots of the matrices E,, Eg are:
0
orı
Br1
E1/2_
E
a- Tn-1
n-2
Fn-2
3 -n-1
gm-2ra 2
(21)
An explicit check yields that these matrices commute, i.c. 1:1/ E3/² E.
1//2 E1/2
Therefore
0
Ca,3 = A
a,3 = -E1/2 x 1/2
C1
1
with Ck = Σ a rar - (22)
1=0
Cn-1
Cn-2
F.2 Proof of Theorem 2
Theorem 2 (Sensitivity for decreasing non-negative Toeplitz matrices). Let M LD'ooep(ma (moe ,,, be a lower triangular Toeplitz matrix with decreasing non-negative entries, i.e.
mo 2 mi m2 2 ... ...m-1 Then its sensitivity 2 in the setting of b-min-separation
2
1/2
4-1 Σ
-11.//
(M
M Mc.14500
Σ
mi-jb
(10)
0=
where M..1+ denotes the (1 +1+ jb)-th column of M.
32

---

Proof. Because all entries of M are positive, so are the entries of M M. Therefore, the condition is fulfilled such that 2 holds with equality, and
( & AEE
nax
Σ
EMM'M)
max [T for [(; ΣE relly
(M M(-1))
RET jEn
(23)
where Igg { C (1, (I,...,, n |-| | c
We now establish that {1, 1 + b,..., 1+ (k ( 11-1 is an optimal index set, which implies the statement of the theorem.
b) v
To see this, let n- be any optimal solution and let 1" E x be a column index that is not as far left as possible, i.e., x²) (i") (i" - 1} would be valid index set in k If such an index i" exists, split n" into the left indices, which are smaller than i", and the remaining. right, ones: x = = Um, R with = = {i E A < )"), T2 = fi E A i f"). Then, we construct new index set in which the left indices are kept but all right ones are shifted by one position to the left: ''' = x U T
with F A {i 1 | T By the condition on i", we know n' IIg.6.
We now prove that ['T'', x') 2 [Tx, r), so n' must also optimal. First, we observe two inequalities: for any i,j 1:
(M 1]. M M = (M,,), Mori) + mm mm-st1mm-5tt 2 (M..)). M M[3))
(24)
and for i 2 1.j > 1:
(Mi.g), Mi-s-11 = Σ Ma.g. MiJ-11 =
my -mm Σ
(25)
1=1-1
Σ
m -11-;1111-+1
My mm-1m0-)
(26)
1=j-1
2 m m
(27)
2 Σ m--Mi--
Ma.g. Mp.1 = (M Maa))
(28)
where the last inequality holds because by assumption mi-i-1 2 mi- for 1 2it1 1.
Now, we split f x²) and f x²) into three terms: the inner products of indices below 1", the ones of terms above i" and the ones between both,
f(=",=²)~ = f( + ( TR) + ( TR). f(7', x') = f( + f( RR FR) + 2f -i) The first term appears identically in both expressions. The second term fulfills
(29)
(30)
f( f( Fi) = Σ
(M M..a) = Σ (Mi.s--), M
(31)
J R
Σ
(M.. Mi-s1) = (
(32)
by Equation 24 The third term fulfills
f(TiL F = Σ Σ (Mr.a) Mi--1) = Σ Σ
(M MM M
(33)
iEN Fe
iETī JEnA
2 Σ
Σ (Ma.)) MiJ11) = ( TT))
(34)
NEW JE
by Equation 28
combination, this establishes [( x') 2 f(x=" m²), and since "t was already
optimal, the same must hold for T'.
33

---

Using the above construction, we can create a new optimal index sets, x², until reaching one that does
not contain any index i" described anymore. Then x* = {1, 1-bb..... ..... 1+ (l-1)6) must hold. If 1 = k. the statement of Theorem i confirmed. Otherwise, n = {I...., + (k 1)6)
1)63 for some E N
superset of x so because of the positivity of entries, ((r',x') 2 f(xx", ==) must hold. Once
again, because x² was optimal, the same must hold for x²' which concludes the proof.
L
F.3 Proof of Corollary L
Corollary 1. The sensitivity of the p-BSR for SGD can be computed using formula (T)p.
Proof. From (T) we know that Ca.3 is a Toeplitz matrix with coefficients (1, C1, .....,.......................... Cn-1) 1). where
CJ o a'r for 0 B S 1, with r; 1(-13) 3t, where BN (
26
is the
i-central binomial coefficient. It suffices to show that cj 2 cj+1 for any j E (1, ,.,...,
First, we show for the r; coefficients:
2i
1
(2
(2i)!
(2i + 2)
Ti Ti+1
4' i
4 (i+1
= 4' i!!!
4 (i+ 1)1 (i+ + 1)!
(35)
1 (2i))
(1
4 (i+ ) =
1 (2i 2) (2i +
(1-
28 + 1
(36)
4' i!
2(i + 1)
. C++1 1
(37)
2(i 1)
where Cj = J
1 (2) is the j-th Catalan number.
J+I
Now, we study the case = 1. If B 0, then c = Gn 1, so monotonicity is fulfilled.
Otherwise, i.e. 0 < 3 < 1, we write
Ck Ck+1
ri
- i
(38)
2=0
1
B,Ck_1gh gh-1
gx+1 BR+ gh+h1
(39)
4k
1h
k
B,Cx-11-+ BK+18
(40)
1
i=0
Using the classic identity between Catalan numbers, 2 2/ Bp BR+1, e.g. Batir et al. 2021 Identity 4.2] we obtain
1b
B
> 0,
(41)
4k+1
where the last inequality follow from the fact that 9 - B > 0 for each O,...,k k and any 3 1. This proves the monotonicity of Ck-
For a < 1, we observe that c = a Σ i=1 T:Th-s 7 for y= 3 < 1. Clearly, the sequence of is decreasing, and by the above argument, the sum is decreasing, too. Consequently, is the product of two decreasing sequences, it also decreasing, which concludes the proof.
F.4 Useful Lemmas
Before the remaining proofs, we establish a number of useful lemmas.
Lemma 2. For any C E Taan with CTC2 it holds for any b {1,...,m) n) that
senss,p(() = ||C'||,00:
(42)
where ||013,00 max ||C 2|2.
34

---

Proof. This follows directly from Theorem 2
max TE LJEE
Σ
max
Cl
max
(43)
Lemma 3. For any C' Rax" with CTC 0 holds for any E {l....,n) and lo E fi..... n)
that
k C < x S &
(44)
Proof.
We first show the upper bound. Observe that for any π C [n]:
Σlc Chu = Σ (9h))) C(a)) - Σ |
(45)
iJER
i,jE
ijET
= ②
([ s | Σ
|| ||]]| |
(46)
KE
iE
Therefore, using Theorem
sens 2
max
Σ IJE
S A| SICTC
(47)
rell,b6
For the lower bound, we introduce some additional notation. Let IIg be the set of 6-separated index 2
sets with exactly elements. Then, from Theorem
we obtain
max T J
ΣICTOh 2
max
ΣICTOI==
max S(x) 2 max S(T),
(48)
DEllibb RET
TETTY
=
with the notation S(1) Eier 11,.1112 for any index set /C (1....,
Now, we prove by backwards induction over = 1....:
max S(T) 2 - I| nlll
(49)
As base case, let 1s = 5. Denote by म, == + b, + 2b. + (n - -b}} for é = 1,... b the uniformly spaced index sets. By construction they all fulfill m, E 1ln/b and U T. [n]]. where the union disjoint. Therefore
max S(T) 2 Σ S(m,) = 1 s((ml) = 1 ICp... -
max S(z) wella/a
(50)
im b
This proves the statement 44 because 1 in this case.
As an induction step, we prove that if 49 holds for some value S n. then it also holds for k -- -121.
Let n* E argmax S(T)) and jr argmingen- S((i)), such that know that S(()) < S(()). Now, set n n"" }i Because n' E &14,1, it follows that
max S(T) 2 S(T') S( S(U)) > - relle-1
(51)
k
where in the last step we used the induction hypothesis. This concludes the proof.
L
Lemma 4. For Ca.9 as in and k: = 1, it holds that
< E(Ca,g, Ca.B) <
(52)
T=0
1=0
35

---

Proof. From Lemmas 3 and we obtain
g( Ca, s |ICa,slle scnst, & (senss./(Ca.s) s ||C12,,0
1 Vn
c
(53)
==0
where the last identify follows from the explicit form of Ca.a. The lower bound follows from
J
|la sensi, (onns,,CC..0) 1 |C1 72
Ea Ca.g)
(54)
Vn
and again the explicit form of ||C||/
Lemma S. For r, = =((-//) = t (3) it holds that:
1
1
1
ro = and
r1 2
and general
2Vj
Srz S
VTj
for j 2
(55)
Proof. The double inequality is a particular case of a more general pair of binomial inequalities when k j and
2j:
(k/mmm < 8k(m- - k)
m
(kH(m/m))
(56)
k 2 - k)
where (k/m) is the binary entropy function, with (1122) = 1. The proof of the general result (5). can be found in Mac Williams and Sloane [977 Chapter 10, Lemma 7, p309].
Lemma 6. Let cx = ai
as in I Then co = 1, and j 2 1:
a 2Vj +1 S S (1 8 VI
a
(57)
Proof. We exploit the upper and lower bounds from Lemma First, we write Ck k Σ}= rytk-j? with y = 2. Then we check immediately that Co = and C1 = [ + B) =
/ 2(1 (1+7) 7) S 1-
For j 2 2 we derive the upper bound by
f-1
9 =r (11++))) +
1 1+75
(58)
2
VTj
πl i i)
1+
V - 1
Vi - -
V +
TVj 73 +
(59)
Vīj
T
1=0
v
1 1
(60)
TVI
TVj - - B)
1)v5
+I
(1+
VJ+I
1
(61)
Vj
TVj
+V
VJ Vj-1 マVj +
11
S1
Sv3g2
3
(62)
πV T (1 (1-
which proves the upper bound aj- The lower bound for j 2 1 follows trivially from
0 2 2'rr, 2 ad
á
N 2
(63)
2
36

---

Lemma 7. For j E {l..... n it holds
log(A
+
1)
1-1
1+ log.
G
<
(64)
4
(1 3)
==
for a and otherwise
3-1 1
1
1
1
2
<
log
(65)
(a
3)2
1-
a-a2
==0
Proof. We first prove the result for a
1. Combining Lemmas
4 and 6 we obtain
j
1
1 + log
E
(66)
9)2
i+1
(1 (1-3) 3)3
(1 3)
1-1
1-1
log(i +1)
Σ
2
(67)
i
1=0
1=0
1=1
For < 1, if follows analogously:
1
j-1
26
1
a
1
log
(68)
1+1
(a
8)2
(a
3)2
a-
1
8
1=0
1=1
j-1
1
1
a
8
2i
8
a
(69)
i i+1
4a
4a
i
2(++1)
log
(70)
4
1 1
(j 1) Q2)
where the last term emerges from 2i=j+1 i
> (25112)
E
25512 1
jtI
J+1
Lemma 8. Let 8 < a S L. Let ... on be the sorted list of singular values of Au, If a < 1, then for j = 1,..., n:
(1+ + 3) So, S (1- (1- a) 3)
(71)
and
n S Aa,88l S (1- a1= - 3)
n
(72)
Ifa 1, then for I,...,n. n,
2 1 n
n
<aj <
(73)
11+8j-
1 ]
and consequently
2 (n + log(n + 1)
S ||A S (n+ + 1) + log n)
(74)
+ 3
1 + 3
Proof. The statements on the singular values follow from the following Lemma 9 because A
E. Es. Because Ea and Eg are diagonalizable and they commute, we have on (Egg), (Ea) S
(E E3) I o1 (EE)), (Ea). For < the lower bound follows from ||A 2 trace Aa, and the upper bound follows from the identity || E Tz.
For 1, the bounds follow from the same identity together with the fact that
log(n (n + 1) 5
< log(n) 1.
(75)
37

---

Lemma 9 (Singular values of Ed). For 0 st let E, = LDToep(1, tt................. t-- E R Then
the singular values 01 (E,) 2 2 on Ex,) fulfill for ,....,
1+t S (E.) - for 0 <t< 1.
and (i
(76)
sin
3
Proof. We follow the steps of SehastienB (017]. and use that the singular values of Fg are the recipro- cals of the singular values of E which themselves are the eigenvalues of ((E))) =: ie., for i = ,,...,7. 7
(Bi)
(77)
Am+1-,;(T)
The E," and can be computed explicitly as
0
0
1
-
0
0
1
0
0
1
1
-1
0
L
1
E
T
(78)
7.
1+ 12
-t
0
-t
0
-t
1+
Lemma 10. All eigenvalues, μ4, of fulfill
(1 -)
S $ (1+
(79)
Proof. By Gershgorin's circle theorem Gersshorri 1-t -tt S 1 t or ii) |
1931
know that μ fulfills i) |1 - /u| < t. i...
S 2f, i.e.
2t + S p S +1 + t2. For t E 0, 1] the
first condition implies the second, 79 must hold.
Case I: For t < 1, the statement (6) follows from Lemma 10 in combination with 77
-1
0
LO
-1
2
0
Case II: Fort 1 the matrix simplifies to T
0
:
:
Note that T is not
0
2
-1
2 exactly Toeplitz, because of the top left entry, so closed-form expressions for the eigenvalues of tridiagonal Toeplitz matrices do not apply Instead, we can compute its eigenvalues explicitly. Matrices of this form have been studied by Elliott [953]] for completeness, we provide full proof
here.
Vn-1). From the eigenvector equation Let be an eigenvalue of T with eigenvector = =( (Po,.... TV D obtain
V V/ V
(80)
/4 T 4/4 = 1 4-1 + 2 Mx41
for k: 1,..., n
(81)
n n-1 Vn- 2T n-
(82)
which yields linear recurrence relation
Ux+1 = (2 1 k- for 1 1.....
(83)
with two boundary conditions
= (1-1)N0
(84)
, = (2- 11) "
(85)
We solve the recurrence relation using the polynomial method Greene and Knuth
1990] The
characteristic polynomial of (83 P(z) = =² (μ - 2) + 1. its roots are
2 1
p2 - i m)
r
10
(86)
2
38

---

for some value 0 E [0, 2T).. Note that the expression under the second square root is positive, because
of Lemma 10 The last equation a consequence of, |r_|² = ((2 14)² - 11)
(2 (2-1)))) = 1.
Consequently.
= 2 28((e")) = 2 cos A.
(87)
From standard results on linear recurrence, it follows that any solution to 88 has the form V aa (r+)² + C2 (r- Y for some constants c1, C2 E C. The fact that V must he real-valued implies that C = C2 =: ae for some values a E RR, Φ E [0, 2x). Dropping the normalization constant (which
we could recover later if needed), we obtain
V ((++)0) + ((++5)) 2 cos( ( + j0).
(88)
Next, we use the boundary conditions to establish values for o and 0. Equation (85) can be rewritten as
cos( Ø + - 2)0) = 2 cos(0) cos((n 1)0)
(89)
which, using 2 cos(a + B) = cas(a + b) + cos(a b), simplifies to
0 = cos( و + n0)
(90)
Consequently, Φ + n0 = उ + kn must hold for some E N.
Equation 84 can be rewritten
cos( 0) (2 cos 1) (o)
(91)
which simplifies
cos(Ø) = cos(B - ②)
(92)
One solution to this would be A = 0, but that would implies μ = 0, which is inconsistent with T being an invertible matrix. So instead, it must hold that e = + ko for some k E N.
Combining both conditions and solving for 0 we obtain
2 + π 2"
0
+ ko for some & E N.
(93)
n+
Each such value 0% for E fo.....n 1} yields an eigenvector with associated eigenvalue f 2 - cos 0g = 4 ssn (0x/2), Now, (76) follows from this in combination with (77).
F5 Proof of Theorem
Theorem 3 (Expected approximation error with single participation). Let Aa, E R be the workload matrix (9) of SGD with momentum 0 < B 3< and weight decay parameter < a where a > B. Assume that each data item can contribute at most once to update vector (e.g. single participation, k = 1). Then, the expected approximation error of the square root factorization, An.8 C² fulfills
< E(Co, Ca.p) I (a 3)=
log
(11)
1- az
for a < 1, and
log(n + +1 - < E(CC,,,, C1,0) <
1 + log(n)
max 1,
(12)
4
B)²
Proof. The proof consists of a combination of Lemmas and k 1, so need just the first column of matrix Ca,
Because in the single participation
1 log
for a 1,
n-1
(1 3)2
E(Ca,g, Ca.s) <
Σ
<
(94)
1
1=0
log
otherwise.
(a 3)
1 o²
39

---

which proves the upper bounds. For the lower bounds, for any a S 1:
CCaa,,,Ca,ol 2
cz 00=
(95)
j=0 i=
j=0
Also, for a = 1:
1 7
log( (n + 1)1)
log(n + +1 - 1
E(C1, 2 π
log(j + +1)
(96)
4n 5=1
4n
4
because log( (n + 1)!) 2 (n + 1) log (n + 1) - n.
F.6 Proof of Theorem
5
Theorem 5. Assume the setting of Theorem and Aa.s = Id Aa, fulfill, for a < 1,
Then, the baseline factorizations Au,3 Aa, IId
V +B
Vni + +olvn)
E(Aa,s,Id) = v(1- aBl -
+ +1) and E(A1,8, -
B2)
v2(1 - B)
(14)
+B+
Vā
E(Id, Aa,s)
+ o[1 and E(Id, A1,s) S 1- 3
+ olVn).
- B2)
(15)
Proof. For a = 1, by Lemma 2 we have:
cens),y(A,,o)) =
n-1
ag
n-1
Bi+1
2
n
2/3
n-1 1
g2
n- 1
325
3
(1 - 3)2 (1 - 3)
(1 3)2
1=0
i=0
i=0
(97)
2
g2n+2
n
+
n
(1 + o()). (1))
(98)
(1 3)2
(1-
(1- (1-1) - B2)
(1- 3)2
For < 1:
a-1
n-1 (a - git222
scnst,
(99)
(a - 3)2
1=0
1=0
n-1
(a0) + 82 n-1
(a²)- - 2aB Σ
> (89))
(100)
a -
B
1=0
i=0
==
a²
20,3
B2
+
(1+
(101)
(a 8)2
- a²
- a3
1- 32
1 aB
(1 (1+00 (1))
(102)
(1-a8)(1- - (1 (2)
Together with
2
1
1 18+1 2
||
/
n j)
8
(n
(103)
1
j=0
1=0
1=0
1
Σ
n-1
(n j)(1 - 28 + 325t2)
(104)
n1l - 3)2
1=0
(n+1)
+ O(L) 2(1 B)3 (1 + o(1))
(105)
2(1 -1)
40

---

3 < 1 and the sum Ej= j8st1 1 is uniformly bounded by
y=== ]3
8²
-1-03z
For < 1:
1
Σ j=0
nn - j) (alt1- - -35t112
n-1 (a--1- -3
+ o(1)
(106)
(a - 8)2
(a 3)2
==
1 1+
(1 o(1)
(107)
(1- (1- a3) B2)
where the second equality due to the fact that we average over the sequence ju which converges to for [ar| 1.
F.7 Proof of Theorem
Theorem 6 (Approximation error of BSR). Let Aa, & Ruxa be the workload matrix
9 of SGD
with momenium 0 < 3 < 1 and weight decay 0 < S 1, with a >> 3. Let A = BB IPp be
a.3 a.
its banded square root factorization as in Definition Then, for any b E {1...., n.. S b and
E (1....; 8 holds:
O3 (v nk logg
+ 03.p(VTi)
for a = 1,
E(Blr) a,8 clel)
(16)
OB
(VT))
for a < 1.
Proof. Consider a Lower Triangular Toeplitz (LTT) matrix multiplication:
a1
0
b1
0
0
C1
0
a22
0
X
62
61
0
=
C1
(108)
On
bn
b
b3
Ca-1
where Cj = El= a 1-i so the LTT structure is preserved with multiplication that allows to work with sequences and their convolutions rather than matrix multiplication. For instance, we would write the previous product in the form:
(a1...... an) (b1..... b,)) (c1 Ci....,,a!) Cn).
(109)
The inverse of the Lower Triangular Toeplitz matrix remains a Lower Triangular Toeplitz (LTT) matrix because we can find a unique sequence (ci,...,C,) such that:
(C1, (ci,..., Cn) x (a....., ....................))))..)..... = (1.0....,0)
(110)
0g = 0 and 61
(111)
111
with the restriction that a₁ * 0; otherwise, the original matrix was not invertible. We consider the banded square root factorization Aa,3 = Blpl C pl which is characterized by the following identity: a, a,
(bu...., Cp-.
0) = (1, cn-1) (1, .t..c.....C.. 1n-1).
(112)
We will bound the Frobenius norm of the LTT matrix (bo,... bn- By the uniqueness of the solution, we obtain that for the first p values we have b, = e, For the next p values we have the following formula:
41

---

bpts +... bpc; + ... + bj+1cp-1 = Cptj + ... ....... cpt1cp-11111 Cpts
(113)
Opts 101 +.........00 bpbj = (cpts +... + Cp(z).
(114)
By induction argument, can that bpt. 2cptj for 0 S values we will prove convergence to constant.
For the remaining
2p
p-
4- = 0g = a++1 35+1
Σ
(115)
==0
We make ansatz for the solution of the linear recurrence in the form:
α (a - B) E qa (a- 3) E c3
+ yo
(116)
where Vj represents the error terms, which will be proven converge to O. The sequence yj satisfies following recurrence formula:
p-1
Fo - yz Vi-
Σ
(117)
i=1
We denote wy c)a") which is decreasing sequence because the values correspond to the C 1/3/a matrix. We rewrite the recurrence in matrix notation:
W1
W2
W3
Wp-2 2
W
yk-1
Ue
0
yk-2
yk-
0
1
0
0
yk-3
yk-2
(118)
0
0
3k--1
yk-3
...
...
...
0
1
yk-b/
W/-b+1/
To show that the error terms Y! goes to as j goes to infinity, we first study the characteristic polynomial of the associate homogeneous relations:
((A) x- wi A0-
Wp-2A Wp-1-
(119)
Wh-1 0, follows from Schur's (relaxed) stability condi- tion INguyen et al. 2007 Theorem 1] that all its (complex) roots lie inside of the open unit circle. Therefore, all solutions the homogeneous relation converge
Because
Wh
W2
rate exponential in j and y = o(1 and Σラ= v/ 3.,3,p (1). Then we can bound the Frobenious norm of the matrix B!nl as:
a.
1
25t2
|
pl a.
ta y
(120)
1=
1=0
(a
3)2
2
We use the following lower bound for the sum of wy:
Σ E. 2 プ -
-
(121)
1=0
1=0
Combining these bounds we can upper bound the Frobenious norm of the matrix Blpl a.3
the following
way:
42

---

a²
(a - 3)2 log
(vp+ - a) 8)2 +0,
for <
[
Paglll//n
1 1- a²
+ log(p)
n n-p
(1 () 1)2 + pp.s(1)
for 1.
3)2
(122)
Simplifying for the leading terms in asymptotics, we have:
(a.
for
[
a. /1/.1////
(123)
OB
pp,s((1) for 1.
Sensitivity of Cl for any p S b:
For the b-min-separrtion participation sensitivity have the following bound
k
log
for a 1
p-1
(a - 3)2 1 + log(p)
1 a²
L n. < k
<
(124)
k: (1 3)²
for 1.
Combining sensitivity with the upper bound for the Frobenious norm we obtain:
for <
((BlP) clPl) a,9:
Ogg (v
nk log
(125)
(vk) for 1.
F.8 Proof of Theorem for Square Root Factorization
Theorem 7 (Approximation eтror of Square Root Factorization). Let Aa,9 Raxx be the workload
matrix of SGD with momentum S < and weight decay < a S with a > B. Let
Aa.8 0 a.
be its square root factorization. Then, for any b E (1....,m) n. and k = it holds:
E(Ca,gs Ca.s) = + (kVlogn + vk Vklogn)
for a = 1,
(17)
(a,o((V) (Vk)
for a <
Proof. We prove the case without weight decay (a = 1) and with weight decay (a < 1) separately. Case 1) no weight decay (a 1).
We start by bounding the b-min separation sensitivity:
k-1k-1 (Cr.))
Σ
((C1,8))..0). (C1)
(126)
2=0
Consider scalar product for general pair of indices, j j>{ i:
n 1-j
((C1, (C1,a) =
Σ
-0g-i4t-
(127)
1=0
43

---

Using the bounds on C 66 for 1 we can lower and upper bound this sum by:
1
dx
c2g 5 (1- 3) t=
(128)
Vt(j i+
(1-
3) 0
Vxlj i+
1 n-j
dæ
Σ
(129)
Vij -
vxlj - (
We can compute the indefinite integral explicitly:
dr
F
(130)
for F(a) = 2oo (v + +1 /2). In combination, we obtain the upper and lower bound for (26):
{1 -0-
j
$
- - (C0)).CCC))) (1 31-
(131)
n
Now we are ready to bound the sensitivity of the matrix C1,3:
1-1
k-1 k-1
Σ
((C1,))ib, (C1,8)ib) +2
Σ Σ ((Ci,)) (C1,s)3b)
(132)
1=0 j=it1
2 k-
$ (1 3)²
(log(n ib) + 1)
(1- (1-
ΣΣ
k
(133)
and, analogously
k sens 2 log(n n ib) +
Σ
(+==))
1
f(blj ((6(-11) m))
(134)
0
Ugi<jsk-1
Firstly, using (5)a kd kk we bound the sum of the logarithms:
log(n -j) log log k llog + log k! k klog + k log k k log n,
(135)
k-1
Σ log(n jbb = k log b + log k! 2 k log b + k log k - h = k loog k.
(136)
3=
To upper bound the last term in sensitivity lower bound 134
we use the auxiliary inequality
f(a) 2 log
(v 1 +1+ <
to derive:
2 k-
k-1
Σ 5(bi 1)) <
Σ
<
<
8 3
3/2
2 Osicisk-1
vb osicisk-1
-5-
Vb
Vt
vb
F=
1=1
(137)
To bound the final we establish the following inequalities for f(a):
f(a) = 2lo V v2 = log { +1) + 2log
(1+
(138)
Va+l
44

---

log (1+1)
< f(a) < log
- (1)
+ 2 log
(139)
Then we can bound the first double sum in [34 as
log
k - i (1:1)
S Σ
k k-
j 1
log
k (k-1) (5-i)
+ 222 log 2. (140)
Osicisk-1
0sicisk-1
Osicisk-1
To bound the term Losiejsh11 log
j-i) we use the following identities:
Σ log
k -
log
k-
(k - -))
(141)
j-i)
(k: 1)!
0giessk-1
i=0
k*
(I: - 1)8-2...2
= loog
= log 1k-1. 2k- ... (k. 1)
(142)
1! 2! ..... (k 1)!
k-1
j ++1 5
k-1
log
Σ
jog(í + - Σ
jlog(k: - j)
(143)
h: j
/=1
j=1
j=1
k-1
= Σ [j 1) log(j) - Σ
(k-
log(i)
(144)
==
==1
k-1
2
Σ
log(j) - log k! + 2k log log k!
(145)
j=1
Now, using that x log x is a monotonically increasing function,
Σ
log
-
s 2
L k
log rdr + k log k + k 122 log ts + 12
(146)
5-
0giejsk-1
k.
= 122 log
2 + k log k - k 1.2 log h: + 222
(147)
3 k²
(148)
As lower bound, we obtain
k: - i
k-1
log
2
log vdv - h: log + klog k k(k - 1) log(( - 1) (149)
5-i
Osicisk-1
11
(k 1)2
= (k 18-1) log - 1)
+ klog - k(k 1) log(A: - 1)
(150)
2
(k-
1)2
= -(k - log 1)
+ k log
(151)
2
22
(152)
Therefore, combining the upper bound (16) and the lower bound 119 yields
Σ
f k j
S (2 log 3/2) < 312
(153)
OStcisk-1
Σ
(
(j-i)
(2 log 1/2)
4
(154)
k
OSicjsk-1
Combining all three terms together we obtain the following bounds for the squared sensitivity 134)
45

---

s - (1-By (log +1) +
(1 - 8)2 k
(155)
k
8 k
2 k²
(C z I (logn - -1 11 3vb
(156)
Now, we recall the bounds for the Frobenius norm of the matrix C1,8
7 and (96
log(n + 1) -1 - ||C1,0ll//m
log n2 +
(157)
4
(1- 3)2
With the auxiliary inequality / + V S vvabbs + b and combining 155 and the bounds
on Frobenius norm and (96) we get that:
8(Ci,8, (13)) S (log + +)+ v5k
Vk
logn+
(158)
(1
B)2
(1 - 3)2
1
k
20
E(C1,8, C1,))
vi(logn - 1)
log(n) 1
(159)
4/2
2/5
3/n
Making the lower bound well-defined requires 2 45, otherwise one can simply take E(Ci, C1,s) 2 As a final step, we combine both inequalities the following asymptotic statement:
E( CI O (Vklogn kvllggn)
(160)
which concludes the proof of the case without weight decay.
Case 2) with weight decay (a < 1).
As above, we first express the 6-min-separation sensitivity of the matrix Ca in terms of inner products.
R- A-1
s = ΣΣ
(Co,,)) ()
(161)
==
and then consider scalar product for general pair indexes j i:
2-1 1-1-3
((Co,,))s (Ca,
Σ
cCj-itt-
(162)
Now, the bounds c from Lemma for < 1, upper and lower bound this sum with the following expression, where 7 8
aj-i n-1-j
2h
a -
((Ca,)),(Ca,d)))
71
V +
-1011)--+1 i+++1
(1-7)?(
a²)
I=
(163)
2
a
j-
(CC.s))
(164)
(1+1)0-i+(+)
4Vj i+
We substitute these bounds into Equation matrix Ca,p:
161 to obtain the following upper bound for sensitivity of
46

---

2
ab--))
sens ( (C, - (CCa, (Ca,plib)
(165)
(1- a²)vb
vJ
j2i
log
+
2
Σ
((5-1)
(166)
(a - 3)2 a
(1- 7)- a²
521
k
1
2ka
(a - 3)2 log
+
(167)
a² (1- (1-
where the second inequality is due to Equation (95).
A lower bound for the sensitivity follows directly from Lemma 3
sensi, (Ca,,a) 2 xlI 2
(168)
The Frohenius norm of the matrix Ca.8 is the same as that for one round participation; thus, we could reuse Inequalities (95
1 log (a: 3-2 1 -22
1 S ||Ca,ell//n
(169)
By merging the bounds for sensitivity and Frobenius norm, we derive the following bounds for error:
1 log (a - 3)2 1 a² (1- a²) (1--)(1-33(1-0)) Vb
2a
E(Ca, Ca.o) S Vi
(170)
E(Co,,, Co.8) 2 Vk.
(171)
The combination of these results yields the following asymptotic statement:
E(Ca, Ca.s) = Da,g(Vk),
(172)
which concludes the proof.
F.9
Proof of Theorem
Theorem 8. Assume the setting of Theorem 6
Then, for any factorization An,9
BC with
CT C 0, the approximation error fulfills
looo
for a =
E(B,C)
(18)
Vk
for a <1
Proof. Let Aa,s = BC be any factorization with CC 2 0. From Lemma
3 it follows that
E(B.C) = Vn _ ( Vk 11
Vk n
AA,,oll.
(173)
where -|. denotes the nuclear norm, and the last inequality follows from its variational form, ||||||] The statement of the Theorem follows by inserting the
||
|+ min (X,YXXXXXXXXXA
corresponding bounds on |A,,, from Lemma 8
F.10 Proof of Theorem
Theorem 9. Assume the setting of Theorem Then, the baseline factorizations Aa, Aa, Id and Aa,s = Id Aa,8 fulfill
k
nk
for a = 1.
E(A,,8, Id) 2
2
for a = 1,
E(Id, An,a) >
V3
(19)
Vk
for I.
V for a 1.
47

---

Proof. Case 1) Aa.s = BC with B = A and C = Id.
It is easy to check that sens
Vk, so E C) V ! An, Because Aa, S A I
a-3 An, componentwise, have for
a=
n(n 1) ||A1,o1f 1| 2 (1- B)2 A A1.01
n(n + 1)
(174)
2
2(1 3)2*
which implies the corresponding statement of the theorem. For 0 < <a < 1, we use that A1,0 > Id componentwise, so c,,0? > ld [| n, which conclude the proof of this case.
Case 2) Aa, BC with B ld and C' Aa.3. We observe that ||d || vn, so
EBB C) (en,s,
(175)
Again, use the fact that Aa,o A a- Aa, Now Aa,o fulfills the conditions of Theorem so from Equation (10) know
2
k-1
(0n8k,6(Aa.))) = Σ 11 j=0
(Aaa)
(176)
We first study 176) for a = 1. Then, from the explicit structure sees that the vectors inside the norm have block structure
11,0 LDToep(11 1, 1,..., 1) one
U1
(110)
with
02
R
(177)
==
Uk
for I,..., ki and v'
E nR=-bk appears only if B.
Now we check
||011 +
||u:11 1 = (n - bk)k² +
12
(178)
k=0
i=1
k:(k: + 1)(2 +
+1)
nk?
bk 3 + b
(179)
6
2 nk.2
b 2 nkk
(180)
3
because bk 2 n. Consequently
seriss senss, (A1,8) 2 l/v
(181)
v3"
which concludes the proof this case. For a < 1, An, 2
componentwise readily implies
senst,6(Aa,o) 2 Vi.
(182)
which implies the statement of the theorem.
48

---

NeurIPS Paper Checklist
1. Claims
Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
Justification: Our abstract and introduction clearly states our claims, the contributions made in the paper, the main assumptions and limitations.
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made in the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made the paper and important assumptions and limitations. A No or NA answer to this question will not be perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalize to other settings.
It is fine to include aspirational goals as motivation as long as is clear that these goals attained by the paper.
2. Limitations
Question:
Does the paper discuss the limitations the work performed by the authors?
Answer: [Yes]
Justification: We discuss the limitation of the approximately optimal matrix factorization methods. We emphasize that the bounds for the banded square root factorization are asymptotic and should be perceived together with the numerical experiments.
Guidelines:
The answer NA means that the paper has limitation while the answer No means that the paper has limitations, but those are not discussed in the paper.
The authors are encouraged to create separate "Limitations" section in their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be.
The authors should reflect the scope of the claims made, e.g., if the approach was only tested on a few datasets or with few runs. In general, empirical results often depend implicit assumptions, which should be articulated.
The authors should reflect the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution low or images are taken in low lighting. Or speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon.
The authors should discuss the computational efficiency the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness.
While the authors might fear that complete honesty about limitations might used by reviewers as grounds for rejection, worse outcome might be that reviewers discover limitations that aren't 't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play impor- tant role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.
3. Theory Assumptions and Proofs
Question: For each theoretical result, does the paper provide the full of assumptions and complete (and correct) proof?
49

---

Answer:
[Yes]
Justification: All lemmas and theorem state the necessary assumptions. We provide proof sketches in the main text and proofs for all results in the appendix.
Guidelines:
The answer NA means that the paper does not include theoretical results. All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated or referenced in the statement of any theorems.
The proof' can either appear in the main paper or the supplemental material, but if they appear in the supplemental material, the authors are encouraged to provide short proof sketch to provide intuition.
Inversell, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix or supplemental material.
Theorems and Lemmas that the proof relies upon should properly referenced.
4. Experimental Result Reproducibility
Does the paper fully disclose all the information needed to reproduce the main ex- perimental results the paper to the extent that affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided or not)?
Question:
Answer: [Yes]
comparison, which we consider be the main experimental results, we have explicit descriptions of the procedures together with the code included the appendix. For the model training, we adapted the public code from prior work (Choquette-Choo et al. 2023c and we describe the changes in the appendix.
Justification: For the factorization
Guidelines:
The answer NA means that the paper does not include experiments. If the paper includes experiments, No answer to this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
If the contribution is dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is novel architecture, describing the architecture fully might suffice, or if the contribution is a specific model and empirical evaluation, it may be necessary to either make possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way to accomplish this, but reproducibility can also be provided via detailed instructions for how to replicate the results, access to a hosted model (e.g.. in the case of large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed.
While NeurlPS does not require releasing code, the conference does require all submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) If the contribution is primarily a new algorithm, the paper should make it clear how to reproduce that algorithm.
(b) the contribution is primarily new model architecture, the paper should describe the architecture clearly and fully.
(c) the contribution is a new model (e.g., a large language model), then there should either be way access this model for reproducing the results or way to reproduce the model (e.g., with an open-source dataset or instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome describe the particular way they provide for reproducibility. In the case of closed-source models, it may be that access to the model is limited in some way (c.g., to registered users), but it should be possible for other researchers to have some path to reproducing verifying the results.
50

---

5. Open access to data and code
Does the paper provide open access to the data and code, with sufficient instruc- tions faithfully reproduce the main experimental results, described in supplemental material?
Question:
Answer: [Yes]
Justification: We provide the code for the approximately optimal matrix factorization. The code for model training is adapted from the publicly accessible code, above.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurlPS code and data submission guidelines (tttts::/nnis..cc)
for more details.
While we encourage the release of code and data, we understand that this might not be possible, so No an acceptable answer. Papers canmot be rejected simply for not including code, unless this is central to the contribution (c.g., for new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurIPS code and data submission guidelines (https for more details.
/nips. htt.lii/gguddeeeddddeaaaaiiisees Coodeeummiseioooooooooo
cumbicc/guided
The authors should provide instructions on data access and preparation, including how access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. If only subset of experiments reproducible, they should state which ones are omitted from the script and why.
At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable).
Providing as much information as possible supplemental material (appended to the paper) is recommended, but including URLs to data and code is permitted.
6. Experimental Setting/Details
Question: Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
Justification: We use the same setting as Choquette-Chon et al. [003cc and report to have the same hyperparameters for the training.
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core of the paper to level of detail that is necessary to appreciate the results and make sense of them.
The full details can be provided either with the code, in appendix, as supplemental material.
7. Experiment Statistical Significance
Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Question:
Answer: [Yes]
Justification: we report error bars and indicate that experiments were performed multiple times to ensure reliability and statistical significance.
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by emror bars, confi- dence intervals, or statistical significance tests, least for the experiments that support the main claims of the paper.
51

---

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.
is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have 96% CI, if the hypothesis of Normality of errors is not verified.
For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out range (e.g. negative еπτог rates).
If error bars reported in tables or plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables text.
8. Experiments Compute Resources
For each experiment, does the paper provide sufficient information the com- puter resources (type of compute workers, memory, time execution) needed to reproduce the experiments?
Question:
Answer:[Yes)
Justification:
Our experiments run on commodity hardware.
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate type of compute workers CPU or GPU, internal cluster, cloud provider, including relevant memory and storage.
The paper should provide the amount of compute required for each of the individual experimental runs as well as estimate the total compute.
The paper should disclose whether the full research project required more compute than the experiments reported the paper (e.g., preliminary or failed experiments that didn't make it into the paper).
9. Code Of Ethics
Question: Does the research conducted the paper conform, in every respect, with the NeurlPS Code of Ethics
Answer: [Yes]
Justification: We have acquainted with the guidelines and preserve the anonymity of the authors in the paper.
Guidelines:
The answer NA means that the authors have not reviewed the NeurlPS Code of Ethics.
If the authors answer No, they should explain the special circumstances that require deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g., there a special consid- eration due to laws or regulations in their jurisdiction).
10. Broader Impacts
Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?
Answer: [Yes]
Justification: Our work is about model training while preserving data privacy. We mention this in the manuscript.
Guidelines:
The answer NA means that there is no societal impact of the work performed.
52

---

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), faimess 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 tied to particular applications, let alone deployments. However, if there is 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, 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 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 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).
11. Safeguards
Question:
Does the paper describe safeguards that have been put in place for responsible release of data models that have high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Answer: [NA]
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 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 best faith effort.
12. Licenses for existing assets
Question: Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?
Answer: [Yes]
Justification:
We have cited all the papers which have results and code used in our paper.
Guidelines:
The answer NA means that the paper does not use existing assets.
The authors should cite the original paper that produced the code package or dataset.
The authors should state which version of the asset is used and, if possible, include URL
The name of the license (c.g., CC-BY 4.0) should be included for each asset.
For scraped data from a particular source (c.g., website), the copyright and terms of service of that source should be provided.
If assets are released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, wwwpppwwthhhooeeeeooee com/ com/datasets
53

---

has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
For existing datasets that are re-packaged. both the original license and the license of the derived asset (if it has changed) should be provided.
If this information is available online, the authors encouraged to reach out the asset's s creators.
13. New Assets
Question: Are new assets introduced in the paper well documented and the documentation provided alongside the assets?
Answer: [NA]
Guidelines:
The answer NA means that the paper does release new assets.
Researchers should communicate the details the datasct/code/model as part of their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent was obtained from people whose asset used.
At submission time, remember to anonymize your assets (if applicable). You can either create anonymized URL or include an anonymized zip file.
14. Crowdsourcing and Research with Human Subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well as details about compensation (if any)?
Answer: [NA]
Guidelines:
The answer NA means that the paper does not involve crowdsourcing research with human subjects.
Including this information in the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then as much detail as possible should be included in the main paper.
According to the NeurlPS Code of Ethics, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector.
15. Institutional Review Board (IRB) Approvals or Equivalent for Research with Human Subjects
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or equivalent approval/review based the requirements of your country or institution) were obtained?
Answer: [NA]
Guidelines:
The answer NA means that the paper does not involve crowdsourcing research with human subjects.
Depending the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
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.
54


Reviews:
[{"review_id": "GA9L1sYGZy", "reviewer": "Reviewer_Cfm2", "summary": "The paper considers the problem of adding correlated noise $C^{-1}z$ instead of independent noise $z$ across iterations in continual counting (equivalently, DP-SGD). Past work gave an algorithm to choose $C^{-1}$ that optimizes some objective on the noise under b-min-separated participation, but this algorithm requires solving an SDP which is infeasible when using a large number of iterations, i.e. optimizing over $C$ that has a large number of rows / columns $n$. The authors propose a choice of $C$, which they call banded square root, which has a succinct representation that is efficiently computable in time independent of $n$. Furthermore, their choice of $C$ is banded, i.e. is non-zero below the $b$-th diagonal for $b \\geq 1$, which gives it some nice properties such as being able to compute $C^{-1}z$ efficiently in a streaming manner. To define the banded square root $C$, they authors take $A$ which represents the workload, effectively a representation of the updates in SGD which can include momentum and weight decay, and take its matrix square root. They then truncate the matrix square root to its first $b$ diagonals. The authors give an explicitly and efficiently calculable formula for the computing the first column of this matrix, which specifies the whole matrix as it is Toeplitz. The authors analyze the asymptotic error guarantees of different choices of $C$, including the banded square root, in the single- and multiple-participation settings. Specifically, their error guarantees show that banded square root improves on $C = I$ or $C = A$, and nearly matches a lower bound on any factorization if using enough bands. The authors conduct experiments training models for classification on CIFAR10 and show that banded square root is competitive with the choice of $C$ given by solving a more expensive SDP.", "questions": "* wrt the first weakness, did the authors consider this alternate approach / if so, do you believe banded square root still offers speedups in this setting?\n* The previous work of Choquette-Choo et al. shows that banded $C$ are compatible with privacy amplification with sampling. Have you considered combining your error analysis with their privacy analysis (e.g., maybe for RDP to simplify) in the setting where batches are sampled in DP-SGD? It might be an interesting question to see how the number of bands affects error in the presence of amplification, since they observed more bands reduce the benefits of amplification.", "limitations": "Yes", "rating": 6, "confidence": 4, "soundness": 4, "presentation": 4, "contribution": 3, "strengths": "* The work makes theoretical progress on a practical problem. Correlated noise/DP-MF is now being used in practice to train models with DP, and especially as model sizes and training runs get larger, more efficient ways to implement DP-MF will lead to better models in practice. In particular, itmakes the results of the DP-MF literature more accessible to those who do not have large amounts of compute to use the techniques they introduced.\n* The theoretical analysis of the error of DP-MF as a function of the number of bands is novel and adds theoretical understanding of banded $C$ that was not present in the previous paper of Choquette-Choo et al. Furthermore, the set of theoretical results is rather extensive and gives a pretty complete theoretical understanding of the authors' methods.\n* While DP-MF is a relatively niche topic, the paper does a good job slowly introducing the problem, their approach, and their theoretical results.", "weaknesses": "* While the previous work did require solving an SDP, to my understanding this is for the case where $C$ can be is only required to be a banded lower-triangular matrix. If $C$ is required to be Toeplitz as well, as the banded square root factorization is, then one only needs to optimize over $b$ variables instead of $\\approx nb$ and it is not clear that the computations in the previous work are expensive, which mitigates the improvements in this paper.\n* It is worth pointing out there is a work of https://arxiv.org/abs/2404.16706, that also gives an efficient-to-compute way to choose $C$, and leads to faster streaming noise generation for than banded matrices. This work I would consider concurrent, so I did not account for any possible overlap between the two papers' results and impact when assigning my score.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for the insightful review. Here, we address the individual concerns:\n\n**\\>  Relation to arXiv:2404.16706.**\n\nIndeed, this interesting preprint (which we cite as [Dvijotham et al., 2024]) is concurrent work. It studies only the case of MF-SGD without momentum or weight decay (“prefix sum”) and is limited to streaming data ($k=1$, “single epoch”).  This setting enables a deeper analytic treatment (in particular the sensitivity computation, which e.g. does not benefit from bandedness), but it comes at the cost of covering fewer real-world scenarios. \n\n**\\> While the previous work did require solving an SDP, to my understanding this is for the case where can be is only required to be a banded lower-triangular matrix. If is required to be Toeplitz as well, as the banded square root factorization is, then one only needs to optimize over variables $b$ instead $nb$ [...]**\n\nThis is an interesting observation that could lead to a hybrid between [Choquette-Choo et al., 2023a] and our work. As a compromise between speed and optimality, one could keep the objective of the problem (4) but restrict it to Toeplitz matrices, thereby parametrizing the problem with $b$ variables. Unfortunately, the objective function would be quite more involved. It would no longer have a standard structure, such as the “matrix fractional function” form of the SDP (4), and it would not be convex in its variables (AOF’s (4) is only convex in the product matrix $S=C^T C$, which enters the objective via $S^{-1}$). Gradient-based optimizers could converge quickly but one would lose the guarantee of being close to a global optimum within a specifiable tolerance. \n\n**\\> Did the authors consider this alternate approach / if so, do you believe banded square root still offers speedups in this setting?**\n\nWe had not considered this approach, but it sounds like promising future work that we would be happy to follow up on (with the permission of the reviewer). That said, we believe it would be unfair to consider it a shortcoming of our work that an alternative path exists which, to our knowledge, has not appeared in any prior work and will come with its own challenges. \n\nRegarding the speed: even if the above idea might accelerate the numerical optimization, it will not surpass the simplicity and efficiency of BSR’s closed-form expressions. Furthermore, even if solved optimally, the quality of such a solution should lie between BSR and AOF, which are close together already. \n\n**\\> Privacy amplification.**\n\nAnother interesting point, thank you for bringing it up. Since the BSR matrix has the same band structure as the AOF matrix, the privacy properties amplification arguments from Choquette-Choo et al. will still be applicable to our results. We will add a note to this in the manuscript. We did not include this aspect in the original submission because the material was already quite dense and we considered it somewhat orthogonal to our main analysis. But indeed, it would be interesting future work.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "khwkr5HQEa", "reviewer": "Reviewer_sFEq", "summary": "This paper addresses optimal matrix factorization mechanism for differential privacy, focussing on stochastic gradient descent (SGD) optimization.\nIt introduces the banded squared root (BSR) factorization and provide formal lower and upper bounds to the approximation error.\nThe BSR achieves competitive formal approximation error compared to state-of-the-art methods and demonstrates practical utility in real-world training tasks. \nThe authors provide theoretical bounds for single and repeated participation training, showing that BSR maintains high training accuracy with provable privacy guarantees while being computationally efficient.", "questions": "./.", "limitations": "./.", "rating": 7, "confidence": 3, "soundness": 4, "presentation": 3, "contribution": 3, "strengths": "- The authors make a significant contribution to the field of differentially private machine learning under matrix factorization mechanisms by leveraging efficient BSR for the SGD workload matrix.\n- The incorporation of the SDG workload matrix in this line of research is elegant.\n- The banded squared root (BSR) factorization is computationally efficient and scalable to large workload matrices, making it a valuable tool.\n- The proof sketches are intuitive and convincing, although I haven't thoroughly reviewed the technical proofs in the appendix.\n- A notable advantage of BSR is its independence from specific training objectives.\n- BSR achieves competitive approximation error, on par with state-of-the-art methods in differentially private training.\n- The approach demonstrates practical utility by maintaining high training accuracy in real-world training tasks.\n- The authors provide theoretical guarantees for both single and repeated participation scenarios, adding to the method's credibility.\n- The theoretical explanations are (for the most part) well-written and easy to follow.\n- The paper includes extensive technical supplementary material.\n- The inclusion of useful, copy-pasteable Python code listings is interesting.\n- The authors show respect for existing work by referencing relevant StackOverflow answers and software packages.", "weaknesses": "- The experimental evaluation is limited to a small set of real-world and synthetic experiments.\n- The range of synthetic datasets used in the experiments is restricted, which may not be representative of diverse scenarios.\n- The data generation process is unclear, and only homogeneous partitions are considered, with no variation in data distributions.\n- The real-world experiments are limited to a single dataset, CIFAR-10, which is insufficient to draw general conclusions.\n- The post-processing step (line 275) lacks clear explanation, leaving it unclear whether it addresses a principled problem or a numerical issue.\n- The experiments do not provide a comparison with non-private approaches, which makes it difficult to understand the practical limitations of the method.\n- The evaluation is limited to a single real-world dataset, CIFAR-10, and a simple ConvNet model.\n- The experiments do not discuss potential practical limitations of BSR, which may impact its applicability.\n- The comparative analysis is limited and does not include a wide range of existing differential privacy methods.\n- The formal proof of Theorem 2 (Appendix D.2) employs unclear and undefined notation (Π), which makes it difficult to follow.\n- Neither the proof sketch nor the formal proof are straightforwardly understandable, which may hinder comprehension.\n\nMinor:\n\n- The tone of \"supposedly optimal\" (line 254) could be perceived as inappropriate and may benefit from rephrasing.\n- The proof sketch of Theorem 2 refers to the wrong appendix (Appendix D.3 instead of Appendix D.2, presumably).", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for the encouraging feedback and detailed comments. \n\n**\\> Experimental Evaluation.**\n\nWe believe there may be a slight misunderstanding. Our main contribution lies in the algorithm and properties of BSR. Specifically, we demonstrate that in the context of MF-SGD, using BSR requires adding less noise compared to the baselines (for identical levels of privacy) and it is essentially on par with the computationally much more costly AOF (e.g., Figure 1). This fact is independent of any subsequent model training steps and requires no data, not even synthetic. We do not introduce a new technique for differential private model training, which would indeed require more extensive experimental validation and comparison to other private training techniques. \n\nThe MF-SGD experiments on CIFAR that we report (Figure 2) are only meant as a proof-of-concept to illustrate that, by keeping all other components fixed, the lower added noise tends to translate into higher accuracy. The exact experimental setup we will be happy to clarify in the manuscript: like prior work, we split the dataset uniformly at random into batches and train in a centralized way. \n\n**\\> Limitations of BSR.**\n\nWe will be happy to extend our discussion in Section 5 to provide a broader view. In particular, we observe that BSR achieves results comparable to AOF, although we currently cannot prove this, as the theoretical properties of AOF are not well understood enough. Regarding real-world usability, the analytic expression for BSR does not apply to variable learning rates, which would be desirable. \n\n**\\> Post-processing step (line 275)**\n\nNote that this step is unrelated to BSR; it was only needed for AOF. We believe the problem is both numerical and principled. By this, we mean that, on the one hand, the problem is indeed numerical: if we could solve (4) with infinite mathematical precision, post-processing would not be needed. On the other hand, the numerical issues arise because of a principled drawback, namely that the AOF procedure (solving (4), computing a Cholesky decomposition, inverting the resulting matrix) is quite sensitive to small numeric deviations. Consequently, errors from solving (4) only to a certain precision, as all numerical solvers do, can lead to disproportionately large undesirable effects. \n\n**\\> Comparison to non-private approaches.**\n \nWe will be happy to include such information in the manuscript. For example, non-private training in the setting of Figure 2 achieves an accuracy of approximately 65%. However, we would like to emphasize that our contribution is orthogonal to the contrast between private and non-private training. Our submission shows that by using a different processing step (BSR), one benefits in terms of speed (over AOF) and accuracy (over the baselines). Outside the context of MF-SGD, such a comparison cannot be made. \n\n**\\> Notation ($\\Pi$), clarity of proof and sketches.**\n\nThank you for bringing this to our attention. We will clarify the proof and sketches and we will remind the reader of the definition of $\\Pi$. \n\n**\\> Minor issues.**\n\nThank you, we will fix these items.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "Kk6HOYhOjj", "reviewer": "Reviewer_1wZc", "summary": "This paper proposes the Banded Square Root Matrix Factorization, an efficient approximation of the optimal banded matrix factorization for differentially private machine learning applications.  The authors give a closed form expression for the BSR C matrix for any SGD + Momentum + Weight Decay workload by exploiting the Toepltiz structure of the workload and factors.  These factorizations can be used in place of the prior work and enjoy some efficiency advantages.", "questions": "* Assuming your method scales well beyond the prior work of n ~ 8192, can you make a convincing case that in modern DP + ML applications we need to scale beyond that point?  I'm not super familiar with how many iterations DP-SGD is usually run for.\n\n* Thm 2 is an interesting result that allows you to compute sensitivity efficiently under this special Toepltiz structure without requiring bandedness, but it's not clear where you are using this Thm.  Is it better to find a banded factorization with bands = minsep or use a non-banded factorization with Thm 2?\n>* Also minor point: might be better to write Eq 10 in terms of m_0, ... m_{n-1} instead of M (or maybe include both expressions).\n>* There is an incorrect reference to Eq 10 in the experiments\n>* Typo: n_{n-1}\n\n* Regarding the statement: \"Apart from the factorization itself, the computational cost of BSR and AOF are nearly identical.\"  Is the factorization time usually the dominating term in the cost for AOF?  Is that still the case for BSR, even when scaling up well beyond n ~ 8192?  What factor(s) other than n (if any) impact the complexity of these mechanisms?\n\n* It is remarkable that the approximation quality of BSR and AOF is so close.  It looks like you evaluated a few settings, one of them being b=100, k=n/100 in Fig 1.  Does this hold more generally across other settings, and can you demonstrate that with an experiment?  There are two extremes where it I might expect *some* gap: (a) 2 bands (b) n bands.  \n>* In Fig 1, the most interesting comparison to me is BSR vs. AOF, the baselines have already been established to be far from optimal in prior work and their presence makes it difficult to perceive difference between the two most interesting lines, would be good to update the plot accordingly.", "limitations": "The authors have discussed some limitations, but perhaps adding some of the things I mentioned in this review (or addressing them directly) would be good.", "rating": 7, "confidence": 5, "soundness": 4, "presentation": 4, "contribution": 3, "strengths": "* The scope, contributions, and key results are clearly written.\n* Near-optimal and efficiency computable banded matrix factorizations could make it feasible to use this mechanism in some new settings where it was previously not possible.\n* Authors demonstrated deep understanding of the area, and made a compelling story with a strong mix of theoretical derivations and empirical observations.", "weaknesses": "* The implementation of AOF seems suboptimal and hence the comparison appears to be biased.  Multiple-day runtimes for n ~ 700 seems high when prior work conducts experiments for n >= 2000 (https://arxiv.org/pdf/2306.08153).  Doing a little digging, I found https://www.kdd.org/kdd2016/papers/files/rpp0224-yuanA.pdf, which ran experiments up to n ~ 8192.  The limits of the BSR approach are not really demonstrated besides the statement \"Even sizes of n = 10,000 or more take at most a few seconds.\"  How much beyond this can you scale?  Would be good to add a scalability experiment.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for the valuable feedback. Below we clarify the raised points:\n\n**\\> Implementation of AOF seems suboptimal.**\n \n(note that Reviewer VxsQ had a similar question, and for easier reading we provide our answer in both our replies)\n\nThere are, of course, many ways to solve AOF (4). We used SCS, because it is a dedicated SPD-solver that allows a simple and reproducible implementation without additional hyperparameters (Algorithm 2). General-purpose gradient-based solvers, such as L-BFGS, were meant by our comment about approximate solvers (lines 122-125). These can be faster for solving (4) up to a certain precision, but they require additional design choices, such as the learning rates, number of update steps, a barrier function to enforce the constraints, and potentially projections to ensure positive definiteness.\n\n**\\> References on larger-scale matrix factorizations.**\n\nRef arXiv:2306.08153 uses an in-house implementation by Google that is not freely available. The manuscript does not state runtimes or hardware. However, the fact that they also only report experiments up to $n=2,052$, we take as confirmation that solving (4) for large $n$ is indeed challenging. Regarding the KDD2016 reference, we do not believe it is comparable. It solves a different optimization problem (in fact, AOF’s formulation (4) was just proposed in 2023), which is also an SDP in the context of DP but with a simpler constraint set. Nevertheless, we’ll be happy to include it in our discussion. \n\nImportantly, however, note that our claim in the work is not that “BSR is x-times faster than AOF”. Actual runtimes always depend on many factors, including hardware, programming language, software packages, etc. We mention runtimes rather as a qualitative guidance to the reader, which is why we do not report a formal table or plot of runtimes. Our core argument is that solving (4) is a challenging optimization problem. *Any* numeric solver will be slower and less convenient to use than our $O(n \\log n)$ closed-form expressions. In real-world problems solving (4) creates a bottleneck, especially because it has to be solved anew for any change of hyperparameters. BSR overcomes this bottleneck.\n\n**\\> Scalability of BSR.**\n\nWe reported experiments to a size of 2000 only due to the computational limits of AOF. For BSR, much larger workload sizes are no problem. We will be happy to add them to the manuscript. The speed advantage of BSR over AOF will increase with $n$, regardless of the solver. We illustrate this in the response PDF: a plain python implementation of Theorem 1 (Alg.3 in PDF) computes the BSR-coefficients for $n=10,000$ in 24ms , for $n=100,000$ in 2.4s and for $n=500,000$ in ~66s (CPU, single-threaded). If a dense matrix is wanted for C, times are approximately 3x higher due to the memory access overhead (Table 1 in PDF).\n\n**\\> Scale of $n$ for modern DP + ML applications?**\n\nIn this context, $n$ is simply the total number of model updates, so large values of $n$ are quite common for standard network training. E.g., training MNIST (60,000 examples) with batchsize 128 for 10 epochs requires $n=(60000/128)*10= 4680$. A standard ResNet-Training on ImageNet1K (1.2M examples, batchsize 512, 90 epochs) yields $n=225,180$. Llama2-like LLM training (2T tokens, batchsize 4M tokens, 1 epoch) would result in $n=500,000$. However, in practice one would not use vanilla MF-SGD for such sizes, because the resulting $C^{-1}$ matrix would require 1 TB of memory to store (in dense form). It remains an open research question how to scale MF-SGD to such sizes. Our contribution only covers a part of that, namely avoiding the need to solve the preparatory (4). \n\n**\\> Theorem 2**\n\nWe use this theorem to compute b-min-sep sensitivity for all cases, including Theorems 6 and 7. The idea of using Theorem 2 to find a better (potentially non-banded) factorization is indeed interesting. However, implementing this is not straightforward, because Theorem 2 requires the matrix $C$ to have decreasing coefficients. It is unclear how one could impose this property during the optimization, where one only has access to $S = C^TC$. \n\n**\\> Is the factorization time the dominating cost?**\n\nFor AOF: yes. The factorization (4) is by far the most costly step, regardless of the solver. The Cholesky-type decomposition (line 118) is also costly but faster (seconds for $n=10000$, minutes for $n=50000$). BSR avoids both steps. Instead, one computes $C$ explicitly, which is efficient due to Theorem 1 (see above), and computing C’s sensitivity is trivial due to Theorem 2. The speed of AOF and BSR are both functions of $n$ and the bandwidth. The runtime of numeric solvers for AOF will also depend on the desired precision and the entries of the workload matrix (which determine the hardness of the optimization problem).\n\n**\\> Notation and typos.**\n\nThank you for the suggestions, we will include them. \n\n**\\> It is remarkable that the approximation quality of BSR and AOF is so close**\n\nWe found BSR to be consistently close to AOF across all settings we tested, not just the one in the manuscript. We’ll be happy to extend the manuscript with more results. \nRegarding bandwidth: Fig.3 in the response PDF shows that for the recommended p=b, BSR is typically close to AOF even for extreme values of $b$ ($b=2$ or $b=n$). For $p\\ll b$ or $p\\gg b$, indeed the distance to (b-banded) AOF grows. Some numeric results for $p=n$ can already be found in Appendix C in column “sqrt” (Tables 10-17 for $b=n$, Tables 2-9 for $b=100$). Unfortunately, their captions were mislabeled, which we will fix.\n\n**\\> Figures**:  the baselines were meant as a reference to help judge the significance of the differences between the methods. However, we could put the plots without baselines in the appendix, or if preferred, we could put the non-banded “sqrt” instead. Please let us know what you would prefer. We hope that also the numeric tables in the appendix allow judging the differences.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "vtxNIyeB1M", "reviewer": "Reviewer_VxsQ", "summary": "The paper proposes the Banded Square Root (BSR) matrix factorization to speed up banded matrix factorization. The authors demonstrate that the workload matrix for SGD with momentum and weight decay can be expressed as a lower triangular Toeplitz matrix. They utilize explicit recursive expressions to compute the square root of this workload matrix, and Theorem 1 provides a closed-form expression for the square root matrix. Theoretical analysis shows that the calculation of the square root matrix is efficient, and the sensitivity of the mechanism is discussed. Detailed analysis of the approximation error is also provided, demonstrating that the square root factorization has asymptotically optimal approximation quality. Experiments show that BSR can achieve performance comparable to previous state-of-the-art methods.", "questions": "1. Enforcing positive semi-definiteness (PSD) for the matrix is crucial for convergence in AOF. In line 274, you mentioned a post-processing step for S to ensure that all its eigenvalues are at least $\\sqrt{1/n}$. Could you provide a more detailed description of this step? I couldn't find it in the source code provided in Appendix A.\n\n2. In equation (6), b should be p, right?", "limitations": "1. If the learning rate is not constant, then the workload matrix will not be a Toeplitz matrix, right? This would significantly impact the efficiency of BSR. During the training process, a proper learning rate schedule can enhance convergence. Therefore, a discussion about varying learning rates should be included.\n\n2. The paper does not provide a general introduction to differential privacy, and Figure 2 presents (epsilon, delta) without sufficient background information. This lack of context makes the paper not self-contained.", "rating": 6, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "1. The discussion regarding SGD with momentum and weight decay is interesting. The workload matrix in this context is a lower triangular Toeplitz matrix, which allows for efficient factorization algorithms. Prior work has not explored this setting; they typically treat the workload matrix with momentum as a general lower triangular matrix. This finding is significant for future research, suggesting potential improvements over current solvers for banded matrix factorization.\n\n2. The theoretical analysis demonstrates that Banded Square Root (BSR) matrix factorization can be computed efficiently even for large problem sizes. The expected approximation error for BSR indicates that it achieves asymptotically optimal approximation quality. These analyses ensure that BSR is both fast and accurate.", "weaknesses": "1. Using CVXPY with SCS to compute AOF is not efficient; L-BFGS can solve this problem more efficiently. For instance, as shown in [1], an implementation using GPUs can solve a problem size of 10,000 in less than an hour. Even with CPUs, the implementation in [2] can solve a problem size of 2,000 in less than an hour, with slight modifications to add the b-min-step constraint.\n\n2. The current experiment only demonstrates results for matrix sizes up to 2000, which is insufficient. A comparison involving matrix sizes of 10,000 or more would be more convincing. How long would it take to obtain a solution for problem sizes of 10,000 and 100,000 using BSR?\n\n[1] https://github.com/apple/pfl-research/blob/2a9d65bd66dc89ef80a24e1be0d28446a0422133/pfl/privacy/ftrl_mechanism.py#L140 \n\n[2] https://github.com/dpcomp-org/hdmm/blob/7a5079a7d4f1a06b0be78019adadf83c538d0514/src/hdmm/templates.py#L391", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for the encouraging feedback. In the following we hope to address all your concerns.\n\n**\\> Using CVXPY with SCS to compute AOF is not efficient.**\n\nThere are, of course, many ways to solve AOF (4). We used SCS, because it is a dedicated SPD-solver which allows for a simple and reproducible implementation without additional hyperparameters (Algorithm 2). General-purpose gradient-based solvers, such as L-BFGS, were meant to be included in our comment about approximate solvers (lines 122-125). These can be faster for solving (4) up to a certain precision, but they require additional design choices, such as the number of update steps, learning rates, a barrier function to enforce the constraints, and potentially projections to ensure positive definiteness.\nImportantly, however, note that our claim in the work is not that “BSR is x-times faster than AOF”. Actual runtimes always depend on many factors, including hardware, programming language, software packages, etc. We mention runtimes rather as a qualitative guidance to the reader, which is why we do not report a formal table or plot of runtimes. Our core argument is that (4) is a challenging optimization problem to solve. *Any* numeric solver will be slower and less convenient to use than our $O(n \\log n)$ closed-form expressions. In real-world problems solving (4) creates a bottleneck, especially because it has to be solved anew for any change of hyperparameters. BSR overcomes this bottleneck.\n\n**\\> References to larger-scale matrix factorizations.**\n \nThank you for the references; we will include them into our manuscript. Ref [1] indeed aims to solve the optimization problem (4). Could you please provide a reference for the timing results you mentioned? We did not spot them in the repository or associated report arXiv:2404.06430, but would be happy to adapt our presentation accordingly. [2] solves a simpler problem (from [McKenna et al., VLDB 2018]). Although adapting the formulation might be possible, the additional constraints could complicate the numerical solution. Nevertheless, we will be happy to include this in our discussion.\n\n**\\> The current experiment only demonstrates results for matrix sizes up to 2000, which is insufficient**\n\nWe report experiments to a size of 2000 due to the computational limits of calculating AOF. For BSR, much larger workload sizes are no problem, and we will be happy to add them to the manuscript. The speed advantage of BSR over AOF will increase with $n$, regardless of the solver. We illustrate this in the response PDF: a plain python implementation of Theorem 1 (Alg.3 in PDF) computes the BSR-coefficients for $n=10,000$ in 24ms , for $n=100,000$ in 2.4s and for $n=500,000$ in ~67s (Table 1 in PDF, single-core CPU, double precision). In practice, one would likely not use vanilla MF-SGD for $n \\ge 500,000$, because the resulting matrices (e.g. $C^{-1}$) would be inconveniently large (1 TB of memory in dense fp32 form). \n\n**\\> Enforcing positive semi-definiteness (PSD).**\n\nWe did not include this step in Algorithm 2, because we considered it a post-processing operation only needed when trying to extract $C$ from the solution matrix $S$ (which despite the use of SCS, is not always perfectly PSD). We will be happy to add it, though. The actual operation is a standard procedure (Alg.4 in the response PDF). Our experiments use $tol=1/\\sqrt{n}$ as a heuristic, because we observed it to give good results across a wide range of problem sizes, in particular $C^{-1}$ remaining stable. Constant values for tol or values proportional to $1/n$ did not work as well. We did not find this issue discussed in the literature. \n\n**\\> Equation (6).**\n\nThanks for spotting this typo. We will fix it. \n\n**\\> Different learning rates.**\n\nIndeed, different learning rates would result in a non-Toeplitz workload matrix, so our analytical expressions no longer apply. However, the banded-square-root method remains applicable, only that the matrix square root has to be computed numerically. This step remains quite efficient because the workload matrix is still lower triangular, for which specialized routines exist, e.g. [Björck& Hammarling. “A Schur method for the square root of a matrix”, 1983]. We will include a discussion of the impact of varying learning rates in the revised version.\n\n**\\> General introduction to differential privacy.**\n\nIndeed, due to space constraints, we were not able to give an introduction in the manuscript, which primarily targets an audience already familiar with DP. We will add a section in the appendix.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}]

Author rebuttal:
Attached is our response PDF; please see the responses to the individual reviews for details.

Meta-review: This paper studies the problem of optimal matrix factorization mechanism for differentially private training. Specifically, the authors develop a banded square root factorization approach and provide formal lower and upper bounds on the approximation error.

All the reviewers agree that the paper presents solid results and makes significant contributions to the field of differentially private training. In light of this, I recommend acceptance and encourage the authors to incorporate the reviewers' suggestions into the revision.
Decision: Accept (poster)