This paper introduces a simple principle for robust statistical inference via appropriate shrinkage on the data. This widens the scope of high-dimensional techniques, reducing the distributional conditions from subexponential or sub-Gaussian to more relaxed bounded second or fourth moment. As an illustration of this principle, we focus on robust estimation of the low-rank matrix ${\mathrm{\Theta }}^{\ast }$ from the trace regression model $\mathit{Y}=Tr({{\mathrm{\Theta }}^{\ast }}^{\top }\mathbf{X})+\mathit{\epsilon }$. It encompasses four popular problems: sparse linear model, compressed sensing, matrix completion and multitask learning. We propose to apply the penalized least-squares approach to the appropriately truncated or shrunk data. Under only bounded $2+\mathit{\delta }$ moment condition on the response, the proposed robust methodology yields an estimator that possesses the same statistical error rates as previous literature with sub-Gaussian errors. For sparse linear model and multitask regression, we further allow the design to have only bounded fourth moment and obtain the same statistical rates. As a byproduct, we give a robust covariance estimator with concentration inequality and optimal rate of convergence in terms of the spectral norm, when the samples only bear bounded fourth moment. This result is of its own interest and importance. We reveal that under high dimensions, the sample covariance matrix is not optimal whereas our proposed robust covariance can achieve optimality. Extensive simulations are carried out to support the theories.

Vector autoregressive (VAR) models aim to capture linear temporal interdependencies among multiple time series. They have been widely used in macroeconomics and financial econometrics and more recently have found novel applications in functional genomics and neuroscience. These applications have also accentuated the need to investigate the behavior of the VAR model in a high-dimensional regime, which will provide novel insights into the role of temporal dependence for regularized estimates of the models parameters. However, hardly anything is known regarding posterior model selection consistency for Bayesian VAR models in such regimes.

In this work, we develop a pseudo-likelihood based Bayesian approach for consistent variable selection in high-dimensional VAR models by considering hierarchical normal priors on the autoregressive coefficients, as well as on the model space. We establish strong selection consistency of the proposed method, namely that the posterior probability assigned to the true underlying VAR model converges to one under high-dimensional scaling where the dimension p of the VAR system grows nearly exponentially with the sample size n.

Further, the result is established under mild regularity conditions on the problem parameters. Finally, as a by-product of these results, we also establish strong selection consistency for the sparse high-dimensional linear regression model with serially correlated regressors and errors.

In this paper, we derive nonasymptotic error bounds for the Lasso estimator when the penalty parameter for the estimator is chosen using K-fold cross-validation. Our bounds imply that the cross-validated Lasso estimator has nearly optimal rates of convergence in the prediction, ${\mathit{L}}^2$, and ${\mathit{L}}^1$ norms. For example, we show that in the model with the Gaussian noise and under fairly general assumptions on the candidate set of values of the penalty parameter, the estimation error of the cross-validated Lasso estimator converges to zero in the prediction norm with the $\sqrt{\mathit{s}log\mathit{p}/\mathit{n}}\times \sqrt{log(\mathit{p}\mathit{n})}$ rate, where n is the sample size of available data, p is the number of covariates and s is the number of nonzero coefficients in the model. Thus, the cross-validated Lasso estimator achieves the fastest possible rate of convergence in the prediction norm up to a small logarithmic factor $\sqrt{log(\mathit{p}\mathit{n})}$, and similar conclusions apply for the convergence rate both in ${\mathit{L}}^2$ and in ${\mathit{L}}^1$ norms. Importantly, our results cover the case when p is (potentially much) larger than n and also allow for the case of non-Gaussian noise. Our paper therefore serves as a justification for the widely spread practice of using cross-validation as a method to choose the penalty parameter for the Lasso estimator.

