Correlation Aware Sparsified Mean Estimation Using Random Projection

In modern machine learning applications, data is naturally distributed across a large number of edge devices or clients. The underlying learning task in such settings is modeled by distributed optimization or the recent paradigm of Federated Learning (FL) konevcny16federated ; fedavg17aistats ; kairouz2021advances ; wang2021field . A crucial subtask in distributed learning is for the server to compute the mean of the vectors sent by the clients. In FL, for example, clients run training steps on their local data and once-in-a-while send their local models (or local gradients) to the server, which averages them to compute the new global model. However, with the ever-increasing size of machine learning models simonyan2014very ; brown2020language , and the limited battery life of the edge clients, communication cost is often the major constraint for the clients. This motivates the problem of (empirical) distributed mean estimation (DME) under communication constraints, as illustrated in Figure 1. Each of the $n$ clients holds a vector ${\mathbf{x}}_{i}\in\mathbb{R}^{d}$, on which there are no distributional assumptions. Given a communication budget, each client sends a compressed version $\widehat{{\mathbf{x}}}_{i}$ of its vector to the server, which utilizes these to compute an estimate of the mean vector $\frac{1}{n}\sum_{i=1}^{n}{\mathbf{x}}_{i}$.

Quantization and sparsification are two major techniques for reducing the communication costs of DME. Quantization gubner1993distributed ; davies2021new_bounds_dme_var_reduction ; vargaftik2022eden ; suresh2017dme_icml involves compressing each coordinate of the client vector to a given precision and aims to reduce the number of bits to represent each coordinate, achieving a constant reduction in the communication cost. However, the communication cost still remains $\Theta(d)$. Sparsification, on the other hand, aims to reduce the number of coordinates each clinet sends and compresses each client vector to only $k\ll d$ of its coordinates (e.g. Rand-$k$ konevcny2018rand_dme ). As a result, sparsification reduces communication costs more aggressively compared to quantization, achieving better communication efficiency at a cost of only $O(k)$. While in practice, one can use a combination of quantization and sparsification techniques for communication cost reduction, in this work, we focus on the more aggressive sparsification techniques. We call $k$, the dimension of the vector each client sends to the server, the per-client communication budget.

Most existing works on sparsification ignore the potential correlation (or similarity) among the client vectors, which often exists in practice. For example, the data of a specific client in federated learning can be similar to that of multiple clients. Hence, it is reasonable to expect their models (or gradients) to be similar as well. To the best of our knowledge, jhunjhunwala2021dme_spatial_temporal is the first work to account for spatial correlation across individual client vectors. They propose the Rand-$k$-Spatial family of unbiased estimators, which generalizes Rand-$k$ and achieves a better estimation error in the presence of cross-client correlation. However, their approach is focused only on the server-side decoding procedure, while the clients do simple Rand-$k$ encoding.

In this work, we consider a more general encoding scheme that directly compresses a vector from $\mathbb{R}^{d}$ to $\mathbb{R}^{k}$ using a (random) linear map. The encoded vector consists of $k$ linear combinations of the original coordinates. Intuitively, this has a higher chance of capturing the large-magnitude coordinates (“heavy hitters”) of the vector than randomly sampling $k$ out of the $d$ coordinates (Rand-$k$), which is crucial for the estimator to recover the true mean vector. For example, consider a vector where only a few coordinates are heavy hitters. For small $k$, Rand-$k$ has a decent chance of missing all the heavy hitters. But with a linear-maps-based general encoding procedure, the large coordinates are more likely to be encoded in the linear measurements, resulting in a more accurate estimator of the mean vector. Guided by this intuition, we ask:

Can we design an improved joint encoding-decoding scheme that utilizes the correlation information and achieves an improved estimation error?

