# Efficient Sparse Least Absolute Deviation Regression with Differential
Privacy

Weidong Liu , Xiaojun Mao , Xiaofei Zhang , and Xin Zhang School of
Mathematical Sciences, MoE Key Lab of Artificial Intelligence, Shanghai Jiao
Tong University, Shanghai, 200240, China (e-mail: weidongl@sjtu.edu.cn).School
of Mathematical Sciences, Ministry of Education Key Laboratory of Scientific
and Engineering Computing, Shanghai Jiao Tong University, Shanghai, 200240,
China (e-mail: maoxj@sjtu.edu.cn).School of Statistics and Mathematics,
Zhongnan University of Economics and Law, Wuhan, 430073, China (e-mail:
zhangxf@zuel.edu.cn).Department of Statistics, Iowa State University, Ames,
50011, USA (e-mail: xinzhang.nac@gmail.com).

###### Abstract

In recent years, privacy-preserving machine learning algorithms have attracted
increasing attention because of their important applications in many
scientific fields. However, in the literature, most privacy-preserving
algorithms demand learning objectives to be strongly convex and Lipschitz
smooth, which thus cannot cover a wide class of robust loss functions (e.g.,
quantile/least absolute loss). In this work, we aim to develop a fast privacy-
preserving learning solution for a sparse robust regression problem. Our
learning loss consists of a robust least absolute loss and an $\ell_{1}$
sparse penalty term. To fast solve the non-smooth loss under a given privacy
budget, we develop a Fast Robust And Privacy-Preserving Estimation (FRAPPE)
algorithm for least absolute deviation regression. Our algorithm achieves a
fast estimation by reformulating the sparse LAD problem as a penalized least
square estimation problem and adopts a three-stage noise injection to
guarantee the $(\epsilon,\delta)$-differential privacy. We show that our
algorithm can achieve better privacy and statistical accuracy trade-off
compared with the state-of-the-art privacy-preserving regression algorithms.
In the end, we conduct experiments to verify the efficiency of our proposed
FRAPPE algorithm.

Keywords: Robust Regression, Sparse Learning, Differential Privacy, Least
Absolute Deviation

##  1 Introduction

Recent years have witnessed that machine learning (ML) technology has rapidly
and fundamentally changed human daily lives. However, as many ML models and
training procedures involve sensitive personal information (e.g. healthcare
data, financial records, web browser history, etc.), data privacy becomes a
crucial concern. To ensure that only population-wise patterns are revealed,
but individual sensitive information is well protected, differential privacy
has been established [1] as a rigorous and standard definition of privacy.
Under the framework of differential privacy, a long line of research on
privacy-preserving ML has been inspired, such as sparse regression [2, 3, 4],
clustering analysis [5, 6], and deep neural network training [7, 8, 9].

Although various differential private ML methods have been developed, few of
them take the model or algorithm robustness into consideration. In real-world
ML applications, data are typically contaminated with heavy-tailed noises and
abnormal values [10]. Naively analyzing the data with classical least square
loss will lead to a fitting bias and deteriorate the estimation performance.
To mitigate the negative impact of outlier data and improve learning
robustness, in the literature, the least absolute loss has become a popular
surrogate for the least square loss in regression analysis. As the sparse
regression problem covers a broad class of supervised ML tasks, there is a
pressing demand to develop an efficient privacy-preserving algorithm for
sparse least absolute deviation (LAD) regression.

Unfortunately, the design of an efficient privacy-preserving sparse LAD
regression algorithm is a non-trivial task. On the one hand, to improve the
privacy of an iterative learning algorithm, a widely adopted strategy is to
inject certain unbiased independent random noise into the intermediate
statistics/variables during iterations [11, 12, 13, 14]. However, such a noise
injection method typically deteriorates the privacy guarantee as the iteration
round $T$ increases. Therefore, with a given privacy budget, there is an upper
limit for the total iteration rounds, which introduces the trade-off between
privacy and accuracy. On the other hand, unlike the classical least square
loss, the least absolute loss is a non-smooth function with no Lipschitz
continuous gradient. Thus, the commonly used Newton or gradient-descent based
algorithms are not applicable to the LAD regression problem. To solve the non-
smooth least absolute loss, a few algorithms were proposed in recent
literature, including the subgradient descent algorithm [15, 16], EM
algorithm[17, 18, 19], ADMM[20, 21] and smoothed approximations [22] But these
algorithms suffer poor convergence performance, which means a large value of
iterations $T$ is needed for the desired estimation accuracy. The limitations
caused by these two conflicting challenges naturally motivate our research
question: _Is it possible to develop an efficient privacy-preserving algorithm
for sparse LAD regression, so that it could achieve excellent performance in
terms of both privacy protection and estimation accuracy?_

In this paper, we will give an _affirmative_ answer to the above question. We
develop a novel algorithm termed FRAPPE for fast robust and privacy-preserving
sparse LAD regression estimation. Our main contributions and their
significance are summarized in the following:

  * •

We study the robust sparse regression problem with privacy constraints. Our
learning loss is a combination of non-smooth least absolute loss and
$\ell_{1}$ sparse penalty. To optimize the loss efficiently, we propose an
efficient algorithm, named FRAPPE. Specifically, our algorithm consists of
three steps: initial estimation, pseudo response transformation, and surrogate
loss optimization. Among them, the key step is the pseudo response
transformation. Our proposed transformation converts the least absolute loss
to the least square loss. Therefore, the problem can be easily solved with a
proximal gradient descent based algorithm.

  * •

To guarantee $(\epsilon,\delta)$-differential privacy, our algorithm is
designed with a three-stage noise injection: The first stage of noise
injection is for the initial estimator. The output perturbation is applied so
that the procedure does not depend on the initialization method. The second
stage is to add random noise to the kernel density estimation, which is
essential for constructing pseudo responses. In the third stage, gradient
perturbation is adopted as our algorithm leverages the proximal gradient
descent for surrogate loss optimization.

  * •

By incorporating the convergence error and injected noises, we establish the
statistical consistency rate with the privacy-accuracy trade-off in Theorem 3.
We show in Theorem 3 that after a constant number of iterations, our
estimators achieve a trade-off between the near-optimal rate of $O(\sqrt{s\log
p/N})$ and the privacy loss
$O(\sqrt{p\log(1/\delta)\log(N\epsilon)}/(N\epsilon))$, where $N,p,s$ are the
size of the total sample, dimension of regression weights, and number of non-
zero weights.

The rest of the paper is organized as follows. Section 2 discusses related
works for robust regression estimation and differential private sparse
regression. In Section 3, we specify the problem formulation and introduce the
proposed FRAPPE algorithm. Section 4 shows the privacy guarantee of the
proposed algorithm and the theoretical guarantee of statistical consistency
and sparsity recovery. Simulation results and real data analysis are reported
in Section 5. The conclusion is given in Section 6, and the technical details
are given in the supplementary materials.

Notation: For any $p$ dimensional vector
$\mathbf{x}=(x_{1},\cdots,x_{p})^{\top}$, we use
$\|\mathbf{x}\|_{1}=\sum_{i=1}^{p}\mbox{$\lvert x_{i}\rvert$}$ to denote the
$\ell_{1}$-norm, and $\|\mathbf{x}\|_{2}=(\sum_{i=1}^{p}x_{i}^{2})^{1/2}$ to
denote the $\ell_{2}$-norm. Given two sequences $a_{n}$ and $b_{n}$, we define
$a_{n}=O(b_{n})$ if there is a constant $M>0$ such that $a_{n}\leq Mb_{n}$. We
say $a_{n}=\Theta(b_{n})$ if and only if $a_{n}=O(b_{n})$ and $b_{n}=O(a_{n})$
both hold. Given a sequence of random variables $X_{n}$ and a sequence of
constants $a_{n}$, we say $X_{n}=O_{\mathbb{P}}(a_{n})$ if for any
$\varepsilon>0$, there exist a finite $M>0$ and a finite $N>0$ such that
$P\left(\left|X_{n}/a_{n}\right|>M\right)<\varepsilon,\forall n>N$.

##  2 Preliminaries and Related Works

###  2.1 Robust Linear Regression

The linear regression model has been one of the most powerful machine learning
tools. Given the independently and identically distributed observations
$\\{\mathbf{x}_{i},y_{i}\\}_{i=1}^{N}$, we assume that the data follow
$y_{i}=\mathbf{x}_{i}^{\top}\bm{\beta}^{*}+e_{i}$, where $\bm{\beta}^{*}$ is
an unknown regression weight and $e_{i}$ is the corresponding random error for
the $i$th observation. Traditionally, the regression weight $\bm{\beta}^{*}$
can be estimated by solving the ordinary least squares (OLS) estimation
problem $\sum_{i=1}^{N}(y_{i}-\mathbf{x}_{i}^{\top}\bm{\beta})^{2}/N$.
However, it is commonly known that the OLS estimation is very sensitive to
outlier data, and thus it is not a robust estimator [23, 24].

In the literature, various methods have been proposed to enhance the
robustness of regression, including M-estimates [23], least median of squares
(LMS) estimates [25], least trimmed squares (LTS) estimates [26] and
S-estimates [27], etc. In this work, we focus on the least absolute deviation
(LAD) regression method to achieve robust regression estimation. LAD
regression method, also known as the median regression method, estimates the
coefficient $\bm{\beta}^{*}$ by minimizing the least absolute loss
$\sum_{i=1}^{N}|y_{i}-\mathbf{x}_{i}^{\top}\bm{\beta}|/N$. Note that the LAD
regression can be viewed as a special case of the well-known quantile
regression by taking the quantile level as $1/2$. Quantile regression has been
widely adopted for robust estimation as its estimator is less likely affected
by the outliers than OLS estimator [28, 29, 30, 31].

