Consider a sequence of real data points $X_{1},\ldots ,X_{n}$ with underlying means $\theta ^{*}_{1},\dots ,\theta ^{*}_{n}$. This paper starts from studying the setting that $\theta ^{*}_{i}$ is both piecewise constant and monotone as a function of the index $i$. For this, we establish the exact minimax rate of estimating such monotone functions, and thus give a nontrivial answer to an open problem in the shape-constrained analysis literature. The minimax rate under the loss of the sum of squared errors involves an interesting iterated logarithmic dependence on the dimension, a phenomenon that is revealed through characterizing the interplay between the isotonic shape constraint and model selection complexity. We then develop a penalized least-squares procedure for estimating the vector $\theta ^{*}=(\theta^{*}_{1},\dots ,\theta ^{*}_{n})^{\mathsf{T}}$. This estimator is shown to achieve the derived minimax rate adaptively. For the proposed estimator, we further allow the model to be misspecified and derive oracle inequalities with the optimal rates, and show there exists a computationally efficient algorithm to compute the exact solution.

We consider the problem of finding anomalies in a $d$-dimensional field of independent random variables $\{Y_{i}\}_{i\in \{1,\ldots,n\}^{d}}$, each distributed according to a one-dimensional natural exponential family $\mathcal{F}=\{F_{\theta }\}_{\theta \in \Theta }$. Given some baseline parameter $\theta _{0}\in \Theta $, the field is scanned using local likelihood ratio tests to detect from a (large) given system of regions $\mathcal{R}$ those regions $R\subset \{1,\ldots,n\}^{d}$ with $\theta _{i}\neq \theta _{0}$ for some $i\in R$. We provide a unified methodology which controls the overall familywise error (FWER) to make a wrong detection at a given error rate.

Fundamental to our method is a Gaussian approximation of the distribution of the underlying multiscale test statistic with explicit rate of convergence. From this, we obtain a weak limit theorem which can be seen as a generalized weak invariance principle to nonidentically distributed data and is of independent interest. Furthermore, we give an asymptotic expansion of the procedures power, which yields minimax optimality in case of Gaussian observations.

In this paper, we introduce new, easily implementable designs for drawing causal inference from randomized experiments on networks with interference. Inspired by the idea of matching in observational studies, we introduce the notion of considering a treatment assignment as a “quasi-coloring” on a graph. Our idea of a perfect quasi-coloring strives to match every treated unit on a given network with a distinct control unit that has identical number of treated and control neighbors. For a wide range of interference functions encountered in applications, we show both by theory and simulations that the classical Neymanian estimator for the direct effect has desirable properties for our designs.

We study the problem of learning a tree Ising model from samples such that subsequent predictions made using the model are accurate. The prediction task considered in this paper is that of predicting the values of a subset of variables given values of some other subset of variables. Virtually all previous work on graphical model learning has focused on recovering the true underlying graph. We define a distance (“small set TV” or ssTV) between distributions $P$ and $Q$ by taking the maximum, over all subsets $\mathcal{S}$ of a given size, of the total variation between the marginals of $P$ and $Q$ on $\mathcal{S}$; this distance captures the accuracy of the prediction task of interest. We derive nonasymptotic bounds on the number of samples needed to get a distribution (from the same class) with small ssTV relative to the one generating the samples. One of the main messages of this paper is that far fewer samples are needed than for recovering the underlying tree, which means that accurate predictions are possible using the wrong tree.

We study the nonparametric maximum likelihood estimator (NPMLE) for estimating Gaussian location mixture densities in $d$-dimensions from independent observations. Unlike usual likelihood-based methods for fitting mixtures, NPMLEs are based on convex optimization. We prove finite sample results on the Hellinger accuracy of every NPMLE. Our results imply, in particular, that every NPMLE achieves near parametric risk (up to logarithmic multiplicative factors) when the true density is a discrete Gaussian mixture without any prior information on the number of mixture components. NPMLEs can naturally be used to yield empirical Bayes estimates of the oracle Bayes estimator in the Gaussian denoising problem. We prove bounds for the accuracy of the empirical Bayes estimate as an approximation to the oracle Bayes estimator. Here our results imply that the empirical Bayes estimator performs at nearly the optimal level (up to logarithmic factors) for denoising in clustering situations without any prior knowledge of the number of clusters.