One naïve solution is to apply the same random rotation matrix ${\bm{G}}\in\mathbb{R}^{d\times d}$ to each client vector, before applying Rand-$k$ or Rand-$k$-Spatial encoding. Indeed, such preprocessing is applied to improve the estimator using quantization techniques on heterogeneous vectors suresh2022correlated_dme_icml ; suresh2017dme_icml . However, as we see in Appendix A.1, for sparsification, we can show that this leads to no improvement. But what happens if every client uses a different random matrix, or applies a random $k\times d$-dimensional linear map? How to design the corresponding decoding procedure to leverage cross-client correlation? As there is no way for one to directly apply the decoding procedure of Rand-$k$-Spatial in such cases. To answer these questions, we propose the Rand-Proj-Spatial family estimator. We propose a flexible encoding procedure in which each client applies its own random linear map to encode the vector. Further, our novel decoding procedure can better leverage cross-client correlation. The resulting mean estimator generalizes and improves over the Rand-$k$-Spatial family estimator.

Next, we discuss some reasonable restrictions we expect our mean estimator to obey. 1) Unbiased. An unbiased mean estimator is theoretically more convenient compared to a biased one horvath2021induced . 2) Non-adaptive. We focus on an encoding procedure that does not depend on the actual client data, as opposed to the adaptive ones, e.g. Rand-$k$ with vector-based sampling probability konevcny2018rand_dme ; wangni2018grad_sparse . Designing a data-adaptive encoding procedure is computationally expensive as this might require using an iterative procedure to find out the sampling probabilities konevcny2018rand_dme . In practice, however, clients often have limited computational power compared to the server. Further, as discussed earlier, mean estimation is often a subroutine in more complicated tasks. For applications with streaming data nokleby2018stochastic , the additional computational overhead of adaptive schemes is challenging to maintain. Note that both Rand-$k$ and Rand-$k$-Spatial family estimator jhunjhunwala2021dme_spatial_temporal are unbiased and non-adaptive.

In this paper, we focus on the severely communication-constrained case $nk\leq d$, when the server receives very limited information about any single client vector. If $nk\gg d$, we see in Appendix A.2 that the cross-client information has no additional advantage in terms of improving the mean estimate under both Rand-$k$-Spatial or Rand-Proj-Spatial, with different choices of random linear maps. Furthermore, when $nk\gg d$, the performance of both the estimators converges to that of Rand-$k$. Intuitively, this means when the server receives sufficient information regarding the client vectors, it does not need to leverage cross-client correlation to improve the mean estimator.

Our contributions can be summarized as follows:

1. 1.

   We propose the Rand-Proj-Spatial family estimator with a more flexible encoding-decoding procedure, which can better leverage the cross-client correlation information to achieve a more general and improved mean estimator compared to existing ones.

2. 2.

   We show the benefit of using Subsampled Randomized Hadamard Transform (SRHT) as the random linear maps in Rand-Proj-Spatial in terms of better mean estimation error (MSE). We theoretically analyze the case when the correlation information is known at the server (see Theorems 4.3, 4.4 and Section 4.3). Further, we propose a practical configuration called Rand-Proj-Spatial(Avg) when the correlation is unknown.

3. 3.

   We conduct experiments on common distributed optimization tasks, and demonstrate the superior performance of Rand-Proj-Spatial compared to existing sparsification techniques.

Quantization and Sparsification. Commonly used techniques to achieve communication efficiency are quantization, sparsification, or more generic compression schemes, which generalize the former two basu2019qsparse . Quantization involves either representing each coordinate of the vector by a small number of bits davies2021new_bounds_dme_var_reduction ; vargaftik2022eden ; suresh2017dme_icml ; alistarh2017qsgd_neurips ; bernstein2018signsgd ; reisizadeh2020fedpaq_aistats , or more involved vector quantization techniques shlezinger2020uveqfed_tsp ; gandikota2021vqsgd_aistats . Sparsification wangni2018grad_sparse ; alistarh2018convergence ; stich2018sparsified ; karimireddy2019error ; sattler2019robust , on the other hand, involves communicating a small number $k<d$ of coordinates, to the server. Common protocols include Rand-$k$ konevcny2018rand_dme , sending $k$ uniformly randomly selected coordinates; Top-$k$ shi2019topk , sending the $k$ largest magnitude coordinates; and a combination of the two barnes2020rtop_jsait . Some recent works, with a focus on distributed learning, further refine these communication-saving mechanisms ozfatura2021time by incorporating temporal correlation or error feedback horvath2021induced ; karimireddy2019error .