Despite the popularity and practical success of total variation (TV) regularization for function estimation, surprisingly little is known about its theoretical performance in a statistical setting. While TV regularization has been known for quite some time to be minimax optimal for denoising one-dimensional signals, for higher dimensions this remains elusive until today. In this paper, we consider frame-constrained TV estimators including many well-known (overcomplete) frames in a white noise regression model, and prove their minimax optimality w.r.t. ${\mathit{L}}^{\mathit{q}}$-risk ($1\le \mathit{q}<\infty$) up to a logarithmic factor in any dimension $\mathit{d}\ge 1$. Overcomplete frames are an established tool in mathematical imaging and signal recovery, and their combination with TV regularization has been shown to give excellent results in practice, which our theory now confirms. Our results rely on a novel connection between frame-constraints and certain Besov norms, and on an interpolation inequality to relate them to the risk functional. Additionally, our results explain a phase transition in the minimax risk for $\mathit{BV}$ functions.

For the sparse vector model, we consider estimation of the target vector, of its ${\ell }_2$-norm and of the noise variance. We construct adaptive estimators and establish the optimal rates of adaptive estimation when adaptation is considered with respect to the triplet “noise level—noise distribution—sparsity.” We consider classes of noise distributions with polynomially and exponentially decreasing tails as well as the case of Gaussian noise. The obtained rates turn out to be different from the minimax nonadaptive rates when the triplet is known. A crucial issue is the ignorance of the noise variance. Moreover, knowing or not knowing the noise distribution can also influence the rate. For example, the rates of estimation of the noise variance can differ depending on whether the noise is Gaussian or sub-Gaussian without a precise knowledge of the distribution. Estimation of noise variance in our setting can be viewed as an adaptive variant of robust estimation of scale in the contamination model, where instead of fixing the “nominal” distribution in advance we assume that it belongs to some class of distributions.

A common goal in statistics and machine learning is to learn models that can perform well against distributional shifts, such as latent heterogeneous subpopulations, unknown covariate shifts or unmodeled temporal effects. We develop and analyze a distributionally robust stochastic optimization (DRO) framework that learns a model providing good performance against perturbations to the data-generating distribution. We give a convex formulation for the problem, providing several convergence guarantees. We prove finite-sample minimax upper and lower bounds, showing that distributional robustness sometimes comes at a cost in convergence rates. We give limit theorems for the learned parameters, where we fully specify the limiting distribution so that confidence intervals can be computed. On real tasks including generalizing to unknown subpopulations, fine-grained recognition and providing good tail performance, the distributionally robust approach often exhibits improved performance.

The aim of the paper is to describe a bootstrap, contrary to the sieve bootstrap, valid under either long memory ($\mathit{LM}$) or short memory ($\mathit{SM}$) dependence. One of the reasons of the failure of the sieve bootstrap in our context is that under $\mathit{LM}$ dependence, the sieve bootstrap may not be able to capture the true covariance structure of the original data. We also describe and examine the validity of the bootstrap scheme for the least squares estimator of the parameter in a regression model and for model specification. The motivation for the latter example comes from the observation that the asymptotic distribution of the test is intractable.

We study exponential families of distributions that are multivariate totally positive of order 2 (${\mathrm{MTP}}_2$), show that these are convex exponential families and derive conditions for existence of the MLE. Quadratic exponential familes of ${\mathrm{MTP}}_2$ distributions contain attractive Gaussian graphical models and ferromagnetic Ising models as special examples. We show that these are defined by intersecting the space of canonical parameters with a polyhedral cone whose faces correspond to conditional independence relations. Hence ${\mathrm{MTP}}_2$ serves as an implicit regularizer for quadratic exponential families and leads to sparsity in the estimated graphical model. We prove that the maximum likelihood estimator (MLE) in an ${\mathrm{MTP}}_2$ binary exponential family exists if and only if both of the sign patterns $(1,-1)$ and $(-1,1)$ are represented in the sample for every pair of variables; in particular, this implies that the MLE may exist with $\mathit{n}=\mathit{d}$ observations, in stark contrast to unrestricted binary exponential families where $2^{\mathit{d}}$ observations are required. Finally, we provide a novel and globally convergent algorithm for computing the MLE for ${\mathrm{MTP}}_2$ Ising models similar to iterative proportional scaling and apply it to the analysis of data from two psychological disorders.

