Information-Based Optimal Subdata Selection for Big Data Linear Regression

• IBOSS deterministically selects the most informative observations, preserving substantial full-data information in a small subdata set without inducing sampling error.The framework borrows the information-maximization motivation of optimal experimental design for data analysis.Source: §1 · p. 3
• IBOSS requires O(np) computation, supports distributed parallel computing, and is significantly faster than existing subdata methods.The algorithm identifies informative points by examining covariates individually.Source: §1 · p. 3
• With fixed subdata size k, IBOSS estimator variances can converge to 0 as full data size n increases, unlike subsampling bounds typically of order 1/k.Subsampling-based estimators also incur sampling error, and their variance need not vanish as n grows.Source: §1 · p. 3
• IBOSS estimators inherit distributional properties of estimators based on the full data, making their sampling distributions easier to assess.The approach is designed for the n ≫ p regime, where ordinary least squares requires O(np^2) time.Source: §1 · pp. 2–3
• The paper develops the IBOSS framework, D-optimal subdata characterization, asymptotic theory, and evaluations through simulations and two real-data applications.It also derives a lower covariance bound for subsampling-based estimators and reports additional numerical results in supplementary material.Source: §1 · pp. 3–4

2 The framework

• For the full-data linear model, the least-squares estimator is BLUE, with covariance equal to the inverse full-data information matrix Mf.The model assumes uncorrelated responses conditional on covariates and errors with E(εi) = 0 and V(εi) = σ2; normality is not required.Source: §2 · p. 4
• The framework selects k subdata points deterministically using covariates alone, then maximizes an optimality criterion of the subdata information matrix M(δ).The criterion includes D-optimality, which maximizes det(M(δ)) and minimizes the expected volume of the joint confidence ellipsoid for β.Source: §2 · pp. 5–6

2.1 Analysis of existing subsampling-based methods

• The information-based framework can assess subsampling without replacement, whose fixed-size properties are more complicated because sampled observations are dependent given the full data.The framework also covers popular UNI, LEV, and LEVUNW methods through their induced information matrices.Source: §2.1 · pp. 6–8
• Common random subsampling estimators have covariance matrices bounded below, in Loewner ordering, by finite quantities of order 1/k that do not approach 0 as n →∞.The bound applies regardless of the weighting choice and whether sampling is with or without replacement, subject to the theorem’s conditions.Source: §2.1 · pp. 7, 9

The D-optimality criterion

• Theorem 2 characterizes D-optimal subdata, with equality when the 2p selected points use each covariate’s minimum and maximum values equally often.When k is much smaller than 2p, exact equality may be unattainable, motivating approximate selection using both small and large covariate values.Source: §3 · pp. 9–10

IBOSS algorithm

• Algorithm 1 selects r smallest and r largest values for each covariate, excluding previously selected points, to construct an IBOSS subdata set.If r = k/(2p) is not an integer, k can be adjusted using floor or ceiling choices, or their combination.Source: §3 · p. 10

Computational properties

• For n > kp, Algorithm 1 obtains the subdata in O(np) time and has full complexity O(np + kp^2 + p^3 + kp) = O(np + kp^2), reducing to O(np).Source: §3 · p. 10

Statistical inference

• Algorithm 1 yields the exact covariance matrix of the estimator, and under normality with known error variance, D has an exact normal distribution.This makes the resulting covariance information directly useful for statistical inference.Source: §3 · p. 10

Implementation and robustness

• The algorithm supports distributed parallel processing by independently finding covariate extremes, while duplicate indexes may reduce subdata size without changing the convergence rate.Selecting extreme covariate values can include outliers, but the subdata follows the same underlying regression model and permits outlier diagnostics.Source: §3 · pp. 10–11

Approximation quality

• Although Algorithm 1 may miss Theorem 2’s unattainable upper bound, under fixed p and reasonable conditions it can achieve the same order.This result links the algorithm’s determinant quality to the theoretical D-optimal benchmark.Source: §3 · p. 11