Estimation of the prediction error of a linear estimation rule is difficult if the data analyst also uses data to select a set of variables and constructs the estimation rule using only the selected variables. In this work, we propose an asymptotically unbiased estimator for the prediction error after model search. Under some additional mild assumptions, we show that our estimator converges to the true prediction error in $L^{2}$ at the rate of $O(n^{-1/2})$, with $n$ being the number of data points. Our estimator applies to general selection procedures, not requiring analytical forms for the selection. The number of variables to select from can grow as an exponential factor of $n$, allowing applications in high-dimensional data. It also allows model misspecifications, not requiring linear underlying models. One application of our method is that it provides an estimator for the degrees of freedom for many discontinuous estimation rules like best subset selection or relaxed Lasso. Connection to Stein’s Unbiased Risk Estimator is discussed. We consider in-sample prediction errors in this work, with some extension to out-of-sample errors in low-dimensional, linear models. Examples such as best subset selection and relaxed Lasso are considered in simulations, where our estimator outperforms both $C_{p}$ and cross validation in various settings.

We study joint estimation of the inverse temperature and magnetization parameters $(\beta ,B)$ of an Ising model with a nonnegative coupling matrix $A_{n}$ of size $n\times n$, given one sample from the Ising model. We give a general bound on the rate of consistency of the bi-variate pseudo-likelihood estimator. Using this, we show that estimation at rate $n^{-1/2}$ is always possible if $A_{n}$ is the adjacency matrix of a bounded degree graph. If $A_{n}$ is the scaled adjacency matrix of a graph whose average degree goes to $+\infty $, the situation is a bit more delicate. In this case, estimation at rate $n^{-1/2}$ is still possible if the graph is not regular (in an asymptotic sense). Finally, we show that consistent estimation of both parameters is impossible if the graph is Erdős–Renyi with parameter $p>0$ independent of $n$, thus confirming that estimation is harder on approximately regular graphs with large degree.

Consider the problem of estimating average treatment effects when a large number of covariates are used to adjust for possible confounding through outcome regression and propensity score models. We develop new methods and theory to obtain not only doubly robust point estimators for average treatment effects, which remain consistent if either the propensity score model or the outcome regression model is correctly specified, but also model-assisted confidence intervals, which are valid when the propensity score model is correctly specified but the outcome model may be misspecified. With a linear outcome model, the confidence intervals are doubly robust, that is, being also valid when the outcome model is correctly specified but the propensity score model may be misspecified. Our methods involve regularized calibrated estimators with Lasso penalties but carefully chosen loss functions, for fitting propensity score and outcome regression models. We provide high-dimensional analysis to establish the desired properties of our methods under comparable sparsity conditions to previous results, which give valid confidence intervals when both the propensity score and outcome models are correctly specified. We present simulation studies and an empirical application which demonstrate advantages of the proposed methods compared with related methods based on regularized maximum likelihood estimation. The methods are implemented in the R package $\mathtt{RCAL}$.

This paper considers a wide range of issues concerning the estimation of the Hurst function of a multifractional Brownian motion when the process is observed on a regular grid. A theoretical lower bound for the minimax risk of this inference problem is established for a wide class of smooth Hurst functions. We also propose a new nonparametric estimator and show that it is rate optimal. Implementation issues of the estimator including how to overcome the presence of a nuisance parameter and choose the tuning parameter from data will be considered. An extensive numerical study is conducted to compare our approach with other approaches.

We study the fundamental limits of detecting the presence of an additive rank-one perturbation, or spike, to a Wigner matrix. When the spike comes from a prior that is i.i.d. across coordinates, we prove that the log-likelihood ratio of the spiked model against the nonspiked one is asymptotically normal below a certain reconstruction threshold which is not necessarily of a “spectral” nature, and that it is degenerate above. This establishes the maximal region of contiguity between the planted and null models. It is known that this threshold also marks a phase transition for estimating the spike: the latter task is possible above the threshold and impossible below. Therefore, both estimation and detection undergo the same transition in this random matrix model. Further information on the performance of the optimal test is also provided. Our proofs are based on Gaussian interpolation methods and a rigorous incarnation of the cavity method, as devised by Guerra and Talagrand in their study of the Sherrington–Kirkpatrick spin-glass model.