The bootstrap, introduced by The Jackknife, the Bootstrap and Other Resampling Plans ((1982), SIAM), has become a very popular method for estimating variances and constructing confidence intervals. A key insight is that one can approximate the properties of estimators by using the empirical distribution function of the sample as an approximation for the true distribution function. This approach views the uncertainty in the estimator as coming exclusively from sampling uncertainty. We argue that for causal estimands the uncertainty arises entirely, or partially, from a different source, corresponding to the stochastic nature of the treatment received. We develop a bootstrap procedure for inference regarding the average treatment effect that accounts for this uncertainty, and compare its properties to that of the classical bootstrap. We consider completely randomized and observational designs as well as designs with imperfect compliance.

We study the principal components of covariance estimators in multivariate mixed-effects linear models. We show that, in high dimensions, the principal eigenvalues and eigenvectors may exhibit bias and aliasing effects that are not present in low-dimensional settings. We derive the first-order limits of the principal eigenvalue locations and eigenvector projections in a high-dimensional asymptotic framework, allowing for general population spectral distributions for the random effects and extending previous results from a more restrictive spiked model. Our analysis uses free probability techniques, and we develop two general tools of independent interest—strong asymptotic freeness of GOE and deterministic matrices and a free deterministic equivalent approximation for bilinear forms of resolvents.

In this paper, we discuss the estimation of a nonparametric component ${\mathit{f}}_1$ of a nonparametric additive model $\mathit{Y}={\mathit{f}}_1({\mathit{X}}_1)+\cdots +{\mathit{f}}_{\mathit{q}}({\mathit{X}}_{\mathit{q}})+\mathit{ϵ}$. We allow the number q of additive components to grow to infinity and we make sparsity assumptions about the number of nonzero additive components. We compare this estimation problem with that of estimating ${\mathit{f}}_1$ in the oracle model $\mathit{Z}={\mathit{f}}_1({\mathit{X}}_1)+\mathit{ϵ}$, for which the additive components ${\mathit{f}}_2,\dots ,{\mathit{f}}_{\mathit{q}}$ are known. We construct a two-step presmoothing-and-resmoothing estimator of ${\mathit{f}}_1$ and state finite-sample bounds for the difference between our estimator and a corresponding smoothing estimator ${\hat{\mathit{f}}}_1^{\text{(oracle)}}$ in the oracle model. In an asymptotic setting, these bounds can be used to show asymptotic equivalence of our estimator and the oracle estimator; the paper thus shows that, asymptotically, under strong enough sparsity conditions, knowledge of ${\mathit{f}}_2,\dots ,{\mathit{f}}_{\mathit{q}}$ has no effect on estimation accuracy. Our first step is to estimate ${\mathit{f}}_1$ with an undersmoothed estimator based on near-orthogonal projections with a group Lasso bias correction. In the second step, we construct pseudo responses $\hat{\mathit{Y}}$ by evaluating this undersmoothed estimator of ${\mathit{f}}_1$ at the design points and then apply the smoothing method of the oracle estimator ${\hat{\mathit{f}}}_1^{\text{(oracle)}}$ to the nonparametric regression problem with “responses” $\hat{\mathit{Y}}$ and covariates ${\mathit{X}}_1$. Our mathematical exposition centers primarily on establishing properties of the presmoothing estimator. We present simulation results demonstrating close-to-oracle performance of our estimator in practical applications.

In this paper, we study high-dimensional sparse Quadratic Discriminant Analysis (QDA) and aim to establish the optimal convergence rates for the classification error. Minimax lower bounds are established to demonstrate the necessity of structural assumptions such as sparsity conditions on the discriminating direction and differential graph for the possible construction of consistent high-dimensional QDA rules.

We then propose a classification algorithm called SDAR using constrained convex optimization under the sparsity assumptions. Both minimax upper and lower bounds are obtained and this classification rule is shown to be simultaneously rate optimal over a collection of parameter spaces, up to a logarithmic factor. Simulation studies demonstrate that SDAR performs well numerically. The algorithm is also illustrated through an analysis of prostate cancer data and colon tissue data. The methodology and theory developed for high-dimensional QDA for two groups in the Gaussian setting are also extended to multigroup classification and to classification under the Gaussian copula model.