• The asymptotic results require covariate distributions in the generalized extreme-value domain and permit fixed r or growing r under specified endpoint and rate conditions.For normal covariates, one condition is log r/log log n →0; for exponential covariates, it is log r/log n →0.Source: §4 · pp. 12–13
• If p/(z(n−r+1)j −z(r)j)2 = oP(1), the slope-variance convergence is faster than the 1/k rate typical of subsampling methods.The bound applies for arbitrary n, r, and p, although the actual convergence may be faster than the stated upper bound.Source: §4 · p. 13
• With fixed subdata size k, the IBOSS slope-estimator variance is governed by the inverse squared full-data covariate range and converges to 0 when that range diverges.The corresponding intercept variance remains bounded below by a term proportional to 1/k.Source: §4 · pp. 12, 14
• Using full-data means to adjust the intercept estimator improves β0 efficiency, with a convergence rate similar to slopes when E(z)≠0 and potentially faster when E(z)=0.The additional computation is O(np).Source: §4 · p. 14
• For special covariance structures, intercept variance remains proportional to 1/k, whereas slope-covariance components converge to 0 at distribution-dependent rates, including 1/log n for normal covariates.For lognormal covariates, the reported component rate is exponential in the stated result.Source: §4 · pp. 14–15
• For tν covariates, D-OPT slope variance converges at n−2/ν; when ν ≤2, this is close to the full-data rate and faster than UNI.Neither D-OPT nor UNI makes the intercept variance converge to 0, while UNI’s slope convergence is much slower.Source: §4 · p. 15

5.1 Simulation studies

• For intercept estimation, D-OPT IBOSS generally outperforms other subdata methods, but its MSE does not decrease with increasing n.In Cases 1 and 3, the intercept estimator performs as well as the full-data approach because E(z) = 0, yielding faster convergence than slope estimators.Source: §5.1 · p. 18
• Increasing k improves every subdata method, while D-OPT IBOSS remains best; at n = 10^6, k = 200 is about 10 times as accurate as LEV with k = 5 × 10^3.The comparison uses MSE for slope estimation under the mixture covariate distribution.Source: §5.1 · pp. 18–19
• IBOSS confidence intervals have empirical coverage close to the nominal 0.95 level, supporting valid inference from IBOSS subdata.The study compares normal and mixture covariate distributions with the full-data analysis, but not subsampling methods because theoretically justified intervals were unavailable.Source: §5.1 · p. 19
• D-OPT IBOSS compares favorably with LEV computationally, and both are more efficient than using the full data for fixed k = 10^3.CPU times are reported across full-data sizes and covariate dimensions, with other cases omitted because results were similar.Source: §5.1 · p. 19
• D-OPT IBOSS uniformly dominates subsampling-based methods, including UNI and LEV, with greater advantage for heavier-tailed covariates; its MSE decreases as n increases despite fixed subdata size k = 103.Source: §5.1 · p. 17
• With mixture covariates, D-OPT IBOSS using k = 10^3 from n = 10^6 has 2.4 times lower MSE than full-data analysis using n = 10^5.This performance is comparable to, and in this comparison better than, the full-data approach at the smaller full-data size.Source: §5.1 · p. 18, Figure 1 · p. 20
• Simulation studies examine intercept and slope MSEs across five covariate distributions, varying n with fixed k = 1000 and varying k with fixed n = 10^6.They also evaluate empirical coverage probabilities and average lengths of 95% confidence intervals for D-OPT IBOSS and full-data methods.Source: Figure 2 · p. 21, §5.1 · pp. 21–23, Figure 3 · p. 22

5.2 Real data

• The real-data evaluation uses a CSFII dietary-intake dataset and a chemical-sensor dataset to compare IBOSS with full-data and subsampling approaches.The CSFII analysis uses k = 10p = 50, while the chemical-sensor analysis uses k = 280.Source: §5.2 · p. 24, §5.2.1 · p. 24, §5.2.2 · p. 27

5.2.1 Example 1: food intakes data

• D-OPT IBOSS slope estimates are not very different from full-data estimates, with consistent signs, while their standard errors are larger but reasonably good given the small subdata size.Source: §5.2.1 · p. 24

5.2.2 Example 2: chemical sensors data

• In chemical-sensor data with n = 4, 188, 261 and p = 14, D-OPT IBOSS performance lies between the full-data and subsampling-based approaches.After log-transforming sensor readings, the log-scale relationships appear linear and similar between the D-OPT subdata and the full-data pattern.Source: §5.2.2 · pp. 27–28

6 Concluding remarks

• IBOSS is developed for big-data linear regression, with D-optimality, asymptotic analysis, and simulated and real-data demonstrations supporting its performance.The authors emphasize that informative subdata selection remains important for extracting useful information from big data.Source: §6 · pp. 28–29
• The analysis derives a covariance lower bound for existing subsampling estimators and develops a computationally efficient algorithm approximating D-optimal IBOSS subdata.The IBOSS subdata are theoretically characterized before the approximate optimization algorithm is introduced.Source: §6 · p. 28
• Future work should develop IBOSS methods under alternative criteria, including A-optimality for average estimator variance and c-optimality for a specified parameter function.These criteria have different inferential purposes and may yield additional efficient IBOSS procedures.Source: §6 · pp. 28–29

A.1 Proof of Theorem 1