We provide statistical guarantees for a family of variational approximations to Bayesian posterior distributions, called $\alpha $-VB, which has close connections with variational approximations of tempered posteriors in the literature. The standard variational approximation is a special case of $\alpha $-VB with $\alpha =1$. When $\alpha \in (0,1]$, a novel class of variational inequalities are developed for linking the Bayes risk under the variational approximation to the objective function in the variational optimization problem, implying that maximizing the evidence lower bound in variational inference has the effect of minimizing the Bayes risk within the variational density family. Operating in a frequentist setup, the variational inequalities imply that point estimates constructed from the $\alpha $-VB procedure converge at an optimal rate to the true parameter in a wide range of problems. We illustrate our general theory with a number of examples, including the mean-field variational approximation to (low)-high-dimensional Bayesian linear regression with spike and slab priors, Gaussian mixture models and latent Dirichlet allocation.

Median-of-means (MOM) based procedures have been recently introduced in learning theory (Lugosi and Mendelson (2019); Lecué and Lerasle (2017)). These estimators outperform classical least-squares estimators when data are heavy-tailed and/or are corrupted. None of these procedures can be implemented, which is the major issue of current MOM procedures (Ann. Statist. 47 (2019) 783–794).

In this paper, we introduce minmax MOM estimators and show that they achieve the same sub-Gaussian deviation bounds as the alternatives (Lugosi and Mendelson (2019); Lecué and Lerasle (2017)), both in small and high-dimensional statistics. In particular, these estimators are efficient under moments assumptions on data that may have been corrupted by a few outliers.

Besides these theoretical guarantees, the definition of minmax MOM estimators suggests simple and systematic modifications of standard algorithms used to approximate least-squares estimators and their regularized versions. As a proof of concept, we perform an extensive simulation study of these algorithms for robust versions of the LASSO.

We present new results for consistency of maximum likelihood estimators with a focus on multivariate mixed models. Our theory builds on the idea of using subsets of the full data to establish consistency of estimators based on the full data. It requires neither that the data consist of independent observations, nor that the observations can be modeled as a stationary stochastic process. Compared to existing asymptotic theory using the idea of subsets, we substantially weaken the assumptions, bringing them closer to what suffices in classical settings. We apply our theory in two multivariate mixed models for which it was unknown whether maximum likelihood estimators are consistent. The models we consider have nonstochastic predictors and multivariate responses which are possibly mixed-type (some discrete and some continuous).

The eigenvector empirical spectral distribution (VESD) is a useful tool in studying the limiting behavior of eigenvalues and eigenvectors of covariance matrices. In this paper, we study the convergence rate of the VESD of sample covariance matrices to the deformed Marčenko–Pastur (MP) distribution. Consider sample covariance matrices of the form $\Sigma ^{1/2}XX^{*}\Sigma ^{1/2}$, where $X=(x_{ij})$ is an $M\times N$ random matrix whose entries are independent random variables with mean zero and variance $N^{-1}$, and $\Sigma $ is a deterministic positive-definite matrix. We prove that the Kolmogorov distance between the expected VESD and the deformed MP distribution is bounded by $N^{-1+\epsilon }$ for any fixed $\epsilon >0$, provided that the entries $\sqrt{N}x_{ij}$ have uniformly bounded 6th moments and $|N/M-1|\ge \tau $ for some constant $\tau >0$. This result improves the previous one obtained in (Ann. Statist. 41 (2013) 2572–2607), which gave the convergence rate $O(N^{-1/2})$ assuming i.i.d. $X$ entries, bounded 10th moment, $\Sigma =I$ and $M<N$. Moreover, we also prove that under the finite $8$th moment assumption, the convergence rate of the VESD is $O(N^{-1/2+\epsilon })$ almost surely for any fixed $\epsilon >0$, which improves the previous bound $N^{-1/4+\epsilon }$ in (Ann. Statist. 41 (2013) 2572–2607).

We consider optimal designs for general multinomial logistic models, which cover baseline-category, cumulative, adjacent-categories and continuation-ratio logit models, with proportional odds, nonproportional odds or partial proportional odds assumption. We derive the corresponding Fisher information matrices in three different forms to facilitate their calculations, determine the conditions for their positive definiteness, and search for optimal designs. We conclude that, unlike the designs for binary responses, a feasible design for a multinomial logistic model may contain less experimental settings than parameters, which is of practical significance. We also conclude that even for a minimally supported design, a uniform allocation, which is typically used in practice, is not optimal in general for a multinomial logistic model. We develop efficient algorithms for searching D-optimal designs. Using examples based on real experiments, we show that the efficiency of an experiment can be significantly improved if our designs are adopted.