This paper is concerned with the limiting spectral behaviors of large dimensional Kendall’s rank correlation matrices generated by samples with independent and continuous components. The statistical setting in this paper covers a wide range of highly skewed and heavy-tailed distributions since we do not require the components to be identically distributed, and do not need any moment conditions. We establish the central limit theorem (CLT) for the linear spectral statistics (LSS) of the Kendall’s rank correlation matrices under the Marchenko–Pastur asymptotic regime, in which the dimension diverges to infinity proportionally with the sample size. We further propose three nonparametric procedures for high dimensional independent test and their limiting null distributions are derived by implementing this CLT. Our numerical comparisons demonstrate the robustness and superiority of our proposed test statistics under various mixed and heavy-tailed cases.

This paper considers both approximate and exact designs for estimating the total effects under one crossover and two interference models. They are different from the traditional block designs in the sense that the assigned treatments also affect their neighboring plots, hence a design is understood as a collection of sequences of treatments. A notable result in literature is that the circular neighbor balanced design (CNBD) is optimal among designs, which do not allow treatments to be neighbors of themselves. However, we find it necessary to allow self-neighboring, and further show that it is the best to allocate each treatment in a subblock of adjacent plots with equal or almost equal numbers of replications. This explains why the efficiency of CNBD drops down to $50\mathrm{\%}$ as the sequence length, say k, increases. Unlike CNBD or the designs for direct effects, our proposed designs do not try to put as many treatments in a sequence as possible. The optimal number of distinct treatments in a sequence is around $\sqrt{2\mathit{k}}$ for crossover designs and $\sqrt{\mathit{k}/0.96}$ for interference models, whenever they are smaller than the total number of treatments under consideration. We systematically study necessary and sufficient conditions for any design to be universally optimal under the approximate design framework, based on which algorithms for deriving optimal or efficient exact designs are proposed. This hybrid nature of cohesively combining theories with algorithms makes our method more flexible than existing ones in the following aspects. (i) Not only symmetric designs are studied, general procedures for producing asymmetric designs are also provided. (ii) Our method applies to any form of within-block covariance matrix instead of specific forms. (iii) We cover all configurations of the numbers of treatments and sequence lengths, especially for large values of them when purely computational methods are not applicable. (iv) On top of the latter, we cover a continuous spectrum of the number of sequences instead of special numbers decided by combinatorial constraints.

We study the sample complexity of learning a high-dimensional simplex from a set of points uniformly sampled from its interior. Learning of simplices is a long studied problem in computer science and has applications in computational biology and remote sensing, mostly under the name of “spectral unmixing.” We theoretically show that a sufficient sample complexity for reliable learning of a K-dimensional simplex up to a total-variation error of ϵ is $\mathit{O}(\frac{{\mathit{K}}^2}{\mathit{ϵ}}log\frac{\mathit{K}}{\mathit{ϵ}})$, which yields a substantial improvement over existing bounds. Based on our new theoretical framework, we also propose a heuristic approach for the inference of simplices. Experimental results on synthetic and real-world datasets demonstrate a comparable performance for our method on noiseless samples, while we outperform the state-of-the-art in noisy cases.

We propose a novel two-regime regression model where regime switching is driven by a vector of possibly unobservable factors. When the factors are latent, we estimate them by the principal component analysis of a panel data set. We show that the optimization problem can be reformulated as mixed integer optimization, and we present two alternative computational algorithms. We derive the asymptotic distribution of the resulting estimator under the scheme that the threshold effect shrinks to zero. In particular, we establish a phase transition that describes the effect of first-stage factor estimation as the cross-sectional dimension of panel data increases relative to the time-series dimension. Moreover, we develop bootstrap inference and illustrate our methods via numerical studies.