Distributed Mean Estimation (DME). DME has wide applications in distributed optimization and FL. Most of the existing literature on DME either considers statistical mean estimation zhang2013lower_bd ; garg2014comm_neurips , assuming that the data across clients is generated i.i.d. according to the same distribution, or empirical mean estimation suresh2017dme_icml ; chen2020breaking ; mayekar2021wyner ; jhunjhunwala2021dme_spatial_temporal ; konevcny2018rand_dme ; vargaftik2021drive_neurips ; vargaftik2022eden_icml , without making any distributional assumptions on the data. A recent line of work on empirical DME considers applying additional information available to the server, to further improve the mean estimate. This side information includes cross-client correlation  jhunjhunwala2021dme_spatial_temporal ; suresh2022correlated_dme_icml , or the memory of the past updates sent by the clients liang2021improved_isit .

Subsampled Randomized Hadamard Transformation (SRHT). SRHT was introduced for random dimensionality reduction using sketching Ailon2006srht_initial ; tropp2011improved ; lacotte2020optimal_iter_sketching_srht . Common applications of SRHT include faster computation of matrix problems, such as low-rank approximation Balabanov2022block_srht_dist_low_rank ; boutsidis2013improved_srht , and machine learning tasks, such as ridge regression lu2013ridge_reg_srht , and least square problems Chehreghani2020graph_reg_srht ; dan2022least_sq_srht ; lacotte2020optimal_first_order_srht . SRHT has also been applied to improve communication efficiency in distributed optimization ivkin2019distSGD_sketch_neurips and FL haddadpour2020fedsketch ; rothchild2020fetchsgd_icml .

Notation. We use bold lowercase (uppercase) letters, e.g. ${\mathbf{x}}$ (${\bm{G}}$) to denote vectors (matrices). ${\mathbf{e}}_{j}\in\mathbb{R}^{d}$, for $j\in[d]$, denotes the $j$-th canonical basis vector. $\|\cdot\|_{2}$ denotes the Euclidean norm. For a vector ${\mathbf{x}}$, ${\mathbf{x}}(j)$ denotes its $j$-th coordinate. Given integer $m$, we denote by $[m]$ the set $\{1,2,\dots,m\}$.

Problem Setup. Consider $n$ geographically separated clients coordinated by a central server. Each client $i\in[n]$ holds a vector ${\mathbf{x}}_{i}\in\mathbb{R}^{d}$, while the server wants to estimate the mean vector $\bar{{\mathbf{x}}}\triangleq\frac{1}{n}\sum_{i=1}^{n}{\mathbf{x}}_{i}$. Given a per-client communication budget of $k\in[d]$, each client $i$ computes $\widehat{{\mathbf{x}}}_{i}$ and sends it to the central server. $\widehat{{\mathbf{x}}}_{i}$ is an approximation of ${\mathbf{x}}_{i}$ that belongs to a random $k$-dimensional subspace. Each client also sends a random seed to the server, which conveys the subspace information, and can usually be communicated using a negligible amount of bits. Having received the encoded vectors $\{\widehat{{\mathbf{x}}}_{i}\}_{i=1}^{n}$, the server then computes $\widehat{{\mathbf{x}}}\in\mathbb{R}^{d}$, an estimator of $\bar{{\mathbf{x}}}$. We consider the severely communication-constrained setting where $nk\leq d$, when only a limited amount of information about the client vectors is seen by the server.