Furthermore, we also consider the problem of identifying the sparse structure
in the regression weight $\bm{\beta}^{*}$. In [29], the authors studied the
Lasso penalized quantile regression problem and proposed a consistent
estimator at a near-oracle rate. The authors of [32] studied the estimator for
the $\ell_{1}$ penalized LAD regression problem and showed that their
estimator could achieve near-oracle risk performance. Later, in [33, 34, 35],
the authors further considered using adaptive weighted Lasso or folded concave
penalty to improve the estimation performance of the sparse LAD regression
problem. Estimating $\bm{\beta}^{*}$ for a sparse LAD regression typically
requires more effort than solving a sparse penalized OLS problem due to the
non-smoothness property of the median loss. The existing works proposed
methods of interior-point algorithm [36], greedy coordinate descent algorithm
[37, 38], subgradient descent algorithm [15], majorization-minimization
algorithm [39] or ADMM algorithm [21], etc, to solve the sparse LAD regression
problem. However, these methods either lack the theoretical convergence
guarantee or suffer a slow convergence rate. Furthermore, none of these works
could provide the guarantee of preserving data privacy. In this work, we aim
to study the efficient algorithm for sparse LAD regression problems and also
to establish the theoretical guarantee of estimation consistency, algorithm
convergence, and data privacy.

###  2.2 Differential Private Regression

Differential privacy (DP) [1] is a canonical metric for data privacy analysis.
Under the DP framework, the privacy is measured by how significant the
distribution of query outcome changes when we only replace one single
sensitive data in the dataset. Mathematically, the definition of DP can be
stated as following:

######  Definition 1 ($(\epsilon,\delta)$-Differential Privacy).

Given any two adjacent data sets
$\mathcal{D},\mathcal{D}^{\prime}\in\mathcal{D}^{n}$ differing by a single
element, a randomized mechanism
$\mathcal{M}:\mathcal{D}^{n}\rightarrow\mathbb{R}^{d}$ provides
$(\epsilon,\delta)$-differential privacy if it produces response in any
measurable set $A$ with probability

$\displaystyle\mathbb{P}[\mathcal{M}(\mathcal{D})\in A]\leq e^{\epsilon}\cdot
P[\mathcal{M}(\mathcal{D}^{\prime})\in A]+\delta.$ (1)

Typically, the literature achieves the $(\epsilon,\delta)$-differential
privacy by perturbing the query mechanism with a certain type of random noise,
such as the Gaussian or Laplacian noise. The scale of injected noise is
jointly determined by the given privacy budget $\epsilon$ the $\ell_{2}$
sensitivity of the query function:

######  Definition 2 ($\ell_{2}$-Sensitivity).

The $\ell_{2}$-Sensitivity is defined as the maximum change in the
$\ell_{2}$-Norm of the function value $f(\cdot)$ on two adjacent datasets
$\mathcal{D}$ and $\mathcal{D}^{\prime}$:

$\displaystyle\Delta_{2}(f)=\max_{\mathcal{D},\mathcal{D}^{\prime}}\|f(\mathcal{D})-f(\mathcal{D}^{\prime})\|_{2},$
(2)

where $\mathcal{D}$ and $\mathcal{D}^{\prime}$ are two adjacent datasets that
differ by only one element.

So far, following the core definition of DP, various privacy concepts and
perturbation methods have been proposed, such as local DP [40], f-DP [41],
output/objective perturbation [42]. Thanks to these data privacy techniques,
privacy-preserving machine learning has achieved enormous success in a large
number of scientific and engineering areas.

As we study the sparse linear regression model, thus, in following, we mainly
review the recent developments of privacy-preserving sparse regression methods
and put our work into a comprehensive comparison. The authors of [3] studied
the trade-off between statistical accuracy and privacy for sparse linear
regression in both low- and high-dimensional cases. They showed that based on
$\ell_{2}$ norm measurement, the optimal private minimax risk is lower bounded
by $\Omega(\sqrt{{s\log p}/{N}}+{s\log p}/{N\epsilon})$ for the high-
dimensional case and $\Omega(\sqrt{p/N}+{p}/{N\epsilon})$ for the low-
dimensional case. In [43], the authors studied the differential private sparse
linear model from the perspective of algorithm design. To achieve the desired
$(\epsilon,\delta)$-DP, they proposed the DP-IGHT algorithm leveraging
recently developed $\ell_{0}$ hard-thresholding and gradient perturbation.
They showed that their algorithm enjoys a linear convergence rate and attains
the lower bound as in [3]. [2] studied the sparse linear model under local
differential privacy. They first showed that the polynomial dependency on $p$
in the lower bound of estimation error can not be avoidable if we need to
protect the whole data privacy, and provided a private algorithm named LDP-IHT
for low-dimensional regression that achieves near-optimal upper bound as
$O({p\log p\log n}/{N\epsilon^{2}})$. Then they established the results that
when only responses’ privacy needs to be protected, their LDP-IHT algorithm
could achieve a better error upper bound as $O(\sqrt{s\log
p}/\sqrt{N}\epsilon)$. However, these works focused on the conventional OLS
loss and did not consider the model robustness. [44] studied differentially
private $l_{1}$-norm linear regression without sparse constraint on the
regression coefficient. Furthermore, they consider a truncated loss, which
introduces an additional approximation error of $O(\sqrt{p/N})$, which is not
negligible when $p$ is at the same order of $N$. [45] extended to enable
sparsity recovery by using a truncated square loss. In our work, we study the
robust sparse linear regression problem. Instead of an OLS loss, we consider
the median loss for linear model fitting. To our best knowledge, differential
private sparse LAD regression is still an open problem and there are few
theoretical analyses on this problem. Table 1 summarizes the related works.
The last column “Oracle rate” indicates whether the work proved that as the
privacy budget $\epsilon$ goes to infinity, the statistical accuracy of the
model weights would achieve the oracle rate $\sqrt{s\log p/N}$.

Table 1: Comparison of different ($\epsilon$, $\delta$)-DP algorithms for
linear regression.

Method | Privacy-accuracy trade-off | Sparsity | Algorithm not | Robustness | Non-smooth | Oracle
---|---|---|---|---|---|---
recovery |  depending on $s$ | loss | rate
[3] | $O\\{s\log p\log N[\log(1/\delta)\log N]^{1/2}/N\epsilon\\}^{*}$ | ✔ | ✗ | ✗ | ✗ | ✔
[43] | $O\\{s[\log(1/\delta)\log p]^{1/2}/N\epsilon\\}^{*}$ | ✔ | ✗ | ✗ | ✗ | ✔
[44] | $O\\{[s^{2}\log N\log(1/\delta)/(N\epsilon)]^{1/2}\\}^{{}^{\prime}}$ | ✗ | – | ✔ | ✗ | ✗
[45] | $O\\{s^{2}\log N\log p^{2}\log(1/\delta)/(N\epsilon)\\}^{{}^{\prime}}$ | ✔ | ✗ | ✔ | ✗ | ✗
Ours | $O\\{[p\log(1/\delta)\log(N\epsilon)]^{1/2}/N\epsilon\\}^{*}$ | ✔ | ✔ | ✔ | ✔ | ✔

  * ∗ Trade-off measured on weights. ${}^{{}^{\prime}}$ Trade-off measured on loss.

##  3 Problem and Algorithm

###  3.1 Problem Formulation

The motivation of the paper is to develop a fast and privacy-preserving
algorithm for the sparse and robust regression problem. More specifically, we
consider the following classical linear model:

$\mathbf{y}=\mathbf{X}\bm{\beta}^{*}+\mathbf{e},$

where $\mathbf{y}=(y_{1},y_{2},\cdots,y_{N})^{\top}$ is the vector of
responses,
$\mathbf{X}=(\mathbf{x}_{1},\mathbf{x}_{2},\cdots,\mathbf{x}_{N})^{\top}$ is
the $N\times p$ feature matrix,
$\mathbf{x}_{i}=(x_{1,i},\cdots,x_{p,i})^{\top}$ is a $p$-dimensional
covariate vector, $\bm{\beta}^{*}=(\beta_{1}^{*},\cdots,\beta_{p}^{*})^{\top}$
is the true sparse regression weight, and
$\mathbf{e}=(e_{1},e_{2},\cdots,e_{N})^{\top}$ is the vector of random errors.
We define $S=\\{1\leq i\leq p:\beta_{i}^{*}\neq 0\\}$ as the support of
$\bm{\beta}^{*}$ and $s=\mbox{$\lvert S\rvert$}$ as the sparsity level. Assume
that the random error has a density function $f$ and is independent of the
covariate vector. It is common to assume that there exists a positive constant
$c_{\bm{\beta}}$ such that for the true regression weight,
$\|\bm{\beta}^{*}\|_{2}\leq c_{\bm{\beta}}$ such as in [3]. Moreover, we
assume that the covariate $\mathbf{x}$ has bounded $\ell_{2}$ norm, which is
common for sparse and robust regressions ([43, 3]).

Note that we have quite a slack assumption for the random error, which allows
heavy-tailed random errors, i.e., abnormal responses such that the variance
could be infinite, such as Cauchy distribution. For such a heavy-tailed error,
the least square based methods are not applicable since the finite second-
moment assumption does not hold. Our error assumption is looser than the
state-of-the-art privacy-preserving algorithms designed for high-dimensional
linear regression such as [43, 3, 45], where finite second or fourth moments
of the error is required.

To tackle the heavy-tailed noise problem, we consider the least absolute
deviation regression, i.e., the regression weight $\bm{\beta}^{*}$ could be
obtained by minimizing the median loss,

$\arg\min_{\bm{\beta}\in\mathbf{R}^{p}}\frac{1}{N}\|\mathbf{y}-\mathbf{X}^{\top}\bm{\beta}\|_{1}.$

Furthermore, to pursue the sparsity structure of $\bm{\beta}^{*}$, we adopt a
regularization term on the weights and obtain the sparse and robust estimated
weights as

$\displaystyle\widehat{\bm{\beta}}=\arg\min_{\bm{\beta}\in\mathbf{R}^{p}}\frac{1}{N}\|\mathbf{y}-\mathbf{X}^{\top}\bm{\beta}\|_{1}+\lambda\|\bm{\beta}\|_{1},$
(3)

