Differentially Private High-dimensional Variable Selection via Integer Programming
Authors: Petros Prastakos
Kayhan Behdin
Rahul Mazumder
Venue: NeurIPS 2025 poster ()
Keywords: Variable Selection
Best Subset Selection
Differential Privacy
Exponential Mechanism
Mixed Integer Programming
Sparse Linear Regression
Sparse Classification
Primary area: social_and_economic_aspects_of_machine_learning
TL;DR: We propose two scalable DP algorithms for high-dimensional sparse variable selection, leveraging modern mixed-integer programming techniques.

Abstract: Sparse variable selection improves interpretability and generalization in high-dimensional learning by selecting a small subset of informative features. Recent advances in Mixed Integer Programming (MIP) have enabled solving large-scale non-private sparse regression—known as Best Subset Selection (BSS)—with millions of variables in minutes. However, extending these algorithmic advances to the setting of Differential Privacy (DP) has remained largely unexplored. In this paper, we introduce two new differentially private estimators for sparse variable selection, levering modern MIP techniques. Our framework is general and applies broadly to problems like sparse regression or classification, and we provide theoretical support recovery guarantees in the case of BSS. Inspired by the exponential mechanism, we develop structured sampling procedures that efficiently explore the non-convex objective landscape, avoiding the exhaustive combinatorial search in the exponential mechanism. We complement our theoretical findings with extensive numerical experiments, using both least squares and hinge loss for our objective function, and demonstrate that our methods achieve state-of-the-art empirical support recovery, outperforming competing algorithms in settings with up to $p=10^4$.