Error Metric. We measure the quality of the decoded vector $\widehat{{\mathbf{x}}}$ using the Mean Squared Error (MSE) $\mathbb{E}\left[\|\widehat{{\mathbf{x}}}-\bar{{\mathbf{x}}}\|_{2}^{2}\right]$, where the expectation is with respect to all the randomness in the encoding-decoding scheme. Our goal is to design an encoding-decoding algorithm to achieve an unbiased estimate $\widehat{{\mathbf{x}}}$ (i.e. $\mathbb{E}[\widehat{{\mathbf{x}}}]=\bar{{\mathbf{x}}}$) that minimizes the MSE, given the per-client communication budget $k$. To consider an example, in rand-$k$ sparsification, each client sends randomly selected $k$ out of its $d$ coordinates to the server. The server then computes the mean estimate as $\widehat{{\mathbf{x}}}^{(\text{Rand-$k$})}=\frac{1}{n}\frac{d}{k}\sum_{i=1}^{n}\widehat{{\mathbf{x}}}_{i}$. By (jhunjhunwala2021dme_spatial_temporal, , Lemma 1), the MSE of Rand-$k$ sparsification is given by

The Rand-$k$-Spatial Family Estimator. For large values of $\frac{d}{k}$, the Rand-$k$ MSE in Eq. 1 can be prohibitive. jhunjhunwala2021dme_spatial_temporal proposed the Rand-$k$-Spatial family estimator, which achieves an improved MSE, by leveraging the knowledge of the correlation between client vectors at the server. The encoded vectors $\{\widehat{{\mathbf{x}}}_{i}\}$ are the same as in Rand-$k$. However, the $j$-th coordinate of the decoded vector is given as

Here, $T:\mathbb{R}\rightarrow\mathbb{R}$ is a pre-defined transformation function of $M_{j}$, the number of clients which sent their $j$-th coordinate, and $\bar{\beta}$ is a normalization constant to ensure $\widehat{{\mathbf{x}}}$ is an unbiased estimator of ${\mathbf{x}}$. The resulting MSE is given by

where $c_{1},c_{2}$ are constants dependent on $n,d,k$ and $T$, but independent of client vectors $\{{\mathbf{x}}_{i}\}_{i=1}^{n}$. When the client vectors are orthogonal, i.e., $\langle{\mathbf{x}}_{i},{\mathbf{x}}_{l}\rangle=0$, for all $i\neq l$, jhunjhunwala2021dme_spatial_temporal show that with appropriately chosen $T$, the MSE in Eq. 3 reduces to Eq. 1. However, if there exists a positive correlation between the vectors, the MSE in Eq. 3 is strictly smaller than that for Rand-$k$ Eq. 1.

While the Rand-$k$-Spatial family estimator proposed in jhunjhunwala2021dme_spatial_temporal focuses only on improving the decoding at the server, we consider a more general encoding-decoding scheme. Rather than simply communicating $k$ out of the $d$ coordinates of its vector ${\mathbf{x}}_{i}$ to the server, client $i$ applies a (random) linear map ${\bm{G}}_{i}\in\mathbb{R}^{k\times d}$ to ${\mathbf{x}}_{i}$ and sends $\widehat{{\mathbf{x}}}_{i}={\bm{G}}_{i}{\mathbf{x}}_{i}\in\mathbb{R}^{k}$ to the server. The decoding process on the server first projects the encoded vectors $\{{\bm{G}}_{i}{\mathbf{x}}_{i}\}_{i=1}^{n}$ back to the $d$-dimensional space and then forms an estimate $\widehat{{\mathbf{x}}}$. We motivate our new decoding procedure with the following regression problem:

To understand the motivation behind Eq. 4, first consider the special case where ${\bm{G}}_{i}={\bm{I}}_{d}$ for all $i\in[n]$, that is, the clients communicate their vectors without compressing. The server can then exactly compute the mean $\bar{\mathbf{x}}=\frac{1}{n}\sum_{i=1}^{n}{\mathbf{x}}_{i}$. Equivalently, $\bar{\mathbf{x}}$ is the solution of $\operatorname*{arg\,min}_{{\mathbf{x}}}\sum_{i=1}^{n}\|{\mathbf{x}}-{\mathbf{x}}_{i}\|_{2}^{2}$. In the more general setting, we require that the mean estimate $\widehat{{\mathbf{x}}}$ when encoded using the map ${\bm{G}}_{i}$, should be “close” to the encoded vector ${\bm{G}}_{i}{\mathbf{x}}_{i}$ originally sent by client $i$, for all clients $i\in[n]$.

We note the above intuition can also be translated into different regression problems to motivate the design of the new decoding procedure. We discuss in Appendix B.2 intuitive alternatives which, unfortunately, either do not enable the usage of cross-client correlation information, or do not use such information effectively. We choose the formulation in Eq. 4 due to its analytical tractability and its direct relevance to our target error metric MSE. We note that it is possible to consider the problem in Eq. 4 in the other norms, such as the sum of $\ell_{2}$ norms (without the squares) or the $\ell_{\infty}$ norm. We leave this as a future direction to explore.

The solution to Eq. 4 is given by $\widehat{{\mathbf{x}}}^{(\text{Rand-Proj})}=(\sum_{i=1}^{n}{\bm{G}}_{i}^{T}{\bm{G}}_{i})^{{\dagger}}\sum_{i=1}^{n}{\bm{G}}_{i}^{T}{\bm{G}}_{i}{\mathbf{x}}_{i}$, where ${\dagger}$ denotes the Moore-Penrose pseudo inverse golub2013matrix_book . However, while $\widehat{{\mathbf{x}}}^{(\text{Rand-Proj})}$ minimizes the error of the regression problem, our goal is to design an unbiased estimator that also improves the MSE. Therefore, we make the following two modifications to $\widehat{{\mathbf{x}}}^{(\text{Rand-Proj})}$: First, to ensure that the mean estimate is unbiased, we scale the solution by a normalization factor $\bar{\beta}$¹¹1We show that it suffices for $\bar{\beta}$ to be a scalar in Appendix B.1. . Second, to incorporate varying degrees of correlation among the clients, we propose to apply a scalar transformation function $T:\mathbb{R}\rightarrow\mathbb{R}$ to each of the eigenvalues of $\sum_{i=1}^{n}{\bm{G}}_{i}^{T}{\bm{G}}_{i}$. The resulting Rand-Proj-Spatial family estimator is given by

Though applying the transformation function $T$ in Rand-Proj-Spatial requires computing the eigendecomposition of $\sum_{i=1}^{n}{\bm{G}}_{i}^{T}{\bm{G}}_{i}$. However, this happens only at the server, which has more computational power than the clients. Next, we observe that for appropriate choice of $\{{\bm{G}}_{i}\}_{i=1}^{n}$, the Rand-Proj-Spatial family estimator reduces to the Rand-$k$-Spatial family estimator jhunjhunwala2021dme_spatial_temporal .

The proof details are in Appendix C.5. We discuss the choice of $T$ and how it compares to Rand-$k$-Spatial in detail in Section 4.3.

To showcase the utility of the Rand-Proj-Spatial family estimator, we propose to set the random linear maps ${\bm{G}}_{i}$ to be scaled Subsampled Randomized Hadamard Transform (SRHT, e.g. tropp2011improved ). Assuming $d$ to be a power of $2$, the linear map ${\bm{G}}_{i}$ is given as