The problem of nonparametric inference on a monotone function has been extensively studied in many particular cases. Estimators considered have often been of so-called Grenander type, being representable as the left derivative of the greatest convex minorant or least concave majorant of an estimator of a primitive function. In this paper, we provide general conditions for consistency and pointwise convergence in distribution of a class of generalized Grenander-type estimators of a monotone function. This broad class allows the minorization or majoratization operation to be performed on a data-dependent transformation of the domain, possibly yielding benefits in practice. Additionally, we provide simpler conditions and more concrete distributional theory in the important case that the primitive estimator and data-dependent transformation function are asymptotically linear. We use our general results in the context of various well-studied problems, and show that we readily recover classical results established separately in each case. More importantly, we show that our results allow us to tackle more challenging problems involving parameters for which the use of flexible learning strategies appears necessary. In particular, we study inference on monotone density and hazard functions using informatively right-censored data, extending the classical work on independent censoring, and on a covariate-marginalized conditional mean function, extending the classical work on monotone regression functions.

Archimax copula models can account for any type of asymptotic dependence between extremes and at the same time capture joint risks at medium levels. An Archimax copula is characterized by two functional parameters: the stable tail dependence function $\ell $, and the Archimedean generator $\psi $ which distorts the extreme-value dependence structure. This article develops semiparametric inference for Archimax copulas: a nonparametric estimator of $\ell $ and a moment-based estimator of $\psi $ assuming the latter belongs to a parametric family. Conditions under which $\psi $ and $\ell $ are identifiable are derived. The asymptotic behavior of the estimators is then established under broad regularity conditions; performance in small samples is assessed through a comprehensive simulation study. The Archimax copula model with the Clayton generator is then used to analyze monthly rainfall maxima at three stations in French Brittany. The model is seen to fit the data very well, both in the lower and in the upper tail. The nonparametric estimator of $\ell $ reveals asymmetric extremal dependence between the stations, which reflects heavy precipitation patterns in the area. Technical proofs, simulation results and $\mathsf{R}$ code are provided in the Online Supplement.

This paper investigates estimation of the mean vector under invariant quadratic loss for a spherically symmetric location family with a residual vector with density of the form $f(x,u)=\eta ^{(p+n)/2}f(\eta \{\|x-\theta \|^{2}+\|u\|^{2}\})$, where $\eta $ is unknown. We show that the natural estimator $x$ is admissible for $p=1,2$. Also, for $p\geq 3$, we find classes of generalized Bayes estimators that are admissible within the class of equivariant estimators of the form $\{1-\xi (x/\|u\|)\}x$. In the Gaussian case, a variant of the James–Stein estimator, $[1-\{(p-2)/(n+2)\}/\{\|x\|^{2}/\|u\|^{2}+(p-2)/(n+2)+1\}]x$, which dominates the natural estimator $x$, is also admissible within this class. We also study the related regression model.

Pairwise comparison data arises in many domains, including tournament rankings, web search and preference elicitation. Given noisy comparisons of a fixed subset of pairs of items, we study the problem of estimating the underlying comparison probabilities under the assumption of strong stochastic transitivity (SST). We also consider the noisy sorting subclass of the SST model. We show that when the assignment of items to the topology is arbitrary, these permutation-based models, unlike their parametric counterparts, do not admit consistent estimation for most comparison topologies used in practice. We then demonstrate that consistent estimation is possible when the assignment of items to the topology is randomized, thus establishing a dichotomy between worst-case and average-case designs. We propose two computationally efficient estimators in the average-case setting and analyze their risk, showing that it depends on the comparison topology only through the degree sequence of the topology. We also provide explicit classes of graphs for which the rates achieved by these estimators are optimal. Our results are corroborated by simulations on multiple comparison topologies.

We analyse the reconstruction error of principal component analysis (PCA) and prove nonasymptotic upper bounds for the corresponding excess risk. These bounds unify and improve existing upper bounds from the literature. In particular, they give oracle inequalities under mild eigenvalue conditions. The bounds reveal that the excess risk differs significantly from usually considered subspace distances based on canonical angles. Our approach relies on the analysis of empirical spectral projectors combined with concentration inequalities for weighted empirical covariance operators and empirical eigenvalues.