• The proof of Theorem 1 establishes unbiasedness and decomposes estimator variance into conditional sampling and selection components, using projection and convexity arguments.The proof also treats sampling with and without replacement through inclusion-probability bounds.Source: §A.1 · pp. 29–30

A.2 Proof of Theorem 2

• Theorem 2’s proof shows that selecting the 2p extreme covariate points equally often yields an identity standardized information matrix and D-optimal subdata.The construction uses each covariate’s minimum or maximum value.Source: §A.2 · p. 31

A.3 Proof of Theorem 3

• Theorem 3’s proof analyzes Algorithm 1 through order statistics, concomitants, subdata means, variances, and the resulting information matrix.Although processing covariates can remove some exact extremes, this does not affect the final result when r remains fixed as n increases.Source: §A.3 · pp. 31–32

A.4 Proof of Theorem 4

• Theorem 4’s proof bounds the relevant spectral quantities using the decomposition R = VΛV^T, eigenvalue inequalities, and bounds on standardized covariates.The argument combines equations (35), (36), (41), and (39).Source: §A.4 · pp. 33–34

A.5 Proof of Theorem 5

• Theorem 5’s equation (22) is established by treating five cases for r, support endpoints, and growth conditions, with equation (21) following directly from (20).The cases include fixed r, r →∞ with r/n →0 and bounded support, finite or infinite support endpoints, and a fifth case combining earlier arguments.Source: §A.5 · pp. 34–35
• For the third case, the proof reduces the required order-statistic ratio to z_(n−r+1)j/z_(n)j = 1 + oP(1) using Hall (1979) and assumption (23).The result follows because r →∞ and the stated condition implies the intermediate order statistic has the required asymptotic ratio.Source: §A.5 · pp. 34–35
• Theorem 5’s fourth case follows by applying the third-case argument to −z, using z_(r)j = −(−z)_(n−r+1)j and the correspondence between conditions (24) and (23).Source: §A.5 · p. 35

A.6 Proof of Theorem 6

• Theorem 6 begins by expressing the selected subdata through coordinate order statistics z_(i)j and their concomitants, then derives the vector component needed for equations (26) and (27).The proof explicitly defines σ_j and ρ_j1j2 and invokes the construction from Algorithm 1 before specializing to two distributions.Source: §A.6 · p. 35

A.6.1 Proof of equation (26) in Theorem 6

• For normal covariates, equations (26) are obtained from Gaussian extreme-value results in Galambos (1987), yielding expansions involving 2 log n and O_P(1) terms.The proof combines results corresponding to the upper order statistics and uses the relation c = k + O_P(1/log n).Source: §A.6.1 · p. 36

A.6.2 Proof of equation (27) in Theorem 6

• For lognormal covariates, the proof writes z_ij = exp(U_ij) with U_i ∼ N(µ, Σ), assumes ρ_lj ≥ 0, and derives exponential extreme-value expansions.The resulting quantities include v_1 = (e^µ1, ..., e^µp)^T and B_5 = diag(e^2µ1, ..., e^2µp).Source: §A.6.2 · p. 37
• The lognormal derivation produces terms proportional to exp(2 log n){e^µj + oP(1)} and exp(2 log n){e^2µj + oP(1)} in the relevant expansions.Source: §A.6.2 · p. 37

A.7 Proof of results in Table 1

• The D-OPT IBOSS variance of the β1 estimator is of the same order as (z(n)1 − z(1)1)−2, whereas UNI’s lower bound is O(1) when ν > 2.Source: §A.7 · p. 38

Supplementary Material

• The supplementary material presents additional numerical results about IBOSS performance.Source: Supplementary Material · p. 41

S.1 Predictive performance

• The predictive evaluation focuses on MSPE because irreducible response variance σ2 cannot be reduced by improving the subdata.The simulations generate new samples of size 5,000, repeat the procedure 1,000 times, and compare log10(MSPE) against log10(n).Source: §S.1 · p. 41
• D-OPT IBOSS uniformly dominates UNI and LEV for prediction, with a larger advantage under heavier-tailed covariate distributions.For Case 6, D-OPT IBOSS is almost identical to the full-data approach, while LEV significantly outperforms UNI.Source: §S.1 · p. 42, Figure S.1 · p. 43

S.2 Column permutation

• IBOSS performance is very similar under the original and randomly permuted covariate-column orders, agreeing with the theoretical results.The experiment uses an unexchangeable correlation structure and unequal column variances.Source: §S.2 · p. 42, Figure S.2 · p. 44

S.3 Interaction model

• The interaction-model comparison fixes k = 1000, varies n, and estimates slope-parameter MSEs from 1,000 simulation iterations.Source: §S.3 · pp. 43, 45
• IBOSS remains the most efficient subdata-based method for estimating slope parameters when the true model includes pairwise interactions omitted from selection.Only main effects guide subdata selection, while interaction terms are included in parameter estimation; covariates follow multivariate normal or lognormal distributions.Source: §S.3 · pp. 42, 45