where $|\cdot|_{1}$ is the $\ell_{1}$-regularization term and $\lambda>0$ is
the regularization parameter. [43, 3, 45] utilized the $\ell_{0}$ hard-
threshold to choose non-zero weights. However, in the simulation, we show that
the $\ell_{1}$-regularization performs more stable than the $\ell_{0}$ hard-
threshold method. In this work, we aim to develop a fast and privacy-
preserving algorithm to solve the sparse LAD problem $\bm{\beta}$ in (3).

###  3.2 Proposed Algorithm

In this section, we propose our FRAPPE algorithm to find the sparse and robust
estimation for $\bm{\beta}^{*}$. Note that the sparse LAD problem has a non-
smooth loss as well as a non-smooth penalty term. To our knowledge, the doubly
non-smooth loss function is typically solved by subgradient descent algorithms
[15, 16]. However, subgradient descent algorithms suffer from a slow
convergence rate. To speed up the computation, we consider transforming the
median loss into a square loss.

Our FRAPPE algorithm has a double-loop structure. In the outer loop, we
construct pseudo responses such that the minimizer of the mean square error
with the pseudo responses can be seen as an update in the Newton-Raphson
iteration to solve the non-smooth problem (3). Thus, the $\ell_{1}$ penalized
LAD loss is transformed as $\ell_{1}$ penalized square loss. In the inner
loop, we solve the transformed $\ell_{1}$ penalized square loss privately.

Given a data-driven initial regression weight, which will be introduced later,
all the parameters are updated iteratively. Without loss of generality, we
introduce the proposed updating steps in the $v$th loop. At the $v$th loop,
given the current estimation $\widehat{\bm{\beta}}_{v}$, we first estimate the
density of the error $\\{e_{i}\\}_{i=1}^{N}$. The density estimation can be
easily achieved by kernel estimation

$\displaystyle\widehat{f}_{v}(0)=\frac{1}{Nh_{v}}\sum_{i=1}^{N}K\Big{(}\frac{y_{i}-\mathbf{x}_{i}^{\top}\widehat{\bm{\beta}}_{v}}{h_{v}}\Big{)},$
(4)

where $\widehat{f}_{v}(0)$ denotes the estimated density evaluated at point
zero, $K(\cdot)$ is a kernel function and $h_{v}$ is the bandwidth at $v$-th
iteration. To guarantee differential privacy, we add Gaussian noise to
$\widehat{f}_{v}(0)$:

$\displaystyle\widehat{f}^{\prime}_{v}(0)=\max\\{\widehat{f}_{v}(0)+u_{\widehat{f}_{v}},c_{f(0)}^{-1}\\},$
(5)

where $u_{\widehat{f}_{v}}\sim N(0,\sigma^{2}_{\widehat{f}_{v}})$ and
$c_{f(0)}$ is some positive constant to ensure that
$\widehat{f}^{\prime}_{v}(0)$ would be bounded away from zero. Then, we
construct the pseudo responses $\widetilde{\mathbf{y}}_{i,v}$ for the $v$th
outer loop

$\displaystyle\widetilde{y}_{i,v}=\mathbf{x}_{i}^{\top}\widehat{\bm{\beta}}_{v}-\big{(}\widehat{f}^{\prime}_{v}(0)\big{)}^{-1}(\mathds{I}[y_{i}\leq\mathbf{x}_{i}^{\top}\widehat{\bm{\beta}}_{v}]-1/2).$
(6)

With the pseudo responses $\widetilde{\mathbf{y}}_{i,v}$, we reformulate the
estimation problem as an $\ell_{1}$ penalized least square estimation problem

$\widetilde{\bm{\beta}}_{v}=\arg\min_{\bm{\beta}\in\mathcal{B}}H_{v}(\bm{\beta})+\lambda_{v}\|\bm{\beta}\|_{1},$

where $\mathcal{B}=\\{\bm{\beta}\in\mathbb{R}^{p}:\|\bm{\beta}\|_{2}\leq
c_{\bm{\beta}}\\}$ and
$H_{v}(\bm{\beta})=\frac{1}{2N}\sum_{i=1}^{N}(\widetilde{y}_{i,v}-\mathbf{x}_{i}^{\top}\bm{\beta})^{2}$.
It has been shown that the solution $\widetilde{\bm{\beta}}_{v}$ can be seen
as the $v$th update in the Newton-Raphson iteration to solve problem (3). We
say that steps (4-6) are in the outer loop to reformulate the median loss as a
square loss.

In the inner loop iterations, we solve the reformulated penalized least square
problem with privacy constraint. We adopt the method of differential private
gradient descent with soft thresholding. At $t$th inner loop, we first update
the regression weight as

$\displaystyle\widehat{\bm{\beta}}_{v,t-\frac{1}{2}}=\text{soft}\big{(}\widehat{\bm{\beta}}_{v,t-1}-\eta(\nabla
H_{v}(\widehat{\bm{\beta}}_{v,t-1})+\mathbf{u}_{\widehat{\bm{\beta}}_{v,t}},\lambda_{v}\eta\big{)},$
(7)

where $\mathbf{u}_{\widehat{\bm{\beta}}_{v,t}}\sim
N(0,\sigma_{\widehat{\bm{\beta}}_{v,t}}^{2}\mathbf{I}_{p}))$ and the soft
operator is defined as $\text{soft}(g,\tau):=\text{sign}(g)(|g|-\tau)_{+}$ and
$(a)_{+}=a\cdot\mathds{I}[a>0]$. Then we do the weight clipping with a
positive constant $c_{\bm{\beta}}$

$\displaystyle\widehat{\bm{\beta}}_{v,t}=\frac{c_{\bm{\beta}}}{\max\\{c_{\bm{\beta}},\|\widehat{\bm{\beta}}_{v,t-\frac{1}{2}}\|_{2}\\}}\cdot\widehat{\bm{\beta}}_{v,t-\frac{1}{2}}.$
(8)

Note that, unlike most of the existing work with gradient clipping, we perform
the clipping on model weight, which projects the weight into a compact set.
According to our problem formulation, the true model weight is bounded by a
constant $c_{\bm{\beta}}$, i.e. $\|\bm{\beta}\|\leq c_{\bm{\beta}}$, so the
above clipping would not make estimation biased. We summarize our algorithm in
Algorithm 1.

Algorithm 1 Fast Robust And Privacy-Preserving Estimation (FRAPPE) algorithm
for least absolute deviation regression

0: Data $\\{y_{i},\mathbf{x}_{i}\\}_{i=1}^{N}$, step size $\eta$, noise scales
$\sigma^{2}_{\widehat{\bm{\beta}}_{0}}$, $\sigma^{2}_{\widehat{f}_{v}}$,
$\sigma^{2}_{\widehat{\bm{\beta}}_{v,t}}$.

1: Initialize regression weight $\widehat{\bm{\beta}}_{0}$, e.g. (9), and
inject the noise
$\widehat{\bm{\beta}}_{1}=\widehat{\bm{\beta}}_{0}+\mathbf{u}_{\widehat{\bm{\beta}}_{0}}$,
where $\mathbf{u}_{\widehat{\bm{\beta}}_{0}}\sim
N(0,\sigma_{\widehat{\bm{\beta}}_{0}}^{2}\mathbf{I}_{p}))$;

2: for $v=1,\cdots,V$ do

3: Estimate the kernel density (4) and inject the noise (5);

4: Construct the pseudo response (6) and set up the weight
$\widehat{\bm{\beta}}_{v,0}=\widehat{\bm{\beta}}_{v}$;

5: for $t=1,\cdots,T$ do

6: Update the weight following (7) and (8);

7: end for

8: Set up the weight $\widehat{\bm{\beta}}_{v+1}=\widehat{\bm{\beta}}_{v,T}$;

9: end for

######  Remark 1.

We would like to preserve data privacy in each step where the data are used,
so it is necessary to add noise to all three steps. Although we add noise to
all three steps, we reduce the total noise addition by performing a loss
transformation. Note that, for the initialization, we only add noise once. For
kernel estimation, we add noise to every outer loop, and the total number of
outer loops is up to a constant, as given in Theorem 3. Thus, most noise is
added in the inner loop gradient descent step. For the inner loop part, we
transform the problem to a least square estimation problem, which can be
solved at a linear speed. Compared to a subgradient-based algorithm, which has
a sublinear convergence rate, our algorithm has fewer inner loop iterations.
Therefore, our algorithm could achieve a better privacy estimation trade-off.

Initialization: In this part, we discuss how to find an initial weight
$\widehat{\bm{\beta}}_{0}$ for Algorithm 1. Instead of the full dataset
$\\{y_{i},\mathbf{x}_{i}\\}_{i=1}^{N}$, we suggest using a subsampled dataset
$\mathcal{S}$ with size $n$. As only the $n$ subsampled data are involved, it
could speed up the computation for initialization and enhance data privacy
with privacy amplification ([46]). We consider estimating
$\widehat{\bm{\beta}}_{0}$ by solving the elastic net penalized LAD regression
problem

$\displaystyle\widehat{\bm{\beta}}_{0}=\arg\min\frac{1}{n}\sum_{(y_{i},\mathbf{x}_{i})\in\mathcal{S}}|y_{i}-\mathbf{x}_{i}^{\top}\bm{\beta}|+\lambda_{01}\|\bm{\beta}\|_{1}+\frac{\lambda_{02}}{2}\|\bm{\beta}\|_{2}^{2},$
(9)

as there exist many algorithms and software packages for it [47, 48]. Because
we do not specify an algorithm for solving the elastic net penalized LAD
regression problem, we use the output perturbation to guarantee DP.

##  4 Main Theoretical Results

In this section, we will first show that our proposed FRAPPE method satisfies
$(\epsilon,\delta)$-DP and then prove the statistical consistency of our
proposed estimator. The theoretical results also provide specific parameter
selections for the proposed FRAPPE algorithm. We start our analysis with the
following assumptions.