Many theoretical results for lasso require the samples to be i.i.d. Recent work has provided guarantees for lasso assuming that the time series is generated by a sparse Vector Autoregressive (VAR) model with Gaussian innovations. Proofs of these results rely critically on the fact that the true data generating mechanism (DGM) is a finite-order Gaussian VAR. This assumption is quite brittle: linear transformations, including selecting a subset of variables, can lead to the violation of this assumption. In order to break free from such assumptions, we derive nonasymptotic inequalities for estimation error and prediction error of lasso estimate of the best linear predictor without assuming any special parametric form of the DGM. Instead, we rely only on (strict) stationarity and geometrically decaying $\beta$-mixing coefficients to establish error bounds for lasso for sub-Weibull random vectors. The class of sub-Weibull random variables that we introduce includes sub-Gaussian and subexponential random variables but also includes random variables with tails heavier than an exponential. We also show that, for Gaussian processes, the $\beta$-mixing condition can be relaxed to summability of the $\alpha$-mixing coefficients. Our work provides an alternative proof of the consistency of lasso for sparse Gaussian VAR models. But the applicability of our results extends to non-Gaussian and nonlinear times series models as the examples we provide demonstrate.

We consider the problem of estimating stochastic volatility for a class of second-order parabolic stochastic PDEs. Assuming that the solution is observed at high temporal frequency, we use limit theorems for multipower variations and related functionals to construct consistent nonparametric estimators and asymptotic confidence bounds for the integrated volatility process. As a byproduct of our analysis, we also obtain feasible estimators for the regularity of the spatial covariance function of the noise.

Functional data analysis is typically conducted within the $L^{2}$-Hilbert space framework. There is by now a fully developed statistical toolbox allowing for the principled application of the functional data machinery to real-world problems, often based on dimension reduction techniques such as functional principal component analysis. At the same time, there have recently been a number of publications that sidestep dimension reduction steps and focus on a fully functional $L^{2}$-methodology. This paper goes one step further and develops data analysis methodology for functional time series in the space of all continuous functions. The work is motivated by the fact that objects with rather different shapes may still have a small $L^{2}$-distance and are therefore identified as similar when using a $L^{2}$-metric. However, in applications it is often desirable to use metrics reflecting the visualization of the curves in the statistical analysis. The methodological contributions are focused on developing two-sample and change-point tests as well as confidence bands, as these procedures appear to be conducive to the proposed setting. Particular interest is put on relevant differences; that is, on not trying to test for exact equality, but rather for prespecified deviations under the null hypothesis.

The procedures are justified through large-sample theory. To ensure practicability, nonstandard bootstrap procedures are developed and investigated addressing particular features that arise in the problem of testing relevant hypotheses. The finite sample properties are explored through a simulation study and an application to annual temperature profiles.

We study polynomial time algorithms for estimating the mean of a heavy-tailed multivariate random vector. We assume only that the random vector $X$ has finite mean and covariance. In this setting, the radius of confidence intervals achieved by the empirical mean are large compared to the case that $X$ is Gaussian or sub-Gaussian.

We offer the first polynomial time algorithm to estimate the mean with sub-Gaussian-size confidence intervals under such mild assumptions. Our algorithm is based on a new semidefinite programming relaxation of a high-dimensional median. Previous estimators which assumed only existence of finitely many moments of $X$ either sacrifice sub-Gaussian performance or are only known to be computable via brute-force search procedures requiring time exponential in the dimension.

In recent years, bootstrap methods have drawn attention for their ability to approximate the laws of “max statistics” in high-dimensional problems. A leading example of such a statistic is the coordinatewise maximum of a sample average of $n$ random vectors in $\mathbb{R}^{p}$. Existing results for this statistic show that the bootstrap can work when $n\ll p$, and rates of approximation (in Kolmogorov distance) have been obtained with only logarithmic dependence in $p$. Nevertheless, one of the challenging aspects of this setting is that established rates tend to scale like $n^{-1/6}$ as a function of $n$.

The main purpose of this paper is to demonstrate that improvement in rate is possible when extra model structure is available. Specifically, we show that if the coordinatewise variances of the observations exhibit decay, then a nearly $n^{-1/2}$ rate can be achieved, independent of $p$. Furthermore, a surprising aspect of this dimension-free rate is that it holds even when the decay is very weak. Lastly, we provide examples showing how these ideas can be applied to inference problems dealing with functional and multinomial data.