Nonlinear relationships

• The data-generating relationships are nonlinear and cannot be linearized, so a finite-dimensional linear model is misspecified.Despite nonlinear full-data models WM1 or WM2, the analysis uses a linear main effects model for the selected subdata.Source: §S.4 · pp. 45–46
• No method, including the full-data approach, dominates others, and larger sample sizes do not necessarily yield more accurate results.The underlying model is incorrect, making it difficult to determine which method will produce satisfactory results.Source: §S.4 · p. 46
• IBOSS does not always produce the worst results, but these studies do not establish that it works better than other methods.The authors explicitly avoid claiming IBOSS superiority in this nonlinear scenario.Source: §S.4 · p. 46
• The study evaluates slope and intercept estimation MSEs and mean-response prediction MSPEs as the full-data size n changes with fixed subdata size k = 1000.Figure S.4 reports log10-transformed metrics for nonlinear models WM1 and WM2, using the linear main effects model for subdata.Source: §S.4 · pp. 46–47, Figure S.4 · p. 47

Accuracy-cost tradeoff

• For the IBOSS method, as required CPU time increases, MSE decreases, indicating a clear tradeoff between computational cost and estimation accuracy.Figure S.5 reports average CPU times and MSEs for different subdata sample sizes k.Source: §S.5 · p. 48, Figure S.5 · p. 48
• For IBOSS at k = 102, CPU time reaches 7.0839 seconds while MSE decreases from 13.57091 to 0.00786855, indicating significantly improved estimation accuracy with little computational-cost increase.Source: §S.5 · p. 48

Comparison with repeated UNI

• At similar CPU times, repeated UNI can produce smaller MSEs than IBOSS, but this difference is not very significant relative to IBOSS’s advantage in Case 4 with heavier-tailed covariates.Full data use n = 5 × 10^5 and p = 500, with subdata sizes k = 10^3, 5 × 10^3, 10^4, 5 × 10^4, and 10^5.Source: §S.5 · p. 49

S.6 Comparison with the divide-and-conquer method

• The divide-and-conquer method divides the full data into S subdata sets, computes ordinary least-squares estimates for each, averages them, and chooses S = floor(n^1/4), following Battey et al. (2015).Average CPU times and MSEs are calculated from 100 repetitions of the simulation.Source: §S.6 · p. 50
• For p = 50, the divide-and-conquer method is better under normal covariates, methods perform similarly under mixtures, and IBOSS D-OPT dominates under a t1 distribution.IBOSS D-OPT is also more computationally efficient, especially for large n, and CPU times for both methods do not depend on covariate distribution.Source: §S.6 · p. 50
• For the setting shown in Figure S.8, IBOSS D-OPT has a more significant computational-cost advantage over DC, while distribution-dependent relative estimation-efficiency patterns are similar to those with p = 50.Figure S.8 reports MSEs and CPU times for estimating the slope parameter, with subdata size k = 1000.Source: §S.6 · p. 51, Figure S.8 · p. 52

Experimental setup

• The IBOSS method selects subdata for elastic-net regularization, implemented with the glmnet R package, with tuning parameters selected by cross-validation; the setup uses p = 60 and k = 1000.Source: §S.7 · pp. 52–53
• Across 100 repetitions, models fitted on each subdata set are evaluated by MSPE on a new sample of size 5,000.Figure S.9 presents the resulting simulation comparisons.Source: §S.7 · p. 53

Results

• D-OPT IBOSS uniformly dominates the subsampling-based UNI and LEV methods, with a larger advantage for heavier-tailed covariate distributions.Figure S.9 reports MSPEs for elastic-net prediction using subdata of size k = 1000.Source: §S.7 · p. 53, Figure S.9 · p. 54
• Ridge regression produces similar results, so its results are omitted.Source: §S.7 · p. 53

Unequal variance

• The heteroscedastic simulation used the main-paper setup, but generated each error variance from a squared exponential random variable.Error-term standard deviations came from an exponential distribution with rate parameter 1.Source: §S.8 · p. 54
• D-OPT IBOSS uniformly dominates the subsampling-based methods UNI and LEV, with a more significant advantage when covariate tails are heavier.Figure S.10 reports slope-estimation MSEs with fixed subdata size k = 1000 and varying full data size n.Source: §S.8 · pp. 54–55, Figure S.10 · p. 55
• Transformations to stabilize variances and weighted least squares were noted as alternatives but left for future work.The paper states that these questions are beyond its scope and will be investigated in another project.Source: §S.8 · p. 55