######  Assumption 1.

There exists a positive constant $c_{\mathbf{x}}$ such that
$\|\mathbf{x}\|_{2}\leq c_{\mathbf{x}}$.

######  Assumption 2.

The error $\mathbf{e}$ follows the density function $f(\cdot)$, which is
bounded and Lipschitz continuous (i.e.,$|f(x)-f(y)|\leq C_{L}|x-y|$ for any
$x,y\in\mathbb{R}$ and some constant $C_{L}>0$). Moreover, we assume that its
derivative $f^{\prime}(\cdot)$ is bounded and $f(0)>c_{f(0)}>0$ for some
constant $c_{f(0)}$.

######  Assumption 3.

Assume that the kernel function $K(\cdot)$ is integrable with
$\int_{-\infty}^{\infty}K(u)\mathrm{d}u=1$. Moreover, assume that $K(\cdot)$
satisfies $K(u)=0$ if $|u|>1$ and $K(u)\leq B$ for a positive constant $B$.
Further, assume $K(\cdot)$ is differentiable and its derivative
$K^{\prime}(\cdot)$ is bounded.

######  Assumption 4.

The total sample size $N$ satisfies $N=\Theta(p^{a+1})$ for some $a>0$. The
subsample size $n$ satisfies $n=\Theta(p^{r})$ for some $0<r<a+1/2$, and the
sparsity level $s$ satisfies $s=O(p^{b})$ for some $0<b<a/3+1/6$.

Assumption 1 is a regular algorithmic assumption. Assumption 2 is a
statistical assumption for the linear model, which assumes the smoothness of
the density function $f(\cdot)$ of the potentially unbounded errors.
Assumptions 3-4 are required for the statistical convergence. The Assumption 3
is a standard assumption for kernel estimation. Assumption 4 states
requirements for the order of dimension $p$, sample size $n$ to calculate the
initial estimator and the sparsity level $s$. The conditions of subsample size
$n$ and sparsity level $s$ are used to ensure that the proposed algorithm will
achieve the near-oracle convergence rate after a finite number of outer loop
iterations, as shown in equation (11).

######  Remark 2.

We considered the sparse LAD problem with $p<N$, as in Assumption 4, the only
requirement for dimension $p$ is that $N=\Theta(p^{a+1})$ for some $a>0$. This
requirement covers the scenario of $p\ll N$ as in the low-dimension regime. In
our framework, we allow $p$ to be in the same order as $N$, as long as the
transformed LSE problem is strongly convex. For example, in the context of
image analysis, our framework requires that the input for the proposed method
should have a pixel count lower than the total number of images, and the pixel
count could be in the same order as the number of images $N$. For the case
with $p>N$, we could add regularizers to make the problem strongly convex [49,
50].

###  4.1 Privacy Guarantee

In this part, we will show that our FRAPPE algorithm satisfies
$(\epsilon,\delta)$-DP with the proposed three-stage noise injection in
Theorem 1. The proof of Theorem 1 utilizes the $\ell_{2}$-sensitivity for each
stage and a standard composition theorem of zero- concentrated differential
privacy (zCDP) [51]. For the first stage, we achieve the privacy guarantee for
the data-driven initial point $\widehat{\bm{\beta}}_{0}$ via the output
perturbation, which does not depend on the initialization algorithm. As we set
up $\widehat{\bm{\beta}}_{0}$ with the elastic net penalized LAD loss, the
output perturbation is differential private by leveraging the
$\lambda_{02}$-strongly convexity of the loss [42]. Then for the second stage,
the DP of kernel density estimation $\widehat{f}_{0}$ follows from the fact
that the kernel function is bounded so that the $\ell_{2}$-sensitivity of the
kernel estimation is at the order of $O(1/N)$. In the last stage, due to the
bounded gradient $\nabla H_{v}(\widehat{\bm{\beta}}_{v,t})$, we adopt the
gradient perturbation method to achieve differential privacy. By properly
selecting the injected noises, each stage achieves $(\epsilon/3,\delta/3)$-DP.
Since we aim to protect data privacy over the three stages without loss of
generality, we use the even split for simplicity. More experiments for the
effect of split of privacy budget $\epsilon$ is given in the supplementary
material. Using the composition property of zCDP, we guarantee that our
proposed algorithm is $(\epsilon,\delta)$-DP.

######  Theorem 1.

Suppose Assumptions 1-3 hold, then Algorithm 1 is $(\epsilon,\delta)$-DP with
$\sigma^{2}_{\widehat{\bm{\beta}}_{0}}=\frac{24c_{\mathbf{x}}^{2}\log(n/(N\delta))}{\epsilon^{2}\lambda_{02}^{2}N^{2}}$,
$\sigma^{2}_{\widehat{f}_{v}}=\frac{24B^{2}\log(1/\delta)V}{\epsilon^{2}N^{2}h_{v}^{2}}$
and
$\sigma^{2}_{\widehat{\bm{\beta}}_{v,t}}=\frac{6G^{2}\log(1/\delta)TV}{\epsilon^{2}N^{2}}$,
where $G=4c_{x}^{2}c_{\beta}+c_{f(0)}$.

We can see that reducing the iteration numbers $V$ and $T$ can improve
privacy. Thus, the fast convergence of the proposed FRAPPE algorithm improves
not only computation efficiency but also data privacy.

###  4.2 Statistical Accuracy

In this section, we provide the statistical convergence rate for the proposed
estimator $\widehat{\bm{\beta}}_{v+1}$. First, we show the rate of
$\widehat{\bm{\beta}}_{2}$ after one iteration of the outer loop.

######  Theorem 2.

Suppose that Assumptions 1-4 hold. Let
$a_{1,N}=\sigma_{\widehat{\bm{\beta}}_{0}}\sqrt{p}\log N+\sqrt{s\log p/n}$ and
take $\lambda_{1}=C_{0}\left(\sqrt{\frac{\log
p}{N}}+a_{1,N}(\sigma_{\widehat{\bm{\beta}}_{0}}\sqrt{p}\log N+\sqrt{s\log
p/n})\right)$ with the constant $C_{0}$. For the estimator
$\widehat{\bm{\beta}}_{2}$ obtained after $T$ iterations of perturbed gradient
descent (7), if the learning rate is $1/2L$, and the number of iterations is
$T=O(\log_{(1-\mu/2L)}\big{(}\frac{pG^{2}\log(1/\delta)}{N^{2}\epsilon^{2}}\big{)})$,
we have

$\displaystyle\|\widehat{\bm{\beta}}_{2}-\bm{\beta}^{*}\|_{2}=O_{\mathbb{P}}$
$\displaystyle\left(a_{1,N}\sqrt{s}\left(\sigma_{\widehat{\bm{\beta}}_{0}}\sqrt{p}\log
N+\sqrt{s\log p/n}\right)+b_{N}+\sqrt{\frac{s\log p}{N}}\right),$ (10)

where $b_{N}=\sqrt{p\log(1/\delta)\log(N\epsilon)}G/(N\epsilon)$.

In the right-hand side of (10), the first term is from the error of the pseudo
response construction $\widetilde{y}_{i,1}$ and the kernel density estimation
$\widehat{f}^{\prime}_{1}(0)$. The second term $b_{N}$ is introduced by the
perturbed gradient descent method. The third term is the error bound of the
sparse least square estimation.

######  Remark 3.

The error bound of the initial estimator is $\sqrt{s\log p/n}$. Note that
$\sigma_{\widehat{\bm{\beta}}_{0}}=O(1/N)$. By Assumption 4, we have
$\sigma_{\widehat{\bm{\beta}}_{0}}\sqrt{p}\log N=O(\log p/p^{a+1/2})$ and
$\sqrt{s\log p/n}=\Theta(\sqrt{\log p/p^{(1-r)b}})$, which implies
$\sigma_{\widehat{\bm{\beta}}_{0}}\sqrt{p}\log N$ is smaller than $\sqrt{s\log
p/n}$. Similarly, we have $a_{1,N}\sqrt{s}=\Theta(\log
p/p^{a+1/2-rb/2}+\sqrt{\log p/p^{(1-2r)b}})$, which is smaller than one.
Therefore, the initial error $\sqrt{s\log p/n}$ gets shrunk. The term
$b_{N}=O(\sqrt{\log p/p^{2a+1}})$ and the last term $\sqrt{s\log
p/N}=O(\sqrt{\log p/p^{a-rb+1}})$ are in an smaller order of $\sqrt{s\log
p/n}$.

By recursively applying Theorem 2, we provide the convergence rate for
$\widehat{\bm{\beta}}_{v+1}$. Let us define, for $1\leq v\leq V$,

$\displaystyle a_{v,N}=$ $\displaystyle\sqrt{\frac{s\log
p}{N}}+b_{N}+\sqrt{s}(\sigma_{\widehat{\bm{\beta}}_{0}}\sqrt{p}\log
N+\sqrt{s\log p/n})^{v}.$

Then, we have the rate for $\widehat{\bm{\beta}}_{v+1}$ shown in the following
theorem.

######  Theorem 3.

Assume that the initial estimator $\widehat{\bm{\beta}}_{0}$ satisfies
$\|\widehat{\bm{\beta}}_{0}-\bm{\beta}^{*}\|_{2}=O_{\mathbb{P}}(\sqrt{s\log
p/n})$. Let $h_{v}=\Theta(a_{v,N})$ for $1\leq v\leq V-1$, and take

$\displaystyle\lambda_{v,N}=C_{0}\left(\sqrt{\frac{\log
p}{N}}+a_{v,N}(\sigma_{\widehat{\bm{\beta}}_{0}}\sqrt{p}\log N+\sqrt{s\log
p/n})\right)$

with $C_{0}$ being a sufficiently large constant. Under Assumptions 1-4, with
inner loop learning rate $1/2L$ and inner loop iteration
$T=O(\log_{(1-\mu/2L)}\big{(}\frac{pG^{2}\log(1/\delta)}{N^{2}\epsilon^{2}}\big{)})$,
we have