where ${\bm{E}}_{i}\in\mathbb{R}^{k\times d}$ is the subsampling matrix, ${\bm{H}}\in\mathbb{R}^{d\times d}$ is the (deterministic) Hadamard matrix and ${\bm{D}}_{i}\in\mathbb{R}^{d\times d}$ is a diagonal matrix with independent Rademacher random variables as its diagonal entries. We choose SRHT due to its superior performance compared to other random matrices. Other possible choices of random matrices for Rand-Proj-Spatial estimator include sketching matrices commonly used for dimensionality reduction, such as Gaussian weinberger2004learning ; tripathy2016gaussian , row-normalized Gaussian, and Count Sketch minton2013improved_bounds_countsketch , as well as error-correction coding matrices, such as Low-Density Parity Check (LDPC) gallager1962LDPC and Fountain Codes Shokrollahi2005fountain_codes . However, in the absence of correlation between client vectors, all these matrices suffer a higher MSE.

In the following, we first compare the MSE of Rand-Proj-Spatial with SRHT against Rand-$k$ and Rand-$k$-Spatial in two extreme cases: when all the client vectors are identical, and when all the client vectors are orthogonal to each other. In both cases, we highlight the transformation function $T$ used in Rand-Proj-Spatial (Eq. 5) to incorporate the knowledge of cross-client correlation. We define

to measure the correlation between the client vectors. Note that ${\mathcal{R}}\in[-1,n-1]$. ${\mathcal{R}}=0$ implies all client vectors are orthogonal, while ${\mathcal{R}}=n-1$ implies identical client vectors.

We consider three practical distributed optimization tasks for evaluation: distributed power iteration, distributed $k$-means and distributed linear regression. We compare Rand-Proj-Spatial(Avg) against Rand-$k$, Rand-$k$-Spatial(Avg), and two more sophisticated but widely used sparsification schemes: non-uniform coordinate-wise gradient sparsification wangni2018grad_sparse (we call it Rand-$k$(Wangni)) and the Induced compressor with Rand-$k$ + Top-$k$ horvath2021induced . The results are presented in Figure 4.

Dataset. For both distributed power iteration and distributed $k$-means, we use the test set of the Fashion-MNIST dataset xiao2017fashion consisting of $10000$ samples. The original images from Fashion-MNIST are $28\times 28$ in size. We preprocess and resize each image to be $32\times 32$. Resizing images to have their dimension as a power of 2 is a common technique used in computer vision to accelerate the convolution operation. We use the UJIndoor dataset ²²2https://archive.ics.uci.edu/ml/datasets/ujiindoorloc for distributed linear regression. We subsample $10000$ data points, and use the first $512$ out of the total $520$ features on signals of phone calls. The task is to predict the longitude of the location of a phone call. In all the experiments in Figure 4, the datasets are split IID across the clients via random shuffling. In Appendix D.1, we have additional results for non-IID data split across the clients.

Setup and Metric. Recall that $n$ denotes the number of clients, $k$ the per-client communication budget, and $d$ the vector dimension. For Rand-Proj-Spatial, we use the first $50$ iterations to estimate $\bar{\beta}$ (see Eq. 5). Note that $\bar{\beta}$ only depends on $n,k,d$, and $T$ (the transformation function in Eq. 5), but is independent of the dataset. We repeat the experiments across 10 independent runs, and report the mean MSE (solid lines) and one standard deviation (shaded regions) for each estimator. For each task, we plot the squared error of the mean estimator $\widehat{{\mathbf{x}}}$, i.e., $\|\widehat{{\mathbf{x}}}-\bar{{\mathbf{x}}}\|_{2}^{2}$, and the values of the task-specific loss function, detailed below.

Tasks and Settings:

1. Distributed power iteration. We estimate the principle eigenvector of the covariance matrix, with the dataset (Fashion-MNIST) distributed across the $n$ clients. In each iteration, each client computes a local principle eigenvector estimate based on a single power iteration and sends an encoded version to the server. The server then computes a global estimate and sends it back to the clients. The task-specific loss here is $\|{\mathbf{v}}_{t}-{\mathbf{v}}_{top}\|_{2}$, where ${\mathbf{v}}_{t}$ is the global estimate of the principal eigenvector at iteration $t$, and ${\mathbf{v}}_{top}$ is the true principle eigenvector.