Differentially Private High-dimensional Variable Selection via Integer Programming
Petros Prastakos Operations Research Center
Kayhan Behdin
Rahul Mazumder
Linkedln
Operations Research Center
MIT Cambridge, MA 02139, USA pprastakêmit. t. edu
Sunnyvale, CA 94085, USA
Sloan School of Management
k com
MIT Cambridge, MA 02139, USA rahulmaz@mit.c edu
Abstract
Sparse variable selection improves interpretability and generalization in high- dimensional learning by selecting a small subset of informative features. Recent ad- vances in Mixed Integer Programming (MIP) have enabled solving large-scale non- private sparse regression-Known as Best Subset Selection (BSS>wwith millions variables in minutes. However, extending these algorithmic advances to the setting Differential Privacy (DP) has remained largely unexplored. In this paper, we in- troduce two pure differentially private estimators for sparse variable selection, levering modern MIP techniques. Our framework is general and applies broadly to problems like sparse regression classification, and provide theoretical support recovery guarantees in the case of BSS. Inspired by the exponential mechanism, we develop structured sampling procedures that efficiently explore the non-convex objective landscape, avoiding the exhaustive combinatorial search in the exponen- tial mechanism. We complement experiments, using both least squares and hinge loss for our objective function, and demonstrate that our methods achieve state-of-the-art empirical support recovery, outperforming competing algorithms in settings with up to = 10 Code is avail- able
theoretical findings with extensive numerical
selection com cot./pperosspraatteee//--aaaaaaa
Introduction
High-dimensional datasets are increasingly common, but extracting meaningful models is challenging due overfitting and lack of interpretability. Statistical regularizations that encourage model simplicity from certain perspectives have been successful in addressing such challenges, becoming staple of high-dimensional statistics and machine learning. One such common regularization sparsity [I6] [5), where one seeks to choose a small subset of features in the data to form the statistical model.
In this paper, we focus on the problem of sparse variable selection. Given the data matrix X E Rnxp and the observations y R", seek to obtain estimator that describes the data well with only a few coordinates of / being nonzero. natural first formulation for this problem is
Σ (11 (3 s.t. || S s, || S2
min BER
where || llo counts the number of nonzero coordinates of a vector. In the case where E(y,, x () = (yi x 2 the objective becomes least squares, and the problem is referred to as Best Subset Selection (BSS, Miller 277 The constraint ||B]|o I enforces sparsity via the sparsity budget
39th Conference on Neural Information Processing Systems (NcurIPS 2025).

---

> 0, and the constraint || S r² for some > 0 serves as an additional (ridge) regularization. A sparse linear estimator can be more interpretable and have better statistical performance [ 16 34
Real-world datasets often contain confidential and personal information, that should be protected. Hence, recent years have seen surge in private learning algorithms, hoping to preserve sensitive in- formation while extracting useful statistical knowledge. In particular, Differential Privacy (DP, Dwork (2] has garnered significant interest in the machine learning and statisties literature. On a high level, DP aims to ensure one cannot obtain too much information from the private dataset, via querying the statistical model in an adversarial way. A significant body of work dedicated to designing DP algorithms for general machine learning tasks [25] 10] [4] [1 [1],, as well as specialized algorithms for specific statistical problems. Particularly, there is long line of work studying the sparse linear
regression problem [31 32, 23 24 29]
In this paper, we develop two scalable pure DP algorithms for variable selection under broad framework where one releases the optimal support in ①T (i.e., the location of nonzero coordinates in the optimal (3). To our knowledge, we are the first to incorporate MIP techniques for this task. While our support recovery results are derived specifically for the Best Subset Selection (BSS) setting, we provide pure DP guarantees for our methods that hold for general loss functions (not just least squares). Specifically, we make the following contributions:
1. Our first method, named top-R, satisfies pure DP under only standard boundedness assumption on the data (achievable via clipping). For BSS, it achieves support recovery with high probability whenever Sinn == minje #0} 183| Z & (log /n, matching the non-private minimax-optimal / /log p)/n threshold in the low-privacy regime.
2. Our second method, named mistakes, also satisfies pure DP, but requires an additional separation assumption on the objective gap for pure-DP guarantees. In BSS, this condition holds with high probability under Bein Z vs log p/n, with the milder condition for support recovery of Smin 2 /m (s log p)/n, which aligns with the condition of [29 in the high-privacy regime.
3. Empirically, our methods outperform the other DP variable selection methods in the literature, including the oftte-of-thh-att approximate DP MCMC approach of [29] for BSS, under wide range of parameter values and up to p = 10,000. We also show strong empirical performance in wider settings, including sparse classification with hinge loss. Our results demonstrate that DP variable selection with provable guarantees and practical scalability is possible by combining optimization and privacy.
1.1 Related Work
DP variable selection. Most of the existing DP literature focuses on non-sparse linear regression, or &g risk excess in sparse regression [31 33 21 7 22). For the specific problem of DP variable selection, previous works have focused on the sparse regression setting. As Lasso tends to promote sparsity, an interesting line of work is based on releasing the variables selected by Lasso in a private fashion 32 23 32] introduce two propose-test-eelasee algorithms for variable selection. However, the failure probability for support recovery for these methods does not approach 0 with growing sample size. 23 propose a computationally efficient
(28) algorithm. which underperforms compared to our methods in practice, and requires a stronger 3.ni condition than in our methods in the case of BSS. Lei et al. [24] propose an algorithm based on the exponential mechanism, requiring to enumerate all feasible supports in (I), limiting the scalability of their method. Recently, Roy and Tewari 29 have proposed a new method based on the notion of Markov chain mixing obtain approximate DP solutions for BSS, resulting in statistically strong estimator. While our 3min conditions are comparable with theirs in the low-privacy regime of top-ft or high-privacy regime of mistakes, we note that we achieve pure-DP guarantees for general loss functions, our algorithms have scope beyond BSS, and our empirical performance is stronger across broad range of parameters.
Modifications to the exponential mechanism. The methods we introduce in this paper involve modifications to the exponential mechanism, fundamental DP algorithm, in order to reduce the outcome set of our sampling distribution. Some other truncations of the exponential mechanism have existed in the DP literature. First, the Restricted Exponential Mechanism (REM) [6] for private

---

mean estimation samples from the exponential mechanism restricted to points of sufficiently large Tukey depth. together with private "safety" check that the restricted set well behaved. Second, the Truncated Exponential Mechanism (TEM) for metric-DP on text [30] restricts selection to a 7-ball around the input and collapses the remainder of the domain into single bucket- equivalenty assigning the outside set a shared score. While the spirit of truncation is analogous to the modifications proposed in this paper, our methods target combinatorial support selection under pure DP, rather than metric-DP over text data mean estimation.
Notation. We let [p] {1.... .P). Data points follow (æ,, 2 X x C RP R, with D = y) z" for a dataset containing observations.
2 Method
Background on Differential Privacy Before continuing with our selection procedure, let us formalize the notion of differential privacy.
Definition (1)). Given the privacy parameters (s, E IR R+ randomized algorithm A(-) is said to satisfy the (e, 8)-DP property
P(A(D) (KKEEPP(PCC K) +
for any measurable event K range(() and for any pair of neighboring datasets D and D'.
We note that in Definition] [ the probability is taken over the randomness of the algorithm A. When ò > the (z, 8)-DP property also commonly referred to as approximate differential privacy, while the special case where 8 = 0 is commonly referred to as pure differential privacy.
Next, let briefly review the exponential mechanism [25. general mechanism achieve pure DP. Consider a general task where the dataset D E 2" is given, and we seek to design a procedure such A: 2n O to choose the outcome of the task, where O is the set of possible outcomes. We also assume we are given an objective function such as R x n + R, where a smaller objective indicates desirable outcome. The global sensitivity of the objective is then defined as
R(o,D) ((o,
(2)
D,D'E D'
D,D' neighbors
Lemma (Exponential Mechanism, McSherry and Talwar (25)). that follows
The exponential mechanism AF(+)
=R(o, D)
IP(Ax(D) o) exp
Vo E O
(3)
2A
ensures OP..)-D
2.1
Selection Procedure
2.1.1
Rop-R Method
The main inspiration for our selection procedure is the exponential mechanism, defined in Lemma In particular, the BSS problem we seek to select subset features with size that are good linear predictor of n observations y. Therefore, natural choice for the outcome set in BSS is the set of subsets of [p) with size 8, = (S C [p) |S| s). Next, a natural choice for the objective in the BSS problem for each S is the least squares loss, when the regression coefficients only be nonzero for features in S. Formally,
R(S,D)
min BERISI
Σ Elv,, (2:) s.t.
|| < 2
(4)
where (xx)s the vector æ; with columns indexed by
Note that, if the elements of the dataset D are unbounded, the global sensitivity of R(S, D) may be unbounded. We thus make the following boundedness assumption on D.
Assumption 1. There cxist positive constants by. by such that supyey 1y| S by. sup ||2/|00 I br.

---

We note that, in practice. one might not know the exact values of bz. by. or such values might not exist. In such cases, one can clip the values of X, y to satisfy the boundedness requirements of Lemma A. In Section A.1 in the appendix, we show that, for the special case where our objective function is least squares, assumption yields < 2 2
Our Proposal With bounded global sensitivity, one can directly apply the exponential mechanism on (S, D) and achieve a (e. O)-DP procedure for problem [ The difficulty in variable selection under DP constraints arises from the need to enumerate all feasible solutions in O. However, one can argue that if a support S is far from the optimal one, the least-squares objective R(S, D) likely to be large, therefore, the probability mass of S in (3) should be small. Therefore, one might ask:
Is it necessary have access to R(S, for all S E O in the exponential mechanism?
Specifically, for the moment, suppose we have access an oracle that for fixed Rt > can return R feasible supports from O that have the smallest objectives. Formally, assume we can access (i ... (D()) where
(D) E min R(S,D) s.t. S S [p]. |S|= = S (,)) 1,... kk 1.
(5)
8
In particular, (D) is the optimal support for BSS in IT. Then, based on discussion above, if R is sufficiently large, the values (S,,, D) for 2 are expected be significantly larger than R(Sx D) for & < R. Therefore, most the probability mass the distribution (3) is concentrated around Sx (D) for k & R. Hence, we might not need to have access to the exact values of (R(S D) for k 2 R, long bound can be taken R(Sx (D), D). To this end, propose the sampling procedure M, shown as Algorithm below, where Po is the probability distribution following
can replace them with suitable lower bound. This lower
exp (E)(
if S R
(o(()) x
(6)
((9) R) exp (-x(
if = R+1.
Algorithm 1 Top- R. method
1: procedure MD, s, be. by R.
2: Clip X, y to b,, by. respectively, in Lemma A.I
Take A in Lemma A.I
Form IP
in (6)
3: Draw (DD) IPo
4: if (D) < R then
5:
return (D)
6: else
7:
return Mo(D, R) 8: procedure Mo(D,
9:
for t < T' do
Draw S O uniformly at random, independent of IPo. if (Sx l R} then
12:
Break
return
Intuitively speaking, M replaces R(Šx(D), D) for k > R with R(Šx(D), D) and then "approxi- mately" samples from the exponential mechanism. To this end, let
R(S,D) R(S,D) R((S{DD), D)
if E (D), ... otherwise
(S).)))... Sn(D)
where we substitute (D), D) for k 2 with R(S)(D), mechanism that uses the objective R. ( < R in Algorithm we return Sa(D) (D). Note that IP ( (D) Po(()(D) in this case, showing M mimics the exponential mechanism
D). Suppose Åg is the exponential

---

Ag. [f( R+ to mimic Ag, we have to sample uniformly from the set S = (Sx(D) R} as IP (AAg) is uniform on S, by the definition of R in 77. However, $ is exponentially large in general. Therefore, we invoke My that the limit of + oc, samples uniformly from S.
Observe that, in the case where R ("), we have that the distribution IP is the same as the exponential mechanism that uses objective R. as in (29), which is (e, O)-DP by Lemma Below, we show this procedure satisfies pure DP for any R E {2, (2) 1} as well. We defer all proofs to the appendix.
Theorem 1 (Privacy for top- R method). Suppose T > 1, < R and that assumption / holds. The procedure M Algorithm is (s',O)-DP where
q² log (1- 60 exp( neby (2
R
E' = log
eº
+
9
dg
In particular, if the procedure M is (e,
Theorem shows that, regardless of the choice of R. as T + oo, we have that + e. However, we note that è increases with R. there is more privacy loss with increasing R, as we have that
--) i (log e + 00 log (1 --))) -T
1
dq
50(e gT/66) 9"
The privacy loss is in contrast to the effect accuracy, we that larger R in Algorithm 1 should intuitively lead to better support recovery. We formalize this intuition Lemma A.5 in the appendix. For k > R, we underestimate R(S,, D) with R(SR, D), consequently increasing the probability mass given to supports S₁ in the procedure ÁF. This reduces the probability mass for the best support $1. Therefore, in practice, like choose larger R explore the objective landscape better, however, a very large can make the computation slower.
The sampling procedure in Algorithm only requires sampling from Ip (which is supported on R + 1 different values), and sampling sparse supports from uniform distribution (in procedure Mo). which can be done efficiently. Therefore, this procedure circumvents the need to sample from non-uniform distribution with exponentially large support. Importantly, Algorithm T satisfies pure (s', 0)-DP, with é' = as T - co. To our knowledge, no such algorithm exists for BSS that can scale to problems with tens of thousands of variables.
2.1.2 Mistakes Method
Our second proposed mechanism assigns probabilities based on the number of mistakes from the optimal solution. Namely, we define (D) S, arg ming R(S, D), and then proceed to partition the (?) 1 supports based the number of mistakes from (D). We denote the partition P2(),, (D)).. Let PD) = (D)). then have that for E [s]
3,(D) = arg min R(S, D). SEP:(D)
(8)
Our mistakes method, denoted M, assigns probabilities according the element of the partition that a support belongs to. Namely, if S E (D) for & {0, 1,.,, 8), we have that
(D)))))
R(L2(D))
PLM(D) S] Σ (P((D) exp( (D)) 2A
exp
2A
exp(
2A
Ei (P)9) exp(
R (.
24
For S E Pz(D), define R(S, D) = R((S(DD)
Below, show this method is (e, DDP under a lower bound assumption the gap in objective value between S, (D) and (D))
Theorem 2 (Privacy for mistakes method). Suppose assumption holds and that R(S²(D), D) R(Š D) > 2A. Then, the mistakes method is (e, -iffferrretlillllyy private.
5

---

We note that, unlike the privacy of our top-R method in Theorem 1 which requires no additional assumptions aside from assumption II Cheorem@rrquiree stronger conditions for the privacy of the mistakes method. In Lemma A.6 in the appendix, we show that, under the sufficient condition that 2 s, where 7 is defined in the following section, and the additional assumptions we have that the inequality R(Š= D) (D)) D) > 2 holds with high probability.
Remark While the privacy of our mistakes method relies on additional assumption that occurs
with high probability (as shown in Lemmal [A.6 in the appendix), providing privacy guarantee with
additional assumptions is not uncommon in the literature. For example, the privacy guarantees of [29]
which is the closest competitor to our method, depends on assumptions that hold with high probability. More specifically, the privacy proof of the Markov Chain Monte Carlo (MCMC) algorithm in [29
relies on the assumption that the mixing of the Markov Chain used for sampling with its stationary distribution has happened. However, this mixing can only be guaranteed with high probability, and
under additional assumptions on the underlying model-see see Theorem 4.3 of [29] for more details. In
contrast, our top-R method always private (assuming by, by are finite), and our mistakes method private under assumptions that are similar to the ones in Theorem 4.3 of [29].
As another example, [32 uses the stability of Lasso, present DP method for support recovery in sparse linear regression. However, the stability of Lasso only holds under certain assumptions on the data, such as the boundedness of the noise and restricted strong convexity. Such assumptions might only hold with high probability in practice, resulting in privacy guarantees that hold with high probability. For more details, we refer Theorem 8 of [32].
Statistical Theory
For the theoretical results in this section, focus the setting of BSS. Consider the model
y 3 +
where {ei)se[n) are i...d. zero-mean sub-Gaussian random variables with parameter and the feature vector 3 is unknown but is assumed to be s-sparse (i.e. its support size J,S* = lé 8 # 0}| =
p). In the remainder of this section, we provide sufficient conditions for our proposed methods to recover S* with high probability. We first state our additional assumptions.
Assumption 2.
There exists positive constant M such that 18 ll2 S M.
Assumption There exist positive constants K-, K4 such that, for all such that |S] = 8, we have
S Amin (X3 X Xs/n) < Amax (X X s/n) S k
Assumption 4. The sparsity level s follows the inequality s S n/ log p. and 2 3.
Assumption tells that the true parameter 3 lies inside an Eg ball. Similar boundedness assumptions are fairly standard in the DP literature [35] 24 7). Assumption the Sparse Riesz Condition (SRC), which is a well-known assumption in the high-dimensional statistics literature (36] 20 26). Finally, Assumption essentially assumes that o(n). i..., sparsity grows slowly relative to sample size.
Define the set of supports that make
[s] mistakes from the true support as
A {S
[p] |S| = s, |S S-1= t
Let = -n X be the sample covariance and Es, be the submatrix of Σ with row indices in Si and column indices in S2. Let Px5 = Xs X X denote the projection to the column space of Xs. Then we have that
y X s. Bg. + Pxs X. B (I, Pxs X Xx Ps- +
and (I, Px5) X Bz. describes the part of the signal that cannot be linearly explained by Xs.
Define also
D(S) Es-1s,s-s5 Es E s,
which the covariance of the residuals of Xs after being regressed on Xs. We have that
(In XX
D (( which we can intuitively consider as the

---

discrimination margin between S and the true support S The larger this quantity, the easier it is for BSS discriminate between S* and any other candidate model S.
We now introduce the central quantity of interest in analyzing support recovery, the identifiability margin, defined as
s (
7 min 1
S S1
We observe that, the more correlated the features, the closer to O, and harder is for BSS to distinguish between the true model and any other candidate support, so exact support recovery is harder. the features are more uncorrelated, 7 increases so exact support recovery is easier.
This intuition is made rigorous in Theorem 2.1 of 13 where it shown that > logg sufficient condition to have
"
{S = arg min Rots (S, where Rots (S) min SEO
lyy BER
X
with high probability, i.e. to have 5 be the unique minimizer for the BSS problem (with unconstrained l2 norm on (3) with high probability. Such a theorem offers us support recovery guarantees for the non-private fp-sparse ordinary least squares problem.
We now transition to the private setting of the fz-constrained version of BSS, and offer support guarantees for our proposed methods under this setting. In the following theorem, we provide sufficient conditions for our private top- R method to recover the true support with high probability.
Theorem 3 (Support recovery for top- R method). Suppose that assumptions
4 hold. Set 2
)M + 4 aab..
Set A = 20, + rb Then, there exists universal constant C > 0 such that,
whenever
2 max{Co? 8A logg
we have that
- 10sp
P(M(D)= =
1+p-
Comparison with previous work: Theorem B shows that, using the appropriate global sensitivity bound and lower bound on т. sufficient condition for recovering the true support with high probability is T Z max f logp s² "ogg], compared to r Z max log 8 #logp] } for the exponential mechanism applied to R) as Theorem 3.5 of [29]. Observe that, in a low privacy regime,
the logg term dominates, aligning with the [13] sufficient condition in the non-private setting. The extra factor of in the second term of our condition is expected, as we not making any additional assumptions on the choice It or on the number of mistakes of the enumerated supports.
In the following theorem, we provide weaker sufficient conditions for our private mistakes method to recover the true support with high probability.
Theorem 4 (Support recovery for mistakes method). Suppose that assumptions hold. Set
r 2 (*+ )M + 40bx, Set A = 2b,, + 2b3r2s. Then, there exists a universal constant C > 0 such
that, whenever
16A logp 6
T max
n
we have that
1 - 18sp
= P(M(D)=
1 2p
Comparison with previous work: Theorem 4 shows that, using the appropriate global sensitivity bound and lower bound on r, a sufficient condition for recovering the true support with high probability is f 2 max { slogg #log2] compared to 7 2 max logg slogg] for the exponential mechanism applied R(S,D D) in Theorem 3.5 of (29 Our condition
ne thus matches that in their paper for high privacy regimes. The strength in our result lies in noting that, unlike the
exponential mechanism in (29], which requires access R(S, D) for all (?) supports S C [p] such that |S] s, method only requires access to the supports that solve minse (P) R(S,D)
for all i E {0, s}.

---

4 Optimization Algorithms
In this section, we discuss how the top R supports, (D), SR(D), which are solutions to the problems ( can be obtained by solving a series of MIPs.
For clarity, we present the case of least squares objective. However, it is worth noting that a key benefit of our MIP approach is its generalizability across different loss functions. Specifically, the only property that we need to have for the loss function f' is that it is convex function of B. We discuss any needed modifications for the case of hinge loss Appendix[C and our method works with Huber or quantile loss as well, and the MIP algorithm would remain unchanged.
To obtain Sx for k E [R] consider:
(z (x) ((k), (6)) E arg min 2,8,0
y - X
s.t. 3,0 6. R2, E 10.1)" 020 0,
Σ
0g
S r?
i=1
i=1
1
s Eie [p].
Σ
24 S 8
j
j [ 1|.
2
1E9,
(D)
In the following proposition, we show that Problems 55) and 9
are equivalent.
Proposition 1. For k 2 1, f 0 = (D).
Problem (9) can be solved to optimality using off-the-shelf solvers like Gurobi for moderately-sized datasets. In order to run our DP methods in even higher dimensions, where p 10,000, we present more tailored algorithm for solving the MIPs in the remainder of this section. This adds the line of work on developing specialized discrete optimization algorithms for solving sparse regression problems and relatives the non-private setting setting-see see for eg [8
We first add a ridge penalty term to the objective, which makes strongly convex and a function of z.
To obtain Sx (D) for k E [R] we define
c(z)
min
11y
x31
>
Σ
2
37
perier 2n
2m
24
and we seek to solve
min c(z)
(10)
subject z E 13',
2, = 5, Σ (ES;(D)
2 2
Vj jk 1].
2
For any 2 E {0, 1),, let & Let
2; if 2,
and 2 2; Unif(a, bj if where 0 and b 1.
A
p
3
B E arg
2n lly - x3112++ ||
min
Σ
2n 1=1
2
By Danskin's 's theorem
3 then have
(Ve(2)), A (8.)
2n 2
Taking 20, 2,, 2n E (0, 1JP, we have by convexity of c that for all (0, 1]2 and for all E (O,.., t).
c(x) 2 c((x)= + Ve(êx)) 2x).
So then the map
(xx) max(c(26) + Vc(2o)"((e (x 20), ...., d(x)) + Vc(&t) -20}
8

---

is a lower bound on the map c. We can now present our outer approximation algorithm for solving
Problem 10 based on g
Algorithm 2 Outer approximation for Sx(DD
1: procedure A(D, &,, tol)
2: Initialize #o E {0,
1)"h s.t. Ei 7 0,t
3: 20 f add add_noise(zo)
4: while ( tol do o(xe)
P
5: 2t41, M1+1 f- arg (0.1)
8.1.
Σ
2 S
Σ iE8y(D)
2, 2 8
fj E [h: 1),
i=1
22 ((2) + Vc(2x) 24))
wk E (4]
6:
2++1 t add noise((2)) ++
7:
return 21
Intuitively, this approach seeks to solve Problem by constructing a sequence of MIP approximations based on cutting planes. At each iteration, the cutting plane 1) 2 ol & + (o(2,)²(( 2ạ) is added, cutting off 2,, the current noisy version of the binary solution 2t, unless 2; happened be optimal as defined by our stopping criterion, which is [--c( < tol. As the algorithm progresses, the c
outer approximation function (² = max
2;) becomes an increasibly better approximation to the loss function making our lower bound ?i: converge the upper bound obtained by evaluating c( Please refer to Appendix BB for additional algorithmic details for the setting of BSS, and refer to Appendix for the necessary modifications in the case of hinge loss.
c(2,) + Vc(2,)
Remark 2.
The approach for obtaining (D), where k E [s], is very analogous to Algorithm 2
except instead of having the constraints Σieš,(D) 2: S 2 Vj - 1]. where E [R], we have (l: 1)
the constraint EI)
5 Numerical Experiments
In our experiments, we draw the data points as y; = =x + for i E [n]. where æ1, 2n iid (0 Σ) E R and the independent noise follows NN x² I,)) where I, is the identity matrix of size Moreover, for i,j [p]. we set Σ, أاج and set nonzero coordinates of 3 take value 1/ /3 at indices (1,3,
2s 1}. We define the Signal to Noise Ratio as SNR || / || In Algorithm[]] we set R = 2 (p 1ip-aa s)s, by = by = 0.5, r = 1.1 and T = co for all our experiments in this paper. In Algorithm 2 we set 0.001, b 0.005, various values of the other parameters.
1.1 and tol 0.005, and consider
In Figure Ta we plot the average proportion of draws from 10 independent trials that recovered the right support for our top-Ri R and mistakes methods using least squares as our objective for p 10, 000. We compare with the MCMC algorithm from [29 and the Samp-Agg algorithm in 23 wherein we use Lasso for the Asupp subroutine. In Figure Tb we show the analogous results for hinge loss, comparing with Lasso Samp-Agg. More experimental results with varying values of SNR, p, 8, P, and e for least squares and hinge loss, as well as results on prediction accuracy. utility loss, and ablation studies corresponding to each algorithm and gathered the proportion of correct supports. For MCMC, we similarly used 50 independent Markov chains from random initialization and gathered the proportion of correct supports after a number of iterations that was chosen to make the runtime comparable to methods.
provided in Appendix D For each trial, we drew 50 times from the distribution
We have that all settings, both of our methods outperform other algorithms for large enough n. The proportion of draws that recover the right support increases with since a larger sample size reduces the threshold required for the identifiability margin 7 (discussed Section 3 have enough separation between the true support and other supports. Furthermore, performance improves much

---

larger values of n when or s is greater or is lower, since the lower bound on T is harder to satisfy those settings. Moreover, we observe that, keeping s, n. p. and fixed, smaller values of SNR and larger values of make recovery harder, as 7 decreases with lower signal strength more correlation between features. Furthermore, the mistakes method numerically outperforms top-R, aligning with which shows it succeeds under milder identifiability condition.
Section
1.0
1.0
correct supports trials
# correct supports total trials
C.8
0.8
C.6
Top
Top R Mistakes
Mistakes
Lasso Subsample
C.4
MCMC
0.4
Lasso Subsample
C.2
0.2
7000
8000
9000
10000
11000
6000
7000
8000
9000
(a) Least squares
(b) Hinge loss
Figure 1: Simulations for p = 10,000, s 5, SNR=5, p 0.1, and E = for least squares and hinge loss. The penalty parameter À in Algorithm 2 was set to 600 and 170 for figures La and lb respectively, and the number MCMC iterations was set to 100,000 for Tal On the æ-axis, we vary the value of and plot the average proportion of draws across 10 independent trials that recovered the right support for each corresponding algorithm. Error bars denote the mean standard error.
Computational resources and license information: All experiments were conducted on computing cluster using 20 cores and 64 GB RAM. The Gurobi Optimizer used under the Gurobi End User License Agreement. CVXPY is distributed under the Apache License, Version 2.0. ABESS package distributed under GNU General Public License, Version 3.
Conclusion
In this paper, introduced two scalable pure DP estimators for variable selection in sparse high- dimensional settings. While we provide utility guarantees specific to the BSS setting, we demonstrate how our methods can be applied more broadly, yielding favorable support recovery in the additional setting of sparse classification with hinge loss. Our contributions enhance privac--rreeeving practices, enabling safer of sensitive datasets in critical areas such as medicine, public health, finance, and personalized recommendation systems.
One limitation of our work is that Theorem requires an additional assumption that holds with high probability for r large enough, and remains an open question whether lighter assumption can be made to yield privacy guarantees for the mistakes method. Furthermore, our theoretical support recovery results yield sufficient conditions for support recovery, but an interesting direction of research may be to find necessary conditions well, to see if our bounds 7 are tight.
7 Acknowledgements
We thank the authors of (29 for sharing their code with us. This research is supported in part by grants from the Office of Naval Research AN00014251200 N0001422126655 A shorter workshop version of the paper appeared in [2). The research started when Kayhan Behdin was a PhD student at MIT.
10

---

References
[1] Kayhan Behdin and Rahul Mazumder. Sparse pea: A new scalable estimator based on integer programming. arXiv preprint arXiv:2109. 11142, 2021.
[2] Kayhan Behdin, Petros Prastakos, and Rahul Mazumder. Differentially private best subset selection via integer programming. International Conference on Learning Representations Privacy in Machine Learning Workshop, 2024.
[3] Dimitri Bertsekas. Nonlinear Programming, volume 4. Athena Scientific, 2016.
[4] Dimitris Bertsimas and Bart Van Parys. Sparse high-dimensional regression: Exact scalable algorithms and phase transitions. 2020.
[5] Dimitris Bertsimas, Angela King, and Rahul Mazumder. Best subset selection via modern optimization lens. 2016.
[6] Gavin Brown, Marco Gaboardi, Adam Smith, Jonathan Ullman, and Lydia Zakynthinou. Covariance-aware private mean estimation without private covariance estimation. Advances in neural information processing systems, 34:7950-7964, 2021.
[7] T T'ny Cai, Yichen Wang, and Linjun Zhang. The cost privacy: Optimal rates of convergence for parameter estimation with differential privacy. The Annals of Statistics, 49(5):2825-288) 2850, 2021.
[8] Marco A Duran and Ignacio Grossmann. An outer-aaproximatiion algorithm for a class of mixed-integer nonlinear programs. Mathematical programming, 36:307-339, 1986.
[9] Cynthia Dwork. Differential privacy. In International colloquium automata, languages, and programming, pages 1-12. Springer, 2006.
[10] Cynthia Dwork, Frank McSherry, Kobbi Nissim, and Adam Smith. Calibrating noise to sensitivity in private data analysis. In Theory of Cryptography: Third Theory of Cryptography Conference, TCC 2006, New York, NY, USA, March 4-7, 2006. Proceedings 3, pages 265-284. Springer, 2006.
[11] Cynthia Dwork, Guy N Rothblum, and Salil Vadhan. Boosting and differential privacy. In 2010 IEEE 51st Annual Symposium on Foundations of Computer Science, pages 51-60. IEEE, 2010.
[12] Cynthia Dwork, Aaron Roth, et The algorithmic foundations differential privacy. Founda- lions and Trends@ in Theoretical Computer Science, 9(3-4)211-4077 2014.
[13] Yongyi Guo, Ziwei Zhu, and Jianqing Fan. Best subset selection is robust against design dependence. arXiv preprint S1V:220070014760 2020.
[14] Moritz Hardt, Katrina Ligett, and Frank McSherry. A simple and practical algorithm for differentially private data release. Advances in neural information processing systems, 25, 2012.
[15] Trevor Hastie, Robert Tibshirani, Jerome H Friedman, and Jerome H Friedman. The elements of statistical learning: data mining, inference, and prediction, volume Springer, 2009.
[16] Trevor Hastie, Robert Tibshirani, and Martin Wainwright. Statistical learning with sparsity: the lasso and generalizations. CRC press, 2015.
[17] Hussein Hazimeh and Rahul Mazumder. Fast best subset selection: Coordinate descent and local combinatorial optimization algorithms. Operations Research, 68(5):1517-1537, 2020.
[18] Hussein Hazimeh, Rahul Mazumder, and Ali Saab. Sparse regression at scale: Branch-and- bound rooted first-order optimization. Mathematical Programming. 196(1):347-38,, 2022.
[19] Hussein Hazimeh, Rahul Mazumder, and Peter Radchenko. Grouped variable selection with discrete optimization: Computational and statistical perspectives. The Annals of Statistics, 51 (1):1-32, 2023.
[20] Jian Huang, Yuling Jiao, Yanyan Liu, and Xiliang Lu. constructive approach to 1_0 penalized regression. Journal of Machine Learning Research, 19(10)1137,, 2018.
11

---

[21] Prateck Jain and Abhradeep Guha Thakurta. (near) dimension independent risk bounds for differentially private learning. In Eric P. Xing and Tony Jebara, editors, Proceedings of the 3 st International Conference on Machine Learning, volume 32 of Proceedings of Machine Learning Research, pages 476-484, Bejing, China, 22-24 Jun 2014. PMLR.
[22] Shiva Prasad Kasiviswanathan and Hongxia Jin. Efficient private empirical risk minimization for high-d -igh-dimensional learning. In Maria Florina Balcan and Kilian Q. Weinberger, editors, Proceedings of The 33rd International Conference on Machine Learning, volume 48 of Pro- ceedings of Machine Learning Research, pages 488-497, New York, New York, USA, 20-22 Jun 2016. PMLR.
[23] Daniel Kifer, Adam Smith, and Abhradeep Thakurta. Private convex empirical risk minimization and high-dimensional regression. In Conference Learning Theory, pages 25-1. JMLR Workshop and Conference Proceedings, 2012.
[24] Jing Lei, Anne-Sophie Charest, Aleksandra Slavkovic, Adam Smith, and Stephen Fienberg. Differentially private model selection with penalized and constrained likelihood. Journal of the Royal Statistical Society Series A: Statistics in Society, 181(3)600 6333, 2018.
[25] Frank McSherry and Kunal Talwar. Mechanism design via differential privacy. In 48th Annual IEEE Symposium on Foundations of Computer Science (FOCS pages 94-103. IEEE, 2007.
[26] Nicolai Meinshausen and Bin Yu. Lasso-type recovery of sparse representations for high- dimensional data. The Annals of Statistics, (1(1:244-2200 2009.
[27] Alan Miller. Subset selection in regression. CRC Press, 2002.
[28] Kobbi Nissim, Sofya Raskhodnikova, and Adam Smith. Smooth sensitivity and sampling in private data analysis. In Proceedings of the thirty-ninth annual ACM symposium Theory of computing, pages 75-84, 2007.
[29] Saptarshi Roy and Ambuj Tewari. On the computational complexity of private high-dimensional model selection via the exponential mechanism. arXiv preprint arXiv:2310.07852, 2023.
[30] Ricardo Silva Carvalho, Theodore Vasiloudis, and Oluwaseyi Feyisetan. Tem: High utility
metric differential privacy on text. arXiv preprint arXiv:2107.07928, 2021.
[31] Kunal Talwar, Abhradeep Guha Thakurta, and Li Zhang. Nearly optimal private lasso. Advances in Neural Information Processing Systems, 28, 2015.
[32] Abhradeep Guha Thakurta and Adam Smith. Differentially private feature selection via stability arguments, and the robustness of the lasso. In Conference on Learning Theory, pages 819-850. PMIR, 2013.
[33] Prateek Varshney, Abhradeep Thakurta, and Prateck Jain. (nearly) optimal private linear regression via adaptive clipping. arXiv preprint arX1v:2207.04686. 2022.
[34] Martin J Wainwright. High-dimensional statistics: non-asymptotic viewpoint, volume 48. Cambridge university press, 2019.
[35] Yu-Xiang Wang. Revisiting differentially private linear regression: optimal and adaptive prediction & estimation in unbounded domain. arXiv preprint orXiv:1803.02554 2018.
[36] Cun-Hui Zhang and Jian Huang. The sparsity and bias of the lasso selection in high-dimensional linear regression. 2008
12

---

A Proofs of Main Results
A.1 Lemma A.1
Lemma A.I. Suppose Assumption 2 holds, with
R(
min BERIS
y X3 s.t. |, r².
(A.1)
Then,
A 26 + 2b2
Proof. Suppose D, D' are two neighboring datasets. Fix support S E O and suppose
E arg min 11y/
1ا s.t. | v r
BER
Then,
R(S, D) - R(S,D')
x - 3
Let us assume without loss of generality that D, D' differ in the n-th observations. Hence, have that
1y - x -
l1 - x'33133 =
Σ [y: (x:),,) - (yi - (2;)
==1
(1/2 (a Y g) < (ya (x (33)
+ (yn- (x 1/,)
< 2y + 33)
$ 26y 262,2²
where the last step uses the Sauchy-Schwartz Schwartz inequality and the fact that |S| s.
A.2 Proof of Theorem 1
First, we prove some technical results that will be used in the proof of Theorem
Define
R(S,D) if (D), (R))))),,,,(()))) otherwise.
R(S
(A.2)
(R(SR(D),))
Lemma A.2. Let A to be taken in 2 Then,
R(Sx(D), -)) s A.
max D.D' D,D' neighbors
Proof. Fix k 2 1 and let consider the following cases:
Case 1: (D)) D') s R(S, Then, by the definition of A,
R()) D) - (((), D') < (D) ( D') +
A.
(A.3)
Case 2: R(Sx D) < R(S) D). Then,
R(S{(D), D) - R($) D') <
+ A
A.
(A.4)
(),,))))) D') > R(S) D'). Trivially, in this case we must have (") 2 2. Then, there must exist So C [p), |aal = such that
Case 3: R(Sx(D),D) > RD')(DD). D) and R (D),
R(S., D) > (D),
D), R(S., D') (D'),
13

---

To this end, define
S (D) ... (D)) 82 (D), ... S = (D', ... (D)))
As (D), > RDS) we have that Sx (D) & so S' U On the other hand, as (D) > R(S((D)) D), (D') E S1 and as SD(D') $ S', we have IS' n1 I l 2. As | 1, we must have 8 n 2 which proves the existence of Sg. Next, note that
R(Sx(D)) D) (D'), S RSo. D) R(So, D)) S..
(A.5)
Lemma A.3. Let A to be taken in 2 Then,
R(S) SA.
max SCIP
mAx D.D'E2
S D,D' are neighbors
Proof. Suppose S = Sk (D) = (D'). Let us consider the following cases:
Case k1 12 2 Then, we have R(S D) RR(Sx(D), D) and RiS, D*) R(ŠS((D)) Therefore,
R(S,D) (OS,D)-R((52)))) - -((Sp()) - (()) SA
(A.6)
by Lemma A.2
Case 2: ln R2 2 R: Then, have R(S,D) R)S,D) R(Šn((D) D) and R(S
RRSS(DD)) D)) Then,
R(S,D) - D)( - R(SR D) R(Sx(D'), D)) A..
(A.7)
= Case 3: kj 2 R2 < R: Then, we have R( = R(S=(D), D) < R(S and R(S, (')
RRS,DD') Then.
( < R(S,D) R(S, D' S..
(A.8)
Case 4: k1, 12 < R: Then, we have R(S, D) RRS,D) and R(S, D') = D'), The result follows.
Lemma A.4. Suppose M is as defined in Algorithm 1 and Ag is an exponential mechanism with the objective R,
&R(S,
P(AA(D)= = 5) x exp
AS
(A.9)
2A
Then, for S EOO
(1- PPM)) = S) < ((A=(D) = +)+
(A.10)
where
R
=
Proof. Fix O and suppose 5 (DD) Moreover, Sn = (Sx(D) k < RJ. Consider the
following cases:
Case 1: Then, based on Algorithm 1
PPM(DD) S) = IP
= = S]).
Therefore,
IP(a(D) -k) PPM() S) = = =) R+,, = S})
(a) 22 P(a(D) = k) + PPM Sn)
N PP(a(D) = k) +
(A.11)
14

---

where (a) is true as S E SR. and (b) is true as Mo return support in SR if it selects some support from nn for all T iterations, showing (P E Sn) q". Note that for k I R. P(((DD) = k) = (A (D) = 5), therefore,
(A (D) = S) = P(M(D) - S) = P(A)((D) - S) +
(A.12)
Case 2: k R: Then, from 66.
(C) R) exp (- R(Sn(D) D))
P(a(DD) R 1)
E exp (- 1-Z3 R(Sk (D): D) + (()) - R) exp
R(Sn D)
(()) - R) exp ( 23 R(
Ek= exp
E (()) ΣxZR-1 exp
(-x))( D))
23
(() (8)
(P)
R
(D)A 5).
(A.13)
Hence, can write
PPM = 5) - (fa(D) = R ()) S/)
g
@ P(a(D) = R+1)
Σ
(P)
1=1
(0) - (3)
() (()
R P ((Ag S) 1- qT
1-9
(1- q²)) )P S)
(A.14)
where (a) is true as
P(Mo(D) = S) = Σ PP (P) = Mo stops after iterations ) =
Σ
1=1
1=1
The proof is complete by (A.12) and A.14
Next, let prove important intermediate result Algorithm [
Theorem A.1. Suppose T > 1, R (?). The procedure M in Algorithm 1 is ()-,8)-DP where
R (e))))) 8 s R// 6=(r/(2 (3)
e' log 1
Proof. From Lemma A.3 and Lemma neighboring datasets. Then, from Lemma
1 we know that Ag is an (e. DD-D procedure. Suppose D, D A.4 A.4
PPM(D 5) s (P =) + 9"
I E (D') = S) +
1
P P(M 5) + -T
(A.15)
where the first and last inequality use Lemma A.4
Proof of Theorem 7 Note that by definition, for S O, we have 0 S R(S,D) < 1 < nt? Then,
exp e ()(S.))
(A(A = =S)
-== exp R (D), D)) + Σ Ex2R+1 exp (- c (() D))
exp( mepp-nebi// (eb))/(A)) (2A))
: 60-
(A.
15

---

Then, from A.15).
PMM(D) - S) ((AA(D) = = +
((An(DD) = (P(A((D)) 5) + -
60
60
E +
q
(((A(DD) SS
on
e
P(M(D') S)
(A.17)
q²
do
where (a) is by (A.16) and (b) by Lemma A.4
A.3 Effect of choice of R on support recovery for top-R method
In this section, we formalize the intuition discussed in Section 2.1.1 regarding the impact that the choice of R has the top- R recovering the optimal BSS support (1 We show that, R increases, the probability of top-R outputting (D) can only improve.
Lemma A.5. Take T do in Algorithm 1 Denote Mi and Ma as two instances of the top- R method, using Rí and R2 enumerated supports, respectively, where Rı R2. Then, we have that
P(M,(D) = Si PPM2)() -
Proof. By Lemma A.4 we have that, when T = 0c.
(PAA $1 PPMM(() SD)
Let AE, and A denote the exponential mechanisms with the objective R. using Ri and Rs enumerated supports, respectively. It then suffices show that
PL (D) (D)) < PI (D) S, (D)).
Define G, = (((((())) - (D)) ) Then, note that
1
(D) S
1 exp (()) R (x))exp((000 R.)esp(-eGa/((2A)
and
(A (D) = (iD)
1 + 2 exp( exp(-eGi/(2A)) + (())
We then have that
Ra
cG₁
cGnn
exp
I RR2
R1) exp
2A
2A
(=n R
and
(?)
cGR,
(?)
eGR,
Rs)
exp
R2) exp
2A
2A
since the supports are sorted in increasing objective value. The result follows.
A.4 Proof of Theorem 2
Proof. First, we claim that if R(S2(D)) D) RR)) 2A, then Si (D') for neighboring datasets D. DD'
Suppose, by way of contradiction, that (D (D'). Then we have that R(Si < R(), D') and R(S1(D'), 2 ((), D). By definition of A. have that
16

---

R(Si(D), D') & R(SI(D), + A and R(S)(D)), D) A - R(Si(D'), D'), Combining we have
(D)) D) + - R(S) (D') D') 2 (D), > (D),
- Thus 2A 2 -(((()))))_R((()))))) This is a contradiction with the assumption that ((S2((): D) - D)
Since (D) = (1(D), i.. Su(D) = (D'),, we have that P.(D) - Pi(D') for all i E {0, 1,..,..... s). Suppose 5 E Ps (D), Pa(D"). Then, by definition of RR we have R(S,D) = (D) DD and R(S, D') - R(S (D'), D"). We have that (D),, - minser. R(S, D). Since S (D') E Pa (D), we have that RD), D) I (D') D). Then, by definition of A, we have that
D) R(S,D') - 7.0.://..................000-0-0
Thus, we have that
R) R(,, SA.
max [s D,D' neighhors
max 2
Since R has bounded global sensitivity A, we have that, by Lemma 1 the mistakes method, which is the exponential mechanism with scoring function R. is (c. -)ddiffrrrettilllyy private.
A.5 Sufficient conditions for privacy of mistakes method
Below, we present sufficient conditions under which R(S2(D). (D). > 24 with high probability, which, by Theorem 2 implies that the mistakes method is (c, DD-DP with high probability.
Lemma A.6. Suppose that assumptions 4 hold. Set 2
gbx, M K_
Then, there exists
universal constant C such that, whenever
5 logp
2 Co
T
n
have that
(P(((S(D))D) R(Š) > 24) 2
1 10s0pp -
Proof. Following the proof of Theorem 2.1 from 13 we have that
n (Roto(())) RS D)) = n 1 T =n +E)
- Pxs
)u T Px..)9)
T( (I Px
+ T Px )E}
( + 2n TT(( PI X sBB.
- -ntle (Pxs -Xx- )E
2 1 (S) + 7 ( + 2n -le T -I- Px X X
+ 7 s DSS) n T - Px..Je.
We then argue that the following two inequalities are true with high probability
1 0 D(S)
2n 1 Pxs) S s.
(A.18)
1 D(S) 35-ss
n -tee T(Pxs - -Pxg- )c < B
(A.19)
4
so that -nnn-'(Raa (Rot D) RS D)) 3 (SS)
Defining n 1/2 - Pxs ) X ss we have that (A.18) is equivalent to
1/2
<
ug 2
|us||2
17

---

Note that since each of the entries of e are ii.d. zero-mean sub- Gaussaan random variables with is sub-Gaussian with parameter so we can apply the Hoeffding to get that, for any 0,
parameter we have that ugel
luslla
bound (Proposition 2.5 in [4] with
just
El
P
ox S 2e
uus
Now, applying union bound over all S Ar, we have that for any 0,
P
3S
Ar,
luge Lus 2
>
o Ets 1 (P=0)
8
2e sta/2 - e 8/8/2 -
Then we have that, whenever
inf sE, ||us|12
{s) 1/2
80
t1/2
n
we have that
jugel |us||2
1/2
P 35 E Ar:
l
S 2pp" -66//2
Regarding AA.19. observe that, as shown in the proof of Theorem 2.1 of 13 exists universal constant e₁ 0 such that, for any O,
we have that there
P le I (Pxs - Px..)el )<g
202x 4e c min(<=/ta/
n
n
Then, we have that for e minfc1, 111
E TPxg (Pxg- Px.. )..cel > 202x
S -a min(a"/(ex)) 4e cmmnmxx ft,z)
n
Noting that 1, we have that for any & 2 1,
20
P T - Px..) )
< -gts
and applying union bound over E Ac, we have
(as E Ar, Le TPxs - PX..)E> 20t
(P= p-8 8) (1)
4e -5ts
< <4p" e
Hence, whenever,
inl SEA,
| 2
8502g
1/2
/1/2
n
we have that
P
SS E Ar, c (Pxx Px.- )>> 4 |lus11
I 4p cotto
n
Combining and taking union bound over all t E [s]. we have that, for any & > 1, whenever
T 2 2
Es
we have that,
1
P V VSSEULIAL (S, D) (St D) 5 21 4e
-6//2 t ccfo).
cEa Now note that, since < 2. we have that 1 2sp 2 - -dt) > 1 8sp -ot= Furthermore, choosing log p. we have that p² 20 log + 0 as p
18

---

Combining the above, define = =/4 b, + ber², where r M1 MM + 4 as assumed, and let
C' max{82co $ Set { co log p. Then, whenever
2 Co logp
n
we have that,
1
P VS E Ui (S, Rolo (5", (D) 2 nt
8sp As
where x = cco
2>> > 2. Note that,
Co 8 logp n
- n (8b) + 86
4
where the first inequality uses that max(1, d2 Jo Σ second inequality uses that logp > by [ssumption and the last inequality follows from Lemma A.1 Hence, we have that
F PVS E Ar, Rnt D) Rots(S" D) > 2A)= 2 8sp 2
Taking union bound using Lemma A.7 we have that
P PVS E 11,,, R(S" D) > 24) 2 1 2p- 8sp 2 22 10sp
and noting
P(R(S2(D), D) R($;(D))) 24) 2 PYS E
- RRSS" DD > 24)
concludes the proof.
A.6 Proof of Theorem a
Before proving the theorems, we will first setup some preliminaries.
Observe that the solution to the unconstrained least squares problem with support restricted to S is given by
X Xs ) Bs = (X X5y
X x-
X
X X
n
n
1
n
442
and the constrained estimator on the same support is given by
Bs.r = arg min l1 x33322 B:11Bllssr
(A.20)
For each support S, define the event Es,r (Bs, Bs,ots} and the intersection of events across all supports as E,- = &
The lemma below shows that if a sufficiently high bound on the l2 norm of / in the constrained optimization in A.20 chosen, the solution to the unconstrained OLS problem is the same as the solution to A.20 for all supports S with high probability.
Lemma A.7. Suppose assumptions hold. Ifr 1- 2p-7
(t M 44 ab. then P Pinsis
Proof. At a high level, we are seeking to bound ||, 2 with high probability. First off, by
assumption 2 and we have that ||u,||s 2 I [X] xs/n)
/ '112 | Xs- /nl//|l || ll s M
Next, note that
X $
Xs
X
X X
X
I Vs
n
n
2
72
2
7
X X
X E
< V8
19

---

Hence, we have |ug l2
-
v
Now, define D.j Xx for all (i.j) [p] [n]. Since
K n
ej is sub- Gabssinn with parameter o, using assumption 1 we have that D.j is uub-Gaussian with parameter σbg. Applying the Hoeffding bound (Proposition 2.5 in (4 with t = 4ob, /n log p,
have that, for all [p].
log n
Da.sl 2 40bz n
V
2p
Observe that
X
max4elel # Digl.
Hence, by union bound, we have that
x
4abs
log
S 2p-7.
By assumption have that
X c
ob2 V
P
2
$ 2p-
n
This yields that [|8s 1|2 S
|uill2 2 + ||22l12
I -+ M + 4 obe with probability at least 2pp-7
Hence, we have that if r 2 ( +++ )
+ 4
then P E, > 2p- as desired.
For the remainder of the paper, define Rots (S, min 1ly x3313 for all S C [p] such that |S] = s. Define the event Eouu == Ni=i (VS Ar. (Rols (S, D) Rots (S", D)) 2 1t7
Lemma A.8. Suppose 2 3. There exists a universal constant C > such that, whenever
Co logp
have that
N sAS E A,, 1 n
(Rnis RSt D)) 2 1 [7)
P[Egap P
2 - 8sp
Proof. By Theorem 2.1 of 133, we have that there exists constant c > such that for any & > 1,
whenever
minsEA, ATT. D(S) S-1S
4
0² logg
t
1 -
we have
P VS E Av, DD) RRos ( > 11 (())) (118)
2 1- (c5--))_ (052
2
Let = 5 = max =), and C' =
4
Then, we have that -(ct - 2)t < 22 and
-(E2 2) 2t, so 4p (re 2)( 2p
(E 2)x
4p
2p > 1 8p 2t. Observe
also that
minSEA, 33 3(S)
minseA. B D(S)
min [E]e
t
Hence, we have that whenever
7 2 Co log
n
we have
P VS VS e A., 1 (R.ss - (S*, D)) > (
2 1 18p 8p-2t
Observe now that B D(S), 2 tr hence we have
P VS E A., A. (Rot D) (S D))
1
tr > 8p-21 21
2
Applying the union bound and the fact that 2 get (P(Epp) proof.
2 1 -8sp- 2, which concludes the
20

---

We now proceed with the proof of the theorems.
Proof of Theorem We will now show that the exponential mechanism with scoring function R. and R 2, denoted as A recovers the true support with high probability.
We first define the event
- = (AS E - D) nlls..
Observe that &, n & C E. By Lemmas A.7 and A.8
we apply union bound, we have that
P[E] 2 PIEE n & 2 1 8sp 2 2p -7
10sp
Furthermore, if we condition on E, we
have that (D) D) = R(S", D).
Then, note that
=
1) exp(- 24 (R(S2(D), D) - R(S)(D).
Then, have that, if we assume
8A s)
logp
> max
we have the following:
(?)
P
1) exp(
(R(S2(D), D) R(S)(D)., D))) s exp(-
cG
24
2A
S p² exp( - 4A
ENIT
S pt exp( -2s log p)
p
where the second inequality uses the fact that t 2 1 for all S E A. Thus, we have
IP(Á,, (D) = () 2 p+
Now consider any exponential mechanism with scoring function R and R 2, denoted (D(D). By the same argument used in Lemma A.5 we have that
IP (A (D) S (PAc(D))
Applying the law of total probability, we have that
(AAE(D) = S*) = ((E))(AE = =(E))=(())))))) = ST189)
2 (AE)) - S118)
10sp-* 2
2 1 1+ +
Now we apply Lemma A.4
Let
R (3)
Then, we have that (1 - (PA =S*) S
PPM = S*), and as T + we have that PPAA((D) = 5*) < P(M((D) S))). and we
conclude that
- 10sp -2
P(M(D)) = 5*) 2
1 p-
21

---

Proof of Theorem 4 Let A be the bounded global sensitivity of R as in Lemma A.1 First, define
the event مم s E , (. D) Rots (S D) 24}. In the proot of Lemma A.6
we show that, given assumptions whenever
T 2 Co2 slogp
n
we have that
P P(YS E U 1 Rota(S, - Roto (D) > 24) 2 1 - 8sp
Now consider the event
8:= ([vy YS E A., R(S,D) R(5", D) N - [
R(S,D) R(S,D') s }
max sSCPl [s|=s D,D' neighhors
max
D.D' E'
By Theorem 2 we have that if R(S,D) - R(S" D) 2A for all S E A which implies that (D), D) - R()) D) > 24, then has the same bound on global sensitivity as R in Lemma A. Hence, we have that &, n Egap n E2 C E. Then, by Lemmas (A.7 and A.8 if we apply
union bound, we have that P[E]= 2 P(E, n Egap n E2a] > 8sp -22p 2p
8sp
1 18sp"
Furthermore, if we condition on E, we have that (D), D) R(S"
Then, note that
PLM(D) 1-1
1+ Ei= (1) exp(
R
24
Then, have that, if we assume
16A logp
T 2 max{Coos. s,
we have the following:
p 8) :) ( exp(
(8)
(R((S D) R(36 D))
2 exp(
-nnT
2A
4A
21 p 11
t=1
21
< 2p 2
p
where we use assumption in the last inequality. Thus, we have that
1
PPM = 2 + 2p-
Thus, we conclude that
= (((8) +
- ((EE))(((()) = =SYEE)
2 ((E)P(M(D) S718)
1- 18sp 2
1+ 2p 2
22

---

A.7 Proof of Proposition 1
Proof. Take k 2 1. Let S and (z, 0) be feasible for Problems and equivalently. Observe that,
since 2 E {0, 1] and Di=1 2; = 8, we have { 2; * 0}} = 8. Furthermore, the constraint that <8 Vj E [k 1] implies that {i 2 + 03 + (D(D) vi E [k 1]. Finally, observe that, from the constraint 3 I 8,2g E [p] it follows that B can be nonzero only on indices i for which 2, 1, and combined with the constraint Σ'= 8 2 r².
(D) S
have that || - r² We then have that Problem 9 solves the problem of minimizing lyy X x among all /3 such that || & r², and such that (3) C. {i: 2, * 03. This then becomes analogous to restricting X to the columns indexed by {i: 24 # 0} and we have that the optimization formulations in Problems and are exactly equivalent, with = {i : 2 # 0}.
B BSS algorithmic details
Experimentally. we make number of adjustments Algorithm 2 to facilitate a faster process of obtaining (D) for & E [R]. We find that, very commonly in simulated experiments, the enumerated supports (z(D), ( (p+()=))s(2)) are the (p-s)s supports that make mistake from $1 (D), and that the largest gaps in objective value across two consecutive enumerated supports often occur when the consecutive supports belong to different elements of the partition PD(D), (D)), P2(D).... (D(). Thus, in order to make the best choice of R to explore the objective landscape better but without significantly increasing computational costs, we pursued the following strategy:
1. We use Algorithm 2 to obtain $1 = (D))
2. We solve e(zx) = p
2n y X x
2n E
for all (p - 8)8 binary (zv).
vectors 2α corresponding to the supports in (D)
3. We use A Loo to obtain the optimal support that makes at least 2 mistakes from So with corresponding binary vector Σ2.
4. We check that [ax c( (zx) S c(ā²).
5. We run [ (a) the optimal support, (b) all the supports with mistake, and (c) the optimal with 2 mistakes, i.c. R 2 + (p- s)s.
Furthermore, to reduce the number of iterations of outer approximation needed in step 3 above, we added additional cuts corresponding to some of the 1-mistake vectors 2ą prior to the start of the while loop in Algorithm 2 We selected cuts by first sorting the values of X y in absolute value, then taking the top m% of the features from this sorting, and using those entries make two swaps the binary vector corresponding to $ (D), thus generating 1-mistake vectors that then used for cuts.
B.1 Solving c(ε) using PGD
At each iteration of outer approximation in Algorithm 2 solve
we use projected gradient descent (PGD) to
32
c(2)
min | s
ll
x
2n
2n
2
9(1)
We have that
Vg(
_ x x (X,3 y)
AB
-x
(X3 (XT(XX 9) + ADiag(
-/221 )3)
1
2 ...
n
n
n
23

---

where for vectors Iº denotes element-wist wise division. Then note that
(X y) + ADiag(
1
5)3)
]3)
(x
=xT(3X3' (x33' - y) + +ADiag
183) 12
21
2p
"
21 ...
3p
1 x X ADiag(
1 (3 -3)1122
n
2p
Amax (X T X + ADiag(
))[[ - 111l
21
:."
2p
Then setting L = 1 Amax (XTX
X + ADiag 1
21 :.. 2,2 )), the PGD update is then
r(3, 1 (g())))
B+11
max(r: 1| Vg(
B.2 Heuristic to kickstart outer approximation
In order to provide good initialization for zo in Algorithm Algorithm in (5). Specifically, we consider the problem
we use the heuristic taken from 2
1
min AERP:
Clyy X31 + \|13113
oS 2n
M(B)
Note that Vh(33) 1/( (X (X y) + x3) and that
1 X (XB - y) + A 1 (x (X3' 4) AB') = || (X x XT) - B')
n2
Amax
(X x XT) 3 ||2
n
Then, letting L =mmxx (X
+ AT), we run the following heuristic:
1. Initialize B1 IR" such that ||1:10
2. For t 2 1:
(a) Sort the entries of B, / Vh(B,) order of decreasing absolute value, let / denote the index set of the largest entries.
(b) Set (Be+1); (3e 1 Vh())) if E I, and (Bz41); otherwise.
C Modifications pertaining to hinge loss
Similarly to BSS, we now consider the following sparse classification problem:
min BERP
max (x s.t. ||,]|o I s, 118112 S r
(C.1)
Our objective function then becomes
R(S,D) = min
max (((()))33)) s.t. || Sr
,EER #1
where D (X, y) as before, and (æ:)s E R is the i-th row of X with columns indexed by S.
We first present a result analogous to Lemma A.1 in order to bound the global sensitivity for the case of hinge loss.
Lemma C.1. Suppose that |X,,] S for i [m].] E pp. Then,
< n
+ rb.
24

---

Proof. Suppose D, D' are two neighboring datasets. Fix support S E O and suppose
arg min S
max v1((z{()))))
71 1=1
Then note that
R(S, - R(,,))
- - ma (x:
Σ max (((((()))))
Let us assume without loss of generality that D, D'' differ in the n-th observations. Then we have that
max (0, 1 - (x(()))))))>
max (0, 1 v(((<=));)))
(max((, (0, - yn (0<(())))))) (<x))
max (0, (((((()))))))))
< 71
max {0. yn ((xn))
(1 rb.
where the last inequality uses the fact that |S|] = and the Cauchy-Schwarr inequality, since - -yn((=))33)) (,x,((z))))) ( s (((x)))l2ll@@]1
C.1
Optimization formulation
As in the BSS case, we consider the penalized form of Problem C.I in order to obtain S₃ (D) for E [R]. We define
3
czz) = min
ma 1 v
+
||
==1
1=1
and we seek solve
min c(z)
subject to z (0,
Σ 21 = 8,
1=1
Σ
2
Vj jk 1].
(E))(D)
Given z {0, 1), define (0, 1]" in Section
Let
p 82
B E arg
min B
max(0, 1 ye(x B)} +
n
2
We then have
> (B.)2 (Ve(²)); = n 8?
and we run Algorithm with these modifications.
C.2 Solving c(2) using Projected Subgradient method
At each iteration of outer approximation, solve
min max - 22 9)3
25

---

using a projected subgradient method. Using the subgradient
(
1 n
Σ
y: + 2ADiag
1
2 13
up
we run the update
r(3: m=9(3:)
B max(r, 113 (B))2))
where 7/2 Vt'
D
Additional experimental results
D.1 BSS
D.1.1 Prediction accuracy and utility loss
In effort to compare prediction accuracy across methods, performed 70/30 random train/test split and implemented Algorithm [23] on the training data, using half the privacy budget (i.e., r/2) for variable selection with the top-R, mistakes, Samp-Agg, MCMC methods, and the remaining half for model optimization via objective perturbation (Algorithm in (3 to obtain the regression coefficients under the privacy budget (Speiv). We experiments with p 100,
= 2, SNR = 5, and 0.1, and conducted 10 independent trials for each value of n. For each trial, we drew 100 samples from the distribution corresponding each algorithm. For MCMC, used 100 independent Markov chains per trial. After running objective perturbation over the selected supports, we obtained 1000 distinct coefficient vectors for each method and each n. and computed the average MSE on the test data. The choices of 120 and MCMC iterations = 1000 were made keep the runtimes comparable, in our support recovery results. We summarize the results below, showing that our top-At R and mistakes methods outperform the competitor algorithms in prediction accuracy for sufficiently large n.
p= 100
p 1000
2.00
Top R
1.8
Mistakes
Mistakes
Lasso Subsample
Lasso Subsample
1.50
MCMC
1.4
MCMC
Test MSE
1.25
1.2
1.00
1.0
0.15
0.50
0.25
0.4
3000
4000
5000
6000
7000
6000
7000
8000
9000
10000
11000
0
(a) = 100
(b) = 1,000
Figure D.1: Numerical experiments for p 100 and p 1,000, with s 5, SNR=5, / 0.1,
and = 2. The penalty parameter À in Algorithm 2 was set to 120 and 250 for figures D and D.Ib respectively. On the --axis, we vary the value of n and plot the average test MSE across 10 independent trials. Error bars denote the mean standard error.
Furthermore, we ran experiments to evaluate the utility loss, defined as the gap between the objective at 3priv and the objective at 3 (the optimal BSS solution), using the same parameters as in the experiments above. with the prediction accuracy results, our top- and mistakes methods outperform the competitor algorithms in terms of utility loss when is large enough.
26

---

p=100
p=1000
lop
lop
0.20
Mistakes
Mistakes
Lasso Subsample
0.05
Lasso Subsample
Obj utility loss
0.15
MCMC
MCMC
0.04
Obj-Pert
0.03
0.02
0.05
0.00
2000
3000
4000
5000
6000
7000
4000
5000
6000
7000
8000
(a) 100
(b) p 1,000
Figure D.2: Numerical experiments for = 100 and p 1,000, with 5, SNR=5, p 0.1, and The penalty parameter À in Algorithm 2 was 120 and 250 for figures [D22 and D.2b respectively. On the r-axis, we vary the value of n and plot the average utility loss across 10 independent trials. Error bars denote the mean standard error.
D.1.2 Ablation studies
In this section, we present several ablation studies in order show the effect changing R, and (be, by) on the fraction of correctly recovered supports and on the Fl score. For the following results, we ran 10 independent trials and drew 100 samples from the distribution corresponding to each of our algorithms.
Ablation study of R. As noted in Appendix B we observe empirically in the results below that the largest gaps in objective value across two consecutive enumerated supports often occur when the consecutive supports belong to different elements of the partition Po(D), Pi (D), Pa(D)....,PP
Moreover, as intuitively expected and shown formally in Lemma A.5 of the appendix, increasing R
has a positive effect on support recovery. We used 120 and (bg. by) = (0.5.0.5)
0.8
0.7
0.6
0.4
0.3
0.1
Fraction correct
0.0
3000
6000
9000
12000
15000
R
Vertical bars denote an objective gap corresponding an increase Hamming distance from the optimal support. On the ar-axis, we vary the value of R and plot the average Fl score and average fraction of correct supports across 10 independent trials for top- R algorithm. Error bars denote the mean standard error.
Figure D.3: Numerical experiments for 1,000, p 100, = 3, P=
27

---

Ablation study of X. We present results below for tuning X. We used R = 2 (p -pss s)s and (b,, by) (0.5, 0.5). We observe that choosing a very large X can have negative effect on support recovery. However, increasing À is beneficial for the runtime of our outer approximation solver. Therefore, we choose a moderate value for À ensure that the runtime of our method remains comparable to that of MCMC.
1.00
0.975
0.96
0.950
0.94
supports
0.92
0.900
0.90
0.675
0.88
Top
Top
Mistakes
Mistakes
0.825
120
140
160
180
200
160
X
(a) F1 score
(b) Proportion correct supports
Figure D.4: Numerical experiments for n 4,000, p 100, 5, SNR=5, p = 0.1, ar-axis, we vary the value of X, and in Figures D.4a and D.4b and respectively plot the average F1 score and average fraction of correct supports across 10 independent trials for the top- R and mistakes algorithm. Error bars denote the mean standard error.
On the
Ablation study of (b., by)- We present results below for tuning (ba, by). For simplicity, we set by = by. We used R + (p - s)s and À 300. Our results indicate that choosing the clipping constants too small too large negatively affect support recovery quality. In practice, these hyperparameters can tuned via cross-validation.
1.0
1.0
0.8
supports total
0.8
0.7
0.7
0.6
0.5
0.4
0.4
0.3
Top
Mistakes
Mistakes
0.2
0.50
0.75
1.00
1.25
1.50
0.50
0.75
1.00
1.25
1.50
6,
(a) score
(b) Proportion correct supports
Figure D.5: Numerical experiments for 4,000. 100. 5, SNR=5, 0.1, 1. On the --axis, we vary the value (bg. by), and Figures D.5a and D.5b and respectively plot the average Fl score and average fraction of correct supports across 10 independent trials for the top-R and mistakes algorithm. Error bars denote the mean standard error.
D.1.3 Support recovery
In this section, we present additional experimental results for support recovery in BSS, for varying values of s. p, and SNR. As noted Section 5 the magnitude of the penalty parameter À in Algorithm and the number of MCMC iterations for the algorithm by [2 were chosen such that the average runtimes of the methods across different trials were comparable. That is, our average runtime for each value of n was most the runtime of the MCMC algorithm.
28

---

To demonstrate the power of our MIP-based estimator, we note that for 1000, s 8, which is one of the settings that we ran, we have (L ~ 10th To put this in perspective, to use the standard exponential mechanism and enumerate all feasible supports, assuming computing each feasible support takes 10 seconds and 16 bits of storage, one would need 764 thousand years and 48 million terabytes of storage. This shows the usefulness our MIP-based estimator in practice, enabling us to solve BSS with DP for problem sizes that otherwise would be prohibitive.
p=1000
1.0
contect supports tottll tr iralls
C.8
0.6
Top R
Mistakes
Lasso Subsamole
0.4
MCMC.
0.0
3000
3500 4000 4500 5000 5500 6000 7000
(a) s = 5, SNR=5, 0.1,
p=1000
cortect curports tatal trtake
0.8
supports total
0.0-
0.85
Top R.
0.4
0.80
assso Subhaa mphs
Top R
0.75
Mistakes Lesso Subsamp
0.70
MCMC
0.65
5000
6000
7030
5000
9000
3000
3500
4000
4500
5300
6000 8500 7000
(b) 3 SNR=2, p= 011,==1
(c) s SNR=5,
0.1,
p=1000
1-1000
1.0
1.0
supports
0.6
0.4
0.4
R
R. Mistaktes Lasso Subsample
0.2
Mistakes
0.2
Lesso Subsample
0.0
MCMC
0.0
MCMC
6030
7000
E000
9000
11000
17000
17500 18000 16500 19000 19500 20000
(d) = 5, SNR=5, = 0.5,
(e) = 8,
p 0.1,
Figure D.6: Numerical experiments for different values of s, p. e and SNR, with p = 000. The penalty parameter X in Algorithm 2 was set to 250 for figures D 6a D.6d and 600 for figure D.6c The number of MCMC iterations was set to 8000 for figures D.6a D.6d and 10000 for figure D.6e On the r-axis, we vary the value of and plot the average proportion of draws across 10 independent trials that recovered the right support for each corresponding algorithm. Error bars denote the mean standard error.
29

---

p=1000
1.0
0.8
sccore
0.6
0.4
Top Mistakes
0.2
Lasso Subsample MCMC
3000
3500
4000
4500 5000 5500 6000 6500 7000
(a) + 5, SNR- 0=
p=1000
p=1000
0.99
0.98
0.6
0.97
Mistakes
0.96
Lasso Subsample
0.4
Top
0.95
0.2
Lasso Subsample
0.94
0.0
MCME
0.93
4000
5000
6000
7000
8000
9000
3000
3500
1000
4500
5000
5500 6000 6500 7000
(b)
5, SNR=2, p c.1.
(c) 5. SNR=5,
0.1.
p-1000
p-1000
1.0
1.0
Mistakes
score
Lasso Subsample
MCMC
0.2
Mistakes Lasso Subsample
0.0
7000
BE0D
17000 17500
19500
11
(d) 8 = 5, SNR=5, p 0.5, e=
(c) 8 = 8, SNR=5, p= 0.1, =1
Figure D.7: Numerical experiments for different values of s, p. c and SNR, with p = 000000 The penalty parameter X in Algorithm 2 was set to 250 for figures D. 7a D.7d and 600 for figure D.7e
The number MCMC iterations was set to 8000 for figures D.7a 7d and 10000 for figure D.7e On the r-axis, we vary the value of and plot the average Fl score across 10 independent trials for each corresponding algorithm. Error bars denote the mean standard error.
30

---

p=100
correct supports /totall trials
0.8
0.6
0.4
Top
0.2
Mistakes
Lasso Subsample
MCMC
0.0
2000
3000
4000
5000
6000
7000
(a) 8 5, p 0.1,c=1
p=100
0.95
0.6
correct sJpports
0.90
Mistakes
Lasco Subsample
0.4
0.85
Top
0.2
Lasso Subsample
0.0
MCMC
0.75
2000
3030
4000
5000
6000
7000
2000
3000
1000
5000
6000
7000
(b)
5, SNR=2, p 0.1.
(c) 5, SNR=S, 0.1,
p-100
o-100
1.0
1.0
Top R
(atal trials
correct supports total bries
0.8
Lasso Subsample
MCMC
0.6
Top
Misakkes
Lasso Subsample
4000
11
(d) == 5, SNR=5, p 0.5,==1
(c) 8 = 8, SNR=5, 0.1, =
Figure D.8: Numerical experiments for different values of s, p. and SNR, with p 100. The penalty parameter X in Algorithm 2 was set to 120 for figures D. D.8d and 350 for figure D.Se The number of MCMC iterations was set 1000 for figures D 8a D.6el On the æ-axis, we vary value of and plot the average proportion of draws across 10 independent trials that recovered the right support for each corresponding algorithm. Error bars denote the mean standard error.
31

---

p=100
1.0
0.8
score
0.4
Top
Mistakes
0.2
Lasso Subsample
MCMC
0.0
2000
3000
4000
5000
6000
7000
(a) 8 5, SNR=5, p 0=
p=100
0.8
0.99
FLscore
0.98
Mistakes
Lasso Subsample
0.97
Top
0.2
Lasso Subsample MCMC
0.95
2000
3030
4000
5000
6000
7000
2000
3000
1000
5000
6000
7000
(b)
5, SNR=2, 0.1.
(c) 5. SNR=5, p 0.1,
p-100
--100
1.0
0.8
0.6
Mistakes
score
0.6
Subsample
0.4
MCNC
0.4
Mistakes Lasso Subsample
0.0
2030
4000
5000
600C
7000
12000
11
(d) == 5, 5 0.5, =1
(c) 8 = 8, SNR=5, = 0.1, =1
Figure D.9: Numerical experiments for different values of s, p. and SNR, with p 100. The penalty parameter À in Algorithm 2 was set to 120 for figures D.9a D.9d and 350 for figure D.9e| The number of MCMC iterations was set to 1000 for figures D.9a D.7e On the r-axis, we vary the value of n and plot the average Fl score across 10 independent trials for each corresponding algorithm. Error bars denote the mean standard error.
32

---

p=10000
1.0
0.8
scooe
0.6
0.4
Top
0.2
Mistakes
Lasso Subsample
MCMC
0.0
7000
8000
9000
10000
11000
Figure 1:10: Numerical experiments for 10, 000, with 5. SNR=5, p 0.1, and c = 1. The penalty parameter Algorithm was set to 600, and the number of MCMC iterations was set to 100, 000. On the z-axis, we vary the value of n and plot the average FI score across 10 independent trials for each corresponding algorithm. Error bars denote the mean standard error.
D.2 Hinge Loss
In this section, present experimental results in the setting of sparse classification, presented in Problem C.I
To generate our data, we first generate 2; = af 3 + for in), where #1.... Fn ~ N( Σ) IRP and the independent noise follows 6 <~ N~ N0, (0,00,,, ²I,) where I, is the identity matrix of size n. We then draw u, ~ Uniform(0, 1] for E [n]. and we set
üd
u, > (24)
yi {-1
where (()) = 1te=*
otherwise
Moreover, for i,j E (p). we Ei, أزكا and nonzero coordinates of 3 to take value 1/v5 indices 31, 2s 1}. We define the Signal to Noise Ratio SNR ||x 1
As with the BSS results, our methods show favorable empirical support recovery in both low and high-dimensional settings, with our mistakes method outperforming our top- R method.
p=
p=1000
trials trials
* correct suppertsitotal
Top Mistkkes
Top
Mistakes
1- Lasso Subsample
Lasso Subsample
0.2
2000 2500 3030 3500 1000 4500 5000 5500 6000
3000
3500
4000
4500
3000
3500
6000
(a) 100
(b)
1,000
Figure D. I: Numerical experiments for p = 100 and p = 1,000, with 5, SNR=5, == 0.1, and 1. The penalty parameter À in Algorithm 2 was set to 90 and 100 for figures) D. and D. respectively. On the z-axis, we vary the value of n and plot the average proportion of draws across 10 independent trials that recovered the right support for each corresponding algorithm. Error bars denote the mean standard error.
33

---

p-100
p-1000
0.8
Top
11 scooeee
Top
Lasto Subsample
Lasso Subsample
0.2
0.2
0.0
2000 2500 3000 3500 1000 4500 5000 5500 6000
3000
3500
1000
1500
5000
5500
6000
(a) p 100
(b) p 1,000
Figure D.12: Numerical experiments for p 100 and p 1,000, with - 5, SNR=5, 0.1, and The penalty parameter À in Algorithm 2 was to 90 and 100 for figures D.12a and DD.22]] respectively. On the w-axis, we vary the value of n and plot the average Fl score across 10 independent trials for each corresponding algorithm. Error bars denote the mean standard error.
p=10000
1.0
0.8
F1 score
0.6
Top Mistakes
0.4
Lasso Subsample
6000
7000
8000
Figure D.13: Numerical experiments for 10,000, with = 5, SNR=5, p 0.1, and e = 1. The penalty parameter À in Algorithm 2 was to 1 170. On the r-axis, vary the value of and plot the average Fl score across 10 independent trials for each corresponding algorithm. Error bars denote the mean standard error.
NeurIPS Paper Checklist
1. Claims
Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Question:
Answer: [Yes]
In the abstract and introduction we claim that present two pure-DP methods that can be broadly applied to loss functions and yield strong theoretical guarantees in the case of BSS. We also claim they achieve state-of-the- thet---rtt empirical performance. This matches with the theoretical and experimental demonstrations in the paper.
Justification:
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made in the paper.
34

---

The abstract and/or introduction should clearly state the claims made, including the contributions made the paper and important assumptions and limitations. A No or NA answer to this question will not be perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalize to other settings.
It is fine to include aspirational goals as motivation as long as it is clear that these goals are not attained by the paper.
2. Limitations
Question: Does the paper discuss the limitations of the work performed by the authors?
Answer: [Yes]
Justification: We point out some of the limitations of our work and potential future research directions in Section G
Guidelines:
The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those are not discussed in the paper.
The authors are encouraged to create separate "Limitations" section in their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be.
The authors should reflect on the scope of the claims made, e.g., if the approach was only tested few datasets with few runs. In general, empirical results often depend on implicit assumptions, which should be articulated.
The authors should reflect on the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution is low or images are taken in low lighting. Or speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails handle technical jargon.
The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness.
While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, worse outcome might be that reviewers discover limitations that aren't It acknowledged judgment and recognize that individual actions in favor of transparency play an impor- tant role in developing norms that preserve the integrity of the community, Reviewers will be specifically instructed to not penalize honesty concerning limitations.
the paper. The authors should use their best
3. Theory assumptions and proofs
Question: For each theoretical result, does the paper provide the full set of assumptions and complete (and correct) proof?
Answer: [Yes]
Justification: We explicitly mention all the assumptions in the statement of the theorems and lemmas. The proofs of these results can be found the supplementary materials sections.
Guidelines:
The answer NA means that the paper does not include theoretical results.
All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated or referenced in the statement of any theorems.
The proofs can either appear in the main paper or the supplemental material, but if they appear in the supplemental material, the authors are encouraged provide short proof sketch provide intuition.
35

---

Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix or supplemental material.
Theorems and Lemmas that the proof relies upon should be properly referenced.
4. Experimental result reproducibility
Question: Does the paper fully disclose all the information needed to reproduce the main ex- perimental results the paper to the extent that affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided or not)?
Answer: [Yes]
Justification: We clearly state the values of all the parameters needed reproduce our experimental results in Section 5 and in the additional experiments section of the appendix.
Guidelines:
The answer NA means that the paper does not include experiments.
If the paper includes experiments, a No answer to this question will not perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
If the contribution is dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is novel architecture, describing the architecture fully might suffice, or if the contribution is specific model and empirical evaluation, it may be necessary to either make possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way to accomplish this, but reproducibility can also be provided via detailed instructions for how to replicate the results, access to a hosted model (e.g.. in the case of large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed.
While NeurIPS does not require releasing code, the conference does require all submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) the contribution is primarily a new algorithm, the paper should make it clear how reproduce that algorithm. (b) If the contribution is primarily the architecture clearly and fully.
model architecture, the paper should describe
(c) the contribution is a new model (e.g., a large language model), then there should either be way to access this model for reproducing the results or way to reproduce the model (e.g., with an open- open-source dataset or instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in authors are welcome to deseribe the particular way they provide for reproducibility. In the case of closed-source models, it may be that access to the model is limited in some way (e.g., to registered users), but it should possible for other researchers have some path reproducing verifying the results.
cases, in which
5. Open access to data and code
Question: Does the paper provide open access to the data and code, with sufficient instruc- tions to faithfully reproduce the main experimental results, as described in supplemental material?
Answer: [Yes]
Justification:
Our code is provided in the supplementary materials.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurlPS code and data submission guidelines (https://nips.cd
for more details.
36

---

While we encourage the release of code and data, we understand that this might not be possible, so "No" an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (e.g., for new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurIPS code and data submission guidelines (https for more details.
7 //nips. http://cpppppccccceeeeeeee/eeeee
The authors should provide instructions data access and preparation, including how to access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. only subset of experiments are reproducible, they should state which ones are omitted from the script and why.
At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable).
Providing as much information as possible in supplemental material (appended the paper) is recommended, but including URLs to data and code is permitted.
6. Experimental setting/details
Question: Does the paper specify all the training and test details (c.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
Justification: Our experiments are run on simulated data. Therefore, technically, we just have the testing step. We clearly define our data generation process in Section
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core of the paper to a level of detail that is necessary to appreciate the results and make sense of them.
The full details can be provided either with the code, in appendix, as supplemental material.
7. Experiment statistical significance
Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Question:
Answer: [Yes]
Justification: Our experiments include error bars reflecting the standard error of the mean across independent trials, capturing variability due to randomness the experimental setup.
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results accompanied by enror bars, confi- dence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper.
The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, overall run with given experimental conditions).
The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.)
The assumptions made should be given (e.g., Normally distributed errors).
It should be clear whether the error bar is the standard deviation the standard error of the mean.
It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have 96% CI, if the hypothesis of Normality of errors is not verified.
37

---

For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out range (e.g. negative enror rates).
If error bars are reported in tables or plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables in the text.
8. Experiments compute resources
For each experiment, does the paper provide sufficient information on the com- puter resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?
Question:
Answer: [Yes]
Justification: We have added information about the computational resources we used for our
experiments in Section 5
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage.
The paper should provide the amount of compute required for each of the individual experimental runs as well as estimate the total compute.
The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary or failed experiments that didn't make into the paper).
9. Code of ethics
Question: Does the research conducted the paper conform, in every respect, with the NeurlPS Code of Ethics http::///haarrssscccpppppcccccct
Answer: [Yes]
Justification: This research work does violate any NeurIPS Code of ethics.
Guidelines:
The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.
If the authors answer No, they should explain the special circumstances that require deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g.. if there is special consid- eration due to laws or regulations their jurisdiction).
10. Broader impacts
Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?
Answer: [Yes]
Justification: We point toward the potential positive impacts of our work in Section To the best of our knowledge, we could not think of any negative societal impact of our work.
Guidelines:
The answer NA means that there is no societal impact of the work performed.
If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact.
Examples of negative societal impacts include potential malicious unintended uses (e.g., disinformation. generating fake profiles, surveillance), fairness considerations (e.g.. deployment of technologies that could make decisions that unfairly impact specific groups). privacy considerations, and security considerations.
The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster.
38

---

The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional unintentional) misuse of the technology.
there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of MI.).
11. Safeguards
Does the paper describe safeguards that have been put in place for responsible release of data models that have high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [NA]
Justification:
The paper poses no such risks.
Guidelines:
The answer NA means that the paper poses such risks.
Released models that have high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere usage guidelines or restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recognize that providing effective safeguards challenging, and many papers do require this, but we encourage authors to take this into account and make best faith effort.
12. Licenses for existing assets
Are the creators or original owners of assets (e.g.. code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?
Question:
Answer: [Yes]
Justification: We clearly cite the relevant papers for any used code in the paper along with the license information.
Guidelines:
The answer NA means that the paper does existing assets.
The authors should cite the original paper that produced the code package or dataset.
The authors should state which version of the asset used and, possible, include a URL.
The name of the license (e.g., CC-BY 4.0) should be included for each
For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should be provided.
If assets released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, wwwppsai thcode. com/ datasets has curated licenses for some datasets. Their licensing guide can help determine the license dataset.
For existing datasets that are re- e-packaged, both the original license and the license of the derived asset (if has changed) should be provided.
If this information is not available online, the authors are encouraged to reach out to the asset's 's creators.
13. New assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [NA]
39

---

Justification:
The paper does not release new assets.
Guidelines:
The answer NA means that the paper does not release new assets.
Researchers should communicate the details the dataset/codeemoodll as part their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent was obtained from people whose asset used.
At submission time, remember anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file.
Crowdsourcing and research with human subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well details about compensation (if any)?
Answer: [NA]
Justification:
The paper does not involve crowdsourcing nor research with human subjects.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing research with human subjects.
Including this information in the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then as much detail as possible should be included in the main paper.
According to the NeurlPS Code of Ethics, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector.
15. Institutional review board (IRB) approvals or equivalent for research with human subjects
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or equivalent approval/review based the requirements of your country or institution) were obtained?
Answer: [NA]
Justification: The paper does not involve crowdsourcing nor research with human subjects.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Depending the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurIPS Code of Ethics and the guidelines for their institution.
For initial submissions, not include any information that would break anonymity (if applicable), such as the institution conducting the review.
16. Declaration of LLM usage
Question: Does the paper describe the usage of L.MMs if is an important, original, or non-standard component of the core methods in this research? Note that if the L.M is used only for writing, editing, formatting purposes and does not impact the core methodology, scientific rigorousness, originality of the research, declaration is not required.
Answer: [NA]
Justification: The core method development this research does not involve L.MMs as any important, original, or non-standard components.
40

---

Guidelines:
The answer NA means that the core method development in this research does not involve LLMs as any important, original, non-standard components.
Please refer to our LLM policy for what should or should not be described.
41


Reviews:
[{"review_id": "CQSDxoueWU", "reviewer": "Reviewer_fXeC", "summary": "Sparse variable selection is important for high-dimensional data, it can improve the interpretability and reduce the overfitting of the models. This paper proposes two new and scalable differentiable privacy estimators for sparse variable selection, “top-R” and “mistake”. And both of them employ the idea of mixed integer programming to generate the candidate supports, avoiding the exhaustive search in the exponential mechanism. They also prove that these two estimators satisfy pure (ε,0)-DP. And the best subset selection matches the non-private minimax rate of top-R in the low privacy regime and those of mistakes in the high privacy regime. Experiments on least squares regression and hinge classification prove that the proposed method is able to recover the correct support better than SOTA, and is still computationally feasible.", "questions": "How to choose the parameters R and λ? An ablation study would be useful.\nIn practice, how should one set the clipping bounds bx and by to balance bias and privacy?", "limitations": "Yes.", "rating": 4, "confidence": 2, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 4, "originality": 3, "strengths_and_weaknesses": "Strength: \n1) This paper gave us an insight into using MIP solvers to approximate the exponential mechanism and solve the problem that effectively solves combinatorial optimization problems under the strong guarantee of pure-DP. \n2) The author provided solid theoretical proof. They provide a detailed proof for pure (ε,0)-DP and analyze the support recovery for BSS. \n3) They gave an empirical comparison with state-of-the-art methods, showing great results.\n\nWeakness: \n1) The proposed methods rely on several assumptions (1-4), but may not be practical for highly correlated features. How will privacy guarantees or support recovery degrade when these assumptions are not satisfied?\n2) Experiments only focus on support recovery (i.e., identifying non-zero coefficients), but do not mention the prediction accuracy. Furthermore, the paper should address how to choose several hyperparameters and discuss the sensitivity of choosing them.", "was_revised": false, "final_justification": "", "rebuttal": "Dear reviewer fXeC,\n\nWe sincerely thank you for your constructive feedback on our paper. We address your concerns below:\n\n**Practicality of assumptions**\n\nAssumption 1 regarding the boundedness of the dataset is a common assumption in the high-dimensional statistics and DP literature (we refer to [1,2] and references therein). We also make a note in the paper that, while in practice one may not know the exact values of $b_x,b_y$, or such values might not exist, one can clip the values of $X,y$ to satisfy the assumption. Importantly, the only assumption required to ensure the privacy of the top-R method is Assumption 1, which can always be satisfied using clipping as discussed.\n\nSimilarly, assumptions 2-4 regarding the boundedness of the true parameter vector, the Sparse Riesz Condition (SRC), and the growth of the sparsity parameter, are standard assumptions often found in the sparse linear regression literature, and we offer numerous citations in section 3 of the paper. \n\n**Experiments on prediction accuracy**\n\nIn an effort to compare prediction accuracy across methods, we did a 70/30 random train test/split and implemented Algorithm 2 in [3] on the train data, wherein we use half the privacy budget (i.e., $\\epsilon/2$) for variable selection using the top-R, mistakes, Samp-Agg, or MCMC methods, and the remaining half of the budget for model optimization using objective perturbation (Algorithm 1 in [3]) to obtain the regression coefficients under the privacy budget. We ran experiments with $p=100, s=5, \\epsilon = 2, SNR= 5, \\rho = 0.1$, and ran 10 independent trials for each value of $n$, where, for each trial, we drew 50 times from the distribution corresponding to each algorithm, and for MCMC, we used 50 independent Markov chains for each trial. After running objective perturbation over the supports, we thus had 500 distinct coefficient vectors for each method and for each $n$, and we computed the average MSE on the test data for every method. Our choices of $\\lambda=120$ and MCMC iterations = 1000 were chosen to make the runtimes comparable, as done in our paper. We summarize our results below, showing that our top-R and mistakes methods outperform the competitor algorithms in terms of prediction accuracy for $n$ large enough.  \n\nAverage Test MSE:\n\n|                             | $n=3000$ | $n=4000$ | $n=5000$ | $n=6000$ |\n|-----------------------------|-----------|-----------|-----------|-----------|\n| top-R                       | 1.894    |  0.650    | 0.290    | 0.276    |\n| mistakes                    | 1.064    | 0.408    | 0.300    | 0.275    |\n| Samp-Agg                    | 1.490   | 0.980    | 0.677    | 0.480    |\n| MCMC                        | 0.532    | 0.386    | 0.329   | 0.317   |\n\n\n**Ablation study to choose hyperparameters**\n\nFor more details on how the choice of R was selected for our experiments, we provide additional information in section B of the appendix. In summary, choosing a larger value of R improves the support recovery quality, but increases the computational cost of the method. Our choice of $R=2+(p-s)s$ is designed to encompass the optimal solution, the solutions with the support differing from the optimal by one, and at least a solution with two differences in the support relative to the optimal. This choice seems to result in reasonable runtime, while leading to good model selection accuracy.\n\nIn addition, we have performed an ablation study to show the effect of changing $R, \\lambda$, and $(b_x, b_y)$ on the fraction of correct supports as well as F1 score, fixing $n=3500, p=100, s=5,  \\epsilon = 1, SNR= 5$, and $\\rho = 0.1$. We ran 10 independent trials and drew 100 times from the distribution corresponding to our algorithms. \n\n- Ablation study of $R$: We present our results below for $R$. Note we used $\\lambda=120$ and $(b_x, b_y) = (0.5,0.5)$. As intuitively expected and shown formally in section A.3 in the appendix, increasing $R$ has a positive effect on support recovery. \n\n| | $R = 177$ | $R = 277$ | $R = 377$ | $R = 477 ~(\\text{i.e.}  ~ 2+(p-s)s)$ |\n|-|------------|-----------|-----------|-----------|\n|Average # correct supports / number of trials|0|0.0025|0.005|0.57\n|Average F1 score|0.047|0.049|0.0535|0.59\n\n- Ablation study of $\\lambda$: We present our results below for tuning $\\lambda$. Note we used $R=2+(p-s)s$ and $(b_x, b_y) = (0.5,0.5)$.\nWe observe that choosing a very large $\\lambda$ can have a negative effect on support recovery. However, increasing $\\lambda$ is beneficial to the runtime of our outer approximation solver. Therefore, we choose a moderate value for $\\lambda$ to ensure the runtime of our method remains comparable to MCMC.\n\nAverage F1 scores:\n|  |$\\lambda=120$ | $\\lambda=140$ | $\\lambda=160$ | $\\lambda=180$|\n|--|---|---|----|---|\n|top-R|0.550|0.534|0.374|0.284|\n|mistakes|0.987|0.988|0.982|0.970|\n\nAverage # correct supports / number of trials:\n|  |$\\lambda=120$ | $\\lambda=140$ | $\\lambda=160$ | $\\lambda=180$|\n|--|---|---|----|---|\n|top-R|0.519|0.501|0.328|0.236|\n|mistakes|0.938|0.941|0.910|0.852|\n\n\n- Ablation study of $(b_x,b_y)$: For simplicity, we set $b_x = b_y$. Note we used $R=2+(p-s)s$ and $\\lambda = 300$. We present the results below. Our results imply that choosing the clipping constants to be too small or large can have a negative effect on the support recovery quality. In practice, one can tune these hyper parameters via cross-validation. \n\n\nAverage F1 scores:\n|  |$b_x=0.5$ | $b_x=1$ | $b_x=1.5$ | $b_x=2$|\n|--|---|---|----|---|\n|top-R|0.111|0.3915|0.0905|0.0505|\n|mistakes|0.952|0.9795|0.9445|0.8065|\n\nAverage # correct supports / number of trials:\n|  |$\\lambda=120$ | $\\lambda=140$ | $\\lambda=160$ | $\\lambda=180$|\n|--|---|---|----|---|\n|top-R|0.0625|0.35|0.04|0|\n|mistakes|0.7675|0.9|0.735|0.2575|\n\n\nReferences:\n\n[1] Roy, S., Wang, Z., & Tewari, A. (2024). On the computational complexity of private high-dimensional model selection via the exponential mechanism. In Advances in Neural Information Processing Systems (NeurIPS).\n\n[2] Chen, L. Y., & Lee, S. (2023). Sparse quantile regression. Journal of Econometrics, 235(2), 2195–2217.\n\n[3] Kifer, D., Smith, A., & Thakurta, A. (2012). Private convex empirical risk minimization and high-dimensional regression. In Proceedings of the Conference on Learning Theory (COLT), JMLR Workshop and Conference Proceedings."}, {"review_id": "KlM7laTG0P", "reviewer": "Reviewer_JE2d", "summary": "The paper focuses on the problem of sparse variable selection with differential privacy. Here we have some per-example loss (such as least squares or the hinge loss) on a model $\\beta$, a dataset $X, y$, and our goal is to find a model $\\beta$ that minimizes the empirical loss on $X$ subject to the model being chosen with differential privacy, and $\\beta$ having at most $s$ non-zero coordinates. Variable selection focuses just on selecting the right $s$ coordinates; in particular, the success of a variable selection algorithm is the fraction of non-zero coordinates in $\\beta^*$, the optimal model, that are also chosen by the algorithm.\n\nThe authors consider an exponential mechanism-based approach to this problem. In their first variant, the loss for choosing a certain set of $s$ coordinates is the minimum loss achieved by any model only supported on those coordinates. They are able to show a bound on the sensitivity of the exponential mechanism under mild assumptions (that can be enforced via clipping), but since there are $\\binom{p}{s}$ different sets of coordinates to choose, running this exponential mechanism is infeasible. They instead propose the following modification: rather than compute the loss for all sets, compute it for the best r sets. For all other sets, set their loss to the r-th best set's loss. They show this is also private, and in fact gets better privacy than the \"default\" value of $r = \\binom{p}{s}$ improving as r decreases. They also show that the top R supports can be computed by solving a series of mixed integer programs which are efficient in practice for p ~= 10000. This eliminates the need to enumerate all solutions. The second variant the authors propose is to set the loss in the exponential mechanism to be the number of wrong coordinates in the set. This is easier to sample from but requires a separation between the loss of the first and second for privacy.\n\nTheoretically, the authors show \"top-r\" recovers the best selection with high probability when the smallest non-zero coordinate of the optimal model is $\\Omega(\\sqrt{\\max\\{1, s^2 / \\epsilon\\} log p / n})$ for a moderate value of r. The non-private requirement is $\\Omega(\\sqrt{\\log p / n})$, so top-R matches the non-private rate up to constants when $\\epsilon \\rightarrow \\infty$ or when the number of coordinates to pick is small. Mistakes requires the non-zero coordinates to be $\\Omega(\\sqrt{(s / \\epsilon) \\cdot \\log p / n})$, i.e. again is optimal when $\\epsilon$ is sufficiently large but requires the additional separation condition for privacy. Empirically, the authors compare top-r and mistakes to an MCMC-based approach of the previous work of Roy et al., and also the Samp-Agg algorithm of Smith et al. Mistakes recovers almost all coordinates in most settings, and top-R outperforms Samp-Agg and outperforms the MCMC based approach when the dataset size is large enough.", "questions": "-Can the separation between the loss of the first and second-best solution be enforced in any clipping-like way? As discussed before, requiring a high-probability condition to hold for the DP guarantee of mistakes is somewhat problematic.\n-If we were to e.g. run objective perturbation or another DP-ERM algorithm over the coordinates selected by top-r or mistakes using the same overall privacy budget as Samp-Agg, would top-r and mistakes still outperform Samp-Agg?", "limitations": "Yes.", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 4, "significance": 3, "originality": 4, "strengths_and_weaknesses": "Strengths:\n* Best coordinate selection is an interesting and well-studied problem.\n* A number of theoretical insights in combination are needed for the results. e.g. top-R requires the combination of a sensitivity bound which is non-trivial, a cleverly-designed MIP to isolate the top-R sets without brute-force enumeration, and the authors also give a privacy proof which shows that imposing a ceiling on the loss function in the exponential mechanism improves its utility.\n* The algorithms are tight in that match the non-private best rates of recoery as $\\epsilon \\rightarrow \\infty$.\n* The presentation makes the flow of ideas appear to be natural which makes for an easy read.\n\nWeaknesses:\n* The mistakes algorithm's privacy guarantee relies on a separation between the loss of the best and second-best datasets, i.e. the privacy guarantee doesn't hold for all pairs of datasets. While the authors argue under certain sampling assumptions for the datasets, this separation holds with high probability, DP is inherently a worst-case statement about all pairs of distributions so \"average-case\" DP definitions are generally problematic in nuanced ways (see e.g. https://differentialprivacy.org/average-case-dp/ for more discussion) and such definitions shouldn't be promoted in the literature. So the mistakes algorithm is probably not desirable (except maybe as a theoretical ideal) unless e.g. the separation condition can be forced using clipping-style arguments.\n* The empirical comparison studies the chance of recovering the right set of coordinates, but the final downstream objective is actually to minimize the model loss and not just find the best coordinates. Samp-Agg in particular does minimize the final loss, so a comparison to Samp-Agg is a bit unfair because (i) Samp-Agg is not optimizing for the highest % of correct coordinates, and (ii) Samp-Agg is spending part of its privacy budget on the model optimization, whereas top-R and mistakes do not spend any privacy budget on this aspect of the problem. \n* The combination of the two above points makes the empirical comparisons a bit unconvincing, since (i) mistakes has a problematic privacy guarantee, (ii) top-R does not uniformly beat MCMC, and (iii) top-R beats Samp-Agg, but that comparison is unfair.", "was_revised": true, "final_justification": "The rebuttal has addressed my concerns. It would be nice if mistakes could be made unconditionally private, or if the authors could include suggestions for dealing with the potential for it to not be private, but the good performance of top-K and the fact that the previous work also included such a conditional privacy guarantee makes this a more minor weakness than I had originally thought. The new experiments make the empirical improvement of the authors' method more convincing and addressed my main concern about the paper.", "rebuttal": "Dear reviewer JE2d,\n\nWe sincerely appreciate your effort in providing detailed feedback on our paper. We address your comments below.\n\n**Privacy of mistakes method**\n\nWe would like to note that while the privacy of our mistakes method relies on an additional assumption that occurs with high probability (as shown in section A.5 in the appendix), providing a privacy guarantee with additional assumptions is not uncommon in the literature. For example, the privacy guarantees of [1], which is the closest competitor to our method, depends on assumptions that hold with high probability. More specifically, the privacy proof of the Markov Chain Monte Carlo (MCMC) algorithm in [1] relies on the assumption that the mixing of the Markov Chain used for sampling with its stationary distribution has happened. However, this mixing can only be guaranteed with high probability, and under additional assumptions on the underlying model—see Theorem 4.3 of [1] for more details. In contrast, our top-R method is always private (assuming $b_x,b_y$ are finite), and our mistakes method is private under assumptions that are similar to the ones in Theorem 4.3 of [1].\n\nAs another example, [2] uses the stability of Lasso, to present a DP method for support recovery in sparse linear regression. However, the stability of Lasso only holds under certain assumptions on the data, such as the boundedness of the noise and restricted strong convexity. Such assumptions might only hold with high probability in practice, resulting in privacy guarantees that hold with high probability. For more details, we refer to Theorem 8 of [2].\n\n\n\n**Comparison with Samp-Agg**\n\nWe would first like to clarify that in our experiments, we do not compare with Algorithm 3 in [3], which allocates half the privacy budget on feature selection and half on model optimization, but rather, we are comparing with Algorithm 5 which exclusively is aimed at feature selection for the same problem and is itself (ε,0)-DP. \n\nNevertheless, as you suggested, in an effort to compare prediction accuracy across methods, we did a 70/30 random train test/split and implemented Algorithm 2 in [3] on the train data, wherein we use half the privacy budget (i.e., $\\epsilon/2$) for variable selection using the top-R, mistakes, Samp-Agg, or MCMC methods, and the remaining half of the budget for model optimization using objective perturbation (Algorithm 1 in [3]) to obtain the regression coefficients under the privacy budget ($\\beta_\\text{priv}$). We ran experiments with $p=100, s=5, \\epsilon = 2, SNR= 5, \\rho = 0.1$, and ran 10 independent trials for each value of $n$, where, for each trial, we drew 50 times from the distribution corresponding to each algorithm, and for MCMC, we used 50 independent Markov chains for each trial. After running objective perturbation over the supports, we thus had 500 distinct coefficient vectors for each method and for each $n$, and we computed the average MSE on the test data for every method. Our choices of $\\lambda=120$ and MCMC iterations = 1000 were chosen to make the runtimes comparable, as done in our paper. We summarize our results below, showing that our top-R and mistakes methods outperform the competitor algorithms in terms of prediction accuracy for $n$ large enough.  \n\nAverage test MSE:\n\n|                             | $n=3000$ | $n=4000$ | $n=5000$ | $n=6000$ |\n|-----------------------------|-----------|-----------|-----------|-----------|\n| top-R                       | 1.894    |  0.650    | 0.290    | 0.276    |\n| mistakes                    | 1.064    | 0.408    | 0.300    | 0.275    |\n| Samp-Agg                    | 1.490   | 0.980    | 0.677    | 0.480    |\n| MCMC                        | 0.532    | 0.386    | 0.329   | 0.317   |\n\nFurthermore, we ran experiments to evaluate the utility loss, or gap between the objective at $\\beta_\\text{priv}$ and the objective at $\\hat{\\beta}$ (the optimal BSS solution), with the same parameters as in the experiments above and we present our results below. As with the prediction accuracy, our top-R and mistakes methods outperform the competitor algorithms in terms of utility loss for $n$ large enough.\n\nAverage utility loss:\n \n|                             | $n=3000$ | $n=4000$ | $n=5000$ | $n=6000$ |\n|-----------------------------|-----------|-----------|-----------|-----------|\n| top-R                       | 0.192   |  0.020    | 0.005   | 0.003    |\n| mistakes                    |0.061   | 0.016   | 0.005    | 0.003   |\n| Samp-Agg                    | 0.082   | 0.037    | 0.018    | 0.011    |\n| MCMC                        | 0.017    | 0.013    | 0.010   | 0.008   |\n\nReferences:\n\n[1] Roy, S., Wang, Z., & Tewari, A. (2024). On the computational complexity of private high-dimensional model selection via the exponential mechanism. In Advances in Neural Information Processing Systems (NeurIPS).\n\n[2] Thakurta, A. G., & Smith, A. (2013). Differentially private feature selection via stability arguments, and the robustness of the lasso. In Proceedings of the Conference on Learning Theory (COLT), PMLR.\n\n[3] Kifer, D., Smith, A., & Thakurta, A. (2012). Private convex empirical risk minimization and high-dimensional regression. In Proceedings of the Conference on Learning Theory (COLT), JMLR Workshop and Conference Proceedings."}, {"review_id": "VpVTLXy7ez", "reviewer": "Reviewer_vS6c", "summary": "The paper introduces recent developments in Mixed Integer Programming (MIP) for the Best Subset Selection (BSS) problem to the differential privacy community. It presents two selection procedures: Top-R, which is based on the exponential mechanism and either selects one of the top R subsets or samples a subset uniformly at random; and the \"Mistakes\" method, which applies the exponential mechanism based on the distance (the number of mistakes) from the optimal solution. To run the Top-R and Mistakes mechanisms, one needs to get the top R subsets or the best subsets with a fixed number of mistakes, which the authors do by sequentially solving modified MIP problems. This approach allows private BSS problems to scale up to $10^4$ parameters.", "questions": "You state that Algorithm 1 satisfies $(\\epsilon, 0)$-DP, but according to Theorem 1, it actually satisfies $(\\epsilon', 0)$-DP unless $T = \\infty$. Is this a typo?\n\nAre there similar modifications of the exponential mechanism in the literature? If so, it would be advisable to mention them in the related works section.", "limitations": "yes", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 4, "clarity": 3, "significance": 3, "originality": 4, "strengths_and_weaknesses": "The authors developed an efficient algorithm for differentially private Best Subset Selection by modifying the exponential mechanism to sample not from the full exponential space of all possible subsets, but only from the best ones. To my knowledge, this adaptation of the exponential mechanism is original, though similar ideas may potentially exist. I consider this paper an important contribution: in numerical experiments, the Mistakes and Top-R methods significantly outperform existing alternatives. Furthermore, the paper provides a non-trivial theoretical analysis of the utility guarantees of the proposed methods in terms of support recovery (Theorems 3 and 4). For these reasons, I support the acceptance of the paper.\n\nWeaknesses:\n\n1) The theoretical guarantees for support recovery appear to be weaker than those provided by the MCMC method [22]. \n\n2) The privacy analysis of the Mistakes method (Theorem 2) relies on extra assumptions that may not hold in practice, potentially limiting the method’s applicability.\n\n3) The presentation of the paper could also be improved. The conclusion section is unnecessarily short, and the font sizes in the plots are too small, making them difficult to read. I would also suggest renaming the \"Mistakes\" method, as the current name might be misleading. Furthermore, the formulas in the Appendix need to be reformatted to fit within the page margins.", "was_revised": true, "final_justification": "The authors developed an efficient algorithm for differentially private Best Subset Selection by modifying the exponential mechanism to sample not from the full exponential space of all possible subsets, but only from the best ones. I consider this paper an important contribution: in numerical experiments, the Mistakes and Top-R methods significantly outperform existing alternatives. Furthermore, the paper provides a non-trivial theoretical analysis of the utility guarantees of the proposed methods in terms of support recovery (Theorems 3 and 4). For these reasons, I support the acceptance of the paper.\n\nDuring the rebuttal, the authors addressed my concerns about the theoretical guarantees and the presentation of the paper.\n\nHowever, the extra assumptions required for the Mistakes methods remain a major weakness.", "rebuttal": "Dear reviewer vS6c,\n\nWe would like to sincerely thank you for your positive feedback on our paper. We would like to address your concerns below.\n\n**Theoretical guarantees for support recovery** \n\nWe would like to clarify that our two methods match the best theoretical guarantees in the literature in complementary privacy regimes, as discussed in the contributions section of our paper and in section 3. Namely, our top-R method matches with the non-private minimax-optimal signal threshold for support recovery in the low privacy regime, and our mistakes method matches the required signal threshold in [1] in the high-privacy regime. \n\n**Privacy of mistakes method**\n\nWe would like to note that while the privacy of our mistakes method relies on an additional assumption that occurs with high probability (as shown in section A.5 in the appendix), providing a privacy guarantee with additional assumptions is not uncommon in the literature. For example, the privacy guarantees of [1], which is the closest competitor to our method, depends on assumptions that hold with high probability. More specifically, the privacy proof of the Markov Chain Monte Carlo (MCMC) algorithm in [1] relies on the assumption that the mixing of the Markov Chain used for sampling with its stationary distribution has happened. However, this mixing can only be guaranteed with high probability, and under additional assumptions on the underlying model—see Theorem 4.3 of [1] for more details. In contrast, our top-R method is always private (assuming $b_x,b_y$ are finite), and our mistakes method is private under assumptions that are similar to the ones in Theorem 4.3 of [1].\n\nAs another example, [2] uses the stability of Lasso, to present a DP method for support recovery in sparse linear regression. However, the stability of Lasso only holds under certain assumptions on the data, such as the boundedness of the noise and restricted strong convexity. Such assumptions might only hold with high probability in practice, resulting in privacy guarantees that hold with high probability. For more details, we refer to Theorem 8 of [2].\n\n**Presentation and typos**\n\nThank you very much for correcting our typo in regards to the privacy of our top-R method. As you correctly point out, top-R is $(\\epsilon’,0)$-DP unless $T=\\infty$. Thank you also for noting how our paper presentation could be improved. We will make sure to edit the formulas in the appendix to fit within the page margins, grow the font sizes in the plots, and expand our conclusion section further. \n\n**Similar modifications of the exponential mechanism**\n\nIn regards to your last question, there are no similar modifications to the exponential mechanism in the literature to our knowledge.\n\nReferences:\n\n[1] Roy, S., Wang, Z., & Tewari, A. (2024). On the computational complexity of private high-dimensional model selection via the exponential mechanism. In Advances in Neural Information Processing Systems (NeurIPS).\n\n[2] Thakurta, A. G., & Smith, A. (2013). Differentially private feature selection via stability arguments, and the robustness of the lasso. In Proceedings of the Conference on Learning Theory (COLT), PMLR."}, {"review_id": "zq7vQPx1WC", "reviewer": "Reviewer_PKV9", "summary": "This paper focuses on the sparse variable selection in the context of Differential Privacy (DP), which aims to improve interpretability and generalization in high-dimensional learning while protecting the confidential and personal information. This paper introduces two new differentially private estimators for sparse variable selection based on modern MIP techniques. Their framework is broadly applicable, and they provide theoretical support recovery guarantees in the case of sparse regression. Finally, extensive numerical experiments, using both least squares and hinge loss, demonstrate the effectiveness and scalability of their methods.", "questions": "Please see weakness above.\nA minor comment is: the notation $\\rho$  is used to represent different concepts.", "limitations": "Some results of this paper are limited to the sparse regression problem, as detailed in the weakness. It would be beneficial to extend these results to broader settings and include numerical results for a wider class of loss functions.", "rating": 3, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 2, "clarity": 3, "significance": 2, "originality": 3, "strengths_and_weaknesses": "This paper is clearly written and well-structured. This paper is the first to apply MIP techniques for sparse variable selection in the context of DP. The work addresses an important problem in the field of variable selection.\n\nHowever, the computational results could be improved. Specifically, the authors do not discuss why the Top-R method performs better as the sample size n increases. Throughout the computational experiments, the sparsity parameter s is set to 5 or 8, much smaller compared to the number of features. This may not be representative. It is necessary to try larger values of s to evaluate the methods. \n\nWhile the two DP algorithms are general, the support recovery results are limited to the sparse regression setting. In addition, they only provide MIP formulations and algorithms for least squares and hinge loss objectives. It remains unclear how MIP techniques extend to general loss function.", "was_revised": true, "final_justification": "While I believe this paper brings forth an interesting idea, there exist technical limitations in their methodology. Thus, I maintain my score.", "rebuttal": "Dear reviewer PKV9, \n\nWe sincerely appreciate your feedback on our paper. We would like to address your concerns below.\n\n**Effect of sample size on Top-R performance**\n\nThank you for your comment! The performance improves with increasing $n$ due to the decrease in the threshold needed for the identifiability margin $\\tau$ to have enough separation between the true support and other supports (Theorem 3 in our paper). We will make sure to update this in our manuscript. \n\n**Empirical performance with larger s**\n\nIn order to evaluate the performance of our methods for larger values of $s$, we ran experiments with $p=100, s=20, \\epsilon = 1, SNR= 5,$ and $\\rho = 0.1$. For each value of $n$, we ran 10 independent trials, where, for each trial, we drew 100 times from the distribution corresponding to each algorithm, and for MCMC, we used 100 independent Markov chains for each trial. Our choice of $\\lambda=80,000$ and MCMC iterations = 1000 resulted in our mistakes method having a considerably smaller average runtime than the MCMC method.\n\nWe present our results below. As expected, we observed that for $n$ large enough, our mistakes method outperformed both the Samp-Agg and MCMC methods. \n\nAverage # correct supports / number of trials:\n\n|               | $n=120,000$ | $n=130,000$ | $n=140,000$ | $n=150,000$ |\n|---------------|-----------|-----------|-----------|-----------|\n| mistakes  | 0.31      | 0.49      | 0.88      | 0.97      |\n| Samp-Agg | 0.71      | 0.79      | 0.82      | 0.83      |\n| MCMC      | 0.13      | 0.12      | 0.06      | 0.10      |\n\nAverage F1 scores:\n\n|               | $n=120,000$ | $n=130,000$ | $n=140,000$ | $n=150,000$ |\n|---------------|-----------|-----------|-----------|-----------|\n| mistakes  | 0.9535    | 0.9735    | 0.994     | 0.998     |\n| Samp-Agg  | 0.985     | 0.9895    | 0.99      | 0.992     |\n| MCMC      | 0.936     | 0.936     | 0.923     | 0.94      |\n\n**MIP extension to general loss functions**\n\nA key benefit of our MIP approach is its generalizability across different loss functions. Specifically, the only property that we need to have for the loss function $\\rho$ is that it is a convex function of $\\beta$. As long as the loss function is convex, our outer approximation approach can be applied to any general loss function. For example, our method works with Huber or quantile loss as well, and the MIP algorithm would remain unchanged.\n\n**Support recovery results are limited to the sparse regression setting**\n\nWe acknowledge that our theoretical support recovery guarantees are currently limited to the sparse regression setting. However, extending these results rigorously to more general loss functions would necessitate a fundamentally different analytical framework and set of assumptions than that used in the least squares setting [1]. Such an extension would be substantial enough to warrant dedicated treatment in future work. Additionally, we emphasize that, to our knowledge, the existing theoretical literature on differentially private sparse variable selection explicitly focuses on linear regression. Therefore, restricting our theoretical support recovery analysis to the sparse regression context aligns naturally with the current state of research and represents a meaningful contribution within this established framework.\n\n**Comment regarding notational conflict** \n\nThank you for noticing our notational conflict with the use of $\\rho$. We will make sure to adjust this in our manuscript.\n\nReferences:\n\n[1] Chen, L. Y., & Lee, S. (2023). Sparse quantile regression. Journal of Econometrics, 235(2), 2195–2217."}]

Author rebuttal:

Meta-review: In this paper the authors propose new algorithms for differentially-private variable selection. Their methods are based on the exponential mechanism but incorporate novel modifications to avoid exhaustive calculation. After discussion, all of the reviewers felt positively about this paper. In particular, the problem is interesting and unique, there are compelling theoretical insights (particularly the technique for saving computation in the exponential mechanism), the methods appear to perform well empirically, and the paper is clear and well-written. The authors are encouraged to incorporate the reviewers' comments in their final revision (along with the new experiments they shared in their rebuttals).
Decision: Accept (poster)