Clustering with Bregman divergences encompasses a wide family of clustering procedures that are well suited to mixtures of distributions from exponential families (J. Mach. Learn. Res. 6 (2005) 1705–1749). However, these techniques are highly sensitive to noise. To address the issue of clustering data with possibly adversarial noise, we introduce a robustified version of Bregman clustering based on a trimming approach. We investigate its theoretical properties, showing for instance that our estimator converges at a sub-Gaussian rate $1/\sqrt{\mathit{n}}$, where n denotes the sample size, under mild tail assumptions. We also show that it is robust to a certain amount of noise, stated in terms of breakdown point. We also derive a Lloyd-type algorithm with a trimming parameter, along with a heuristic to select this parameter and the number of clusters from sample. Some numerical experiments assess the performance of our method on simulated and real datasets.

We consider the estimation and inference in a system of high-dimensional regression equations allowing for temporal and cross-sectional dependency in covariates and error processes, covering rather general forms of weak temporal dependence. A sequence of regressions with many regressors using LASSO (Least Absolute Shrinkage and Selection Operator) is applied for variable selection purpose, and an overall penalty level is carefully chosen by a block multiplier bootstrap procedure to account for multiplicity of the equations and dependencies in the data. Correspondingly, oracle properties with a jointly selected tuning parameter are derived. We further provide high-quality de-biased simultaneous inference on the many target parameters of the system. We provide bootstrap consistency results of the test procedure, which are based on a general Bahadur representation for the Z-estimators with dependent data. Simulations demonstrate good performance of the proposed inference procedure. Finally, we apply the method to quantify spillover effects of textual sentiment indices in a financial market and to test the connectedness among sectors.

Multiple testing of a single hypothesis and testing multiple hypotheses are usually done in terms of p-values. In this paper, we replace p-values with their natural competitor, e-values, which are closely related to betting, Bayes factors and likelihood ratios. We demonstrate that e-values are often mathematically more tractable; in particular, in multiple testing of a single hypothesis, e-values can be merged simply by averaging them. This allows us to develop efficient procedures using e-values for testing multiple hypotheses.

Causal questions are omnipresent in many scientific problems. While much progress has been made in the analysis of causal relationships between random variables, these methods are not well suited if the causal mechanisms only manifest themselves in extremes. This work aims to connect the two fields of causal inference and extreme value theory. We define the causal tail coefficient that captures asymmetries in the extremal dependence of two random variables. In the population case, the causal tail coefficient is shown to reveal the causal structure if the distribution follows a linear structural causal model. This holds even in the presence of latent common causes that have the same tail index as the observed variables. Based on a consistent estimator of the causal tail coefficient, we propose a computationally highly efficient algorithm that estimates the causal structure. We prove that our method consistently recovers the causal order and we compare it to other well-established and nonextremal approaches in causal discovery on synthetic and real data. The code is available as an open-access $\mathtt{R}$ package.

This paper investigates the statistical estimation of a discrete mixing measure ${\mathit{\mu }}^0$ involved in a kernel mixture model. Using some recent advances in ${\ell }_1$-regularization over the space of measures, we introduce a “data fitting and regularization” convex program for estimating ${\mathit{\mu }}^0$ in a grid-less manner from a sample of mixture law, this method is referred to as Beurling-LASSO.

Our contribution is two-fold: we derive a lower bound on the bandwidth of our data fitting term depending only on the support of ${\mathit{\mu }}^0$ and its so-called “minimum separation” to ensure quantitative support localization error bounds; and under a so-called “nondegenerate source condition” we derive a nonasymptotic support stability property. This latter shows that for a sufficiently large sample size n, our estimator has exactly as many weighted Dirac masses as the target ${\mathit{\mu }}^0$, converging in amplitude and localization towards the true ones. Finally, we also introduce some tractable algorithms for solving this convex program based on “Sliding Frank–Wolfe” or “Conic Particle Gradient Descent”.

Statistical performances of this estimator are investigated designing a so-called “dual certificate”, which is appropriate to our setting. Some classical situations as, for example, mixtures of super-smooth distributions (see, e.g., Gaussian distributions) or ordinary-smooth distributions (see, e.g., Laplace distributions), are discussed at the end of the paper.