2. Distributed $k$-means. We perform $k$-means clustering balcan2013distributed with the data distributed across $n$ clients (Fashion-MNIST, 10 classes) using Lloyd’s algorithm. At each iteration, each client performs a single iteration of $k$-means to find its local centroids and sends the encoded version to the server. The server then computes an estimate of the global centroids and sends them back to the clients. We report the average squared mean estimation error across 10 clusters, and the $k$-means loss, i.e., the sum of the squared distances of the data points to the centroids.

For both distributed power iterations and distributed $k$-means, we run the experiments for $30$ iterations and consider two different settings: $n=10,k=102$ and $n=50,k=20$.

3. Distributed linear regression. We perform linear regression on the UJIndoor dataset distributed across $n$ clients using SGD. At each iteration, each client computes a local gradient and sends an encoded version to the server. The server computes a global estimate of the gradient, performs an SGD step, and sends the updated parameter to the clients. We run the experiments for $50$ iterations with learning rate $0.001$. The task-specific loss is the linear regression loss, i.e. empirical mean squared error. To have a proper scale that better showcases the difference in performance of different estimators, we plot the results starting from the 10th iteration.

Results. It is evident from Figure 4 that Rand-Proj-Spatial(Avg), our estimator with the practical configuration $T$ (see Section 4.3) that does not require the knowledge of the actual degree of correlation among clients, consistently outperforms the other estimators in all three tasks. Additional experiments for the three tasks are included in Appendix D.1. Furthermore, we present the wall-clock time to encode and decode client vectors using different sparsification schemes in Figure 5. Though Rand-Proj-Spatial(Avg) has the longest decoding time, the encoding time of Rand-Proj-Spatial(Avg) is less than that of the adaptive Rand-$k$(Wangni) sparsifier. In practice, the server has more computational power than the clients and hence can afford a longer decoding time. Therefore, it is more important to have efficient encoding procedures.

We note two practical limitations of the proposed Rand-Proj-Spatial.
1) Computation Time of Rand-Proj-Spatial. The encoding time of Rand-Proj-Spatial is $O(kd)$, while the decoding time is $O(d^{2}\cdot nk)$. The computation bottleneck in decoding is computing the eigendecomposition of the $d\times d$ matrix $\sum_{i=1}^{n}{\bm{G}}_{i}^{T}{\bm{G}}_{i}$ of rank at most $nk$. Improving the computation time for both the encoding and decoding schemes is an important direction for future work.
2) Perfect Shared Randomness. It is common to assume perfect shared randomness between the server and the clients in distributed settings zhou2022ldp_sparse_vec_agg . However, to perfectly simulate randomness using Pseudo Random Number Generator (PRNG), at least $\log_{2}d$ bits of the seed need to be exchanged in practice. We acknowledge this gap between theory and practice.

In this paper, we propose the Rand-Proj-Spatial estimator, a novel encoding-decoding scheme, for communication-efficient distributed mean estimation. The proposed client-side encoding generalizes and improves the commonly used Rand-$k$ sparsification, by utilizing projections onto general $k$-dimensional subspaces. On the server side, cross-client correlation is leveraged to improve the approximation error. Compared to existing methods, the proposed scheme consistently achieves better mean estimation error across a variety of tasks. Potential future directions include improving the computation time of Rand-Proj-Spatial and exploring whether the proposed Rand-Proj-Spatial achieves the optimal estimation error among the class of non-adaptive estimators, given correlation information. Furthermore, combining sparsification and quantization techniques and deriving such algorithms with the optimal communication cost-estimation error trade-offs would be interesting.

We would like to thank the anonymous reviewer for providing valuable feedback on the title of this work, interesting open problems, alternative motivating regression problems and practical limitations of shared randomness. This work was supported in part by NSF grants CCF 2045694, CCF 2107085, CNS-2112471, and ONR N00014-23-1- 2149.