$\displaystyle\|\widehat{\bm{\beta}}_{v+1}-\bm{\beta}^{*}\|_{2}=O_{\mathbb{P}}$
$\displaystyle\left(\sqrt{\frac{s\log
p}{N}}+b_{N}+\sqrt{s}(\sigma_{\widehat{\bm{\beta}}_{0}}\sqrt{p}\log
N+\sqrt{s\log p/n})^{v+1}\right),$

where $b_{N}=\sqrt{p\log(1/\delta)\log(N\epsilon)}G/(N\epsilon)$.

######  Remark 4.

We give a fast way to find an initial estimator that satisfies the assumption
of initialization in Theorem 3. The error bound of initialization is
$\sqrt{s\log p/n}$, where $n$ could be much smaller than $N$. Then, Theorem 3
shows that under Assumption 4, after a constant number of iterations

$\displaystyle V\geq\frac{2\log(N/\log p)}{\log(n/s\log p)},$ (11)

the estimation accuracy will be improved and eventually achieve a rate
$\sqrt{s\log p/N}+b_{N}$, and $\sqrt{s\log p/N}$ matches the oracle
convergence rate up to a logarithmic factor and the term $b_{N}$ is due to the
privacy accuracy trade-off. Compared with the initialization, the error bound
is almost improved by a rate of $\sqrt{n/N}$.

######  Remark 5.

The privacy cost of our proposed estimator is of the order of $O(\sqrt{p})$,
which polynomially depends on the dimensionality $p$. This order is
intrinsically caused by the soft thresholding for the $\ell_{1}$ penalty.
Under the worse scenario, the soft thresholding operator would keep all the
elements of $\hat{\bm{\beta}}_{v,t-\frac{1}{2}}$ as non-zero, so that the
injected $p$-dimensional noises would be accumulated during the estimation.
Some previous works achieved the privacy cost at the order of $O(s\log p)$ by
using peeling algorithm [3] or hard thresholding operation [43], which only
remain a fixed number of the largest elements in regression weight during
optimization. However, these methods require prior knowledge of the sparsity
level and extra convexity condition on the loss, which cannot easily be
extended to our case. As it is a non-trivial task to improve the privacy cost,
we would like to leave it as a future study.

##  5 Numerical Evaluation

In this section, we empirically evaluate the statistical and computational
efficiency of the proposed FRAPPE algorithm on both synthetic data and real
datasets.

###  5.1 Synthetic Data Experiments

Data genaration: We generate data from the linear model:
$y_{i}=\mathbf{x}_{i}^{\top}\bm{\beta}^{*}+e_{i}$, $i=1,\cdots,N$, where
$\mathbf{x}_{i}$ is a $p$-dimensional vector,
$\bm{\beta}^{*}=(\beta_{1}^{*},\cdots,\beta_{p}^{*})^{\top}$ is the true
regression weight, and $e_{i}$ is the random noise. The vectors
$\\{\mathbf{x}_{i}\\}$ are drawn i.i.d. from multivariate normal distribution
with mean zeros and covariance matrix $\Sigma$, which is constructed as
$\Sigma_{ij}=0.1^{|i-j|}$, for $1\leq i,j\leq p$. Let $s$ be the sparsity
level, and then the true regression weight is set as
$\bm{\beta}^{*}=s^{-1}10\cdot(1,2,\cdots,s,0,,\cdots,0)$. The random noise
$e_{i}$ is i.i.d. sample from some given distributions.

Compared algorithms: We compare the proposed FRAPPE algorithm with five other
state-of-art algorithms designed for private sparse linear regression:

  * •

SgpLAD: using subgradient perturbation descent (SGP) algorithm to solve the
$\ell_{1}$ penalized least absolute deviation loss (LAD). The subgradient
descent method is obtained from [15]. To achieve
$(\epsilon,\delta)$-differential privacy, we add random Gaussian noise to the
subgradient in each iteration of Algorithm 1 in [15].

  * •

GpLASSO: using gradient perturbation (GP) to minimize $\ell_{1}$ penalized
square loss. The minimizer of $\ell_{1}$ penalized square loss is known as the
Lasso estimator. To keep privacy, random Gaussian noise is added to the
gradient in each iteration when solving the Lasso problem.

  * •

DPSLR: Differentially Private Sparse Linear Regression given in Algorithm 4.2
of [3]. This estimator is designed to minimize the square loss. To recover the
sparsity of the regression weights, the “peeling” algorithm, proposed by [52],
is used to select $s$ weights with the top-$s$ largest absolute values in a
differentially private way.

  * •

HtPSLR: Heavy-tailed Private Sparse Linear Regression given in Algorithm 3 of
[45]. First, they truncate the original data at a given level to avoid heavy-
tailed data. Then, gradient descent and ”peeling” are performed to obtain a
sparse regression weight.

  * •

DPIGHT: Differentially Private Iterative Gradient Hard Thresholding proposed
by [43] and the specific algorithm for the sparse linear regression problem is
given in section 5.2 of [43]. This estimator aims to minimize the square loss.
In each iteration, the weights are updated by perturbed gradient and then
perform the hard thresholding operator.

Simulation setup: The initial estimator is obtained by solving (9) with a
randomly selected sample with size $n=200$. We use a biweight kernel function

