Conditional Sparse $\ell_p$-norm Regression With Optimal Probability

We consider the following conditional linear regression problem: the task is to identify both (i) a $k$-DNF condition $c$ and (ii) a linear rule $f$ such that the probability of $c$ is (approximately) at least some given bound $\mu$, and $f$ minimizes the $\ell_p$ loss of predicting the target $z$ in the distribution of examples conditioned on $c$. Thus, the task is to identify a portion of the distribution on which a linear rule can provide a good fit. Algorithms for this task are useful in cases where simple, learnable rules only accurately model portions of the distribution. The prior state-of-the-art for such algorithms could only guarantee finding a condition of probability $\Omega(\mu/n^k)$ when a condition of probability $\mu$ exists, and achieved an $O(n^k)$-approximation to the target loss, where $n$ is the number of Boolean attributes. Here, we give efficient algorithms for solving this task with a condition $c$ that nearly matches the probability of the ideal condition, while also improving the approximation to the target loss. We also give an algorithm for finding a $k$-DNF reference class for prediction at a given query point, that obtains a sparse regression fit that has loss within $O(n^k)$ of optimal among all sparse regression parameters and sufficiently large $k$-DNF reference classes containing the query point.


Introduction
In areas such as advertising, it is common to break a population into segments, on account of a belief that the population as a whole is heterogeneous, and thus better modeled as separate subpopulations. A natural question is how we can identify such subpopulations. A related question arises in personalized medicine. Namely, we need existing cases to apply data-driven methods, but how can we identify cases that should be grouped together in order to obtain more accurate models?
In this work, we consider a problem of this sort recently formalized by Juba 2017: we are given a data set where each example has both a real-valued vector together with a Boolean-valued vector associated with it. Let Y denote the real-valued vector for some example, and let X denote the vector of Boolean attributes. We are also given some target prediction Z of interest, and we wish to identify a k-DNF c such that c( X) specifies a subpopulation on which we can find a linear predictor for Z, a, Y , that achieves small loss. In particular, it may be that neither a, nor any other linear predictor achieves small loss over the entire data distribution, and we simply wish to find the subset of the distribution, determined by our k-DNF c evaluated on the Boolean attributes of our data, on which accurate regression is possible. We refer to this subset of the distribution as a "segment". Juba proposed two algorithms for different families of loss functions, one for the sup norm (that only applies to O(1)-sparse regression) and one for the usual ℓ p norms. But, the algorithm for the ℓ p norm has a major weakness, namely that it can only promise to identify a condition c ′ ( X) that has a probability polynomially smaller than optimal. This is especially problematic if the segment had relatively small probability to begin with: it may be that there is no such condition with enough examples in the data for adequate generalization. By contrast, the algorithm for the sup norm recovers a condition c ′ with probability essentially the same as the optimal condition, but of course the sup norm is an undesirable norm to use for regression, as it is maximally sensitive to outliers. In general, ℓ p regression is more tolerant to outliers as p decreases, and the sup norm is roughly "p = ∞." In this work, in Theorem 10, we show how to give an algorithm for O(1)-sparse ℓ p regression with a two-fold improvement over the previous ℓ p regression algorithm (mentioned above) when the coefficient vector is sparse: 1. The new algorithm recovers a condition c ′ with probability essentially matching the optimal condition. 2. We also obtain a smaller blow-up of the loss relative to the optimal regression parameters a and optimal condition c, reducing the degree of this polynomial factor. More concretely, Juba's previous algorithm identifies k-DNF conditions, as does ours. This algorithm only identified a condition of probability Ω(µ/n k ) when a condition with probability µ exists, and only achieved loss bounded by O(ǫn k ) when it is possible to achieve loss ǫ. We improve the probability of the condition to (1 − η)µ for any desired η > 0, and also reduce the bound on the loss toÕ(ǫn k/2 ), or further toÕ(ǫg log log(n k )) if the condition has only g "terms" (ANDs of literals), in Corollary 12. This latter algorithm furthermore features a smaller sample complexity, that only depends logarithmically on the number n of Boolean attributes. We include a synthetic data experiment demonstrating that the latter algorithm can successfully recover more terms of a planted solution than Juba's sup norm algorithm.
We also give an algorithm for the closely related reference class regression problem in Theorem 14. In this problem, we are given a query point x * (along with values for the real attributes y * ), and we wish to estimate the corresponding z * . In order to compute this estimate for a single point, we find some k-DNF c ′ such that c ′ ( x * ) = 1 so x * is in the support of the distribution conditioned on c ′ ( X), and such that there are (sparse) regression parameters a ′ for c ′ that (nearly) match the optimal loss for sparse regression under any condition c with c( x * ) = 1. In this case, the condition c ′ is a reference class containing x * , obtaining the tightest possible fit among k-DNF classes containing x * . We think of this as describing a collection of "similar cases" for use in estimating the target value z * from the attributes y * . For this problem, we guarantee that we find a condition c ′ and sparse regression parameters a ′ that achieve loss within O(n k ) of optimal.
Our algorithms are based on sparsifiers for linear systems, for example as obtained for the ℓ 2 norm by Batson et al. 2012, and as obtained for non-Euclidean ℓ p norms by Cohen and Peng 2015 based on Lewis weights [Lewis, 1978]. The strategy is similar to Juba's sup norm algorithm: there, by enumerating the vertices of the polytopes obtained on various subsets of the data, he is able to obtain a list of all possible candidates for the optimal estimates of the sup-norm regression coefficients. Here, we can similarly obtain such a list of estimates for the ℓ p norms by enumerating the possible sparsified linear systems. Sparsifiers are often used to accelerate algorithms, e.g., to run in time in terms of the number of nonzero entries rather than the size of the overall input, though here we use that the size of a sparsified representation is small in a rather different way. Namely, the small representation allows for a feasible enumeration of possible candidates for the regression coefficients.
We note that this means that we also obtain an algorithm for O(1)-sparse ℓ p -norm regression in the list-learning model [Balcan et al., 2008, Charikar et al., 2017. Although Charikar et al. in particular are able to solve a large family of problems in the list-learning model, they observe that their technique can only obtain trivial estimates for regression problems. We show that once we have such a list, algorithms for the conditional distribution search problem [Juba, 2016, Zhang et al., 2017, Juba et al., 2018 can be generically used to extract a description of a good event c for the given candidate linear rule. Similar algorithms can also then be used to find a good reference class.

Relationship to other work
The conditional linear regression problem is most closely related to the problem of prediction with a "reject" option. In this model, a prediction algorithm has the option to abstain from making a prediction. In particular, El-Yaniv and Weiner 2012 considered linear regression in such a model. However, like most of the work in this area, their approach is based on scoring the "confidence" of a prediction function, and only making a prediction when the confidence is sufficiently high. This does not necessarily yield a nice description of the region in which the predictor will actually make predictions. On the other hand, Cortes et al. 2016 consider an alternative variant of learning with rejection that does obtain such nice descriptions; but, in their model, they assume that abstention comes at a fixed, known cost, and they simply optimize the overall loss when their space of possible prediction values includes this fixed-cost "abstain" option. By contrast, we can explicitly consider various rates of coverage or loss.
Our work is also related to the field of robust statistics [Huber, 1981, Rousseeuw andLeroy, 1987]: in this area, one seeks to mitigate the effect of outliers. That is, roughly, we suppose that a small fraction of the data has been corrupted, and we wish for our models and inferences to be relatively unaffected by this corrupted data. The difference between robust statistics and the conditional regression problem is that in conditional regression we are willing to ignore most of the data. We only wish to find a small fraction, described by some rule, on which we can obtain a good regression fit.
Another work that similarly focuses on learning for small subsets of the data is the the work by Charikar et al. 2017 in the list-learning model of Balcan et al. 2008. In that work, one seeks to find a setting of some parameters that (nearly) minimizes a given loss function on an unknown small subset of the data. Since it is generally impossible to produce a single set of parameters to solve this task -many different small subsets could have different choices of optimal parameters -the objective is to produce a small list of parameters that contains a near-optimal solution for the unknown subset somewhere in it. As we mentioned above, our approach actually gives a list-learning algorithm for O(1)-sparse ℓ p regression; but, we show furthermore how to identify a conditional distribution such that some regression fit in the list obtains a good fit, thus also solving the conditional linear regression task. Indeed, in general, we find that algorithms for solving the list-learning task for regression yield algorithms for conditional linear regression. The actual technique of Charikar et al. does not solve our problem since their approximation guarantee essentially always admits the zero vector as a valid approximation. Thus, they also do not obtain an algorithm for list-learning (sparse) regression either.
Another task similar to list-learning of linear regression is the problem of finding dense linear relationships, as solved by RANSAC [Fischler and Bolles, 1981]. But, these techniques only work in constant dimension. By contrast, although we are seeking sparse regression rules, this is a sparse fit in high dimension. As with list-learning, in contrast to the conditional regression problem, such algorithms do not provide a description of the points fitting the dense linear relationship.
Finally, yet another task similar to list-learning for regression is fitting linear mixed models [McCulloch andSearle, 2001, Jiang, 2007]. In this approach, one seeks to explain all (or almost all) of the data as a mixture of several linear rules. The guarantee here is incomparable to ours: in contrast to list-learning, the linear mixed model simply needs a list of linear rules that accounts for nearly all of the data; it does not need to find a list that accounts for all possible sufficiently large subsets of the data. So, there is no guarantee that any of the mixture components represent an approximation to the regression parameters corresponding to an event of interest. On the other hand, in linear mixed models, one does not need to give any description at all of which points should lie in which of the mixture components. In applications, one usually assigns points to the linear rule that gives it the smallest residual, but this may be less useful for predicting the values for new points.

Preliminaries
We now formally define the problems we consider in this work, and recall the relevant background.

The conditional regression and distribution search tasks
Formally, we focus on the following task: Definition 1 (Conditional ℓ p -norm Regression) Conditional ℓ p -norm linear regression is the following task. Suppose we are given access to i.i.d. examples drawn from a joint distribution over If a * has at most s nonzero components, and we requireˆ a likewise has at most s nonzero components, then this is the conditional s-sparse ℓ p -norm regression task.
As stressed by Juba 2017, the restriction to k-DNF conditions is not arbitrary. If we could identify conditions that capture arbitrary conjunctions, even for one-dimensional regression, this would yield PAC-learning algorithms for general DNFs. In addition to this being an unexpected breakthrough, recent work by Daniely and Shalev-Shwartz 2016 shows that such algorithms would imply new algorithms for random k-SAT, falsifying a slight strengthening of Feige's hypothesis [Feige, 2002]. We thus regard it as unlikely that any algorithm can hope to find conditions of this kind. Of the classes of Boolean functions that do not contain arbitrary conjunctions, k-DNFs are the most natural large class, and hence are the focus for this model.
The difficulty of the problem lies in the fact that initially we are given neither the condition nor the linear predictor. Naturally, if we are told what the relevant subset is, we can just use standard methods for linear regression to obtain the linear predictor; conversely, if we are given the linear predictor, then we can use algorithms for the conditional distribution search (learning abduction) task (introduced by Juba 2016, recalled below) to identify a condition on which the linear predictor has small error. Our final algorithm actually uses this connection, by considering a list of candidates for the linear predictors to use for labeling the data, and choosing the linear predictor that yields a condition that selects a large subset. In particular, we will use algorithms for the following weighted variant of the conditional distribution search task: Definition 2 (Conditional Distribution Search) Weighted conditional distribution search is the following problem. Suppose we are given access to examples drawn i.i.d. from a distribution over This task is closely related to agnostic conditional distribution search, which is the special case where the weights only take values 0 or 1. The current state of the art for agnostic conditional distribution search is an algorithm given by Zhang et al. 2017, achieving anÕ( √ n k )-approximation to the optimal error. That work built on an earlier algorithm due to Peleg 2007. Peleg's work already showed how to extend his original algorithm to a weighted variant of the problem, and we observe that an analogous modification of the algorithm used by Zhang et al. will obtain an algorithm for the more general weighted conditional distribution search problem we are considering here: Theorem 3 (Peleg 2007, Zhang et al. 2017) There is a polynomial-time algorithm for weighted conditional distribution search achieving anÕ(n k/2 (log b+ log 1/η + log log 1/δ))-approximation with probability 1 − δ given m = Θ b 3 µǫη 2 (n k + log 1 δ ) examples.
We note that it is possible to obtain much stronger guarantees when we are seeking a small formula. Juba et al. 2018 present an algorithm that, when there is a k-DNF with g terms achieving error ǫ, uses only m =Õ( gk µǫη 2 log n δ ) examples and obtains a k-DNF with probability at least (1−η)µ and errorÕ(ǫg log m). Although this is stated for the unweighted case (i.e., ǫ is a probability), it is easy to verify that since our loss is nonnegative and bounded, by rescaling the losses to lie in [0, 1], we can obtain an analogous guarantee for the weighted case: Theorem 4 (Juba et al. 2018) There is a polynomial-time algorithm for weighted conditional distribution search when the condition has g terms using m =Õ( bgk µǫη 2 log n δ ) examples achieving añ O(g log m)-approximation using a k-DNF withÕ(g log m) terms with probability 1 − δ.

Reference class regression
Using a similar approach, we will also solve a related problem, selecting a best k-DNF "reference class" for regression. In this task, we are not merely seeking some k-DNF event of probability µ on which the conditional loss is small. Rather, we are given some specific observed Boolean attribute values x * , and we wish to find a k-DNF condition c that is satisfied by x * solving the previous task. That is, c should have probability at least µ and our sparse regression fit has small conditional loss, conditioned on c. Naturally, the motivation here is that we have some specific point ( x * , y * ) for which we are seeking to predict z * , and so we are looking for a "reference class" c such that we can get the tightest possible regression estimate of z * from y * ; to do so, we need to take µ large enough that we have enough data to get a high-confidence estimate, and we need x * to lie in the support of the conditional distribution for which we are computing this estimate.
Definition 5 (Reference Class ℓ p -Regression) Reference class ℓ p -norm regression is the following task. We are given a query point x * ∈ {0, 1} n , target density µ ∈ (0, 1), ideal loss bound ǫ 0 > 0 approximation parameter η > 0, confidence parameter δ > 0, and access to i.i.d. examples drawn from a joint distribution over We wish to findˆ a ∈ R d with ˆ a 2 ≤ b and a reference class k-DNFĉ such that with probability and k-DNFs c * such that c * ( x * ) = 1 and Pr[c * ( X)] ≥ µ. If we also require bothˆ a and a * to have at most s nonzero components, then this is the reference class s-sparse ℓ p -norm regression task.
The selection and use of such reference classes for estimation goes back to work by Reichenbach 1949. Various refinements of this approach were proposed by Kyburg 1974 andPollock 1990, e.g., to choose the estimate provided by the highest-accuracy reference class that is consistent with the most specific reference class containing the point of interest x * . Our approach is not compatible with these proposals, as they essentially disallow the use of the kind of disjunctive classes that are our exclusive focus. Along the lines we noted earlier, it is unlikely that there exist efficient algorithms for selecting reference classes that capture arbitrary conjunctions, so k-DNFs are essentially the most expressive class for which we can hope to solve this task. Bacchus et al. 1996 give a nice discussion of other unintended shortcomings of disallowing disjunctions. A concrete example discussed by Bacchus et al. is the genetic disease Tay-Sachs. Tay-Sachs only occurs in two very specific, distinct populations: Eastern European Jews and French Canadians. Thus, a study of Tay-Sachs should consider a reference class at least partially defined by a disjunction over membership in these two populations.

ℓ p sparsification
Our approach is based on techniques for extracting low-dimensional sketches of small subspaces in high dimensions. The usual ℓ 2 norm uses much simpler underlying techniques, and we describe it first. The extension to ℓ p norms for p = 2 is obtained via Lewis weights [Lewis, 1978].

Euclidean sparsifiers
The kind of sketches we need originate in the work of Batson et al. 2012. Specifically, it will be convenient to start from the following variant due to Boutsidis et al. 2014: Lemma 6 (BSS weights [Boutsidis et al., 2014] where λ j denotes the jth largest eigenvalue. In particular, taking r = t/γ 2 for some γ ∈ (0, 1), we obtain that for the Proof: Observe that since the each ith row of [u] has unit norm, it must have an entry u i,j * that is at least 1/ √ d in magnitude. By the above argument, where notice in particular, the ith row of [s] ⊤ contributes at least ( √ s i u i,j * ) 2 ≥ s i /d to the norm.

Sparsifiers for non-Euclidean norms
It is possible to obtain an analogue of the BSS weights for p = 2 using techniques based on Lewis weights [Lewis, 1978]. Lewis weights are a general way to reduce problems involving ℓ p norms to analogous ℓ 2 computations. Cohen and Peng 2015 applied this to sparsification to obtain the following family of sparsifiers:

Theorem 8 (ℓ p weights [Cohen and Peng, 2015]) Given a d × t matrix [u]
there exists a set of r(p, t, γ) weights s 1 , . . . , s r such that for the d × r matrix [s] which has as its ith column s i e i , where r(p, t, γ) is asymptotically bounded as in Table 1. p Required dimension r p = 1 t log t γ 2 1 < p < 2 1 γ 2 t log(t/γ) log 2 log(t/γ) p = 2 t/γ 2 p > 2 log 1/γ γ 5 t p/2 log t Cohen and Peng also show how to construct the sparsifiers for a given matrix efficiently, but we won't be able to make use of this, since we will be searching for the sparsifier for an unknown subset of the rows.
We furthermore obtain an analogue of Lemma 7 for the ℓ p weights, using essentially the same argument: Lemma 9 Suppose the rows of [u] are orthonormal. Then the matrix [s] obtained by Theorem 8 has entries of magnitude at most (1 + γ) √ d.

The weighting algorithms for conditional and reference class regression
The results from sparsification of linear systems tell us that we can estimate the loss on the subset of the data selected by the unknown condition by computing the loss on an appropriate "sketch," a weighted average of the losses on a small subset of the data. The dimensions in Table 1 give the size of these sparsified systems, i.e., the number of examples from the unknown subset we need to use to estimate the loss over the whole subset. The key point is that since the predictors are sparse, the dimension is small; and, since we are willing to accept a constant-factor approximation to the loss, a small number of points suffice. Therefore, it is feasible to enumerate these small tuples of points to obtain a list of candidate sets of points for use in the sketch. We also need to enumerate the weights for these points, but since we have also argued that the weights are bounded and we are (again) willing to tolerate a constant-factor approximation to the overall expression, there is also a small list of possible approximate weights for the points. The collection of points, together with the weights, gives a candidate for what might be an appropriate sketch for the empirical loss on the unknown subset. We can use each such candidate approximation for the loss to recover a candidate for the linear predictor. Thus, we obtain a list of candidate linear predictors that we can use to label our data as described above. More precisely, our algorithm is as shown in Algorithm 1. In the following, let Π d 1 ,...,ds denote the projection to coordinates d 1 , d 2 , . . . , d s .
In our proof of this theorem, we will find it convenient to use the Rademacher generalization bounds for linear predictors (note that x → |x| p is pb p−1 -Lipschitz on [−b, b]): Theorem 11 ( Bartlett andMendelson 2002, Kakade et al. 2009 b], and any δ ∈ (0, 1), let L p ( a) denote E[| a, Y − Z| p ], and for an an i.i.d. sample of size m letL p ( a) be the empirical loss 1 m m j=1 | a, y (j) − z (j) | p . We then have that with probability 1 − δ for all a with a 2 ≤ b, Note that although this bound is stated in terms of the ℓ 2 norm of the attribute and parameter vectors y and a, we can obtain a bound in terms of the dimension s of the sparse rule if we are given a bound B on the magnitude of the entries: b ≤ √ sB.
Proof of Theorem 10: Given that we are directly checking the empirical ℓ p loss before returning a and c, for the quoted number of examples m it is immediate by a union bound over the iterations that any a and c we return are satisfactory with probability 1 − δ. All that needs to be shown is that the algorithm will find a pair that passes this final check.
By Theorem 11, we note that it suffices to have b 2p γ 2p ǫ 2p (2pb + 2 ln(12/δ)) 2 examples from the distribution conditioned on the unknown k-DNF event c * to obtain that the ℓ p loss of each candidate for a is estimated to within an additive γǫ with probability 1 − δ/3. By Hoeffding's inequality, therefore when we draw m 0 examples, there is a sufficiently large subset satisfying c * with probability 1 − δ/3.
We let [u] be an orthonormal basis for span{(Π d 1 ,...,ds y (j) , z (j) ) : c * ( x (j) ) = 1, j ≤ m 0 } and invoke Lemma 6 for ℓ 2 or Theorem 8 for p = 2. In either case, there is some set of weights s 1 , . . . , s r 0 for a subset of r 0 coordinates j 1 , . . . , j r 0 such that for any v in the column span of [u], [s] ⊤ [u] v has ℓ p norm that is a 1 ± γ-approximation to the ℓ p norm of [u] v. In particular, for any a, observing v = [y] a − z is in the column span of [u] by construction, we obtain Now, we observe that we can discard weights (and dimensions) from s 1 , . . . , s r 0 of magnitude smaller than γ/r 1/p 0 , since for any unit vector v, the contribution of such entries to [s] ⊤ [u] v p p (recalling there are at most r 0 nonzero entries) is at most γ p . So we may assume the r ≤ r 0 remaining weights all have magnitude at least γ/r 1/p . Furthermore, if we round each weight to the nearest power of (1 + γ), this only changes [s] ⊤ [u] v| p p by an additional (1 ± γ) factor. Finally, we note that since (Π d 1 ,...,ds y (j) , z (j) ) has dimension s+1, Lemmas 7 and 9 guarantee that the magnitude is also at most (1 + γ) √ s + 1. Thus it indeed suffices to find the powers (q 1 , . . . , q r ) for our r examples j 1 , . . . , j r such that (1 + γ) q ℓ is within (1 + γ) of s ℓ , and the resulting set of weights will approximate the ℓ p -norm of every v in the column span to within a 1 + 3γ-factor. Now, when the loop in Algorithm 1 considers (i) the dimensions d * 1 , . . . , d * s contained in the optimal s-sparse regression rule a * (ii) the set of examples j * 1 , . . . , j * r used for the sparse approximation for these coordinates and (iii) the appropriate weights (1 + γ) q * 1 , . . . , (1 + γ) q * r , the algorithm will obtain a vector a that achieves a (1 + 3γ)-approximation to the empirical ℓ p -loss of a * on the same s coordinates. It then follows from Theorem 3 that with probability at least 1 − δ/3 over the data, WtCond will in turn return to us a k-DNF c with probability (1 − η)µ that selects a subset of the data on which a achieves an αǫ =Õ((1 + γ) √ n k (log b + log 1/η + log log 1/δ)ǫ) approximation to the empirical ℓ p loss of a * on c * . This choice of a and c passes the final check and is thus sufficient.
By simply plugging in the algorithm from Theorem 4 for WtCond, we can obtain the following improvement when the desired k-DNF condition is small: Corollary 12 For any constant s and γ > 0, r as given in Table 1 for t = (s + 1), and m = Θ (1+γ)bgk µǫη 2 log n δ + s log d + r log m 0 log(γ 1/p s/r) γ + m 0 examples, if there is a g-term k-DNF solution to the conditional s-sparse ℓ p -regression task, then the modified Algorithm 1 runs in polynomial time and solves the task with aÕ(g log m)-term k-DNF with α =Õ((1 + γ)ǫg log m).
Note that this guarantee is particularly strong in the case where the k-DNF would be small enough to be reasonably interpretable by a human user.
The extension to reference class ℓ p -norm regression proceeds by replacing the weighted condition search algorithm with a variant of the tolerant elimination algorithm from Juba 2016, given in Algorithm 2.
Proof: For convenience, let N ≤ (1 + log 1+η b/ǫ 0 )(1 + log 1+η 1/µ 0 ) denote the total number of iterations. Consider first what happens when the loop considers the largest µ ≤ Pr[c * ( X)]/(1 + η) and the smallest ǫ that is at least (1 + η) 2 ǫ * . On this iteration, for each term T of c * , we observe that So, since W is bounded by b, by a Chernoff bound, 1 . Since this is in turn at most ǫµ, T will be included in c on this iteration.
But similarly, for T not in c with E[W · T ( X)] > (1 + η)µǫ, the Chernoff bound also yields that j:T (x (j) )=1 w (j) ≥ ǫµm with probability 1 − δ/2 n ≤k N . By a union bound over all T ∈ c * and T not in c * with such large error, we see that with probability at least 1 − δ/2N , all of the terms of c * are included in c and only terms with E[W · T ( X)] ≤ (1 + η)µǫ are included in c. So, Furthermore, by yet another application of a Chernoff bound, c * is true of at least µm examples with probability at least 1 − δ/2N . Thus, with probability 1 − δ/N , after this iterationĉ is set to some k-DNF andǫ ≤ (1 + η) 2 max{ǫ * , ǫ 0 }. Now, furthermore, on every iteration, we see more generally that with probability 1 − δ/N , only . Thus, with probability 1 − δ overall, since we found above thatǫ ≤ (1 + η) 2 max{ǫ * , ǫ 0 }, we return a k-DNFĉ as claimed.

Experimental evaluation
To evaluate our algorithm's performance in practice, we two kinds of experiments: one using synthetic data with a planted solution, and another using some standard benchmark data sets from the LIBSVM repository [Chang and Lin, 2011].
Synthetic data. To generate the synthetic data sets, we first chose a random 2-DNF over Boolean attributes by sampling terms uniformly at random. We also fixed two out of d coordinates at random and sampled parameters a from a mean 0, variance σ 2 Gaussian for our target regression fit. For the actual data, we then sampled Boolean attributes X that, with 25% probability are a uniformly random satisfying assignment to the DNF and with 75% probability are a uniformly random falsifying assignment. We sampled Y from a standard Gaussian in R d . Finally, for the examples where X was a satisfying assignment, we let Z be given by the linear rule a, Y + ν where ν is a mean 0, variance σ 2 Gaussian; otherwise, if X was a falsifying assignment, Z was set to ν * where ν * is a standard mean 0, variance 1 Gaussian. (Parameters in Table 2.) Thus, conditioned on the planted formula, there is a small-error regression fit a-the expected error of a is |ν|. Off of this planted formula, the expected error for the optimal prediction z = 0 is |ν * | = 1. For these simple data sets, we made several modifications to simplify and accelerate the algorithm. First, we simply fixed all of the weights (of the form (1 + γ) q i ) to 1, i.e., q i = 0 for all i, since the leverage scores for Gaussian data (and indeed, most natural data sets) are usually small. Second, we set m 0 = 200 for the 1000-example data set and m 0 = 500 for the 5000-example data set. Third, we used the version of the algorithm for small DNFs considered in Corollary 12, that uses algorithm of Juba et al. 2018 (recalled in Theorem 4) for WtCond.
With these modifications, we ran Algorithm 1 for ℓ 2 regression with s = 2, γ = 1, and µ = 0.22 on the 1000-example data set. The algorithm was able to identify the planted 2-DNF together with the two components used in the regression rule in this case, and thus also we approximately identify the linear predictor. The mean squared error of the approximate predictor returned by our algorithm was 0.0149 on the planted 2-DNF (c.f. the ideal linear fit has variance-0.01 Gaussian noise added to it). However, we note that once the planted 2-DNF has been identified, a good regression fit can be found by any number of methods. We also ran Algorithm 1 for ℓ 2 regression with s = 2, γ = 1, and µ = 0.2465 on the 5000-example data set. The algorithm was able to identify a subset of 11 out of the 16 planted k-DNF terms, picking up about 98% of the planted k-DNF data. The mean squared error of the approximate predictor was 0.4221.
As a baseline, we also used the algorithm from Juba 2017 for conditional sparse sup-norm regression on our synthetic data. In order to obtain reasonable performance, we needed to make some modifications: we modified the algorithm along the lines of Algorithm 1, to only search over tuples of a random m 0 = 200 examples for the 1000-example data set and m 0 = 500 for the 5000-example data set in producing its candidate regression parameters. We also chose to take the parameters that achieved the smallest residuals under the sup norm, among those that obtained a condition satisfied by the desired fraction of the data. Again, we used s = 2, and for the norm bound we used a threshold of ǫ = 0.24 for the 1000-example data set and 1.1 for the 5000-example data set. (We considered a range of possible thresholds; when ǫ ≈ 1.5 on the σ 2 = 0.1 data, the sup norm algorithm starts adding terms outside the planted DNF, significantly increasing its error.) For the 1000-example data set, this algorithm was able to identify the planted 2-DNF. However, for the 5000-example data set, it only found 10 out 16 terms of the planted 2-DNF, picking up about 96% of the planted k-DNF data.
Real-world benchmarks. We also compared Algorithm 1 against selective regressors [El-Yaniv and Wiener, 2012] on some of the LIBSVM [Chang and Lin, 2011] regression data sets that were used to evaluate that work. We split each of those datasets into a training set (1/3 of the data) and test set (2/3 of the data). We generated Boolean attributes by choosing different binary splits (median, quartile) on the numerical features, excluding the target attribute. We then run our algorithm 1 on the training set to obtain a 2-DNF. Next, we filtered both the training and test data satisfying the 2-DNF, using the first one to train a new regression fit on the selected subset and the other as a holdout to estimate the error of this resulting regression fit. We compared this to the test error for selective regression. Risk-Coverage (RC) curves of the results are shown in Figure 1. We outperform the baseline on the Boston housing dataset, and we generally achieved lower error than the baseline for lower coverage (less than 0.5) and somewhat higher error for higher coverage for the other three data sets (Body fat, Space, and Cpusmall). Note that in addition, we obtained a 2-DNF that describes the subpopulation on all four, which is the main advantage of our method. Since selective regression first tries to fit the entire data set, and then chooses a subset where that fixed predictor does well, it is to be expected that it may miss a small subset where a different predictor can do better, but that its freedom to abstain on a somewhat arbitrary subset may give it an advantage at high coverage.

Directions for future work
There are several natural open problems. First, although sparsity is desirable, our exponential dependence of the running time (or list size) on the sparsity is problematic. (The sup norm regression algorithm [Juba, 2017] also suffered this deficiency.) Is it possible to avoid this? Second, our algorithm for reference class regression has an O(n k ) blow-up of the loss, as compared toÕ(n k/2 ) for conditional regression. Can we achieve a similar approximation factor for reference class regression? Finally, we still do not know how close to optimal this blow-up of the loss is; in particular, we do not have any lower bounds. Note that this is a computational and not a statistical issue, since we can obtain uniform convergence over all k-DNF conditions.