$K(x)=\left\\{\begin{array}[]{ll}0,&x\leq-1,\\\
-\frac{315}{64}x^{6}+\frac{735}{64}x^{4}-\frac{525}{64}x^{2}+\frac{105}{64},&-1\leq
x\leq 1,\\\ 0,&x\geq 1\end{array}\right.$

It is easy to verify that the biweight kernel function satisfies Assumption 3.
Other common kernel functions, such as uniform, triweight, and Epanechnikov,
also satisfy Assumption 3 and share similar experimental results, which are
shown in the experiment for the effect of kernel in the supplement material.
For the choice of bandwidth, as shown in Theorem 3, it should have the same
order as $a_{v,N}$, and thus we set the bandwidth as $h_{v}=\sqrt{\frac{s\log
N}{N}}+s^{-1/2}\left(0.9\right)^{(v+1)/2}$. For the
$(\epsilon,\delta)$-differential privacy, we set $\delta=10^{-3}$.
Additionally, we set the inner loop iteration $T=50$ and the outer loop
iteration $V=10$. To give a fair comparison, we keep the same number of total
iterations for other methods, i.e., the total number of iterations for other
methods is 500.

We specify 20 values for the regularization parameter $\lambda$ according to
the input data for the soft-thresholding based methods FRAPPE, SgpLAD, and
GpLASSO. Also, 20 values for the sparsity level $s$ are specified for the
hard-thresholding based methods DPSLR, HtPSLR, and DPIGHT. The $\lambda$ and
$s$ with the smallest Bayesian information criterion (BIC) value are chosen
for each algorithm, respectively.

Performance measurements: In the experiment with the synthetic dataset, we
evaluate the performance of the six methods in terms of the following two
measurements:

  1. 1.

Mean square error (MSE) of the estimated regression weights to examine the
model accuracy.

  2. 2.

$F_{1}$-score to evaluate the sparsity recovery. The $F_{1}$-score is defined
as

$F_{1}=\left(\frac{\text{recall}^{-1}+\text{precision}^{-1}}{2}\right)^{\\!-1},$

where recall is the fraction of correctly identified non-zero weights over
true sparsity and precision is the fraction of correctly identified non-zero
weights over total identified non-zero ones. Note that the $F_{1}$ score
ranges between 0 and 1, and the closer to 1, the better the sparsity recovery.

Effect of Heavy-Tailed Noises:

The proposed FRAPPE algorithm is designed for data contaminated with heavy-
tailed noises and extreme values, which are pervasive in the real world. In
this set of experiments, we aim to compare the practical performance of the
proposed FRAPPE algorithm with the five other state-of-art differentially
private sparse linear regression algorithms under different choices of noise
distributions. We fix the dimension of $\mathbf{x}_{i}$ at $p=100$, sparsity
level $s=10$, and privacy budget $\epsilon=0.5$. We consider the following
three noise distributions:

  * •

Normal: the noise $e_{i}\sim N(0,1)$. Normal distribution has a light tail and
finite moments.

  * •

Student’s t: the noise $e_{i}\sim t(2)$. Student’s t distribution has a tail
between Normal and $\mathrm{Cauchy}$ distributions, since $t(1)$ is equivalent
to $\mathrm{Cauchy}$ distribution and $t(\infty)$ is Normal distribution.

  * •

$\mathrm{Cauchy}$: the noise $e_{i}\sim\mathrm{Cauchy}(0,1)$.
$\mathrm{Cauchy}$ distribution has a heavy tail and infinite variance.

Figure 1 shows the probability densities for the three noise distributions. We
can easily observe that $\mathrm{Cauchy}$ distribution has a heavier tail than
t(2), and N(0,1) has the lightest tail behavior.

Figure 1: The probability densities for the three noise distributions.
$\mathrm{Cauchy}$ distribution has a heaviest tail and N(0,1) has the lightest
tail.

Table 2 shows the MSE and $F_{1}$-score of all methods under three different
noise distributions. Remember that FRAPPE and SgpLAD are based on LAD loss,
and the other four methods are based on square loss. For Normal noise, the
best methods are FRAPPE and DPIGHT, which indicates that under light-tailed
noise, the LAD loss and square loss have similar performance. SgpLAD is worse
than FRAPPE due to its slow convergence, especially when the data size
decreases. DPIGHT outperforms GpLASSO, DPSLR, and HtPSLR because its iterative
gradient hard thresholding procedure is more efficient than soft thresholding
and “peeling” procedure for sparse learning.

For the Student’s t(2) noise, FRAPPE, HtPSLR, and DPIGHT have similar
performance. Because HtPSLR and DPIGHT employ data trimming to adopt abnormal
observations and this method works for noise with finite second moment.
GpLASSO and DPSLR perform worst as they are not designed for heavy-tailed
noise.

For heavy-tailed $\mathrm{Cauchy}$ noise, the proposed FRAPPE method has the
best performance among all methods. The heavy-tailed noise deteriorates the
estimation performance of the methods based on square loss. Because Heavy-
tailed noise breaks the finite second-moment assumption for least square
methods. HtPSLR shrunk the original data to avoid extreme values, which leads
to a better performance than the no-shrinkage algorithm DPSLR, but HtPSLR
still needs a finite second moment assumption that does not hold for
$\mathrm{Cauchy}$ noise. More precisely, FRAPPE reduces at least about 63.6%
of MSE of the square loss based methods.

Note that as the total data size $N$ increases, the MSEs of all methods
decrease, and as the dimension of the covariates increases, the MSEs of all
methods increase, which matches the oracle convergence error $\sqrt{s\log
p/N}$. Moreover, when the data size $N$ gets larger or the dimension $p$ gets
smaller, the MSE of SgpLAD increases much faster than that of FRAPPE. Because
SgpLAD converges much slower than FRAPPE, and we will show that in the later
section.

Thus, the consistent performance of the FRAPPE algorithm confirms that the
proposed FRAPPE is robust to both light-tailed and heavy-tailed noises.

Table 2: MSE and $F_{1}$ score for method comparison. We report results under three different noise distributions and sample sizes $N$. We fix privacy budget $\epsilon$ as 0.5, sparsity level $s$ as 10, and dimension $p$ as 100. The bold values under the same setting (row) are the two smallest MSEs. Noise | $N$ | FRAPPE | SgpLAD | GpLASSO | DPSLR | HtPSLR | DPIGHT
---|---|---|---|---|---|---|---
MSE | $F_{1}$ | MSE | $F_{1}$ | MSE | $F_{1}$ | MSE | $F_{1}$ | MSE | $F_{1}$ | MSE | $F_{1}$
Normal | 2000 | 0.05 | 0.91 | 6.15 | 0.64 | 1.76 | 0.20 | 2.25 | 0.75 | 0.10 | 0.93 | 0.05 | 0.98
5000 | 0.01 | 0.92 | 0.08 | 0.91 | 0.32 | 0.21 | 0.87 | 0.85 | 0.03 | 0.96 | 0.01 | 0.99
10000 | 0.01 | 0.95 | 0.01 | 0.96 | 0.09 | 0.21 | 0.90 | 0.83 | 0.02 | 0.96 | 0.01 | 1
Student’s t | 2000 | 0.31 | 0.96 | 5.11 | 0.63 | 1.88 | 0.22 | 2.31 | 0.75 | 0.37 | 1.00 | 0.22 | 0.99
5000 | 0.18 | 0.96 | 0.24 | 0.96 | 0.99 | 0.63 | 0.88 | 0.92 | 0.22 | 1.00 | 0.13 | 1.00
10000 | 0.12 | 0.96 | 0.12 | 0.97 | 0.44 | 0.91 | 0.93 | 0.96 | 0.17 | 1.00 | 0.10 | 1.00
$\mathrm{Cauchy}$ | 2000 | 0.44 | 0.99 | 5.91 | 0.66 | 5.50 | 0.55 | 2.83 | 0.82 | 1.45 | 0.96 | 1.19 | 0.91
5000 | 0.23 | 0.98 | 0.32 | 0.99 | 2.39 | 0.96 | 1.08 | 0.99 | 0.95 | 0.99 | 0.85 | 0.92
10000 | 0.15 | 0.98 | 0.14 | 0.98 | 1.08 | 0.98 | 0.96 | 1.00 | 0.62 | 1.00 | 0.51 | 0.93
Table 3: MSE and $F_{1}$ score for method comparison. We report results under three different noise distributions and dimensions $p$. We fix privacy budget $\epsilon$ as 0.5, sparsity level $s$ as 10, and data size $N$ as 5000. The bold values under the same setting (row) are the two smallest MSEs. Noise | $p$ | FRAPPE | SgpLAD | GpLASSO | DPSLR | HtPSLR | DPIGHT
---|---|---|---|---|---|---|---
MSE | $F_{1}$ | MSE | $F_{1}$ | MSE | $F_{1}$ | MSE | $F_{1}$ | MSE | $F_{1}$ | MSE | $F_{1}$
Normal | 50 | 0.01 | 0.91 | 0.01 | 0.96 | 0.16 | 0.34 | 0.85 | 0.83 | 0.03 | 0.93 | 0.01 | 0.99
100 | 0.01 | 0.91 | 0.09 | 0.92 | 0.32 | 0.20 | 0.83 | 0.82 | 0.02 | 0.95 | 0.01 | 1.00
200 | 0.02 | 0.90 | 0.25 | 0.96 | 0.65 | 0.12 | 0.85 | 0.81 | 0.03 | 0.94 | 0.01 | 0.99
Student’s t | 50 | 0.16 | 0.97 | 0.16 | 0.97 | 0.46 | 0.61 | 0.91 | 0.91 | 0.24 | 1.00 | 0.13 | 1.00
100 | 0.18 | 0.96 | 0.24 | 0.96 | 0.99 | 0.63 | 0.88 | 0.92 | 0.22 | 1.00 | 0.13 | 1.00
200 | 0.21 | 0.97 | 0.69 | 0.98 | 1.26 | 0.30 | 0.89 | 0.88 | 0.24 | 1.00 | 0.14 | 0.99
$\mathrm{Cauchy}$ | 50 | 0.19 | 0.98 | 0.21 | 0.98 | 1.19 | 0.81 | 1.07 | 0.99 | 0.98 | 0.99 | 0.71 | 0.95
100 | 0.22 | 0.99 | 0.27 | 0.98 | 2.32 | 0.96 | 1.06 | 0.99 | 0.95 | 0.99 | 0.81 | 0.92
200 | 0.25 | 0.99 | 1.05 | 0.95 | 2.49 | 0.69 | 1.24 | 0.96 | 1.04 | 0.98 | 0.88 | 0.91

Effect of the Sparsity Level: In this set of experiments, we evaluate the
performance of all algorithms under different levels of sparsity. The oracle
convergence error is $\sqrt{s\log p/N}$, and thus the MSE should increase with
the sparsity level.

We fix sample size $N=5000$, dimension $p=100$, and privacy budget
$\epsilon=0.5$. The results of MSE versus sparsity level $s$ under different
noise distributions are reported in Figure 2. We observe that the MSE
increases as the sparsity level $s$ increases for all methods except for
GpLASSO under heavy-tailed loss. Because the heavy-tailed data destroy the
statistical convergence theory for GpLASSO, the GpLASSO does not follow the
oracle property anymore. The MSE of SGPLAD has a faster increase rate than the
oracle rate due to its slow convergence when $s$ gets larger. We can see that
when the sparsity level gets larger than 20, the MSE of DPSLR, HtPSLR and
DPIGHT have a faster increase rate than the oracle. This finding indicates the
unstable performance of hard thresholding when the number of non-zero
regression weights is large.

Note that the proposed FRAPPE has the best performance among all the methods
under all noise types and sparsity levels. Also note that the outer loop
iteration numbers should increase with the sparsity $s$ to achieve the same
convergence rate as shown in Theorem 3. The results are based on 10 outer loop
iterations, and hence, FRAPPE only requires a small number of outer loop
refinements to achieve its theoretical statistical error. Thus, overall,
FRAPPE performs consistently well in spite of light-tailed or heavy-tailed
noise and different sparsity level.

|  |
---|---|---
(a) Normal noise | (b) Student’s t noise | (c) $\mathrm{Cauchy}$ noise
Figure 2: The MSE vs the sparsity level $s$ ranging from 1 to 30. The three figures are corresponding to different noise distributions under sample size $N=5000$, dimension $p=100$, and privacy budget $\epsilon=0.5$. Compared with the existing algorithms, our method achieved the best performance. |  |
---|---|---
(a) Normal noise | (b) Student’s t noise | (c) $\mathrm{Cauchy}$ noise
Figure 3: The MSE vs the privacy budget $\epsilon$ ranging from 0 to 1. The
three figures are corresponding to different noise distributions under sample
size $N=5000$, dimension $p=100$, and sparsity level $s=10$. Compared with the
existing algorithms, our method achieved the best performance.

Effect of the Privacy Budget: In this set of experiments, we investigate the
performance of all algorithms under different privacy budgets $\epsilon$. The
differential private budget $\epsilon$ controls the privacy level, and the
smaller value indicates a more stringent privacy constraint. We also compare
the results obtained from the non-private version of the proposed FRAPPE
algorithm (Non-private).

We fix sample size $N=5000$, dimension $p=100$, and sparsity level $s=10$. The
results of MSE versus privacy budget $\epsilon$ under different noise
distributions are reported in Figure 3. We can observe the privacy and
statistical accuracy trade-off for all methods. For FRAPPE, the decrease speed
of MSE becomes slower as increasing $\epsilon$, which is consistent with our
theoretical result. As shown in Theorem 3, for FRAPPE, the order of MSE is
$O(\sqrt{s\log p/N}+\sqrt{p\log(N\epsilon)}/(N\epsilon))$. Note that the MSE
of HtPSLR and DPIGHT changes slightly as the privacy budget $\epsilon$
changes. There is no statistical accuracy of the regression weights for
HtPSLR. For DPIGHT, [43] showed that the order of MSE is $O(\sqrt{s\log
p/N}+s\sqrt{\log p}/(N\epsilon))$. The second term $s\sqrt{\log
p}/(N\epsilon)$ shows the privacy accuracy trade-off and is smaller than
FRAPPE one at the given setting in the experiment. However, except for the
smallest $\epsilon$, FRAPPE performs the best among all privacy budget
$\epsilon$ and different noise types.

Computation Time Comparison: In this section, we further study the computation
efficiency of our proposed estimator. The motivation of the proposed method is
to shorten the computation time for the LAD problem. Hence, we compare the
proposed FRAPPE algorithm with SgpLAD algorithm. We fix the privacy budget
$\epsilon=0.5$, sparsity level $s=10$, data size $N=5000$, and dimension
$p=100$. In Figure 4, we report the MSE versus computation time in seconds of
FRAPPE and SgpLAD with $\mathrm{Cauchy}$ noise. For the Normal and Student’s t
noise, the results are similar. We can observe that to achieve the same MSE,
the computation time for FRAPPE is only about half of that for SgpLAD.

|
---|---
(a) $N=2000$ | (b) $N=5000$
Figure 4: The MSE vs computation time. We report the estimation MSE as
computation time increases. The two figures correspond to different sample
sizes $N=2000$ and $N=5000$ under $\mathrm{Cauchy}$ noise, sparsity level
$s=10$, dimension $p=100$, and privacy budget $\epsilon=0.5$. Compared with
the SgpLAD algorithm, our method performed better.

Some more experiments evaluated the effect of inner loop and outer loop
iterations, the effect of different splits of the privacy budget, the effect
of the initial sample size, and the effect of different choices of the kernel
functions for density estimation are shown in the supplementary material.

###  5.2 Real Data Analysis

In this section, we study the performance of the FRAPPE algorithm on two real
datasets: the Communities and Crime dataset from the UCI Machine Learning
Repository [53] and the Ames Housing dataset, which is publicly available at
Kaggle.

Figure 5: The boxplots for the normalized responses in the training data for
Ames Housing and Communities and Crime dataset.

Dataset Description: For the Communities and Crime dataset, we analyze the
linear relationship between the demographic variables and the community crime
number and identify the significant demographic variables. In our analysis,
the response variable is the total number of violent crimes per 100K
population (ViolentCrimesPerPop). The final dataset includes 99 features and
1993 communities. For the Ames Housing dataset, we analyze the linear
relationship between the sale price of a house and the house’s attributes,
such as building class, number of rooms, and original and construction year.
The dataset includes 1121 observations, where each observation contains a sale
price and 37 attributes of the corresponding house. Figure 5 shows boxplots
for the normalized responses in the training data for the Ames Housing and
Communities and Crime dataset. We can observe that both data have skewed
distribution and that the Ames housing data is more spread than the
Communities and Crime data. So we may say that both data have heavy-tailed
data,

Estimation Results: We randomly select 80% of data points as the training
dataset, and we use the other 20% as the test dataset. We apply all six
methods to the training dataset and obtain the estimated model weights. The
parameter settings are the same as in the synthetic data experiment. For
evaluation, we calculate the mean square error (MSE) and mean absolute error
(MAE) between the predicted responses and the true responses in the test
dataset. We also report the number of non-zero estimated model weights as
Sparsity. Tables 4 and 5 give the MSE, MAE, and Sparsity for all six
estimators under different privacy budgets $\epsilon=0.10$, 0.15, 0.20, 0.25
and 0.30.

For the Ames Housing dataset, FRAPPE and SgpLAD outperform other square loss
based methods since the data has a heavy tail as shown in Figure 5. FRAPPE and
SgpLAD have similar results, however, FRAPPE took a shorter computation time
as shown in Figure 4. More specifically, the FRAPPE estimator reduces MSE by
at least about 22%, and MAE by at least about 31% compared with square loss
based methods.

For the Communities and Crime dataset, FRAPPE outperforms other methods as it
is robust and fast. For this dataset, the results are quite sensitive to the
initialization point, and hence the MSEs and MAEs for HtPSLR are not stable.

Hard-threshold methods (DPSLR, HtPSLR, DPIGHT) select few covariates and
suffer from a larger MSE. FRAPPE selects a reasonable number of covariates and
has the best model performance in terms of both MSE and MAE. The results show
that the proposed FRAPPE has the best performance in estimating as well as in
sparsity recovery.

Table 4: The testing results for the Ames Housing dataset. We report MSE, MAE, and sparsity under different privacy budgets $\epsilon$ ranging from 0.1 to 0.3. The bold values under the same setting (row) are the two smallest MSEs or MAEs. $\epsilon$ |  | FRAPPE | SgpLAD | GpLASSO | DPSLR | HtPSLR | DPIGHT
---|---|---|---|---|---|---|---
0.10 | MSE | 0.32 | 0.29 | 1.12 | 32.74 | 1.02 | 0.77
MAE | 0.32 | 0.30 | 0.71 | 3.98 | 0.62 | 0.60
Sparsity | 33 | 28 | 14 | 4 | 4 | 4
0.15 | MSE | 0.30 | 0.29 | 0.87 | 22.72 | 0.69 | 0.58
MAE | 0.30 | 0.30 | 0.62 | 3.15 | 0.52 | 0.52
Sparsity | 31 | 25 | 22 | 4 | 4 | 4
0.20 | MSE | 0.28 | 0.30 | 0.70 | 10.46 | 0.66 | 0.48
MAE | 0.29 | 0.30 | 0.54 | 2.21 | 0.49 | 0.46
Sparsity | 28 | 23 | 30 | 4 | 4 | 5
0.25 | MSE | 0.30 | 0.29 | 0.57 | 6.59 | 0.78 | 0.44
MAE | 0.30 | 0.30 | 0.48 | 1.79 | 0.51 | 0.44
Sparsity | 23 | 22 | 34 | 4 | 3 | 5
0.30 | MSE | 0.32 | 0.28 | 0.47 | 5.35 | 0.76 | 0.41
MAE | 0.32 | 0.29 | 0.43 | 1.58 | 0.52 | 0.42
Sparsity | 20 | 23 | 35 | 4 | 4 | 5
Table 5: The testing results for the Communities and Crime dataset. We report MSE, MAE, and sparsity under different privacy budgets $\epsilon$ ranging from 0.1 to 0.3. The bold values under the same setting (row) are the two smallest MSEs or MAEs. $\epsilon$ |  | FRAPPE | SgpLAD | GpLASSO | DPSLR | HtPSLR | DPIGHT
---|---|---|---|---|---|---|---
0.10 | MSE | 0.49 | 1.74 | 1.08 | 8.83 | 19.27 | 4.85
MAE | 0.44 | 0.54 | 0.71 | 2.65 | 2.16 | 0.78
Sparsity | 11 | 31 | 8 | 3 | 3 | 3
0.15 | MSE | 0.48 | 0.77 | 0.73 | 1.18 | 36.64 | 4.30
MAE | 0.44 | 0.55 | 0.57 | 1.55 | 3.02 | 0.82
Sparsity | 11 | 21 | 60 | 7 | 3 | 4
0.20 | MSE | 0.49 | 0.67 | 0.90 | 1.15 | 47.14 | 4.26
MAE | 0.45 | 0.52 | 0.63 | 0.75 | 3.46 | 0.81
Sparsity | 10 | 18 | 35 | 3 | 3 | 3
0.25 | MSE | 0.51 | 0.76 | 0.73 | 1.17 | 25.98 | 4.04
MAE | 0.46 | 0.55 | 0.57 | 0.76 | 2.37 | 0.76
Sparsity | 10 | 0 | 41 | 3 | 3 | 4
0.30 | MSE | 0.54 | 0.62 | 0.50 | 0.99 | 33.25 | 4.29
MAE | 0.47 | 0.50 | 0.46 | 0.67 | 2.91 | 0.83
Sparsity | 19 | 30 | 64 | 5 | 3 | 4

##  6 Conclusion

In this paper, we studied the robust LAD regression problem with sparsity and
privacy constraints. We proposed an efficient learning method named FRAPPE,
which includes surrogate loss reformulation and three-stage noise injection
for fast computation and privacy-preserving. We showed that with properly
selected parameters, our estimators achieve a trade-off between the near-
optimal rate of $\sqrt{s\log p/N}$ and
$\sqrt{p\log(1/\delta)\log(N\epsilon)}G/(N\epsilon)$. Extensive numerical
studies were conducted and verified that our proposed FRAPPE method
outperforms the state-of-the-art privacy-preserving learning algorithms for
robust sparse regression analysis. An interesting extension of our work is to
generalize the proposed algorithm to classification problems, such as logistic
regression.

## References

  * [1] C. Dwork, A. Roth _et al._ , “The algorithmic foundations of differential privacy.” _Foundations and Trends in Theoretical Computer Science_ , vol. 9, no. 3-4, pp. 211–407, 2014.
  * [2] D. Wang and J. Xu, “On sparse linear regression in the local differential privacy model,” in _Proceedings of the 36th International Conference on Machine Learning_ , vol. 97, 2019, pp. 6628–6637.
  * [3] T. T. Cai, Y. Wang, and L. Zhang, “The cost of privacy: Optimal rates of convergence for parameter estimation with differential privacy,” _The The Annals of Statistics_ , vol. 49, no. 5, pp. 2825–2850, 2021.
  * [4] D. Wang, A. Smith, and J. Xu, “High dimensional sparse linear regression under local differential privacy: Power and limitations,” in _2018 NIPS workshop in Privacy-Preserving Machine Learning_ , vol. 235, 2018.
  * [5] J. Ren, J. Xiong, Z. Yao, R. Ma, and M. Lin, “Dplk-means: A novel differential privacy k-means mechanism,” in _2017 IEEE Second International Conference on Data Science in Cyberspace (DSC)_ , 2017, pp. 133–139.
  * [6] Y. Mülle, C. Clifton, and K. Böhm, “Privacy-integrated graph clustering through differential privacy.” in _EDBT/ICDT Workshops_ , vol. 157, 2015.
  * [7] M. Abadi, A. Chu, I. Goodfellow, H. B. McMahan, I. Mironov, K. Talwar, and L. Zhang, “Deep learning with differential privacy,” in _Proceedings of the 2016 ACM SIGSAC Conference on Computer and Communications Security_ , 2016, pp. 308–318.
  * [8] T. Ha, T. K. Dang, T. T. Dang, T. A. Truong, and M. T. Nguyen, “Differential privacy in deep learning: An overview,” in _2019 International Conference on Advanced Computing and Applications (ACOMP)_ , 2019, pp. 97–102.
  * [9] Y. Zhu, X. Yu, M. Chandraker, and Y.-X. Wang, “Private-knn: Practical differential privacy for computer vision,” in _Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition_ , 2020, pp. 11 854–11 862.
  * [10] J. L. Horowitz, “Bootstrap methods for median regression models,” _Econometrica_ , vol. 66, no. 6, pp. 1327–1351, 1998.
  * [11] B. Jayaraman, L. Wang, D. Evans, and Q. Gu, “Distributed learning without distress: Privacy-preserving empirical risk minimization,” in _Advances in Neural Information Processing Systems_ , vol. 31, 2018.
  * [12] J. Zhang, K. Zheng, W. Mou, and L. Wang, “Efficient private erm for smooth objectives,” in _Proceedings of the 26th International Joint Conference on Artificial Intelligence_ , 2017, pp. 3922–3928.
  * [13] L. Wang, B. Jayaraman, D. Evans, and Q. Gu, “Efficient privacy-preserving stochastic nonconvex optimization,” in _Proceedings of the Thirty-Ninth Conference on Uncertainty in Artificial Intelligence_ , vol. 216, 2023, pp. 2203–2213.
  * [14] X. Zhang, M. Fang, J. Liu, and Z. Zhu, “Private and communication-efficient edge learning: a sparse differential gaussian-masking distributed sgd approach,” in _Proceedings of the Twenty-First International Symposium on Theory, Algorithmic Foundations, and Protocol Design for Mobile Networks and Mobile Computing_ , 2020, pp. 261–270.
  * [15] H. Wang and C. Li, “Distributed quantile regression over sensor networks,” _IEEE Transactions on Signal and Information Processing over Networks_ , vol. 4, no. 2, pp. 338–348, 2017.
  * [16] H. Wang, L. Xia, and C. Li, “Distributed online quantile regression over networks with quantized communication,” _Signal Processing_ , vol. 157, pp. 141–150, 2019.
  * [17] Y. Tian, M. Tian, and Q. Zhu, “Linear quantile regression based on em algorithm,” _Communications in Statistics-Theory and Methods_ , vol. 43, no. 16, pp. 3464–3484, 2014.
  * [18] Y. Tian, Q. Zhu, and M. Tian, “Estimation of linear composite quantile regression using em algorithm,” _Statistics & Probability Letters_, vol. 117, pp. 183–191, 2016.
  * [19] Y.-h. Zhou, Z.-x. Ni, and Y. Li, “Quantile regression via the em algorithm,” _Communications in Statistics-Simulation and Computation_ , vol. 43, no. 10, pp. 2162–2172, 2014.
  * [20] L. Yu and N. Lin, “Admm for penalized quantile regression in big data,” _International Statistical Review_ , vol. 85, no. 3, pp. 494–518, 2017.
  * [21] Y. Gu, J. Fan, L. Kong, S. Ma, and H. Zou, “Admm for high-dimensional sparse penalized quantile regression,” _Technometrics_ , vol. 60, no. 3, pp. 319–331, 2018.
  * [22] X. He, X. Pan, K. M. Tan, and W.-X. Zhou, “Smoothed quantile regression with large-scale inference,” _Journal of Econometrics_ , vol. 232, no. 2, pp. 367–388, 2023.
  * [23] P. J. Huber, “Robust statistics,” in _International Encyclopedia of Statistical Science_. Springer, 2011, pp. 1248–1251.
  * [24] C. Yu and W. Yao, “Robust linear regression: A review and comparison,” _Communications in Statistics-Simulation and Computation_ , vol. 46, no. 8, pp. 6261–6282, 2017.
  * [25] A. F. Siegel, “Robust regression using repeated medians,” _Biometrika_ , vol. 69, no. 1, pp. 242–244, 1982.
  * [26] P. J. Rousseeuw and C. Croux, “Alternatives to the median absolute deviation,” _Journal of the American Statistical Association_ , vol. 88, no. 424, pp. 1273–1283, 1993.
  * [27] P. Rousseeuw and V. Yohai, “Robust regression by means of s-estimators,” in _Robust and Nonlinear Time Series Analysis_. Springer, 1984, pp. 256–272.
  * [28] Y. Li and J. Zhu, “$\ell_{1}$-norm quantile regression,” _Journal of Computational and Graphical Statistics_ , vol. 17, no. 1, pp. 163–185, 2008.
  * [29] A. Belloni and V. Chernozhukov, “$\ell_{1}$-penalized quantile regression in high-dimensional sparse models,” _The Annals of Statistics_ , vol. 39, no. 1, pp. 82–130, 2011.
  * [30] L. Wang, Y. Wu, and R. Li, “Quantile regression for analyzing heterogeneity in ultra-high dimension,” _Journal of the American Statistical Association_ , vol. 107, no. 497, pp. 214–222, 2012.
  * [31] Q. Zheng, L. Peng, and X. He, “High dimensional censored quantile regression,” _The Annals of Statistics_ , vol. 46, no. 1, pp. 308–343, 2018\.
  * [32] L. Wang, “The $\ell_{1}$ penalized lad estimator for high dimensional linear regression,” _Journal of Multivariate Analysis_ , vol. 120, pp. 135–151, 2013.
  * [33] Y. Wu and Y. Liu, “Variable selection in quantile regression,” _Statistica Sinica_ , vol. 19, no. 2, pp. 801–817, 2009.
  * [34] J. Fan, Y. Fan, and E. Barut, “Adaptive robust variable selection,” _The Annals of Statistics_ , vol. 42, no. 1, pp. 324–351, 2014.
  * [35] J. Fan, L. Xue, and H. Zou, “Strong oracle optimality of folded concave penalized estimation,” _The Annals of Statistics_ , vol. 42, no. 3, pp. 819–849, 2014.
  * [36] R. Koenker and P. Ng, “A frisch-newton algorithm for sparse quantile regression,” _Acta Mathematicae Applicatae Sinica_ , vol. 21, no. 2, pp. 225–236, 2005.
  * [37] T. T. Wu and K. Lange, “Coordinate descent algorithms for lasso penalized regression,” _The Annals of Applied Statistics_ , vol. 2, no. 1, pp. 224–244, 2008.
  * [38] B. Peng and L. Wang, “An iterative coordinate descent algorithm for high-dimensional nonconvex penalized quantile regression,” _Journal of Computational and Graphical Statistics_ , vol. 24, no. 3, pp. 676–694, 2015.
  * [39] D. R. Hunter and K. Lange, “Quantile regression via an mm algorithm,” _Journal of Computational and Graphical Statistics_ , vol. 9, no. 1, pp. 60–77, 2000.
  * [40] J. C. Duchi, M. I. Jordan, and M. J. Wainwright, “Minimax optimal procedures for locally private estimation,” _Journal of the American Statistical Association_ , vol. 113, no. 521, pp. 182–201, 2018.
  * [41] Q. Zheng, S. Chen, Q. Long, and W. Su, “Federated f-differential privacy,” in _Proceedings of The 24th International Conference on Artificial Intelligence and Statistics_ , vol. 130. PMLR, 2021, pp. 2251–2259.
  * [42] K. Chaudhuri, C. Monteleoni, and A. D. Sarwate, “Differentially private empirical risk minimization.” _Journal of Machine Learning Research_ , vol. 12, no. 3, pp. 1069–1109, 2011.
  * [43] L. Wang and Q. Gu, “Differentially private iterative gradient hard thresholding for sparse learning,” in _Proceedings of the 28th International Joint Conference on Artificial Intelligence_ , 2019, pp. 3740–3747.
  * [44] D. Wang and J. Xu, “Differentially private $\ell_{1}$-norm linear regression with heavy-tailed data,” in _2022 IEEE International Symposium on Information Theory (ISIT)_ , 2022, pp. 1856–1861.
  * [45] L. Hu, S. Ni, H. Xiao, and D. Wang, “High dimensional differentially private stochastic optimization with heavy-tailed data,” in _Proceedings of the 41st ACM SIGMOD-SIGACT-SIGAI Symposium on Principles of Database Systems_ , 2022, pp. 227–236.
  * [46] B. Balle, G. Barthe, and M. Gaboardi, “Privacy amplification by subsampling: Tight analyses via couplings and divergences,” _Advances in Neural Information Processing Systems_ , vol. 31, 2018.
  * [47] C. Yi and J. Huang, “Semismooth newton coordinate descent algorithm for elastic-net penalized huber loss regression and quantile regression,” _Journal of Computational and Graphical Statistics_ , vol. 26, no. 3, pp. 547–557, 2017.
  * [48] M. Su and W. Wang, “Elastic net penalized quantile regression model,” _Journal of Computational and Applied Mathematics_ , vol. 392, p. 113462, 2021\.
  * [49] S. Shalev-Shwartz, O. Shamir, N. Srebro, and K. Sridharan, “Stochastic convex optimization,” in _22nd Annual Conference on Learning Theory (COLT)_ , vol. 2, no. 4, 2009, pp. 5–15.
  * [50] S. Neel, A. Roth, and S. Sharifi-Malvajerdi, “Descent-to-delete: Gradient-based methods for machine unlearning,” in _Proceedings of the 32nd International Conference on Algorithmic Learning Theory_ , vol. 132, 2021, pp. 931–962.
  * [51] M. Bun and T. Steinke, “Concentrated differential privacy: Simplifications, extensions, and lower bounds,” in _Theory of Cryptography Conference_ , vol. 9985. Springer, 2016, pp. 635–658.
  * [52] C. Dwork, W. Su, and L. Zhang, “Differentially private false discovery rate control,” _Journal of Privacy and Confidentiality_ , vol. 11, no. 2, 2021\.
  * [53] M. Redmond, “Communities and Crime Unnormalized,” UCI Machine Learning Repository, 2011, DOI: https://doi.org/10.24432/C5PC8X.