Papers by Topic

Abstract

Transfer learning, which leverages knowledge from an auxiliary source dataset to improve performance in a primary target domain, has emerged as a pivotal machine learning technique. In this paper, we consider minimax and adaptive estimation of large bandable covariance matrices within the transfer learning framework. We first establish the minimax rate of convergence under the spectral norm and propose a rate-optimal estimation procedure. Our findings reveal intriguing phase transition phenomena that highlight the effectiveness of transfer learning and the use of source samples. We then address the problem of adaptation, establishing the adaptive rate of convergence up to a logarithmic factor. Our results demonstrate that, in sharp contrast to conventional settings, the cost of adaptation in transfer learning can be substantial in certain cases. We propose a novel data-driven algorithm that dynamically adapts to unknown model parameters. These theoretical insights are further validated by a simulation study, demonstrating the practicality and efficiency of the proposed adaptive algorithm.

Abstract

This paper considers minimax and adaptive transfer learning for nonparametric classification under the posterior drift model with distributed differential privacy constraints. Our study is conducted within a heterogeneous framework, encompassing diverse sample sizes, varying privacy parameters, and data heterogeneity across different servers. We first establish the minimax misclassification rate, precisely characterizing the effects of privacy constraints, source samples, and target samples on classification accuracy. The results reveal interesting phase transition phenomena and highlight the intricate trade-offs between preserving privacy and achieving classification accuracy. We then develop a data-driven adaptive classifier that achieves the optimal rate within a logarithmic factor across a large collection of parameter spaces while satisfying the same set of differential privacy constraints. Simulation studies and real-world data applications further elucidate the theoretical analysis with numerical results.

Abstract

Transfer learning provides a powerful tool for incorporating related data into a target study of interest. In epidemiology and medical studies, the classification of a target disease could borrow information across diseases and populations. In this work, we consider transfer learning for high-dimensional generalized linear models (GLMs). A novel algorithm, TransGLM, that incorporates data from the target study as well as the auxiliary studies is proposed. Minimax rate of convergence for estimation is established and the proposed estimator is shown to be rate-optimal. Statistical inference for the target regression coefficients is also studied. Asymptotic normality for a debiased estimator is established and confidence intervals are constructed. Numerical studies show significant improvements in estimation and inference accuracy. Proposed methods are applied to a real data study concerning the classification of colorectal cancer using microbiomes, and are shown to enhance the classification accuracy in comparison to the single-task methods.

Abstract

This paper studies transfer learning for estimating the mean of random functions based on discretely sampled data, where, in addition to observations from the target distribution, auxiliary samples from similar but distinct source distributions are available. The paper considers both common and independent designs and establishes the minimax rates of convergence for both designs. The results reveal an interesting phase transition phenomenon under the two designs and demonstrate the benefits of utilizing the source samples in the low sampling frequency regime. For practical applications, this paper proposes novel data-driven adaptive algorithms that attain the optimal rates of convergence within a logarithmic factor simultaneously over a large collection of parameter spaces. The theoretical findings are complemented by a simulation study that further supports the effectiveness of the proposed algorithms.

Abstract

Motivated by a range of applications, we study in this paper the problem of transfer learning for nonparametric contextual multi-armed bandits under the covariate shift model, where we have data collected on source bandits before the start of the target bandit learning. The minimax rate of convergence for the cumulative regret is established and a novel transfer learning algorithm that attains the minimax regret is proposed. The results quantify the contribution of the data from the source domains for learning in the target domain in the context of nonparametric contextual multi-armed bandits. In view of the general impossibility of adaptation to unknown smoothness, we develop a data-driven algorithm that achieves near-optimal statistical guarantees (up to a logarithmic factor) while automatically adapting to the unknown parameters over a large collection of parameter spaces under an additional self-similarity assumption. A simulation study is carried out to illustrate the benefits of utilizing the data from the auxiliary source domains for learning in the target domain.

Abstract

Transfer learning for high-dimensional Gaussian graphical models (GGMs) is studied with the goal of estimating the target GGM by utilizing the data from similar and related auxiliary studies. The similarity between the target graph and each auxiliary graph is characterized by the sparsity of a divergence matrix. An estimation algorithm, Trans-CLIME, is proposed and shown to attain a faster convergence rate than the minimax rate in the single-task setting. Furthermore, we introduce a universal debiasing method which can be coupled with many initial graph estimators and can be analytically computed in one step. Applying such a debiasing method, a debiased Trans-CLIME estimator is obtained and is shown to be element-wise asymptotically normal. We use the latter fact to construct a multiple testing procedure for edge detection with false discovery rate control. The proposed estimation and multiple testing procedures demonstrate superior numerical performance in simulations and are applied to infer the gene networks in a target brain tissue by leveraging the gene expressions from multiple other brain tissues. A significant decrease in prediction errors and a significant increase in power for link detection are observed.

Abstract

This paper considers estimation and prediction of a high-dimensional linear regression in the setting of transfer learning where, in addition to observations from the target model, auxiliary samples from different but possibly related regression models are available. When the set of informative auxiliary studies is known, an estimator and a predictor are proposed and their optimality is established. The optimal rates of convergence for prediction and estimation are faster than the corresponding rates without using the auxiliary samples. This implies that knowledge from the informative auxiliary samples can be transferred to improve the learning performance of the target problem. When the set of informative auxiliary samples is unknown, we propose a data-driven procedure for transfer learning, called Trans-Lasso, and show its robustness to non-informative auxiliary samples and its efficiency in knowledge transfer. The proposed procedures are demonstrated in numerical studies and are applied to a dataset concerning the associations among gene expressions. It is shown that Trans-Lasso leads to improved performance in gene expression prediction in a target tissue by incorporating data from multiple different tissues as auxiliary samples.

Abstract

Human learners have the natural ability to use knowledge gained in one setting for learning in a different but related setting. This ability to transfer knowledge from one task to another is essential for effective learning. In this paper, we study transfer learning in the context of nonparametric classification based on observations from different distributions under the posterior drift model, which is a general framework and arises in many practical problems. We first establish the minimax rate of convergence and construct a rate-optimal two-sample weighted \(K\)-NN classifier. The results characterize precisely the contribution of the observations from the source distribution to the classification task under the target distribution. A data-driven adaptive classifier is then proposed and is shown to simultaneously attain within a logarithmic factor of the optimal rate over a large collection of parameter spaces. Simulation studies are carried out and the numerical results further illustrate the theoretical analysis. Extensions to the case of multiple source distributions are also considered.

Abstract

Estimating covariance matrices is fundamental to a wide range of statistical applications. This paper studies minimax and adaptive estimation of high-dimensional covariance matrices under \(\rho\)-zero-concentrated differential privacy (\(\rho\)-zCDP) over three nested classes: the pointwise-decay class \(H_\alpha\), the row-tail class \(G_\alpha\), and the separated-block class \(F_\alpha\). We consider both squared operator norm loss and normalized squared Frobenius norm loss. For \(H_\alpha\) and \(G_\alpha\), we develop center--outer dyadic estimators tailored to the refined geometry of the two classes, while for \(F_\alpha\), we develop a blockwise tridiagonal estimator. The resulting minimax-optimal rates reveal a nontrivial interplay among the smoothness \(\alpha\), the loss, the geometry of the covariance class, and the privacy constraint. In contrast to the non-private setting, privacy distinguishes covariance classes that share the same leading non-private rate and induces a polynomial dependence on the ambient dimension. We further develop procedures that adapt to the unknown decay parameter over all three covariance classes under both losses, at the cost of at most polylogarithmic factors. To establish minimax lower bounds, we develop a novel differentially private van Trees inequality that connects Fisher information with the \(\rho\)-zCDP constraint and may be useful for other private estimation problems. We also construct carefully designed prior distributions to obtain matching minimax lower bounds.

Abstract

Privacy-preserving data analysis has become a central challenge in modern statistics. At the same time, a long-standing goal in statistics is the development of adaptive procedures-- methods that achieve near-optimal performance across diverse function classes without prior knowledge of underlying smoothness or complexity. While adaptation is often achievable at no extra cost in the classical non-private setting, this naturally raises a fundamental question: to what extent is adaptation still possible under privacy constraints? We address this question in the context of density estimation under federated differential privacy (FDP), a framework that encompasses both central and local DP models. We establish sharp results that characterize the cost of adaptation under FDP for both global and pointwise estimation, revealing fundamental differences from the non-private case. We then propose an adaptive FDP estimator that achieves explicit performance guarantees by introducing a new noise mechanism, enabling one-shot adaptation via post-processing. This approach strictly improves upon existing adaptive DP methods. Finally, we develop new lower bound techniques that capture the limits of adaptive inference under privacy and may be of independent interest beyond this problem. Our findings reveal a sharp contrast between private and non-private settings. For global estimation, where adaptation can be achieved for free in the classical non-private setting, we prove that under FDP an intrinsic adaptation cost is unavoidable. For pointwise estimation, where a logarithmic penalty is already known to arise in the non-private setting, we show that FDP introduces an additional logarithmic factor, thereby compounding the cost of adaptation. Taken together, these results provide the first rigorous characterization of the adaptive privacy-accuracy trade-off.

Abstract

We study functional linear regression under differential privacy in both fully observed and discretely observed settings. For the fully observed case, we propose a private algorithm that integrates private covariance estimation with noisy gradient descent, and establish its convergence rate. We further derive a matching minimax lower bound, demonstrating that the proposed estimator is minimax-optimal up to a logarithmic factor. In the discrete setting, we analyze the trade-offs among privacy, statistical accuracy, and discretization error, revealing a notable phase transition phenomenon. Simulations and real-data applications confirm the effectiveness of our approach in preserving privacy while maintaining statistical accuracy. Differential privacy in functional linear regression presents unique challenges due to the infinite-dimensional structure of functional data, which complicates sensitivity analysis, noise calibration, and algorithm design. Standard privacy mechanisms developed for finite-dimensional settings are inadequate in this context. To address these difficulties, we use basis expansions to reduce the problem to high-dimensional linear regression and introduce a preconditioning strategy that ensures robustness, even under ill-conditioned covariance structures. This framework broadens the applicability of differential privacy to complex functional models and ill-posed settings, advancing privacy-preserving methods in infinite-dimensional data analysis.

Abstract

Federated Learning (FL) has gained significant recent attention in machine learning for its enhanced privacy and data security, making it indispensable in fields such as healthcare, finance, and personalized services. This paper investigates federated PCA and estimation for spiked covariance matrices under distributed differential privacy constraints. We establish minimax rates of convergence, with a key finding that the central server's optimal rate is the harmonic mean of the local clients' minimax rates. This guarantees consistent estimation at the central server as long as at least one local client provides consistent results. Notably, consistency is maintained even if some local estimators are inconsistent, provided there are enough clients. These findings highlight the robustness and scalability of FL for reliable statistical inference under privacy constraints. To establish minimax lower bounds, we derive a matrix version of van Trees' inequality, which is of independent interest. Furthermore, we propose an efficient algorithm that preserves differential privacy while achieving near-optimal rates at the central server, up to a logarithmic factor. We address significant technical challenges in analyzing this algorithm, which involves a three-layer spectral decomposition. Numerical performance of the proposed algorithm is investigated using both simulated and real data.

Abstract

Nonparametric goodness-of-fit testing in the white-noise-with-drift model under federated differential privacy (FDP) constraints is studied. In this framework, data is distributed across multiple locations, with each submitting a differentially private summary to a central server--encompassing both local and central differential privacy. We quantify the cost of privacy in the federated setting by establishing matching lower and upper bounds, up to a logarithmic factor, on the minimax separation rate. This optimal rate benchmarks the difficulty of the testing problem, factoring in model characteristics such as the number of observations, noise level, and regularity of the signal class, along with the strictness of the \((\varepsilon,\delta)\)-DP requirement and the degree to which the data is distributed. Our results demonstrate interesting and novel phase transition phenomena: where the cost of DP is minimal for testing in more centralized settings, it can significantly affect distributed scenarios. Furthermore, it is revealed that distributed one-shot protocols with access to shared randomness outperform those without access to shared randomness. We also construct a data-driven testing procedure that can adapt to an unknown regularity parameter over a large collection of function classes with minimal additional cost, while adhering to the same set of DP constraints.

Abstract

Estimating a covariance matrix and its associated principal components is a fundamental problem in contemporary statistics. While optimal estimation procedures have been developed with well-understood properties, the increasing demand for privacy preservation introduces new complexities to this classical problem. In this paper, we study optimal differentially private Principal Component Analysis (PCA) and covariance estimation within the spiked covariance model. We precisely characterize the sensitivity of eigenvalues and eigenvectors under this model and establish the minimax rates of convergence for estimating both the principal components and covariance matrix. These rates hold up to logarithmic factors and encompass general Schatten norms, including spectral norm, Frobenius norm, and nuclear norm as special cases. We introduce computationally efficient differentially private estimators and prove their minimax optimality, up to logarithmic factors. Additionally, matching minimax lower bounds are established. Notably, in comparison with existing literature, our results accommodate a diverging rank, necessitate no eigengap condition between distinct principal components, and remain valid even if the sample size is much smaller than the dimension.

Abstract

Achieving optimal statistical performance while ensuring the privacy of personal data is a challenging yet crucial objective in modern data analysis. However, characterizing the optimality, particularly the minimax lower bound, under privacy constraints is technically difficult. To address this issue, we propose a novel approach called the score attack, which provides a lower bound on the differential-privacy-constrained minimax risk of parameter estimation. The score attack method is based on the tracing attack concept in differential privacy and can be applied to any statistical model with a well-defined score statistic. It can optimally lower bound the minimax risk of estimating unknown model parameters, up to a logarithmic factor, while ensuring differential privacy for a range of statistical problems. We demonstrate the effectiveness and optimality of this general method in various examples, such as the generalized linear model in both classical and high-dimensional sparse settings, the Bradley-Terry-Luce model for pairwise comparisons, and nonparametric regression over the Sobolev class.

Abstract

This paper studies federated learning for nonparametric regression in the context of distributed samples across different servers, each adhering to distinct differential privacy constraints. The setting we consider is heterogeneous, encompassing both varying sample sizes and differential privacy constraints across servers. Within this framework, both global and pointwise estimation are considered, and optimal rates of convergence over the Besov spaces are established. Distributed privacy-preserving estimators are proposed and their risk properties are investigated. Matching minimax lower bounds, up to a logarithmic factor, are established for both global and pointwise estimation. Together, these findings shed light on the tradeoff between statistical accuracy and privacy preservation. In particular, we characterize the compromise not only in terms of the privacy budget but also concerning the loss incurred by distributing data within the privacy framework as a whole. This insight captures the folklore wisdom that it is easier to retain privacy in larger samples, and explores the differences between pointwise and global estimation under distributed privacy constraints.

Abstract

Data privacy is a central concern in many applications involving ranking from incomplete and noisy pairwise comparisons, such as recommendation systems, educational assessments, and opinion surveys on sensitive topics. In this work, we propose differentially private algorithms for ranking based on pairwise comparisons. Specifically, we develop and analyze ranking methods under two privacy notions: edge differential privacy, which protects the confidentiality of individual comparison outcomes, and individual differential privacy, which safeguards potentially many comparisons contributed by a single individual. Our algorithms-– including a perturbed maximum likelihood estimator and a noisy count-based method–-are shown to achieve minimax optimal rates of convergence under the respective privacy constraints. We further demonstrate the practical effectiveness of our methods through experiments on both simulated and real-world data.

Abstract

This paper considers minimax and adaptive transfer learning for nonparametric classification under the posterior drift model with distributed differential privacy constraints. Our study is conducted within a heterogeneous framework, encompassing diverse sample sizes, varying privacy parameters, and data heterogeneity across different servers. We first establish the minimax misclassification rate, precisely characterizing the effects of privacy constraints, source samples, and target samples on classification accuracy. The results reveal interesting phase transition phenomena and highlight the intricate trade-offs between preserving privacy and achieving classification accuracy. We then develop a data-driven adaptive classifier that achieves the optimal rate within a logarithmic factor across a large collection of parameter spaces while satisfying the same set of differential privacy constraints. Simulation studies and real-world data applications further elucidate the theoretical analysis with numerical results.

Abstract

Privacy-preserving data analysis is a rising challenge in contemporary statistics, as the privacy guarantees of statistical methods are often achieved at the expense of accuracy. In this paper, we investigate the tradeoff between statistical accuracy and privacy in mean estimation and linear regression, under both the classical low- dimensional and modern high-dimensional settings. A primary focus is to establish minimax optimality for statistical estimation with the \((\epsilon,\delta)\)-differential privacy constraint. By refining the “tracing adversary” technique for lower bounds in the theoretical computer science literature, we improve existing minimax lower bound for low-dimensional mean estimation and establish new lower bounds for high-dimensional mean estimation and linear regression problems. We also design differentially private algorithms that attain the minimax lower bounds up to logarithmic factors. In particular, for high-dimensional linear regression, a novel private iterative hard thresholding algorithm is proposed. The numerical performance of differentially private algorithms is demonstrated by simulation studies and applications to real data sets.

Abstract

Privacy-preserving data analysis has become a central challenge in modern statistics. At the same time, a long-standing goal in statistics is the development of adaptive procedures-- methods that achieve near-optimal performance across diverse function classes without prior knowledge of underlying smoothness or complexity. While adaptation is often achievable at no extra cost in the classical non-private setting, this naturally raises a fundamental question: to what extent is adaptation still possible under privacy constraints? We address this question in the context of density estimation under federated differential privacy (FDP), a framework that encompasses both central and local DP models. We establish sharp results that characterize the cost of adaptation under FDP for both global and pointwise estimation, revealing fundamental differences from the non-private case. We then propose an adaptive FDP estimator that achieves explicit performance guarantees by introducing a new noise mechanism, enabling one-shot adaptation via post-processing. This approach strictly improves upon existing adaptive DP methods. Finally, we develop new lower bound techniques that capture the limits of adaptive inference under privacy and may be of independent interest beyond this problem. Our findings reveal a sharp contrast between private and non-private settings. For global estimation, where adaptation can be achieved for free in the classical non-private setting, we prove that under FDP an intrinsic adaptation cost is unavoidable. For pointwise estimation, where a logarithmic penalty is already known to arise in the non-private setting, we show that FDP introduces an additional logarithmic factor, thereby compounding the cost of adaptation. Taken together, these results provide the first rigorous characterization of the adaptive privacy-accuracy trade-off.

Abstract

Federated Learning (FL) has gained significant recent attention in machine learning for its enhanced privacy and data security, making it indispensable in fields such as healthcare, finance, and personalized services. This paper investigates federated PCA and estimation for spiked covariance matrices under distributed differential privacy constraints. We establish minimax rates of convergence, with a key finding that the central server's optimal rate is the harmonic mean of the local clients' minimax rates. This guarantees consistent estimation at the central server as long as at least one local client provides consistent results. Notably, consistency is maintained even if some local estimators are inconsistent, provided there are enough clients. These findings highlight the robustness and scalability of FL for reliable statistical inference under privacy constraints. To establish minimax lower bounds, we derive a matrix version of van Trees' inequality, which is of independent interest. Furthermore, we propose an efficient algorithm that preserves differential privacy while achieving near-optimal rates at the central server, up to a logarithmic factor. We address significant technical challenges in analyzing this algorithm, which involves a three-layer spectral decomposition. Numerical performance of the proposed algorithm is investigated using both simulated and real data.

Abstract

Nonparametric goodness-of-fit testing in the white-noise-with-drift model under federated differential privacy (FDP) constraints is studied. In this framework, data is distributed across multiple locations, with each submitting a differentially private summary to a central server--encompassing both local and central differential privacy. We quantify the cost of privacy in the federated setting by establishing matching lower and upper bounds, up to a logarithmic factor, on the minimax separation rate. This optimal rate benchmarks the difficulty of the testing problem, factoring in model characteristics such as the number of observations, noise level, and regularity of the signal class, along with the strictness of the \((\varepsilon,\delta)\)-DP requirement and the degree to which the data is distributed. Our results demonstrate interesting and novel phase transition phenomena: where the cost of DP is minimal for testing in more centralized settings, it can significantly affect distributed scenarios. Furthermore, it is revealed that distributed one-shot protocols with access to shared randomness outperform those without access to shared randomness. We also construct a data-driven testing procedure that can adapt to an unknown regularity parameter over a large collection of function classes with minimal additional cost, while adhering to the same set of DP constraints.

Abstract

Federated learning (FL) is a distributed machine learning technique designed to preserve data privacy and security, and it has gained significant importance due to its broad range of applications. This paper addresses the problem of optimal functional mean estimation from discretely sampled data in a federated setting. We consider a heterogeneous framework where the number of individuals, measurements per individual, and privacy parameters vary across one or more servers, under both common and independent design settings. In the common design setting, the same design points are measured for each individual, whereas in the independent design setting, each individual has their own random collection of design points. Within this framework, we establish minimax upper and lower bounds for the estimation error of the underlying mean function, highlighting the differences between common and independent designs under distributed privacy constraints. We propose algorithms that achieve the optimal trade-off between privacy and accuracy and provide optimality results that quantify the fundamental limits of private functional mean estimation. We further support the theory with simulations and a real-data illustration using BMI trajectories from the Health and Retirement Study.

Abstract

This paper studies federated learning for nonparametric regression in the context of distributed samples across different servers, each adhering to distinct differential privacy constraints. The setting we consider is heterogeneous, encompassing both varying sample sizes and differential privacy constraints across servers. Within this framework, both global and pointwise estimation are considered, and optimal rates of convergence over the Besov spaces are established. Distributed privacy-preserving estimators are proposed and their risk properties are investigated. Matching minimax lower bounds, up to a logarithmic factor, are established for both global and pointwise estimation. Together, these findings shed light on the tradeoff between statistical accuracy and privacy preservation. In particular, we characterize the compromise not only in terms of the privacy budget but also concerning the loss incurred by distributing data within the privacy framework as a whole. This insight captures the folklore wisdom that it is easier to retain privacy in larger samples, and explores the differences between pointwise and global estimation under distributed privacy constraints.

Abstract

We consider the problem of estimating the correlation of two random variables X and Y, where the pairs (X, Y) are not observed together, but are instead separated coordinate-wise at two servers: server 1 contains all the X observations, and server 2 contains the corresponding Y observations. In this vertically distributed setting, we assume that each server has its own privacy constraints, owing to which they can only share suitably privatized statistics of their own component observations. We consider differing privacy budgets \((\epsilon_1,\delta_1)\) and \((\epsilon_2,\delta_2)\) for the two servers and determine the minimax optimal rates for correlation estimation allowing for both non-interactive and interactive mechanisms. We also provide correlation estimators that achieve these rates and further develop inference procedures, namely, confidence intervals, for the estimated correlations. Our results are characterized by an interesting rate in terms of the sample size \(n\), \(\epsilon_1\), and \(\epsilon_2\), which is strictly slower than the usual central privacy estimation rates. More interestingly, we find that the interactive mechanism is always better than its non-interactive counterpart whenever the two privacy budgets are different. Results from extensive numerical experiments support our theoretical findings.

Abstract

Distributed estimation of a Gaussian mean under communication constraints is studied in a decision theoretical framework. Minimax rates of convergence, which characterize the tradeoff between communication costs and statistical accuracy, are established under the independent protocols. Communication-efficient and statistically optimal procedures are developed. In the univariate case, the optimal rate depends only on the total communication budget, so long as each local machine has at least one bit. However, in the multivariate case, the minimax rate depends on the specific allocations of the communication budgets among the local machines. Although optimal estimation of a Gaussian mean is relatively simple in the conventional setting, it is quite involved under communication constraints, both in terms of the optimal procedure design and the lower bound argument. An essential step is the decomposition of the minimax estimation problem into two stages, localization and refinement. This critical decomposition provides a framework for both the lower bound analysis and optimal procedure design. The optimality results and techniques developed in the present paper can be useful for solving other problems such as distributed nonparametric function estimation and sparse signal recovery.

Abstract

Distributed estimation of a Gaussian mean with unknown variance under communication constraints is studied. Necessary and sufficient communication costs under different types of distributed protocols are derived for any estimator that is adaptively rate-optimal over a range of possible values for the variance. Communication-efficient and statistically optimal procedures are developed. The analysis reveals an interesting and important distinction among different types of distributed protocols: compared to the independent protocols, interactive protocols such as the sequential and blackboard protocols require less communication costs for rate-optimal adaptive Gaussian mean estimation. The lower bound techniques developed in the present paper are novel and can be of independent interest.

Abstract

Distributed minimax estimation and distributed adaptive estimation under communication constraints for Gaussian sequence model and white noise model are studied. The minimax rate of convergence for distributed estimation over a given Besov class, which serves as a benchmark for the cost of adaptation, is established. We then quantify the exact communication cost for adaptation and construct an optimally adaptive procedure for distributed estimation over a range of Besov classes. The results demonstrate significant differences between nonparametric function estimation in the distributed setting and the conventional centralized setting. For global estimation, adaptation in general cannot be achieved for free in the distributed setting. The new technical tools to obtain the exact characterization for the cost of adaptation can be of independent interest.

Abstract

This paper investigates the fundamental limits for detecting a high-dimensional sparse matrix contaminated by white Gaussian noise from both the statistical and computational perspectives. We consider \(p\times p\) matrices whose rows and columns are individually \(k\)-sparse. We provide a tight characterization of the statistical and computational limits for sparse matrix detection, which precisely describe when achieving optimal detection is easy, hard, or impossible, respectively. Although the sparse matrices considered in this paper have no apparent submatrix structure and the corresponding estimation problem has no computational issue at all, the detection problem has a surprising computational barrier when the sparsity level \(k\) exceeds the cubic root of the matrix size \(p\): attaining the optimal detection boundary is computationally at least as hard as solving the planted clique problem. The same statistical and computational limits also hold in the sparse covariance matrix model, where each variable is correlated with at most \(k\) others. A key step in the construction of the statistically optimal test is a structural property for sparse matrices, which can be of independent interest.

Abstract

The interplay between computational efficiency and statistical accuracy in high-dimensional inference has drawn increasing attention in the literature. In this paper, we study computational and statistical boundaries for submatrix localization. Given one observation of a signal submatrix (of magnitude \(\lambda\) and size \(k_m\times k_n\)) contaminated with a noise matrix (of size \(m\times n\)), we establish two transition thresholds for the signal to noise ratio \(\lambda/\sigma\) in terms of \(m,n,k_m\), and \(k_n\). The first threshold, \(\mathrm{SNR}_c\), corresponds to the computational boundary. Below this threshold, it is shown that no polynomial time algorithm can succeed in identifying the submatrix, under the hidden clique hypothesis. We introduce an adaptive linear time algorithm that identifies the submatrix with high probability when the signal strength is above the threshold \(\mathrm{SNR}_c\). The second threshold, \(\mathrm{SNR}_s\), captures the statistical boundary, below which no method can succeed with probability going to one in the minimax sense. The exhaustive search method successfully finds the submatrix above this threshold. The results show an interesting phenomenon that \(\mathrm{SNR}_c\) is always significantly larger than \(\mathrm{SNR}_s\), which implies an essential gap between statistical optimality and computational efficiency for submatrix localization.

Abstract

Federated Learning (FL) has gained significant recent attention in machine learning for its enhanced privacy and data security, making it indispensable in fields such as healthcare, finance, and personalized services. This paper investigates federated PCA and estimation for spiked covariance matrices under distributed differential privacy constraints. We establish minimax rates of convergence, with a key finding that the central server's optimal rate is the harmonic mean of the local clients' minimax rates. This guarantees consistent estimation at the central server as long as at least one local client provides consistent results. Notably, consistency is maintained even if some local estimators are inconsistent, provided there are enough clients. These findings highlight the robustness and scalability of FL for reliable statistical inference under privacy constraints. To establish minimax lower bounds, we derive a matrix version of van Trees' inequality, which is of independent interest. Furthermore, we propose an efficient algorithm that preserves differential privacy while achieving near-optimal rates at the central server, up to a logarithmic factor. We address significant technical challenges in analyzing this algorithm, which involves a three-layer spectral decomposition. Numerical performance of the proposed algorithm is investigated using both simulated and real data.

Abstract

Estimating a covariance matrix and its associated principal components is a fundamental problem in contemporary statistics. While optimal estimation procedures have been developed with well-understood properties, the increasing demand for privacy preservation introduces new complexities to this classical problem. In this paper, we study optimal differentially private Principal Component Analysis (PCA) and covariance estimation within the spiked covariance model. We precisely characterize the sensitivity of eigenvalues and eigenvectors under this model and establish the minimax rates of convergence for estimating both the principal components and covariance matrix. These rates hold up to logarithmic factors and encompass general Schatten norms, including spectral norm, Frobenius norm, and nuclear norm as special cases. We introduce computationally efficient differentially private estimators and prove their minimax optimality, up to logarithmic factors. Additionally, matching minimax lower bounds are established. Notably, in comparison with existing literature, our results accommodate a diverging rank, necessitate no eigengap condition between distinct principal components, and remain valid even if the sample size is much smaller than the dimension.

Abstract

A general framework for principal component analysis (PCA) in the presence of heteroskedastic noise is introduced. We propose an algorithm called HeteroPCA, which involves iteratively imputing the diagonal entries of the sample covariance matrix to remove estimation bias due to heteroskedasticity. This procedure is computationally efficient and provably optimal under the generalized spiked covariance model. A key technical step is a deterministic robust perturbation analysis on singular subspaces, which can be of independent interest. The effectiveness of the proposed algorithm is demonstrated in a suite of problems in high-dimensional statistics, including singular value decomposition (SVD) under heteroskedastic noise, Poisson PCA, and SVD for heteroskedastic and incomplete data.

Abstract

Driven by a wide range of applications, many principal subspace estimation problems have been studied individually under different structural constraints. This paper presents a unified framework for the statistical analysis of a general structured principal subspace estimation problem which includes as special cases non-negative PCA/SVD, sparse PCA/SVD, subspace constrained PCA/SVD, and spectral clustering. General minimax lower and upper bounds are established to characterize the interplay between the information-geometric complexity of the structural set for the principal subspaces, the signal-to-noise ratio (SNR), and the dimensionality. The results yield interesting phase transition phenomena concerning the rates of convergence as a function of the SNRs and the fundamental limit for consistent estimation. Applying the general results to the specific settings yields the minimax rates of convergence for those problems, including the previous unknown optimal rates for non-negative PCA/SVD, sparse SVD and subspace constrained PCA/SVD.

Abstract

We study the asymptotic distributions of the spiked eigenvalues and the largest nonspiked eigenvalue of the sample covariance matrix under a general covariance matrix model with divergent spiked eigenvalues, while the other eigenvalues are bounded but otherwise arbitrary. The limiting normal distribution for the spiked sample eigenvalues is established. It has distinct features that the asymptotic mean relies on not only the population spikes but also the nonspikes and that the asymptotic variance in general depends on the population eigenvectors. In addition, the limiting Tracy-Widom law for the largest nonspiked sample eigenvalue is obtained. Estimation of the number of spikes and the convergence of the leading eigenvectors are also considered. The results hold even when the number of the spikes diverges. As a key technical tool, we develop a Central Limit Theorem for a type of random quadratic forms where the random vectors and random matrices involved are dependent. This result can be of independent interest.

Abstract

Perturbation bounds for singular spaces, in particular Wedin's \(\sin\Theta\) theorem, are a fundamental tool in many fields including high-dimensional statistics, machine learning, and applied mathematics. In this paper, we establish separate perturbation bounds, measured in both spectral and Frobenius \(\sin\Theta\) distances, for the left and right singular subspaces. Lower bounds, which show that the individual perturbation bounds are rate-optimal, are also given. The new perturbation bounds are applicable to a wide range of problems. In this paper, we consider in detail applications to low-rank matrix denoising and singular space estimation, high-dimensional clustering, and canonical correlation analysis (CCA). In particular, separate matching upper and lower bounds are obtained for estimating the left and right singular spaces. To the best of our knowledge, this is the first result that gives different optimal rates for the left and right singular spaces under the same perturbation. In addition to these problems, applications to other high-dimensional problems such as community detection in bipartite networks, multidimensional scaling, and cross-covariance matrix estimation are also discussed.

Abstract

This paper considers the noisy sparse phase retrieval problem: recovering a sparse signal \(\mathbf{x}\in\mathbb{R}^p\) from noisy quadratic measurements \(y_j=(\mathbf{a}_j^\top\mathbf{x})^2+\epsilon_j,\ j=1,\ldots,m\), with independent sub-exponential noise \(\epsilon_j\). The goals are to understand the effect of the sparsity of \(\mathbf{x}\) on the estimation precision and to construct a computationally feasible estimator to achieve the optimal rates. Inspired by the Wirtinger Flow method for noiseless phase retrieval proposed and analyzed in [12], a novel thresholded gradient descent algorithm is proposed and it is shown to adaptively achieve the minimax rates of convergence over a wide range of sparsity levels when the \(\mathbf{a}_j\)'s are independent standard Gaussian random vectors, provided that the sample size is sufficiently large compared to the sparsity of \(\mathbf{x}\).

Abstract

Matrix completion has attracted significant recent attention in many fields including statistics, applied mathematics and electrical engineering. Current literature on matrix completion focuses primarily on independent sampling models under which the individual observed entries are sampled independently. Motivated by applications in genomic data integration, we propose a new framework of structured matrix completion (SMC) to treat structured missingness by design. Specifically, our proposed method aims at efficient matrix recovery when a subset of the rows and columns of an approximately low-rank matrix are observed. We provide theoretical justification for the proposed SMC method and derive lower bound for the estimation errors, which together establish the optimal rate of recovery over certain classes of approximately low-rank matrices. Simulation studies show that the method performs well in finite sample under a variety of configurations. The method is applied to integrate several ovarian cancer genomic studies with different extent of genomic measurements, which enables us to construct more accurate prediction rules for ovarian cancer survival.

Abstract

This paper presents a unified geometric framework for the statistical analysis of a general ill-posed linear inverse model which includes as special cases noisy compressed sensing, sign vector recovery, trace regression, orthogonal matrix estimation, and noisy matrix completion. We propose computationally feasible convex programs for statistical inference including estimation, confidence intervals and hypothesis testing. A theoretical framework is developed to characterize the local estimation rate of convergence and to provide statistical inference guarantees. Our results are built based on the local conic geometry and duality. The difficulty of statistical inference is captured by the geometric characterization of the local tangent cone through the Gaussian width and Sudakov estimate.

Abstract

This paper studies matrix completion under a general sampling model using the max-norm as a convex relaxation for the rank of the matrix. The optimal rate of convergence is established for the Frobenius norm loss. It is shown that the max-norm constrained minimization method is rate-optimal and it yields a more stable approximate recovery guarantee, with respect to the sampling distributions, than previously used trace-norm based approaches. The computational effectiveness of this method is also studied, based on a first-order algorithm for solving convex programs involving a max-norm constraint.

Abstract

This paper considers sparse spiked covariance matrix models in the high-dimensional setting and studies the minimax estimation of the covariance matrix and the principal subspace as well as the minimax rank detection. The optimal rate of convergence for estimating the spiked covariance matrix under the spectral norm is established, which requires significantly different techniques from those for estimating other structured covariance matrices such as bandable or sparse covariance matrices. We also establish the minimax rate under the spectral norm for estimating the principal subspace, the primary object of interest in principal component analysis. In addition, the optimal rate for the rank detection boundary is obtained. This result also resolves a gap between the upper and lower bounds in a recent paper by Berthet and Rigollet [1] where the special case of rank one is considered.

Abstract

Estimation of low-rank matrices is of significant interest in a range of contemporary applications. In this paper, we introduce a rank-one projection model for low-rank matrix recovery and propose a constrained nuclear norm minimization method for stable recovery of low-rank matrices in the noisy case. The procedure is adaptive to the rank and robust against small low-rank perturbations. Both upper and lower bounds for the estimation accuracy under the Frobenius norm loss are obtained. The proposed estimator is shown to be rate-optimal under certain conditions. The estimator is easy to implement via convex programming and performs well numerically. The main results obtained in the paper also have implications to other related statistical problems. An application to estimation of spike covariance matrices from one-dimensional random projections is considered. The results demonstrate that it is possible to accurately estimate the covariance matrix of a high-dimensional distribution based only on one-dimensional projections.

Abstract

This paper considers compressed sensing and affine rank minimization in both noiseless and noisy cases and establishes sharp restricted isometry conditions for sparse signal and low-rank matrix recovery. The analysis relies on a key technical tool which represents points in a polytope by convex combinations of sparse vectors. The technique is elementary while leads to sharp results. It is shown that for any given constant \(t\ge4/3\), in compressed sensing \(\delta^A_{tk}<\sqrt{(t-1)/t}\) guarantees the exact recovery of all \(k\)-sparse signals in the noiseless case through the constrained \(\ell_1\) minimization, and similarly in affine rank minimization \(\delta^M_{tk}<\sqrt{(t-1)/t}\) ensures the exact reconstruction of all matrices with rank at most \(r\) in the noiseless case via the constrained nuclear norm minimization. Moreover, for any \(\epsilon>0\), \(\delta^A_{tk}<\sqrt{(t-1)/t}+\epsilon\) is not sufficient to guarantee the exact recovery of all \(k\)-sparse signals for large \(k\). Similar result also holds for matrix recovery. In addition, the conditions \(\delta^A_{tk}<\sqrt{(t-1)/t}\) and \(\delta^M_{tr}<\sqrt{(t-1)/t}\) are also shown to be sufficient respectively for stable recovery of approximately sparse signals and low-rank matrices in the noisy case.

Abstract

Principal component analysis (PCA) is one of the most commonly used statistical procedures with a wide range of applications.This paper considers both minimax and adaptive estimation of the principal subspace in the high dimensional setting. Under mild technical conditions, we first establish the optimal rates of convergence for estimating the principal subspace which are sharp with respect to all the parameters, thus providing a complete characterization of the difficulty of the estimation problem in term of the convergence rate. The lower bound is obtained by calculating the local metric entropy and an application of Fano's Lemma. The rate optimal estimator is constructed using aggregation, which, however, might not be computationally feasible. We then introduce an adaptive procedure for estimating the principal subspace which is fully data driven and can be computed efficiently. It is shown that the estimator attains the optimal rates of convergence simultaneously over a large collection of the parameter spaces. A key idea in our construction is a reduction scheme which reduces the sparse PCA problem to a high-dimensional multivariate regression problem. This method is potentially also useful for other related problems.

Abstract

We consider in this paper the problem of noisy 1-bit matrix completion under a general non-uniform sampling distribution using the max-norm as a convex relaxation for the rank. A max-norm constrained maximum likelihood estimate is introduced and studied. The rate of convergence for the estimate is obtained. Information-theoretical methods are used to establish a minimax lower bound under the general sampling model. The minimax upper and lower bounds together yield the optimal rate of convergence for the Frobenius norm loss. Computational algorithms and numerical performance are also discussed.

Abstract

This paper establish new restricted isometry conditions for compressed sensing and affine rank minimization. It is shown for compressed sensing that \(\delta^A_k+\theta^A_{k,k}<1\) guarantees the exact recovery of all \(k\)-sparse signals in the noiseless case through the constrained \(\ell_1\) minimization. Furthermore, the upper bound 1 is sharp in the sense that for any \(\epsilon>0\), the condition \(\delta^A_k+\theta^A_{k,k}<1+\epsilon\) is not sufficient to guarantee such exact recovery using any recovery method. Similarly, for affine rank minimization, if \(\delta^A_k+\theta^A_{k,k}<1\) then all matrices with rank at most \(r\) can be reconstructed exactly in the noiseless case via the constrained nuclear norm minimization; and for any \(\epsilon>0\), \(\delta^A_k+\theta^A_{k,k}<1+\epsilon\) does not ensure such exact recovery using any method. Moreover, in the noisy case the conditions \(\delta^A_k+\theta^A_{k,k}<1\) and \(\delta^A_k+\theta^A_{k,k}<1\) are also sufficient for the stable recovery of sparse signals and low-rank matrices respectively. Applications and extensions are also discussed.

Abstract

This paper establishes a sharp condition on the restrict isometry property (RIP) for both the sparse signal recovery and low-rank matrix recovery. It is shown that if the measurement matrix \(A\) satisfies the RIP condition \(\delta^A_k<1/3\), then all \(k\)-sparse signals \(\beta\) can be recovered exactly via the constrained \(\ell_1\) minimization based on \(y=A\beta\). Similarly, if the linear map \(M\) satisfies the RIP condition \(\delta^M_r<1/3\), then all matrices \(X\) of rank at most \(r\) can be recovered exactly via the constrained nuclear norm minimization based on \(b=M(X)\). Furthermore, in both cases it is not possible to do so in general when the condition does not hold. In addition, noisy cases are considered and oracle inequalities are given under the sharp RIP condition.

Abstract

Estimating covariance matrices is fundamental to a wide range of statistical applications. This paper studies minimax and adaptive estimation of high-dimensional covariance matrices under \(\rho\)-zero-concentrated differential privacy (\(\rho\)-zCDP) over three nested classes: the pointwise-decay class \(H_\alpha\), the row-tail class \(G_\alpha\), and the separated-block class \(F_\alpha\). We consider both squared operator norm loss and normalized squared Frobenius norm loss. For \(H_\alpha\) and \(G_\alpha\), we develop center--outer dyadic estimators tailored to the refined geometry of the two classes, while for \(F_\alpha\), we develop a blockwise tridiagonal estimator. The resulting minimax-optimal rates reveal a nontrivial interplay among the smoothness \(\alpha\), the loss, the geometry of the covariance class, and the privacy constraint. In contrast to the non-private setting, privacy distinguishes covariance classes that share the same leading non-private rate and induces a polynomial dependence on the ambient dimension. We further develop procedures that adapt to the unknown decay parameter over all three covariance classes under both losses, at the cost of at most polylogarithmic factors. To establish minimax lower bounds, we develop a novel differentially private van Trees inequality that connects Fisher information with the \(\rho\)-zCDP constraint and may be useful for other private estimation problems. We also construct carefully designed prior distributions to obtain matching minimax lower bounds.

Abstract

The covariance function is a central object in functional and longitudinal data analysis, yet its optimal estimation remains challenging because of its bivariate structure and its sensitivity to the sampling design. We study covariance function estimation for functional data under two distinct designs: the common design, in which all curves are observed on a shared grid, and the independent design, in which the sampling locations vary across curves. By exploiting the tensor-product structure of the covariance function, we propose a unified smoothing-spline estimator that applies to both settings. We establish minimax optimal convergence rates for the proposed estimator under both sampling designs and identify several phase transition phenomena. In particular, under the independent design with sparse sampling, the estimator attains a nearly “univariate” minimax rate of \((nm)^{-2\alpha/(2\alpha+1)}\), up to a logarithmic factor, despite the intrinsically bivariate nature of the covariance function. By contrast, under the common design, the sparse-regime rate is dominated by the discretization error \(m^{-2\alpha}\). We further develop an adaptive implementation based on cross-validation for selecting the smoothing parameter, which does not require prior knowledge of the smoothness parameter \(\alpha\) or the sampling-frequency regime. Numerical experiments on both simulated and real datasets support the theoretical findings and demonstrate the practical effectiveness of the proposed method.

Abstract

Transfer learning, which leverages knowledge from an auxiliary source dataset to improve performance in a primary target domain, has emerged as a pivotal machine learning technique. In this paper, we consider minimax and adaptive estimation of large bandable covariance matrices within the transfer learning framework. We first establish the minimax rate of convergence under the spectral norm and propose a rate-optimal estimation procedure. Our findings reveal intriguing phase transition phenomena that highlight the effectiveness of transfer learning and the use of source samples. We then address the problem of adaptation, establishing the adaptive rate of convergence up to a logarithmic factor. Our results demonstrate that, in sharp contrast to conventional settings, the cost of adaptation in transfer learning can be substantial in certain cases. We propose a novel data-driven algorithm that dynamically adapts to unknown model parameters. These theoretical insights are further validated by a simulation study, demonstrating the practicality and efficiency of the proposed adaptive algorithm.

Abstract

A general framework for principal component analysis (PCA) in the presence of heteroskedastic noise is introduced. We propose an algorithm called HeteroPCA, which involves iteratively imputing the diagonal entries of the sample covariance matrix to remove estimation bias due to heteroskedasticity. This procedure is computationally efficient and provably optimal under the generalized spiked covariance model. A key technical step is a deterministic robust perturbation analysis on singular subspaces, which can be of independent interest. The effectiveness of the proposed algorithm is demonstrated in a suite of problems in high-dimensional statistics, including singular value decomposition (SVD) under heteroskedastic noise, Poisson PCA, and SVD for heteroskedastic and incomplete data.

Abstract

Driven by a wide range of applications, many principal subspace estimation problems have been studied individually under different structural constraints. This paper presents a unified framework for the statistical analysis of a general structured principal subspace estimation problem which includes as special cases non-negative PCA/SVD, sparse PCA/SVD, subspace constrained PCA/SVD, and spectral clustering. General minimax lower and upper bounds are established to characterize the interplay between the information-geometric complexity of the structural set for the principal subspaces, the signal-to-noise ratio (SNR), and the dimensionality. The results yield interesting phase transition phenomena concerning the rates of convergence as a function of the SNRs and the fundamental limit for consistent estimation. Applying the general results to the specific settings yields the minimax rates of convergence for those problems, including the previous unknown optimal rates for non-negative PCA/SVD, sparse SVD and subspace constrained PCA/SVD.

Abstract

We study the asymptotic distributions of the spiked eigenvalues and the largest nonspiked eigenvalue of the sample covariance matrix under a general covariance matrix model with divergent spiked eigenvalues, while the other eigenvalues are bounded but otherwise arbitrary. The limiting normal distribution for the spiked sample eigenvalues is established. It has distinct features that the asymptotic mean relies on not only the population spikes but also the nonspikes and that the asymptotic variance in general depends on the population eigenvectors. In addition, the limiting Tracy-Widom law for the largest nonspiked sample eigenvalue is obtained. Estimation of the number of spikes and the convergence of the leading eigenvectors are also considered. The results hold even when the number of the spikes diverges. As a key technical tool, we develop a Central Limit Theorem for a type of random quadratic forms where the random vectors and random matrices involved are dependent. This result can be of independent interest.

Abstract

This paper studies the estimation of high-dimensional minimum variance portfolio (MVP) based on high frequency returns which can exhibit heteroskedasticity and possibly be contaminated by microstructure noise. Under certain sparsity assumptions on the precision matrix, we propose an estimator of MVP and prove that our portfolio asymptotically achieves the minimum variance in a sharp sense. In addition, we introduce consistent estimators of the minimum variance, which provide reference targets. Simulation and empirical studies demonstrate that our proposed portfolio performs favorably

Abstract

Making accurate inference for gene regulatory networks, including inferring about pathway by pathway interactions, is an important and difficult task. Motivated by such genomic applications, we consider multiple testing for conditional dependence between subgroups of variables. Under a Gaussian graphical model framework, the problem is translated into simultaneous testing for a collection of submatrices of a high-dimensional precision matrix with each submatrix summarizing the dependence structure between two subgroups of variables. A novel multiple testing procedure is proposed and both theoretical and numerical properties of the procedure are investigated. Asymptotic null distribution of the test statistic for an individual hypothesis is established and the proposed multiple testing procedure is shown to asymptotically control the false discovery rate (FDR) and false discovery proportion (FDP) at the pre-specified level under regularity conditions. Simulations show that the procedure works well in controlling the FDR and has good power in detecting the true interactions. The procedure is applied to a breast cancer gene expression study to identify between pathway interactions.

Abstract

Perturbation bounds for singular spaces, in particular Wedin's \(\sin\Theta\) theorem, are a fundamental tool in many fields including high-dimensional statistics, machine learning, and applied mathematics. In this paper, we establish separate perturbation bounds, measured in both spectral and Frobenius \(\sin\Theta\) distances, for the left and right singular subspaces. Lower bounds, which show that the individual perturbation bounds are rate-optimal, are also given. The new perturbation bounds are applicable to a wide range of problems. In this paper, we consider in detail applications to low-rank matrix denoising and singular space estimation, high-dimensional clustering, and canonical correlation analysis (CCA). In particular, separate matching upper and lower bounds are obtained for estimating the left and right singular spaces. To the best of our knowledge, this is the first result that gives different optimal rates for the left and right singular spaces under the same perturbation. In addition to these problems, applications to other high-dimensional problems such as community detection in bipartite networks, multidimensional scaling, and cross-covariance matrix estimation are also discussed.

Abstract

Missing data occur frequently in a wide range of applications. In this paper, we consider estimation of high-dimensional covariance matrices in the presence of missing observations under a general missing completely at random model in the sense that the missingness is not dependent on the values of the data. Based on incomplete data, estimators for bandable and sparse covariance matrices are proposed and their theoretical and numerical properties are investigated. Minimax rates of convergence are established under the spectral norm loss and the proposed estimators are shown to be rate-optimal under mild regularity conditions. Simulation studies demonstrate that the estimators perform well numerically. The methods are also illustrated through an application to data from four ovarian cancer studies. The key technical tools developed in this paper are of independent interest and potentially useful for a range of related problems in high-dimensional statistical inference with missing data.

Abstract

Covariance structure plays an important role in high dimensional statistical inference. In a range of applications including imaging analysis and fMRI studies, random variables are observed on a lattice graph. In such a setting it is important to account for the lattice structure when estimating the covariance operator. In this paper we consider both minimax and adaptive estimation of the covariance operator over collections of polynomially decaying and exponentially decaying parameter spaces. We first establish the minimax rates of convergence for estimating the covariance operator under the operator norm. The results show that the dimension of the lattice graph significantly affects the optimal rates convergence, often much more so than the dimension of the random variables. We then consider adaptive estimation of the covariance operator. A fully data driven block thresholding procedure is proposed and is shown to be adaptively rate optimal simultaneously over a wide range of polynomially decaying and exponentially decaying parameter spaces. The adaptive block thresholding procedure is easy to implement and numerical experiments are carried out to illustrate the merit of the procedure.

Abstract

Multiple testing of correlations arises in many applications including gene coexpression network analysis and brain connectivity analysis. In this paper, we consider large scale simultaneous testing for correlations in both the one-sample and two-sample settings. New multiple testing procedures are proposed and a bootstrap method is introduced for estimating the proportion of the nulls falsely rejected among all the true nulls. The properties of the proposed procedures are investigated both theoretically and numerically. It is shown that the procedures asymptotically control the overall false discovery rate and false discovery proportion at the nominal level. Simulation results show that the methods perform well numerically in terms of both the size and power of the test and it significantly outperforms two alternative methods. The two-sample procedure is also illustrated by an analysis of a prostate cancer dataset for the detection of changes in coexpression patterns between gene expression levels.

Abstract

Precision matrix is of significant importance in a wide range of applications in multivariate analysis. This paper considers adaptive minimax estimation of sparse precision matrices in the high dimensional setting. Optimal rates of convergence are established for a range of matrix norm losses. A fully data driven estimator based on adaptive constrained \(\ell_1\) minimization is proposed and its rate of convergence is obtained over a collection of parameter spaces. The estimator, called ACLIME, is easy to implement and performs well numerically. A major step in establishing the minimax rate of convergence is the derivation of a rate-sharp lower bound. A "two-directional" lower bound technique is applied to obtain the minimax lower bound. The upper and lower bounds together yield the optimal rates of convergence for sparse precision matrix estimation and show that the ACLIME estimator is adaptively minimax rate optimal for a collection of parameter spaces and a range of matrix norm losses simultaneously.

Abstract

Estimation of high dimensional precision matrix has many important applications. Motivated by the analysis of gene expression data measured in different tissues or disease states, we consider the problem of jointly estimating multiple precision matrices in order to effectively utilize the partially shared graphical structures of the corresponding graphs. The procedure is based on a weighted constrained \(\ell_\infty/\ell_1\) minimization, which can be effectively solved by a second-order cone programming. It is shown that under the entry-wise \(\ell_\infty\) norm loss such a joint estimation of precision matrices leads to a faster convergence rate than estimating the precision matrices separately. After thresholding, the joint estimation can have a faster rate than separate estimation under the matrix \(\ell_1\) norm. We also show that when the collection of the precision matrices have a common support, the proposed procedure leads to the exact graph structure recovery with probability tending to 1. The method is illustrated through an analysis of a human myocardium gene expression data in normal and failing human hearts. The results indicate the genes in the hypertrophic cardiomyopathy pathway are less connected in the failing heart tissues than normal hearts.

Abstract

This is an expository paper that reviews recent developments on optimal estimation of structured high-dimensional covariance and precision matrices. Minimax rates of convergence for estimating several classes of structured covariance and precision matrices, including bandable, Toeplitz, and sparse covariance matrices as well as sparse precision matrices, are given under the spectral norm loss. Data-driven adaptive procedures for estimating various classes of matrices are presented. Some key technical tools including large deviation results and minimax lower bound arguments that are used in the theoretical analyses are discussed. In addition, estimation under other losses and a few related problems such as Gaussian graphical models, sparse principal component analysis, and hypothesis testing on the covariance structure are considered. Some open problems on estimating high-dimensional covariance and precision matrices and their functionals are also discussed.

Abstract

Motivated by differential co-expression analysis in genomics, we consider in this paper estimation and testing of high-dimensional differential correlation matrices. An adaptive thresholding procedure is introduced and theoretical guarantees are given. Minimax rate of convergence is established and the proposed estimator is shown to be adaptively rate-optimal over collections of paired correlation matrices with approximately sparse differences. Simulation results show that the procedure significantly outperforms two other natural methods that are based on separate estimation of the individual correlation matrices. The procedure is also illustrated through an analysis of a breast cancer dataset, which provides evidence at the gene co-expression level that several genes, of which a subset has been previously verified, are associated with the breast cancer. Hypothesis testing on the differential correlation matrices is also considered. A test, which is particularly well suited for testing against sparse alternatives, is introduced. In addition, other related problems, including estimation of a single sparse correlation matrix, estimation of the differential covariance matrices, and estimation of the differential cross-correlation matrices, are also discussed.

Abstract

This paper considers sparse spiked covariance matrix models in the high-dimensional setting and studies the minimax estimation of the covariance matrix and the principal subspace as well as the minimax rank detection. The optimal rate of convergence for estimating the spiked covariance matrix under the spectral norm is established, which requires significantly different techniques from those for estimating other structured covariance matrices such as bandable or sparse covariance matrices. We also establish the minimax rate under the spectral norm for estimating the principal subspace, the primary object of interest in principal component analysis. In addition, the optimal rate for the rank detection boundary is obtained. This result also resolves a gap between the upper and lower bounds in a recent paper by Berthet and Rigollet [1] where the special case of rank one is considered.

Abstract

Model organisms and human studies have led to increasing empirical evidence that interactions among genes contribute broadly to genetic variation of complex traits. In the presence of gene-by-gene interactions, the dimensionality of the feature space becomes extremely high relative to the sample size. This imposes a significant methodological challenge in identifying gene-by-gene interactions. In the present paper, through a Gaussian graphical model framework, we translate the problem of identifying gene-by-gene interactions associated with a binary trait \(D\) into an inference problem on the difference of two high-dimensional precision matrices, which summarize the conditional dependence network structures of the genes. We propose a procedure for testing the differential network globally that is particularly powerful against sparse alternatives. In addition, a multiple testing procedure with false discovery rate control is developed to infer the specific structure of the differential network. Theoretical justification is provided to ensure the validity of the proposed tests and optimality results are derived under sparsity assumptions. A simulation study demonstrates that the proposed tests maintain the desired error rates under the null and have good power under the alternative. The methods are applied to a breast cancer gene expression study.

Abstract

Differential entropy and log determinant of the covariance matrix of a multivariate Gaussian distribution have many applications in coding, communications, signal processing and statistical inference. In this paper we consider in the high dimensional setting optimal estimation of the differential entropy and the log-determinant of the covariance matrix. We first establish a central limit theorem for the log determinant of the sample covariance matrix in the high dimensional setting where the dimension \(p(n)\) can grow with the sample size \(n\). An estimator of the differential entropy and the log determinant is then considered. Optimal rate of convergence is obtained. It is shown that in the case \(p(n)/n\to0\) the estimator is asymptotically sharp minimax. The ultra-high dimensional setting where \(p(n)>n\) is also discussed.

Abstract

It is often of interest to understand how the structure of a genetic network differs between two conditions. In this paper, each condition-specific network is modeled using the precision matrix of a multivariate normal random vector, and a method is proposed to directly estimate the difference of the precision matrices. In contrast to other approaches, such as separate or joint estimation of the individual matrices, direct estimation does not require those matrices to be sparse, and thus can allow the individual networks to contain hub nodes. Under the assumption that the true differential network is sparse, the direct estimator is shown to be consistent in support recovery and estimation. It is also shown to outperform existing methods in simulations, and its properties are illustrated on gene expression data from late-stage ovarian cancer patients.

Abstract

Principal component analysis (PCA) is one of the most commonly used statistical procedures with a wide range of applications.This paper considers both minimax and adaptive estimation of the principal subspace in the high dimensional setting. Under mild technical conditions, we first establish the optimal rates of convergence for estimating the principal subspace which are sharp with respect to all the parameters, thus providing a complete characterization of the difficulty of the estimation problem in term of the convergence rate. The lower bound is obtained by calculating the local metric entropy and an application of Fano's Lemma. The rate optimal estimator is constructed using aggregation, which, however, might not be computationally feasible. We then introduce an adaptive procedure for estimating the principal subspace which is fully data driven and can be computed efficiently. It is shown that the estimator attains the optimal rates of convergence simultaneously over a large collection of the parameter spaces. A key idea in our construction is a reduction scheme which reduces the sparse PCA problem to a high-dimensional multivariate regression problem. This method is potentially also useful for other related problems.

Abstract

This paper considers testing a covariance matrix \(\Sigma\) in the high dimensional setting where the dimension \(p\) can be comparable or much larger than the sample size \(n\). The problem of testing the hypothesis \(H_0:\Sigma=\Sigma_0\) for a given covariance matrix \(\Sigma_0\) is studied from a minimax point of view. We first characterize the boundary that separates the testable region from the non-testable region by the Frobenius norm when the ratio between the dimension \(p\) over the sample size \(n\) is bounded. A test based on a U-statistic is introduced and is shown to be rate optimal over this asymptotic regime. Furthermore, it is shown that the power of this test uniformly dominates that of the corrected likelihood ratio test (CLRT) over the entire asymptotic regime under which the CLRT is applicable. The power of the U-statistic based test is also analyzed when \(p/n\) is unbounded.

Abstract

Toeplitz covariance matrices are used in the analysis of stationary stochastic processes and a wide range of applications including radar imaging, target detection, speech recognition, and communications systems. In this paper, we consider optimal estimation of large Toeplitz covariance matrices and establish the minimax rate of convergence for two commonly used parameter spaces under the spectral norm. The properties of the tapering and banding estimators are studied in detail and are used to obtain the minimax upper bound. The results also reveal a fundamental difference between the tapering and banding estimators over certain parameter spaces. The minimax lower bound is derived through a novel construction of a more informative experiment for which the minimax lower bound is obtained through an equivalent Gaussian scale model and through a careful selection of a finite collection of least favorable parameters. In addition, optimal rate of convergence for estimating the inverse of a Toeplitz covariance matrix is also established.

Abstract

Motivated by the analysis of eQTL data, we introduce a sparse high dimensional multivariate regression model for studying the conditional independent relationships among a set of genes adjusting for possible genetic effects, as well as the genetic architecture that influences the gene expressions. The precision matrix in the model specifies a covariate-adjusted Gaussian graph, which presents the conditional dependency structure of gene expression after the confounding genetic effects on gene expressions are taken into account. We present a covariate-adjusted precision matrix estimation (CAPME) method using constrained \(\ell_1\) minimization, which can be easily implemented by linear programming. Asymptotic convergence rates and sign consistency are established for the estimators of the regression coefficients and the precision matrix, allowing both the number of genes and the number of the genetic variants to diverge. Numerical performance of the estimators is investigated using both simulated and real data sets. Simulation results show that CAPME results in significant improvements in both precision matrix estimation and in graphical structure selection when compared to the standard Gaussian graphical model assuming constant means. CAPME is also applied to analyze a yeast eQTL data for the identification of the gene network among a set of genes in the Mitogen-activated protein kinase (MAPK) pathway.

Abstract

This paper considers in the high-dimensional setting three inter-related problems: (a) testing the equality of two covariance matrices \(\Sigma_1\) and \(\Sigma_2\); (b) recovering the support of \(\Sigma_1-\Sigma_2\) ; and (c) testing the equality of \(\Sigma_1\) and \(\Sigma_2\) row by row. We propose a new test for testing the hypothesis \(H_0:\Sigma_1=\Sigma_2\) and investigate its theoretical and numerical properties. The limiting null distribution of the test statistic is derived and the power of the test is studied. The test is shown to enjoy certain optimality and to be especially powerful against sparse alternatives. The simulation results show that the test significantly outperforms the existing methods both in terms of size and power. Analysis of a prostate cancer dataset is carried out to demonstrate the application of the testing procedures. When the null hypothesis of equal covariance matrices is rejected, it is often of significant interest to further investigate how they differ from each other. Motivated by applications in genomics, we also consider recovering the support of \(\Sigma_1-\Sigma_2\) and testing the equality of the two covariance matrices row by row. New procedures are introduced and their properties are studied. Applications to gene selection is also discussed.

Abstract

This paper considers estimation of sparse covariance matrices and establishes the optimal rate of convergence under a range of matrix operator norm and Bregman divergence losses. A major focus is on the derivation of a rate sharp minimax lower bound. The problem exhibits new features that are significantly different from those that occur in the conventional nonparametric function estimation problems. Standard techniques fail to yield good results and new tools are thus needed. We first develop a lower bound technique that is particularly well suited for treating "two-directional" problems such as estimating sparse covariance matrices. The result can be viewed as a generalization of Le Cam's method in one direction and Assouad's Lemma in another. This lower bound technique is of independent interest and can be used for other matrix estimation problems. We then establish a rate sharp minimax lower bound for estimating sparse covariance matrices under the spectral norm by applying the general lower bound technique. A thresholding estimator is shown to attain the optimal rate of convergence under the spectral norm. The results are then extended to the general matrix \(\ell_w\) operator norms for 1 ≤ w ≤ ∞. In addition, we give a unified result on the minimax rate of convergence for sparse covariance matrix estimation under a class of Bregman divergence losses.

Abstract

Driven by a wide range of applications in high-dimensional data analysis, there has been significant recent interest in the estimation of large covariance matrices. In this paper, we consider optimal estimation of a covariance matrix as well as its inverse over several commonly used parameter spaces under the matrix \(\ell_1\) norm. Both minimax lower and upper bounds are derived. The lower bounds are established by using hypothesis testing arguments, where at the core are a novel construction of collections of least favorable multivariate normal distributions and bounding the affinities between pairs of distributions. The lower bound analysis also provides insight into where the difficulties of the covariance matrix estimation problem come from. Specific thresholding estimator and tapering estimator are constructed and shown to be minimax rate optimal. The optimal rates of convergence established in the paper can serve as a benchmark for the performance of covariance matrix estimation methods.

Abstract

Estimation of large covariance matrices has drawn considerable recent attention and the theoretical focus so far is mainly on developing a minimax theory over a fixed parameter space. In this paper, we consider adaptive covariance matrix estimation where the goal is to construct a single procedure which is minimax rate optimal simultaneously over each parameter space in a large collection. A fully data-driven block thresholding estimator is proposed. The estimator is constructed by carefully dividing the sample covariance matrix into blocks and then simultaneously estimating the entries in a block by thresholding. The estimator is shown to be optimally rate adaptive over a wide range of bandable covariance matrices. A simulation study is carried out and shows that the block thresholding estimator performs well numerically. Some of the technical tools developed in this paper can also be of independent interest.

Abstract

The coherence of a random matrix, which is defined to be the largest magnitude of the Pearson correlation coefficients between the columns of the random matrix, is an important quantity for a wide range of applications including high-dimensional statistics and signal processing. Inspired by these applications, this paper studies the limiting laws of the coherence of \(n\times p\) random matrices for a full range of the dimension \(p\) with a special focus on the ultra high-dimensional setting. Assuming the columns of the random matrix are independent random vectors with a common spherical distribution, we give a complete characterization of the behavior of the limiting distributions of the coherence. More specifically, the limiting distributions of the coherence are derived separately for three regimes: \(\log(p)/n\to0\), \(\log(p)/n\to\beta\in(0,\infty)\), and \(\log(p)/n\to\infty\). The results show that the limiting behavior of the coherence differs significantly in different regimes and exhibits interesting phase transition phenomena as the dimension \(p\) grows as a function of \(n\). Applications to statistics and compressed sensing in the ultra high-dimensional setting are also discussed.

Abstract

Testing covariance structure is of significant interest in many areas of statistical analysis and construction of compressed sensing matrices is an important problem in signal processing. Motivated by these applications, we study in this paper the limiting laws of the coherence of an \(n\times p\) random matrix in the high-dimensional setting where \(p\) can be much larger than \(n\). Both law of large numbers and limiting distribution are derived. We then consider testing the bandedness of the covariance matrix of a high dimensional Gaussian distribution which includes testing for independence as a special case. The limiting laws of the coherence of the data matrix play a critical role in the construction of the test. The asymptotic results is also applied to the construction of compressed sensing matrices.

Abstract

In this paper we consider estimation of sparse covariance matrices and propose a thresholding procedure which is adaptive to the variability of individual entries. The estimators are fully data driven and enjoy excellent performance both theoretically and numerically. It is shown that the estimators adaptively achieve the optimal rate of convergence over a large class of sparse covariance matrices under the spectral norm. In contrast, the commonly used universal thresholding estimators are shown to be sub-optimal over the same parameter spaces. Support recovery is also discussed. The adaptive thresholding estimators are easy to implement. Numerical performance of the estimators is studied using both simulated and real data. Simulation results show that the adaptive thresholding estimators uniformly outperform the universal thresholding estimators. The method is also illustrated in an analysis on a dataset from a small round blue-cell tumors microarray experiment.

Abstract

A constrained \(\ell_1\) minimization method is proposed for estimating a sparse inverse covariance matrix based on a sample of \(n\) iid \(p\)-variate random variables. The resulting estimator is shown to enjoy a number of desirable properties. In particular, it is shown that the rate of convergence between the estimator and the true \(s\)-sparse precision matrix under the spectral norm is \(s(\log p/n)^{1/2}\) when the population distribution has either exponential-type tails or polynomial-type tails. Convergence rates under the elementwise \(\ell_1\) norm and Frobenius norm are also presented. In addition, graphical model selection is considered. The procedure is easily implementable by linear programming. Numerical performance of the estimator is investigated using both simulated and real data. In particular, the procedure is applied to analyze a breast cancer dataset. The procedure performs favorably in comparison to existing methods.

Abstract

Covariance matrix plays a central role in multivariate statistical analysis. Significant advances have been made recently on developing both theory and methodology for estimating large covariance matrices. However, a minimax theory has yet been developed. In this paper we establish the optimal rates of convergence for estimating the covariance matrix under both the operator norm and Frobenius norm. It is shown that optimal procedures under the two norms are different and consequently matrix estimation under the operator norm is fundamentally different from vector estimation. The minimax upper bound is obtained by constructing a special class of tapering estimators and by studying their risk properties. A key step in obtaining the optimal rate of convergence is the derivation of the minimax lower bound. The technical analysis requires new ideas that are quite different from those used in the more conventional function/sequence estimation problems.

Abstract

We study scalar-on-function linear regression when each covariate curve is observed only through finitely many noisy point evaluations. Our goal is to characterize the minimax estimation and prediction risks as joint functions of the number of trajectories \(n\) and the within-trajectory resolution \(m\). Working in a fixed trigonometric eigenbasis, with covariance eigenvalues decaying at rate \(\alpha\) and slope function of Sobolev smoothness \(s\), we derive matching minimax upper and lower bounds under two canonical sampling schemes. Under an independent random design, the minimax prediction rate is \[ n^{-\frac{2\alpha+2s}{2\alpha+2s+1}} + (nm)^{-\frac{2\alpha+2s}{4\alpha+2s+1}}. \] The first term is the fully observed functional linear regression benchmark, while the second term captures the cost of noisy point evaluations after amplification by the inverse covariance operator. Under a common design on an equally spaced grid, the shared sampling geometry introduces additional obstructions, and the minimax prediction rate becomes \[ n^{-\frac{2\alpha+2s}{2\alpha+2s+1}} + (nm)^{-\frac{2\alpha+2s}{4\alpha+2s+1}} + m^{-(2\alpha+2s)} + m^{-4\alpha}. \] Here the third term represents fixed-grid discretization error, whereas the fourth reflects the cost of identifying unknown eigenvalues from observations on a common grid. We further construct data-driven adaptive estimators, based on covariance-scale screening and blockwise prediction-energy thresholding, that attain these rates without prior knowledge of the eigenvalue sequence or the smoothness indices. The results reveal a sharp resolution-dependent phase transition under independent design and a richer phase diagram under common design. Numerical simulations and a real-data example illustrate the theoretical findings.

Abstract

We study functional linear regression under differential privacy in both fully observed and discretely observed settings. For the fully observed case, we propose a private algorithm that integrates private covariance estimation with noisy gradient descent, and establish its convergence rate. We further derive a matching minimax lower bound, demonstrating that the proposed estimator is minimax-optimal up to a logarithmic factor. In the discrete setting, we analyze the trade-offs among privacy, statistical accuracy, and discretization error, revealing a notable phase transition phenomenon. Simulations and real-data applications confirm the effectiveness of our approach in preserving privacy while maintaining statistical accuracy. Differential privacy in functional linear regression presents unique challenges due to the infinite-dimensional structure of functional data, which complicates sensitivity analysis, noise calibration, and algorithm design. Standard privacy mechanisms developed for finite-dimensional settings are inadequate in this context. To address these difficulties, we use basis expansions to reduce the problem to high-dimensional linear regression and introduce a preconditioning strategy that ensures robustness, even under ill-conditioned covariance structures. This framework broadens the applicability of differential privacy to complex functional models and ill-posed settings, advancing privacy-preserving methods in infinite-dimensional data analysis.

Abstract

This paper studies the high-dimensional mixed linear regression (MLR) where the output variable comes from one of the two linear regression models with an unknown mixing proportion and an unknown covariance structure of the random covariates. Building upon a high-dimensional EM algorithm, we propose an iterative procedure for estimating the two regression vectors and establish their rates of convergence. We further construct debiased estimators and establish their asymptotic normality. Confidence intervals centered at the debiased estimators are constructed for individual coordinates. Furthermore, a large-scale multiple testing procedure is proposed for testing the regression coefficients and is shown to control the false discovery rate (FDR) asymptotically. Simulation studies are carried out to examine the numerical performance of the proposed methods and their superiority over existing methods. The proposed methods are further illustrated through an analysis of a dataset of multiplex image cytometry, which investigates the interaction networks among the cellular phenotypes that include the expression levels of 20 epitopes or combinations of markers.

Abstract

Achieving optimal statistical performance while ensuring the privacy of personal data is a challenging yet crucial objective in modern data analysis. However, characterizing the optimality, particularly the minimax lower bound, under privacy constraints is technically difficult. To address this issue, we propose a novel approach called the score attack, which provides a lower bound on the differential-privacy-constrained minimax risk of parameter estimation. The score attack method is based on the tracing attack concept in differential privacy and can be applied to any statistical model with a well-defined score statistic. It can optimally lower bound the minimax risk of estimating unknown model parameters, up to a logarithmic factor, while ensuring differential privacy for a range of statistical problems. We demonstrate the effectiveness and optimality of this general method in various examples, such as the generalized linear model in both classical and high-dimensional sparse settings, the Bradley-Terry-Luce model for pairwise comparisons, and nonparametric regression over the Sobolev class.

Abstract

We consider high dimensional M-estimation in settings where the response \(Y\) is possibly missing at random and the covariates \(X\in\mathbb{R}^p\) can be high dimensional compared to the sample size \(n\). The parameter of interest \(\theta_0\in\mathbb{R}^d\) is defined as the minimizer of the risk of a convex loss, under a fully non-parametric model, and \(\theta_0\) itself is high dimensional which is a key distinction from existing works. Standard high dimensional regression and series estimation with possibly mis-specified models and missing \(Y\) are included as special cases, as well as their counterparts in causal inference using ‘potential outcomes’. Assuming \(\theta_0\) is \(s\)-sparse \((s\ll n)\), we propose an \(L_1\)-regularized debiased and doubly robust (DDR) estimator of \(\theta_0\) based on a high dimensional adaptation of the traditional double robust (DR) estimator’s construction. Under mild tail assumptions and arbitrarily chosen (working) models for the propensity score (PS) and the outcome regression (OR) estimators, satisfying only some high-level conditions, we establish finite sample performance bounds for the DDR estimator showing its (optimal) \(L_2\) error rate to be \(\sqrt{s(\log d)/n}\) when both models are correct, and its consistency and DR properties when only one of them is correct. Further, when both the models are correct, we propose a desparsified version of our DDR estimator that satisfies an asymptotic linear expansion and facilitates inference on low dimensional components of \(\theta_0\). Finally, we discuss various of choices of high dimensional parametric/semi-parametric working models for the PS and OR estimators. All results are validated via detailed simulations.

Abstract

Existing high-dimensional statistical methods are largely developed for analyzing individual-level data. In this work, we study estimation and inference for high-dimensional linear models when only "proxy data" is available. These proxies encompass marginal statistics and sample covariance matrices computed from distinct sets of individuals. We develop a rate optimal method for estimation and inference for the regression coefficient vector and its linear functionals based on the proxy data. We show the intrinsic limitations in the proxy-data based inference: the minimax optimal rate for estimation is slower than that in the conventional case where individual data are observed. These interesting findings are illustrated through simulation studies and an analysis of a dataset concerning the genetic associations of hindlimb muscle weights in a mouse population.

Abstract

We introduce the R package SIHR for statistical inference in high-dimensional generalized linear models with continuous and binary outcomes. The package provides functionalities for constructing confidence intervals and performing hypothesis tests for low-dimensional objectives in both one-sample and two-sample regression settings. We illustrate the usage of SIHR through numerical examples and present real data applications to demonstrate the package’s performance and practicality.

Abstract

This paper studies the problem of statistical inference for genetic relatedness between binary traits based on individual-level genome-wide association data. Specifically, under the high-dimensional logistic regression models, we define parameters characterizing the cross-trait genetic correlation, the genetic covariance and the trait-specific genetic variance. A novel weighted debiasing method is developed for the logistic Lasso estimator and computationally efficient debiased estimators are proposed. The rates of convergence for these estimators are studied and their asymptotic normality is established under mild conditions. Moreover, we construct confidence intervals and statistical tests for these parameters, and provide theoretical justifications for the methods, including the coverage probability and expected length of the confidence intervals, as well as the size and power of the proposed tests. Numerical studies are conducted under both model generated data and simulated genetic data to show the superiority of the proposed methods and their applicability to the analysis of real genetic data. Finally, by analyzing a real data set on autoimmune diseases, we demonstrate the ability to obtain novel insights about the shared genetic architecture between ten pediatric autoimmune diseases.

Abstract

This paper presents a selective survey of recent developments in statistical inference and multiple testing for high-dimensional regression models, including linear and logistic regression. We examine the construction of confidence intervals and hypothesis tests for various low-dimensional objectives such as regression coefficients and linear and quadratic functionals. The key technique is to generate debiased and desparsified estimators for the targeted low-dimensional objectives and estimate their uncertainty. In addition to covering the motivations for and intuitions behind these statistical methods, we also discuss their optimality and adaptivity in the context of high-dimensional inference. In addition, we review the development of statistical inference based on multiple regression models and the advancement of large-scale multiple testing for high-dimensional regression. The R package SIHR has implemented the high-dimensional inference methods reviewed in this paper.

Abstract

This paper develops a unified statistical inference framework for high-dimensional binary generalized linear models (GLMs) with general link functions. Both unknown and known design distribution settings are considered. A two-step weighted bias-correction method is proposed for constructing confidence intervals and simultaneous hypothesis tests for individual components of the regression vector. Minimax lower bound for the expected length is established and the proposed confidence intervals are shown to be rate-optimal up to a logarithmic factor. The numerical performance of the proposed procedure is demonstrated through simulation studies and an analysis of a single cell RNA-seq data set, which yields interesting biological insights that integrate well into the current literature on the cellular immune response mechanisms as characterized by single-cell transcriptomics. The theoretical analysis provides important insights on the adaptivity of optimal confidence intervals with respect to the sparsity of the regression vector. New lower bound techniques are introduced and they can be of independent interest to solve other inference problems in high-dimensional binary GLMs.

Abstract

Linear mixed-effects models are widely used in analyzing clustered or repeated measures data. We propose a quasi-likelihood approach for estimation and inference of the unknown parameters in linear mixed-effects models with high-dimensional fixed effects. The proposed method is applicable to general settings where the dimension of the random effects and the cluster sizes are possibly large. Regarding the fixed effects, we provide rate optimal estimators and valid inference procedures that do not rely on the structural information of the variance components. We also study the estimation of variance components with high-dimensional fixed effects in general settings. The algorithms are easy to implement and computationally fast. The proposed methods are assessed in various simulation settings and are applied to a real study regarding the associations between body mass index and genetic polymorphic markers in a heterogeneous stock mice population.

Abstract

In this paper, we study sparse group Lasso for high-dimensional double sparse linear regression, where the parameter of interest is simultaneously element-wise and group-wise sparse. This problem is an important instance of the simultaneously structured model – an actively studied topic in statistics and machine learning. In the noiseless case, we provide matching upper and lower bounds on sample complexity for the exact recovery of sparse vectors and for stable estimation of approximately sparse vectors, respectively. In the noisy case, we develop upper and matching minimax lower bounds for estimation error. We also consider the debiased sparse group Lasso and investigate its asymptotic property for the purpose of statistical inference. Finally, numerical studies are provided to support the theoretical results.

Abstract

Instrumental variables have been widely used to estimate the causal effect of a treatment on an outcome.Existing confidence intervals for causal effects based on instrumental variables assume that all of the putative instrumental variables are valid; a valid instrumental variable is a variable that affects the outcome only by affecting the treatment and is not related to unmeasured confounders. However, in practice, some of the putative instrumental variables are likely to be invalid. This paper presents two tools to conduct valid inference and tests in the presence of invalid instruments. First, we propose a simple and general approach to construct confidence intervals based on taking unions of well-known confidence intervals. Second, we propose a novel test for the null causal effect based on a collider bias. Our two proposals outperform traditional instrumental variable confidence intervals when invalid instruments are present and can also be used as a sensitivity analysis when there is concern that instrumental variables assumptions are violated. The new approach is applied to a Mendelian randomization study on the causal effect of low-density lipoprotein on globulin levels.

Abstract

Group inference has been a long-standing question in statistics and the development of high-dimensional group inference is an essential part of statistical methods for analyzing complex data sets, including hierarchical testing, tests of interaction, detection of heterogeneous treatment effects and local heritability. Group inference in regression models can be measured with respect to a weighted quadratic functional of the regression sub-vector corresponding to the group. Asymptotically unbiased estimators of these weighted quadratic functionals are constructed and a procedure using these estimator for inference is proposed. We derive its asymptotic Gaussian distribution which allows to construct asymptotically valid confidence intervals and tests which perform well in terms of length or power. The results simultaneously address four challenges encountered in the literature: controlling coverage or type I error even when the variables inside the group are highly correlated, achieving a good power when there are many small coefficients inside the group, computational efficiency even for a large group, and no requirements on the group size. We apply the methodology to several interesting statistical problems and demonstrate its strength and usefulness on simulated and real data.

Abstract

Privacy-preserving data analysis is a rising challenge in contemporary statistics, as the privacy guarantees of statistical methods are often achieved at the expense of accuracy. In this paper, we investigate the tradeoff between statistical accuracy and privacy in mean estimation and linear regression, under both the classical low- dimensional and modern high-dimensional settings. A primary focus is to establish minimax optimality for statistical estimation with the \((\epsilon,\delta)\)-differential privacy constraint. By refining the “tracing adversary” technique for lower bounds in the theoretical computer science literature, we improve existing minimax lower bound for low-dimensional mean estimation and establish new lower bounds for high-dimensional mean estimation and linear regression problems. We also design differentially private algorithms that attain the minimax lower bounds up to logarithmic factors. In particular, for high-dimensional linear regression, a novel private iterative hard thresholding algorithm is proposed. The numerical performance of differentially private algorithms is demonstrated by simulation studies and applications to real data sets.

Abstract

The ability to predict individualized treatment effects (ITEs) based on a given patient's profile is essential for personalized medicine. We propose a hypothesis testing approach to choosing between two potential treatments for a given individual in the framework of high-dimensional linear models. The methodological novelty lies in the construction of a debiased estimator of the ITE and establishment of its asymptotic normality uniformly for an arbitrary future high-dimensional observation, while the existing methods can only handle certain specific forms of observations. We introduce a testing procedure with the type-I error controlled and establish its asymptotic power. The proposed method can be extended to making inference for general linear contrasts, including both the average treatment effect and outcome prediction. We introduce the optimality framework for hypothesis testing from both the minimaxity and adaptivity perspectives and establish the optimality of the proposed procedure. An extension to high-dimensional approximate linear models is also considered. The finite sample performance of the procedure is demonstrated in simulation studies and further illustrated through an analysis of electronic health records data from patients with rheumatoid arthritis.

Abstract

We consider statistical inference for the explained variance \(\beta^\top\Sigma\beta\) under the high-dimensional linear model \(Y=X\beta+\epsilon\) in the semi-supervised setting, where \(\beta\) is the regression vector and \(\Sigma\) is the design covariance matrix. A calibrated estimator, which efficiently integrates both labelled and unlabelled data, is proposed. It is shown that the estimator achieves the minimax optimal rate of convergence in the general semi-supervised framework. The optimality result characterizes how the unlabelled data affects the minimax optimal rate. Moreover, the limiting distribution for the proposed estimator is established and data-driven confidence intervals for the explained variance are constructed. We further develop a randomized calibration technique for statistical inference in the presence of weak signals and apply the obtained inference results to a range of important statistical problems, including signal detection and global testing, prediction accuracy evaluation, and confidence ball construction. The numerical performance of the proposed methodology is demonstrated in simulation studies and an analysis of estimating heritability for a yeast segregant data set with multiple traits.

Abstract

Estimating the genetic relatedness between two traits based on the genome-wide association data is an important problem in genetics research. In the framework of high-dimensional linear models, we introduce two measures of genetic relatedness and develop optimal estimators for them. One is genetic covariance, which is defined to be the inner product of the two regression vectors, and another is genetic correlation, which is a normalized inner product by their lengths. We propose functional de-biased estimators (FDEs), which consist of an initial estimation step with the plug-in scaled Lasso estimator, and a further bias correction step. We also develop estimators of the quadratic functionals of the regression vectors, which can be used to estimate the heritability of each trait. The estimators are shown to be minimax rate-optimal and can be efficiently implemented. Simulation results show that FDEs provide better estimates of the genetic relatedness than simple plug-in estimates. FDE is also applied to an analysis of a yeast segregant data set with multiple traits to estimate the genetic relatedness among these traits.

Abstract

Modern, high dimensional data has renewed investigation on instrumental variables (IV) analysis, primary focusing on estimation of the included endogenous variable under sparsity and little attention towards specification tests. This paper studies in high dimension the Durbin-Wu-Hausman (DWH) test, a popular specification test for endogeneity in IV regression. We show, surprisingly, that the DWH test maintains its size in high dimension, but at an expense in power. We propose a new test that remedies this issue and has better power than the DWH test. Simulation studies reveal that our test achieves near-oracle performance to detect endogeneity.

Abstract

The instrumental variable (IV) method is commonly used to estimate the causal effect of a treatment on an outcome by using IVs that satisfy the assumptions of association with treatment, no direct effect on the outcome and ignorability. A major challenge in IV analysis is to find said IVs, but typically one is unsure of whether all of the putative IVs are in fact valid (i.e. satisfy the assumptions). We propose a general inference procedure that provides honest inference in the presence of invalid IVs, even after controlling for a large number of covariates. The key step of our method is a novel selection procedure, which we call Two-Stage Hard Thresholding (TSHT), where we use hard thresholding to select the set of non-redundant instruments in the first stage and subsequently use hard thresholding to select the valid instruments in the second stage using the thresholding from the first stage. TSHT allows us to not only select invalid IVs, but also provides honest confidence intervals of the treatment effect at \(\sqrt n\) rate. We establish asymptotic properties of our procedure and demonstrate that our procedure performs well in simulation studies compared to traditional IV methods, especially when the instruments are invalid.

Abstract

This paper considers point and interval estimation of the \(\ell_q\) loss of an estimator in high-dimensional linear regression with random design. Both the setting of known identity design covariance matrix and known noise level and the setting of unknown design covariance matrix and noise level are studied. We establish the minimax convergence rate for estimating the \(\ell_q\) loss and the minimax expected length of confidence intervals for the \(\ell_q\) loss of a broad collection of estimators of the regression vector. We also investigate the adaptivity of the confidence intervals for the \(\ell_q\) loss. The results reveal interesting and significant differences between estimating the \(\ell_2\) loss and \(\ell_q\) loss with \(1\le q<2\) as well as the differences between the two settings. A major step in our analysis is to establish rate sharp lower bounds for the minimax estimation error and the expected length of minimax and adaptive confidence intervals for the \(\ell_q\) loss, which requires the development of new technical tools. A significant difference between loss estimation and the traditional parameter estimation is that for loss estimation the constraint is on the performance of the estimator of the regression vector, but the lower bounds are on the difficulty of estimating its \(\ell_q\) loss. The technical tools developed in this paper can also be of independent interest.

Abstract

Motivated by important applications in genomics, we consider in this paper global and coordinatewise tests for the comparisons of two high-dimensional linear regression models. A procedure for testing the equality of the two regression vectors globally is proposed and shown to be particularly powerful against sparse alternatives. In addition, we introduce a multiple testing procedure for identifying unequal coordinates while controlling the false discovery rate (FDR) and false discovery proportion (FDP). Theoretical justifications are provided to guarantee the validity of the proposed tests and optimality results are established under sparsity assumptions on the regression coefficients. The proposed testing procedures are easy to implement. Numerical properties of the procedures are investigated through simulation and real data analysis. Simulation results show that the proposed tests maintain the desired error rates under the null and have good power under the alternative at moderate sample sizes. The procedures are applied to the Framingham Offspring study to investigate the interactions between smoking and cardiovascular related genetic mutations important for an inflammation marker.

Abstract

Confidence sets play a fundamental role in statistical inference. In this paper, we consider confidence intervals for high dimensional linear regression with random design. We first establish the convergence rates of the minimax expected length for confidence intervals in the oracle setting where the sparsity parameter is given. The focus is then on the problem of adaptation to sparsity for the construction of confidence intervals. Ideally, an adaptive confidence interval should have its length automatically adjusted to the sparsity of the unknown regression vector, while maintaining a prespecified coverage probability. It is shown that such a goal is in general not attainable, except when the sparsity parameter is restricted to a small region over which the confidence intervals have the optimal length of the usual parametric rate. It is further demonstrated that the lack of adaptivity is not due to the conservativeness of the minimax framework, but is fundamentally caused by the difficulty of learning the bias accurately.

Abstract

This paper presents a unified geometric framework for the statistical analysis of a general ill-posed linear inverse model which includes as special cases noisy compressed sensing, sign vector recovery, trace regression, orthogonal matrix estimation, and noisy matrix completion. We propose computationally feasible convex programs for statistical inference including estimation, confidence intervals and hypothesis testing. A theoretical framework is developed to characterize the local estimation rate of convergence and to provide statistical inference guarantees. Our results are built based on the local conic geometry and duality. The difficulty of statistical inference is captured by the geometric characterization of the local tangent cone through the Gaussian width and Sudakov estimate.

Abstract

Instrumental variables have been widely used for estimating the causal effect be- tween exposure and outcome. Conventional estimation methods require complete knowledge about all the instruments' validity; a valid instrument must not have a direct effect on the outcome and not be related to unmeasured confounders. Often, this is impractical as highlighted by Mendelian randomization studies where genetic markers are used as instruments and complete knowledge about instruments' validity is equivalent to complete knowledge about the involved genes' functions. In this paper, we propose a method for estimation of causal effects when this complete knowledge is absent. It is shown that identification and estimation is possible under a weaker requirement that more than 50% of instruments are valid, without precisely knowing which of the 50%+ instruments are valid. Sharp identification limits with invalid instruments are given. A fast penalized \(\ell_1\) estimation method, called sisVIVE, is introduced for estimating the causal effect without knowing which instruments are valid, with theoretical guarantees on its performance. The proposed method is demonstrated on simulated data and a real Mendelian randomization study concerning the effect of body mass index on health-related quality of life index. An R package sisVive is available on CRAN. Supplementary materials for this article are available online.

Abstract

Integrative genomics offers a promising approach to more powerful genetic association studies. The hope is that combining outcome and genotype data with other types of genomic information can lead to more powerful SNP detection. We present a new association test based on a statistical model that explicitly assumes that genetic variations affect the outcome through perturbing gene expression levels. It is shown analytically that the proposed approach can have more power to detect SNPs that are associated with the outcome through transcriptional regulation, compared to tests using the outcome and genotype data alone, and simulations show that our method is relatively robust to misspecification. We also provide a strategy for applying our approach to high-dimensional genomic data. We use this strategy to identify a potentially new association between a SNP and a yeast cell's response to the natural product tomatidine, which standard association analysis did not detect.

Abstract

This paper establish new restricted isometry conditions for compressed sensing and affine rank minimization. It is shown for compressed sensing that \(\delta^A_k+\theta^A_{k,k}<1\) guarantees the exact recovery of all \(k\)-sparse signals in the noiseless case through the constrained \(\ell_1\) minimization. Furthermore, the upper bound 1 is sharp in the sense that for any \(\epsilon>0\), the condition \(\delta^A_k+\theta^A_{k,k}<1+\epsilon\) is not sufficient to guarantee such exact recovery using any recovery method. Similarly, for affine rank minimization, if \(\delta^A_k+\theta^A_{k,k}<1\) then all matrices with rank at most \(r\) can be reconstructed exactly in the noiseless case via the constrained nuclear norm minimization; and for any \(\epsilon>0\), \(\delta^A_k+\theta^A_{k,k}<1+\epsilon\) does not ensure such exact recovery using any method. Moreover, in the noisy case the conditions \(\delta^A_k+\theta^A_{k,k}<1\) and \(\delta^A_k+\theta^A_{k,k}<1\) are also sufficient for the stable recovery of sparse signals and low-rank matrices respectively. Applications and extensions are also discussed.

Abstract

This paper establishes a sharp condition on the restrict isometry property (RIP) for both the sparse signal recovery and low-rank matrix recovery. It is shown that if the measurement matrix \(A\) satisfies the RIP condition \(\delta^A_k<1/3\), then all \(k\)-sparse signals \(\beta\) can be recovered exactly via the constrained \(\ell_1\) minimization based on \(y=A\beta\). Similarly, if the linear map \(M\) satisfies the RIP condition \(\delta^M_r<1/3\), then all matrices \(X\) of rank at most \(r\) can be recovered exactly via the constrained nuclear norm minimization based on \(b=M(X)\). Furthermore, in both cases it is not possible to do so in general when the condition does not hold. In addition, noisy cases are considered and oracle inequalities are given under the sharp RIP condition.

Abstract

Testing covariance structure is of significant interest in many areas of statistical analysis and construction of compressed sensing matrices is an important problem in signal processing. Motivated by these applications, we study in this paper the limiting laws of the coherence of an \(n\times p\) random matrix in the high-dimensional setting where \(p\) can be much larger than \(n\). Both law of large numbers and limiting distribution are derived. We then consider testing the bandedness of the covariance matrix of a high dimensional Gaussian distribution which includes testing for independence as a special case. The limiting laws of the coherence of the data matrix play a critical role in the construction of the test. The asymptotic results is also applied to the construction of compressed sensing matrices.

Abstract

We consider the orthogonal matching pursuit (OMP) algorithm for the recovery of a high-dimensional sparse signal based on a small number of noisy measurements. The OMP is an iterative greedy algorithm that selects at each step the variable which is most correlated with the current residuals. Stopping rules are given and the algorithm is fully data driven. It is shown that under conditions on the mutual incoherence and the minimum magnitude of the nonzero components of the signal, the support of the signal can be recovered exactly by the OMP algorithm with high probability. In addition, we also consider the problem of identifying significant variables in the case where some of the nonzero components are possibly small. It is shown that in this case the OMP algorithm will still select all the significant variables before possibly selecting incorrect ones. Moreover, with modified stopping rules, the OMP algorithm can ensure that no incorrect variables are selected.

Abstract

This article considers sparse signal recovery in the presence of noise. A mutual incoherence condition which was previously used for exact recovery in the noiseless case is shown to be sufficient for stable recovery in the noisy case. Both bounded noise and Gaussian noise settings are considered. Furthermore, the condition is shown to be sharp. In addition, an oracle inequality is given under the mutual incoherence condition.

Abstract

In this paper we present a concise and coherent analysis of the constrained \(\ell_1\) minimization method for stable recovering of high-dimensional sparse signals both in the noiseless case and noisy case. The analysis is surprisingly simple and elementary, while leads to strong results. In particular, it is shown that the sparse recovery problem can be solved via \(\ell_1\) minimization under weaker conditions than what is known in the literature. A key technical tool is an elementary inequality, called Shifting Inequality, which, for a given nonnegative decreasing sequence, bounds the \(\ell_2\) norm of a subsequence in terms of the \(\ell_1\) norm of another subsequence by shifting the elements to the upper end.

Abstract

This article considers constrained \(\ell_1\) minimization methods for the recovery of high dimensional sparse signals in three settings: noiseless, bounded error and Gaussian noise. A unified and elementary treatment is given in these noise settings for two \(\ell_1\) minimization methods: the Dantzig selector and \(\ell_1\) minimization with an \(\ell_2\) constraint. The results of this paper improve the existing results in the literature by weakening the conditions and tightening the error bounds. The improvement on the conditions shows that signals with larger support can be recovered accurately. This paper also establishes connections between restricted isometry property and the mutual incoherence property. Some results of Candes, Romberg and Tao (2006) and Donoho, Elad, and Temlyakov (2006) are extended.

Abstract

The rapid growth of digital textual data from a variety of fields has made it increasingly important to develop statistical methods for the analysis of such data with rigorous theoretical guarantees. In this paper, we focus on supervised topic modeling within the framework of generalized linear models (GLMs) and probabilistic latent semantic indexing (pLSI) models. One of the major challenges of the analysis is that the covariates are unobservable. We propose a novel bias-adjusted estimator of the covariates and use it to estimate the regression vector. We establish minimax optimal rates of convergence and show that the proposed estimator is rate-optimal up to a logarithmic factor. In addition, we consider statistical inference for individual regression coefficients and construct confidence intervals based on an asymptotically unbiased and normally distributed estimator. The effectiveness of our proposed algorithms is demonstrated through simulation studies and applications to the analysis of a movie review dataset and a gut microbiome dataset.

Abstract

Sparse topic modeling under the probabilistic latent semantic indexing (pLSI) model is studied. Novel and computationally fast algorithms for estimation and inference of both the word-topic matrix and the topic-document matrix are proposed and their theoretical properties are investigated. Both minimax upper and lower bounds are established and the results show that the proposed algorithms are rate- optimal, up to a logarithmic factor. Moreover, a refitting algorithm is proposed to establish asymptotic normality and construct valid confidence intervals for the individual entries of the word-topic and topic-document matrices. Simulation studies are carried out to investigate the numerical performance of the proposed algorithms. The results show that the proposed algorithms perform well numerically and are more accurate in a range of simulation settings comparing to the existing literature. In addition, the methods are illustrated through an analysis of the COVID-19 Open Research Dataset (CORD-19).

Abstract

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 multi-group classification and to classification under the Gaussian copula model.

Abstract

This paper aims to develop an optimality theory for linear discriminant analysis in the high-dimensional setting. A data-driven and tuning free classification rule, which is based on an adaptive constrained \(\ell_1\) minimization approach, is proposed and analyzed. Minimax lower bounds are obtained and this classification rule is shown to be simultaneously rate optimal over a collection of parameter spaces. In addition, we consider classification with incomplete data under the missing completely at random (MCR) model. An adaptive classifier with theoretical guarantees is introduced and optimal rate of convergence for high-dimensional linear discriminant analysis under the MCR model is established. The technical analysis for the case of missing data is much more challenging than that for the complete data. We establish a large deviation result for the generalized sample covariance matrix, which serves as a key technical tool and can be of independent interest. An application to lung cancer and leukemia studies is also discussed.

Abstract

We propose a general semi-supervised inference framework focused on the estimation of the population mean. We consider both the ideal semi-supervised setting where infinitely many unlabeled samples are available, as well as the ordinary semi-supervised setting in which only a finite number of unlabeled samples is available. As usual in semi-supervised settings, there exists an unlabeled sample of covariate vectors and a labeled sample consisting of covariate vectors along with real-valued responses (“labels"). Otherwise the formulation is 'assumption-lean" in that no major conditions are imposed on the statistical or functional form of the data. Estimators are proposed along with corresponding confidence intervals for the population mean. Theoretical analysis on both the asymptotic behavior and \(\ell_2\)-risk for the proposed procedures are given. Surprisingly, the proposed estimators, based on a simple form of the least squares method, outperform the ordinary sample mean. The method is further extended to a nonparametric setting, in which the oracle rate can be achieved asymptotically. The proposed estimators are further illustrated by simulation studies and a real data example involving estimation of the homeless population.

Abstract

Unsupervised learning is an important problem in statistics and machine learning with a wide range of applications. In this paper, we study clustering of high-dimensional Gaussian mixtures and propose a procedure, called CHIME, that is based on the EM algorithm and a direct estimation method for the sparse discriminant vector. Both theoretical and numerical properties of CHIME are investigated. We establish the optimal rate of convergence for the excess mis-clustering error and show that CHIME is minimax rate optimal. In addition, the optimality of the proposed estimator of the discriminant vector is also established. Simulation studies show that CHIME outperforms the existing methods under a variety of settings. The proposed CHIME procedure is also illustrated in an analysis of a glioblastoma gene expression data set and shown to have superior performance. Clustering of Gaussian mixtures in the conventional low-dimensional setting is also considered. The technical tools developed for the high-dimensional setting are used to establish the optimality of the clustering procedure that is based on the classical EM algorithm.

Abstract

This paper considers sparse linear discriminant analysis of high-dimensional data. In contrast to the existing methods which are based on separate estimation of the precision matrix Ω and the difference δ of the mean vectors, we introduce a simple and effective classifier by estimating the product Ωδ directly through constrained \(\ell_1\) minimization. The estimator can be implemented efficiently using linear programming and the resulting classifier is called the linear programming discriminant (LPD) rule. The LPD rule is shown to have desirable theoretical and numerical properties. It exploits the approximate sparsity of Ωδ and as a consequence allows cases where it can still perform well even when Ω and/or δ cannot be estimated consistently. Asymptotic properties of the LPD rule are investigated and consistency and rate of convergence results are given. The LPD classifier has superior finite sample performance and significant computational advantages over the existing methods that require separate estimation of Ω and δ. The LPD rule is also applied to analyze real datasets from lung cancer and leukemia studies. The classifier performs favorably in comparison to existing methods.

Abstract

We study the misclassification error for community detection in general heterogeneous stochastic block models (SBM) with noisy or partial label information. We establish a connection between the misclassification rate and the notion of minimum energy on the local neighborhood of the SBM. We develop an optimally weighted message passing algorithm to reconstruct labels for SBM based on the minimum energy flow and the eigenvectors of a certain Markov transition matrix. The general SBM considered in this paper allows for unequal-size communities, degree heterogeneity, and different connection probabilities among blocks. We focus on how to optimally weigh the message passing to improve misclassification.

Abstract

In this paper, we study detection and fast reconstruction of the celebrated Watts-Strogatz (WS) small-world random graph model (Watts and Strogatz, 1998) which aims to describe real-world complex networks that exhibit both high clustering and short average length properties. The WS model with neighborhood size \(k\) and rewiring probability \(\beta\) can be viewed as a continuous interpolation between a deterministic ring lattice graph and the Erdös-Rényi random graph. We study both the computational and statistical aspects of detecting the deterministic ring lattice structure (or local geographical links, strong ties) in the presence of random connections (or long range links, weak ties), and for its recovery. The phase diagram in terms of \((k,\beta)\) is partitioned into several regions according to the difficulty of the problem. We propose distinct methods for the various regions.

Abstract

Community detection, which aims to cluster \(N\) nodes in a given graph into \(r\) distinct groups based on the observed undirected edges, is an important problem in network data analysis. In this paper, the popular stochastic block model (SBM) is extended to the generalized stochastic block model (GSBM) that allows for adversarial outlier nodes, which are connected with other nodes in the graph in an arbitrary way. Under this model, we introduce a procedure using convex optimization followed by \(k\)-means algorithm with \(k=r\). Both theoretical and numerical properties of the method are analyzed. A theoretical guarantee is given for our methodology to accurately detect the communities with small misclassification rate under the setting where the number of clusters can grow with \(N\). This theoretical result admits to the best known result in the literature of computationally feasible community detection. Numerical results show that our method is both computationally fast and robust to different kinds of outliers, while some popular computationally fast community detection algorithms, such as spectral clustering applied to adjacency matrices or graph Laplacians, may fail due to a very small portion of outliers. We apply a slight modification of our method to a political blogs data set, showing that our method is competent in practice, and comparable to existing computationally feasible methods in the literature. To the best of the authors' knowledge, our result is the first in the literature in terms of clustering communities with fast growing numbers under the generalized stochastic block model where a portion of arbitrary outlier nodes exist.

Abstract

Motivated by a range of applications in genomics and engineering, we consider in this paper detection and identification of sparse signal segments in three settings: signals with known shape, arbitrary signals, and smooth signals. Optimal rates of detection are established for the three cases and rate-optimal detectors are constructed. The detectors are easily implementable and are based on scanning with linear and quadratic statistics. Identification of the location of signal segments is also considered. Our analysis reveals both similarities and differences in the strategy and fundamental difficulty of detection and identification among these three settings.

Abstract

This paper investigates the fundamental limits for detecting a high-dimensional sparse matrix contaminated by white Gaussian noise from both the statistical and computational perspectives. We consider \(p\times p\) matrices whose rows and columns are individually \(k\)-sparse. We provide a tight characterization of the statistical and computational limits for sparse matrix detection, which precisely describe when achieving optimal detection is easy, hard, or impossible, respectively. Although the sparse matrices considered in this paper have no apparent submatrix structure and the corresponding estimation problem has no computational issue at all, the detection problem has a surprising computational barrier when the sparsity level \(k\) exceeds the cubic root of the matrix size \(p\): attaining the optimal detection boundary is computationally at least as hard as solving the planted clique problem. The same statistical and computational limits also hold in the sparse covariance matrix model, where each variable is correlated with at most \(k\) others. A key step in the construction of the statistically optimal test is a structural property for sparse matrices, which can be of independent interest.

Abstract

This paper studies the problem of detecting dependence between two mixture distributions, motivated by questions arising from statistical genomics. The fundamental limits of detecting weak positive dependence are derived and an oracle test statistic is proposed. It is shown that for mixture distributions whose components are stochastically ordered, the oracle test statistic is asymptotically optimal. Connections are drawn between dependency detection and signal detection, where the goal of the latter is to detect the presence of non-null components in a single mixture distribution. It is shown that the oracle test for dependency can also be used as a signal detection procedure in the two-sample setting, and there can achieve detection even when detection using each sample separately is provably impossible. A nonparametric data-adaptive test statistic is then proposed, and its closed-form asymptotic distribution under the null hypothesis of independence is established. Simulations show that the adaptive procedure performs as well as the oracle test statistic, and that both can be more powerful than existing methods. In an application to the analysis of the shared genetic basis of psychiatric disorders, the adaptive test is able to detect genetic relationships not detected by other procedures.

Abstract

Genome-wide association studies (GWAS) and differential expression analyses have had limited success in finding genes that cause complex diseases such as heart failure (HF), a leading cause of death in the United States. This paper proposes a new statistical approach that integrates GWAS and expression quantitative trait loci (eQTL) data to identify important HF genes. For such genes, genetic variations that perturb its expression are also likely to influence disease risk. The proposed method thus tests for the presence of simultaneous signals: SNPs that are associated with the gene's expression as well as with disease. An analytic expression for the \(p\)-value is obtained, and the method is shown to be asymptotically adaptively optimal under certain conditions. It also allows the GWAS and eQTL data to be collected from different groups of subjects, enabling investigators to integrate public resources with their own data. Simulation experiments show that it can be more powerful than standard approaches and also robust to linkage disequilibrium between variants. The method is applied to an extensive analysis of HF genomics and identifies several genes with biological evidence for being functionally relevant in the etiology of HF. It is implemented in the R package ssa.

Abstract

The interplay between computational efficiency and statistical accuracy in high-dimensional inference has drawn increasing attention in the literature. In this paper, we study computational and statistical boundaries for submatrix localization. Given one observation of a signal submatrix (of magnitude \(\lambda\) and size \(k_m\times k_n\)) contaminated with a noise matrix (of size \(m\times n\)), we establish two transition thresholds for the signal to noise ratio \(\lambda/\sigma\) in terms of \(m,n,k_m\), and \(k_n\). The first threshold, \(\mathrm{SNR}_c\), corresponds to the computational boundary. Below this threshold, it is shown that no polynomial time algorithm can succeed in identifying the submatrix, under the hidden clique hypothesis. We introduce an adaptive linear time algorithm that identifies the submatrix with high probability when the signal strength is above the threshold \(\mathrm{SNR}_c\). The second threshold, \(\mathrm{SNR}_s\), captures the statistical boundary, below which no method can succeed with probability going to one in the minimax sense. The exhaustive search method successfully finds the submatrix above this threshold. The results show an interesting phenomenon that \(\mathrm{SNR}_c\) is always significantly larger than \(\mathrm{SNR}_s\), which implies an essential gap between statistical optimality and computational efficiency for submatrix localization.

Abstract

A common feature in large-scale scientific studies is that signals are sparse and it is desirable to significantly narrow down the focus to a much smaller subset in a sequential manner. In this paper, we consider two related data screening problems: One is to find the smallest subset such that it virtually contains all signals and another is to find the largest subset such that it essentially contains only signals. These screening problems are closely connected to but distinct from the more conventional signal detection or multiple testing problems. We develop data-driven screening procedures which control the error rates with near optimality properties and study how to design the experiments efficiently to achieve the goals in data screening. An application to multistage high-throughput studies is given to illustrate the merits of the proposed screening methods.

Abstract

Copy number variations (CNVs) are alterations of the DNA of a genome that results in the cell having an abnormal number of copies of one or more sections of the DNA. Germline CNVs have been shown to be associated with many complex diseases. Detecting and identifying all the CNVs in a given sample or in multiple population-based samples is an important first step in many CNV analyses. In this chapter, we review statistical methods for CNV identification, focusing on latest developed methods for sparse segment identifications in various settings. We review methods for optimal CNV identification for a single sample based on SNP allele intensity data, methods for robust CNV identification based on the next generation sequence (NGS) data. and methods for detection of recurrent CNVs in a population when a large set of samples are available. Our review focuses on problem formulations and optimal statistical properties of the procedures. We illustrate these methods using data from the 1000 Genomes Project and data from a large genome-wide association study of neuroblastoma. Areas that need further research are also presented.

Abstract

This paper considers testing the equality of multiple high-dimensional mean vectors under dependency. We propose a test that is based on a linear transformation of the data by the precision matrix which incorporates the dependence structure of the variables. The limiting null distribution of the test statistic is derived and is shown to be the extreme value distribution of type I. The convergence to the limiting distribution is, however, slow when the number of groups is relatively large. An intermediate correction factor is introduced which significantly improves the accuracy of the test. It is shown that the test is particularly powerful against sparse alternatives and enjoys certain optimality. A simulation study is carried out to examine the numerical performance of the test and compare with other tests given in the literature. The numerical results show that the proposed test significantly outperforms those tests against sparse alternatives.

Abstract

This paper considers in the high dimensional setting a canonical testing problem in multivariate analysis, namely testing the equality of two mean vectors. We introduce a new test statistic that is based on a linear transformation of the data by the precision matrix which incorporates the correlations among the variables. Limiting null distribution of the test statistic and the power of the test are analyzed. It is shown that the test is particularly powerful against sparse alternatives and enjoys certain optimality. A simulation study is carried out to examine the numerical performance of the test and compare with other tests given in the literature. The results show that the proposed test significantly outperforms those tests in a range of settings.

Abstract

Detection of sparse signals arises in a wide range of modern scientific studies. The focus so far has been mainly on Gaussian mixture models. In this paper, we consider the detection problem under a general sparse mixture model and obtain an explicit expression for the detection boundary. It is shown that the fundamental limits of detection is governed by the behavior of the log-likelihood ratio evaluated at an appropriate quantile of the null distribution. We also establish the adaptive optimality of the higher criticism procedure across all sparse mixtures satisfying certain mild regularity conditions. In particular, the general results obtained in this paper recover and extend in a unified manner the previously known results on sparse detection far beyond the conventional Gaussian model and other exponential families.

Abstract

Copy number variant is an important type of genetic structural variation appearing in germline DNA, ranging from common to rare in a population. Both rare and common copy number variants have been reported to be associated with complex diseases, and it is therefore important to simultaneously identify both based on a large set of population samples. We develop a proportion adaptive segment selection procedure that automatically adjusts to the unknown proportions of the carriers of the segment variants. We characterize the detection boundary that separates the region where a segment variant is detectable by some method from the region where it cannot be detected. Although the detection boundaries are very different for the rare and common segment variants, it is shown that the proposed procedure can reliably identify both whenever they are detectable. Compared with methods for single sample analysis, this procedure gains power by pooling information from multiple samples. The method is applied to analyze neuroblastoma samples and identifies a large number of copy number variants that are missed by single-sample methods.

Abstract

Motivated by copy number variation analysis based on next generation sequencing data, we consider the problem of detecting and identifying sparse short segments hidden in an ultra long linear sequence of data with unspecified noise distribution. Based on a local median transformation, we propose a computationally efficient method called robust segment identifier (RSI), which provides a robust and optimal solution for segment identification over a wide range of noise distributions. We theoretically quantify the conditions for detecting the segment signals and show that the RSI consistently estimates the signal segments whenever it is possible to detect their existence. We present simulations to demonstrate the effect of data transformation and the efficiency of our method under different noise distributions. We also present results from an application to copy number variant analysis using next generation sequencing data of the HapMap Yoruban sample NA19240 to further illustrate the theory and the methods.

Abstract

The problem of detecting heterogeneous and heteroscedastic Gaussian mixtures is considered. The focus is on how the parameters of heterogeneity, heteroscedasticity, and proportion of non-null component influence the difficulty of the problem. We establish an explicit detection boundary which separates the detectable region where the likelihood ratio test is shown to reliably detect the presence of non-null effect, from the undetectable region where no method can do so. In particular, the results show that the detection boundary changes dramatically when the proportion of non-null component shifts from the sparse regime to the dense regime. Furthermore, it is shown that the Higher Criticism test, which does not require the specific information of model parameters, is optimally adaptive to the unknown degrees of heterogeneity and heteroscedasticity in both the sparse and dense cases.

Abstract

Genome-wide association studies have successfully identified hundreds of novel genetic variants associated with many complex human diseases. However, there is a lack of rigorous work on evaluating the statistical power for identifying these variants. In this paper, we consider the problem of sparse signal identification in genome-wide association studies and present two analytical frameworks for detailed analysis of the statistical power for detecting and identifying the disease-associated variants. We present an explicit sample size formula for achieving a given false non-discovery rate while controlling the false discovery rate based on an optimal false discovery rate procedure. The problem of sparse genetic variants recovery is also considered and a boundary condition is established in terms of sparsity and signal strength for almost exact recovery of disease-associated variants as well as nondisease-associated variants. A data-adaptive procedure is proposed to achieve this bound. These results provide important tools for sample size calculation and power analysis for large-scale multiple testing problems. The analytical results are illustrated with a genome-wide association study of neuroblastoma.

Abstract

Motivated by DNA copy number analysis in genetics, we consider the problem of detecting and identifying sparse short segments in a long one-dimensional sequence of data with additive Gaussian white noise, where the number, length and location of the segments are unknown. We present a statistical characterization of the identifiable region of a segment where it is possible to reliably separate the segment from Gaussian noise. An efficient likelihood ratio selection (LRS) procedure for identifying the segments is developed and the asymptotic optimality of this method is presented in the sense that the LRS can separate the signal segments from the noise as long as the signal segments can be estimated. The proposed method is demonstrated with simulations and analysis of a real data set on identification of copy number variants based on high-density single nucleotide polymorphisms (SNP) data. The results show that the LRS procedure can yield greater gain in power of detecting the true segments than some standard signal identification methods.

Abstract

For high dimensional statistical models, researchers have begun to focus on situations which can be described as having relatively few moderately large coefficients. Such situations lead to some very subtle statistical problems. In particular, Ingster and Donoho and Jin have considered a sparse normal means testing problem, in which they described the precise demarcation, or the detection boundary. Meinshausen and Rice have shown that it is even possible to estimate consistently the fraction of nonzero coordinates on a subset of the detectable region, but leave unanswered the question of exactly which parts of the detectable region that consistent estimation is possible. In the present paper we develop a new approach for estimating the fraction of nonzero means for problems where the nonzero means are moderately large. We show that the detection region described by Ingster and Donoho and Jin turns out to be the region where it is possible to consistently estimate the expected fraction of nonzero coordinates. This theory is developed further and minimax rates of convergence are derived. A procedure is constructed which attains the optimal rate of convergence in this setting. Furthermore, the procedure also provides an honest lower bound for confidence intervals while minimizing the expected length of such an interval. Simulations are used to enable comparison with the work of Meinshausen and Rice, where a procedure is given but where rates of convergence have not been discussed. Extensions to more general Gaussian mixture models are also given.

Abstract

Population-level single-cell gene expression data captures the gene expressions of thousands of cells for each individual within a sizable cohort. This data enables the construction of cell-type- and individual-specific gene co-expression network by estimating the covariance matrices. It is important to understand how such co-expression networks are associated with individual-level covariates. This paper considers Fréchet regression with multivariate Gaussian distribution as an outcome and vector covariates, where the Wasserstein distance between distributions is used as a replacement for the Euclidean distance. A test statistic is defined based on Fréchet mean and covariate weighted Fréchet mean. The asymptotic distribution of the test statistic is derived under the assumption of simultaneously diagonalizable covariance matrices. Although the proposed test statistic is motivated by considering the multivariate normal distribution as the outcome, it can be applied for testing the association between covariance matrices and covariates, where permutation can be used for assessing its statistical significance. Simulations show that the proposed test has correct type 1 error and adequate power. Results from an analysis of large-scale single-cell data reveal an association between the gene co-expression network of genes in the nutrient sensing pathway and age, indicating the perturbed gene co-expression network as people age.

Abstract

Phylogenetic association analysis is an essential and powerful tool for studying the association between microbial compositions and the outcome of interest in microbiome studies. However, existing methods for testing such associations are more sensitive to a linear association in a high-dimensional setting and the assumptions of confounding effects. Methods that are capable of characterizing complex association, including non-monotonic association, are therefore needed. This paper proposes a new phylogenetic association analysis framework to address these challenges. The new framework introduces conditional rank correlation as a measure of association to detect a wide range of dependencies, which is robust to the outlier, accounts for confounders in a fully nonparametric way. The new framework aggregates conditional rank correlations for subtrees as the weighted sum and maximum to capture dense and sparse signals. To determine the significance level, we calibrate the test statistics by a nearest neighbor bootstrapping method, which is easy to use and can incorporate extra data sets when available. The practical merits of the new framework are demonstrated by numerical experiments using both simulated and real microbiome data sets.

Abstract

Motivated by applications in text mining and discrete distribution inference, we investigate the testing for equality of probability mass functions of K groups of high-dimensional multinomial distributions. Special cases of this problem include global testing for topic models, two-sample testing in authorship attribution, and closeness testing for discrete distributions. A test statistic, which is shown to have an asymptotic standard normal distribution under the null, is proposed. This parameter-free null distribution holds true without requiring identical multinomial parameters within each group or equal group sizes. The optimal detection boundary for this testing problem is established, and the proposed test is shown to achieve this optimal detection boundary across the entire parameter space of interest. The proposed method is demonstrated in simulation studies and applied to analyze two real-world datasets to examine, respectively, variation among consumer reviews of Amazon movies and the diversity of statistical paper abstracts.

Abstract

Motivated by applications in single-cell biology and metagenomics, we investigate the problem of matrix reordering based on a noisy disordered monotone Toeplitz matrix model. We establish the fundamental statistical limit for this problem in a decision-theoretic framework and demonstrate that a constrained least squares estimator achieves the optimal rate. However, due to its computational complexity, we analyze a popular polynomial-time algorithm, spectral seriation, and show that it is suboptimal. To address this, we propose a novel polynomial-time adaptive sorting algorithm with guaranteed performance improvement. Simulations and analyses of two real single-cell RNA sequencing datasets demonstrate the superiority of our algorithm over existing methods.

Abstract

This paper investigates the theoretical foundations of the t-distributed stochastic neighbor embedding (t-SNE) algorithm, a popular nonlinear dimension reduction and data visualization method. A novel theoretical framework for the analysis of t-SNE based on the gradient descent approach is presented. For the early exaggeration stage of t-SNE, we show its asymptotic equivalence to power iterations based on the underlying graph Laplacian, characterize its limiting behavior, and uncover its deep connection to Laplacian spectral clustering, and fundamental principles including early stopping as implicit regularization. The results explain the intrinsic mechanism and the empirical benefits of such a computational strategy. For the embedding stage of t-SNE, we characterize the kinematics of the low-dimensional map throughout the iterations, and identify an amplification phase, featuring the intercluster repulsion and the expansive behavior of the low-dimensional map, and a stabilization phase. The general theory explains the fast convergence rate and the exceptional empirical performance of t-SNE for visualizing clustered data, brings forth the interpretations of the t-SNE visualizations, and provides theoretical guidance for applying t-SNE and selecting its tuning parameters in various applications.

Abstract

Integrating the summary statistics from a genome-wide association study and expression quantitative trait loci data provides a powerful way of identifying genes with expression levels that are potentially associated with complex diseases. We introduce a parameter called T-score that quantifies the genetic overlap between a gene and the disease phenotype based on the summary statistics, based on the mean values of two Gaussian sequences. Specifically, given two independent samples \(X_n\sim N(\theta,\Sigma_1)\) and \(Y_n\sim N(\mu,\Sigma_2)\), the T-score is defined as \(\sum_{i=1}^n |\theta_i\mu_i|\), a nonsmooth functional, that characterizes the number of shared signals between two absolute normal mean vectors \(|\theta|\) and \(|\mu|\). Using approximation theory, estimators are constructed and shown to be minimax rate-optimal and adaptive over various parameter spaces. Simulation studies demonstrate the superiority of the proposed estimators over existing methods. Lastly, the method is applied to an integrative analysis of heart failure genomics data sets and we identify several genes and biological pathways that are potentially causal to human heart failure.

Abstract

This paper focuses on the non-asymptotic concentration of the heteroskedastic Wishart-type matrices. Suppose $Z$ is a $p_1$-by-$p_2$ random matrix and $Z_{ij} \sim N(0,\sigma_{ij}^2)$ independently, we prove the expected spectral norm of Wishart matrix deviations (i.e., $\mathbb{E} \left\|ZZ^\top - \mathbb{E} ZZ^\top\right\|$) is upper bounded by \[ (1+\epsilon)\left\{2\sigma_C\sigma_R + \sigma_C^2 + C\sigma_R\sigma_*\sqrt{\log(p_1 \wedge p_2)} + C\sigma_*^2\log(p_1 \wedge p_2)\right\}, \] where $\sigma_C^2 := \max_j \sum_{i=1}^{p_1}\sigma_{ij}^2$, $\sigma_R^2 := \max_i \sum_{j=1}^{p_2}\sigma_{ij}^2$ and $\sigma_*^2 := \max_{i,j}\sigma_{ij}^2$. A minimax lower bound is developed that matches this upper bound. Then, we derive the concentration inequalities, moments, and tail bounds for the heteroskedastic Wishart-type matrix under more general distributions, such as sub-Gaussian and heavy-tailed distributions. Next, we consider the cases where $Z$ has homoskedastic columns or rows (i.e., $\sigma_{ij} \approx \sigma_i$ or $\sigma_{ij} \approx \sigma_j$) and derive the rate-optimal Wishart-type concentration bounds. Finally, we apply the developed tools to identify the sharp signal-to-noise ratio threshold for consistent clustering in the heteroskedastic clustering problem.

Abstract

Consider a standard white Wishart matrix with parameters \(n\) and \(p\). Motivated by applications in high-dimensional statistics and signal processing, we perform asymptotic analysis on the maxima and minima of the eigenvalues of all the \(m\) by \(m\) principal minors, under the asymptotic regime that \(n,p,m\) go to infinity. Asymptotic results concerning extreme eigenvalues of principal minors of real Wigner matrices are also obtained. In addition, we discuss an application of the theoretical results to the construction of compressed sensing matrices, which provides insights in high dimensional linear regression.

Abstract

Motivated by the problem of estimating the bacterial growth rates for genome assemblies from shotgun metagenomic data, we consider the permuted monotone matrix model \(Y=\Theta\Pi+Z\), where \(Y\in\mathbb{R}^{n\times p}\) is observed, \(\Theta\in\mathbb{R}^{n\times p}\) is an unknown approximately rank-one signal matrix with monotone rows, \(\Pi\in\mathbb{R}^{p\times p}\) is an unknown permutation matrix, and \(Z\in\mathbb{R}^{n\times p}\) is the noise matrix. This paper studies the estimation of the extreme values associated to the signal matrix Θ, including its first and last columns, as well as their difference (the range vector). Treating these estimation problems as compound decision problems, minimax rate-optimal estimators are constructed using the spectral column sorting method. Numerical experiments through simulated and synthetic microbiome metagenomic data are presented, showing the superiority of the proposed methods over the alternatives. The methods are illustrated by comparing the growth rates of gut bacteria between inflammatory bowel disease patients and normal controls.

Abstract

Motivated by recent research on quantifying bacterial growth dynamics based on genome assemblies, we consider a permuted monotone matrix model \(Y=\Theta\Pi+Z\), where the rows represent different samples, the columns represent contigs in genome assemblies and the elements represent log-read counts after preprocessing steps and Guanine-Cytosine (GC) adjustment. In this model, Θ is an unknown mean matrix with monotone entries for each row, Π is a permutation matrix that permutes the columns of Θ, and Z is a noise matrix. This paper studies the problem of estimation/recovery of \(\Pi\) given the observed noisy matrix \(Y\). We propose an estimator based on the best linear projection, which is shown to be minimax rate-optimal for both exact recovery, as measured by the 0-1 loss, and partial recovery, as quantified by the normalized Kendall’s tau distance. Simulation studies demonstrate the superior empirical performance of the proposed estimator over alternative methods. We demonstrate the methods using a synthetic metagenomics dataset of 45 closely related bacterial species and a real metagenomic dataset to compare the bacterial growth dynamics between the responders and the non-responders of the IBD patients after 8 weeks of treatment.

Abstract

Weighted UniFrac distance, a plug-in estimator of the Wasserstein distance of read counts on a tree, has been widely used to measure the microbial community difference in microbiome studies. Our investigation however shows that such a plug-in estimator, although intuitive and commonly used in practice, suffers from potential bias. Motivated by this finding, we study the problem of optimally estimating the Wasserstein distance between two distributions on a tree from the sampled data in the high-dimensional setting. Minimax rate of convergence is established. To overcome the bias problem, we introduce a new estimator, referred to as moment-screening estimator on a tree (MET), by conducting implicit best polynomial approximation that incorporates the tree structure. The new estimator is computationally efficient and is shown to be minimax rate-optimal. Numerical studies using both simulated and real biological datasets demonstrate the practical merits of MET, including reduced biases and statistically more significant differences in microbiome between inactive Crohn's disease patients and the normal controls.

Abstract

Quantitative comparison of microbial composition from different populations is a fundamental task in various microbiome studies. We consider two-sample testing for microbial compositional data by leveraging the phylogenetic tree information. Motivated by existing phylogenetic distances, we take a minimum-cost flow perspective to study such testing problems. Our investigation shows that multivariate analysis of variance with permutation using phylogenetic distances, one of the most commonly used methods in practice, is essentially a sum-of-squares type test and has better power for dense alternatives. However, empirical evidence from real data sets suggests that the phylogenetic microbial composition difference between two populations is usually sparse. Motivated by this observation, we propose a new maximum type test, Detector of Active Flow on a Tree, and investigate its properties. It is shown that the proposed method is particularly powerful against sparse phylogenetic composition difference and enjoys certain optimality. The practical merit of the proposed method is demonstrated by simulation studies and an application to a human intestinal biopsy microbiome data set for patients with ulcerative colitis.

Abstract

Quantum state tomography, which aims to reconstruct quantum states described by density matrices, is becoming increasingly important in many scientific studies involving quantum systems. This paper considers the reconstruction of high-dimensional low-rank density matrices based on Pauli measurements. In particular it focuses on estimation of eigenspace for a large low-rank density matrix. Both ordinary principal component analysis (PCA) and iterative thresholding sparse PCA (ITSPCA) are studied and optimal rates of convergence are established. Minimax lower bounds for eigenspace estimation under the spectral and Frobenius norms are derived. It is shown that the convergence rates of the ITSPCA algorithm matches the minimax lower bounds and the procedure is thus rate-optimal. With these PCA estimators, we reconstruct the large low-rank density matrix and obtain the optimal convergence rate. A simulation study is carried out to investigate the finite sample performance of the proposed estimators of the density matrices.

Abstract

We develop adaptive estimation and inference methods for high-dimensional Gaussian copula regression that achieve the same performance without the knowledge of the marginal transformations as that for high-dimensional linear regression. Using a Kendall's tau based covariance matrix estimator, an \(\ell_1\)-regularized estimator is proposed and a corresponding de-biased estimator is developed for the construction of the confidence intervals and hypothesis tests. Theoretical properties of the procedures are studied and the proposed estimation and inference methods are shown to be adaptive to the unknown monotone marginal transformations. Prediction of the response for a given value of the covariates is also considered. The procedures are easy to implement and perform well numerically. The methods are also applied to analyze the Communities and Crime Unnormalized Data from the UCI Machine Learning Repository.

Abstract

Perturbation bounds for singular spaces, in particular Wedin's \(\sin\Theta\) theorem, are a fundamental tool in many fields including high-dimensional statistics, machine learning, and applied mathematics. In this paper, we establish separate perturbation bounds, measured in both spectral and Frobenius \(\sin\Theta\) distances, for the left and right singular subspaces. Lower bounds, which show that the individual perturbation bounds are rate-optimal, are also given. The new perturbation bounds are applicable to a wide range of problems. In this paper, we consider in detail applications to low-rank matrix denoising and singular space estimation, high-dimensional clustering, and canonical correlation analysis (CCA). In particular, separate matching upper and lower bounds are obtained for estimating the left and right singular spaces. To the best of our knowledge, this is the first result that gives different optimal rates for the left and right singular spaces under the same perturbation. In addition to these problems, applications to other high-dimensional problems such as community detection in bipartite networks, multidimensional scaling, and cross-covariance matrix estimation are also discussed.

Abstract

The interplay between computational efficiency and statistical accuracy in high-dimensional inference has drawn increasing attention in the literature. In this paper, we study computational and statistical boundaries for submatrix localization. Given one observation of a signal submatrix (of magnitude \(\lambda\) and size \(k_m\times k_n\)) contaminated with a noise matrix (of size \(m\times n\)), we establish two transition thresholds for the signal to noise ratio \(\lambda/\sigma\) in terms of \(m,n,k_m\), and \(k_n\). The first threshold, \(\mathrm{SNR}_c\), corresponds to the computational boundary. Below this threshold, it is shown that no polynomial time algorithm can succeed in identifying the submatrix, under the hidden clique hypothesis. We introduce an adaptive linear time algorithm that identifies the submatrix with high probability when the signal strength is above the threshold \(\mathrm{SNR}_c\). The second threshold, \(\mathrm{SNR}_s\), captures the statistical boundary, below which no method can succeed with probability going to one in the minimax sense. The exhaustive search method successfully finds the submatrix above this threshold. The results show an interesting phenomenon that \(\mathrm{SNR}_c\) is always significantly larger than \(\mathrm{SNR}_s\), which implies an essential gap between statistical optimality and computational efficiency for submatrix localization.

Abstract

This paper considers the noisy sparse phase retrieval problem: recovering a sparse signal \(\mathbf{x}\in\mathbb{R}^p\) from noisy quadratic measurements \(y_j=(\mathbf{a}_j^\top\mathbf{x})^2+\epsilon_j,\ j=1,\ldots,m\), with independent sub-exponential noise \(\epsilon_j\). The goals are to understand the effect of the sparsity of \(\mathbf{x}\) on the estimation precision and to construct a computationally feasible estimator to achieve the optimal rates. Inspired by the Wirtinger Flow method for noiseless phase retrieval proposed and analyzed in [12], a novel thresholded gradient descent algorithm is proposed and it is shown to adaptively achieve the minimax rates of convergence over a wide range of sparsity levels when the \(\mathbf{a}_j\)'s are independent standard Gaussian random vectors, provided that the sample size is sufficiently large compared to the sparsity of \(\mathbf{x}\).

Abstract

This paper presents a unified geometric framework for the statistical analysis of a general ill-posed linear inverse model which includes as special cases noisy compressed sensing, sign vector recovery, trace regression, orthogonal matrix estimation, and noisy matrix completion. We propose computationally feasible convex programs for statistical inference including estimation, confidence intervals and hypothesis testing. A theoretical framework is developed to characterize the local estimation rate of convergence and to provide statistical inference guarantees. Our results are built based on the local conic geometry and duality. The difficulty of statistical inference is captured by the geometric characterization of the local tangent cone through the Gaussian width and Sudakov estimate.

Abstract

In this paper, an over-sampled periodogram higher criticism (OPHC) test is proposed for the global detection of sparse periodic effects in a complex-valued time series. An explicit minimax detection boundary is established between the number and magnitude of the complex sinusoids hidden in the series. The OPHC test is shown to be asymptotically powerful in the detectable region. Numerical simulations illustrate and verify the effectiveness of the proposed test. Furthermore, the periodogram over-sampled by \(O(\log N)\) is proven universally optimal in global testing for periodicities under a mild minimum separation condition. Connections to the problem of detecting a stream of pulses from frequency measurements in signal processing is also discussed.

Abstract

Quantum state tomography aims to determine the state of a quantum system as represented by a density matrix. It is a fundamental task in modern scientific studies involving quantum systems. In this paper, we study estimation of high-dimensional density matrices based on a relatively small number of Pauli measurements. In particular, under appropriate notion of sparsity, we establish the minimax optimal rates of convergence for estimation of the density matrix under both the spectral and Frobenius norm losses; and show how these rates can be achieved by a common thresholding approach. Numerical performance of the proposed estimator is also investigated.

Abstract

This paper considers testing the equality of multiple high-dimensional mean vectors under dependency. We propose a test that is based on a linear transformation of the data by the precision matrix which incorporates the dependence structure of the variables. The limiting null distribution of the test statistic is derived and is shown to be the extreme value distribution of type I. The convergence to the limiting distribution is, however, slow when the number of groups is relatively large. An intermediate correction factor is introduced which significantly improves the accuracy of the test. It is shown that the test is particularly powerful against sparse alternatives and enjoys certain optimality. A simulation study is carried out to examine the numerical performance of the test and compare with other tests given in the literature. The numerical results show that the proposed test significantly outperforms those tests against sparse alternatives.

Abstract

This paper considers in the high dimensional setting a canonical testing problem in multivariate analysis, namely testing the equality of two mean vectors. We introduce a new test statistic that is based on a linear transformation of the data by the precision matrix which incorporates the correlations among the variables. Limiting null distribution of the test statistic and the power of the test are analyzed. It is shown that the test is particularly powerful against sparse alternatives and enjoys certain optimality. A simulation study is carried out to examine the numerical performance of the test and compare with other tests given in the literature. The results show that the proposed test significantly outperforms those tests in a range of settings.

Abstract

This paper studies the asymptotic behaviors of the pairwise angles among \(n\) randomly and uniformly distributed unit vectors in \(\mathbb{R}^p\) as the number of points \(n\to\infty\), while the dimension \(p\) is either fixed or growing with \(n\). For both settings, we derive the limiting empirical distribution of the random angles and the limiting distributions of the extreme angles. The results reveal interesting differences in the two settings and provide a precise characterization of the folklore that "all high-dimensional random vectors are almost always nearly orthogonal to each other". Applications to statistics and machine learning and connections with some open problems in physics and mathematics are also discussed.

Abstract

The coherence of a random matrix, which is defined to be the largest magnitude of the Pearson correlation coefficients between the columns of the random matrix, is an important quantity for a wide range of applications including high-dimensional statistics and signal processing. Inspired by these applications, this paper studies the limiting laws of the coherence of \(n\times p\) random matrices for a full range of the dimension \(p\) with a special focus on the ultra high-dimensional setting. Assuming the columns of the random matrix are independent random vectors with a common spherical distribution, we give a complete characterization of the behavior of the limiting distributions of the coherence. More specifically, the limiting distributions of the coherence are derived separately for three regimes: \(\log(p)/n\to0\), \(\log(p)/n\to\beta\in(0,\infty)\), and \(\log(p)/n\to\infty\). The results show that the limiting behavior of the coherence differs significantly in different regimes and exhibits interesting phase transition phenomena as the dimension \(p\) grows as a function of \(n\). Applications to statistics and compressed sensing in the ultra high-dimensional setting are also discussed.

Abstract

Testing covariance structure is of significant interest in many areas of statistical analysis and construction of compressed sensing matrices is an important problem in signal processing. Motivated by these applications, we study in this paper the limiting laws of the coherence of an \(n\times p\) random matrix in the high-dimensional setting where \(p\) can be much larger than \(n\). Both law of large numbers and limiting distribution are derived. We then consider testing the bandedness of the covariance matrix of a high dimensional Gaussian distribution which includes testing for independence as a special case. The limiting laws of the coherence of the data matrix play a critical role in the construction of the test. The asymptotic results is also applied to the construction of compressed sensing matrices.

Abstract

Theory and methodology for nonparametric regression have been particularly well developed in the case of additive homoscedastic Gaussian noise. Inspired by the asymptotic equivalence theory, there have been ongoing efforts in recent years to construct explicit procedures that turn other function estimation problems into a standard nonparametric regression with Gaussian noise. Then in principle any good Gaussian nonparametric regression method can be used to solve those more complicated nonparametric models. In particular, Brown, Cai and Zhou (2010) considered nonparametric regression in the natural exponential families with a quadratic variance function. In this paper we extend the scope of Brown, Cai and Zhou (2010) to general natural exponential families by introducing a new explicit procedure that is based on the variance stabilizing transformation. The new approach significantly reduces the bias of the inverse transformation and as a consequence it enables the method to be applicable to a wider class of exponential families. Combining this procedure with a wavelet block thresholding estimator for Gaussian nonparametric regression, we show that the resulting estimator enjoys a high degree of adaptivity and spatial adaptivity with near-optimal asymptotic performance over a broad range of Besov spaces.

Abstract

Most results in nonparametric regression theory are developed only for the case of additive noise. In such a setting many smoothing techniques including wavelet thresholding methods have been developed and shown to be highly adaptive. In this paper we consider nonparametric regression in exponential families which include, for example, Poisson regression, binomial regression, and gamma regression. We propose a unified approach of using a mean-matching variance stabilizing transformation to turn the relatively complicated problem of nonparametric regression in exponential families into a standard homoscedastic Gaussian regression problem. Then in principle any good nonparametric Gaussian regression procedure can be applied to the transformed data. In this paper we use a wavelet block thresholding rule to construct the final estimator of the regression function. The procedure is easily implementable. Both theoretical and numerical properties of the estimator are investigated. The estimator is shown to enjoy a high degree of adaptivity and spatial adaptivity. It simultaneously attains the optimal rates of convergence under integrated squared error over a wide range of Besov spaces and achieves adaptive local minimax rate for estimating functions at a point. The estimator also performs well numerically.

Abstract

Density estimation has traditionally been treated separately from nonparametric regression. In this paper we propose and implement a density estimation procedure which begins by turning density estimation into a regression problem. This regression problem is created by binning the original observations into many small size bins, and by then applying a suitable form of root transformation to the binned data counts. A nonparametric regression estimator is then applied to the transformed data. Finally, the estimated regression function is un-rooted by squaring and normalizing. In principle many common nonparametric regression estimators could be used in the implementation of this algorithm. We propose use of a wavelet block thresholding estimator. The entire algorithm is then convenient to implement. We show that the resulting density estimator enjoys a high degree of adaptivity. A numerical example and a practical data example are discussed to illustrate and explain the use of this density estimation procedure.

Abstract

We consider a blockwise James-Stein estimator for nonparametric function estimation in suitable wavelet or Fourier bases. The estimator can be readily explained and implemented. We show that the estimator is asymptotically sharp-adaptive in minimax risk over any Sobolev ball containing the true function. Further, for a moderately broad range of bounded sets in Besov space our estimator is asymptotically nearly sharp adaptive in the sense that it comes within the Donoho-Liu constant, 1.24, of being exactly sharp adaptive. Other parameter spaces are also considered. The paper concludes with a Monte-Carlo study comparing the performance of our estimator to that of three other popular wavelet estimators. Our procedure generally (but not always) outperforms two of these and is overall comparable, or perhaps slightly superior, to the third.

Abstract

Asymptotic equivalence results developed in the literature so far are only for bounded loss functions. This limits the potential applications of the theory because many commonly used loss functions in statistical inference are unbounded. In this paper we develop asymptotic equivalence results for robust nonparametric regression with unbounded loss functions. The results imply that all the Gaussian nonparametric regression procedures can be robustified in a unified way. A key step in our equivalence argument is to bin the data and then take the median of each bin. Through binning and taking the medians of the binned data, the general nonparametric regression model is turned into a standard Gaussian regression model, and then in principle any procedures for Gaussian nonparametric regression can be applied. The asymptotic equivalence results have significant practical implications. To illustrate the general principles of the equivalence argument we consider two important nonparametric inference problems: robust estimation of the regression function and the estimation of a quadratic functional. In both cases easily implementable procedures are constructed and are shown to enjoy a high degree of robustness and adaptivity. Other problems such as construction of confidence sets and nonparametric hypothesis testing can be handled in a similar fashion.

Abstract

A data-driven block thresholding procedure for wavelet regression is proposed and its theoretical and numerical properties are investigated. The procedure empirically chooses the block size and threshold level at each resolution level by minimizing Stein's unbiased risk estimate. The estimator is sharp adaptive over a class of Besov bodies and achieves simultaneously within a small constant factor of the minimax risk over a wide collection of Besov Bodies including both the dense and sparse cases. The procedure is easy to implement. Numerical results show that it has superior finite sample performance in comparison to the other leading wavelet thresholding estimators.

Abstract

We consider a wavelet thresholding approach to adaptive variance function estimation in heteroscedastic nonparametric regression. A data-driven estimator is constructed by applying wavelet thresholding to the squared first-order differences of the observations. We show that the variance function estimator is nearly optimally adaptive to the smoothness of both the mean and variance functions. The estimator is shown to achieve the optimal adaptive rate of convergence under the pointwise squared error simultaneously over a range of smoothness classes. The estimator is also adaptively within a logarithmic factor of the minimax risk under the global mean integrated squared error over a collection of spatially inhomogeneous function classes. Numerical implementation and simulation results are also discussed.

Abstract

In this paper we develop a nonparametric regression method that is simultaneously adaptive over a wide range of function classes for the regression function and robust over a large collection of error distributions, including those that are heavy-tailed, and may not even possess variances or means. Our approach is to first use local medians to turn the problem of nonparametric regression with unknown noise distribution into a standard Gaussian regression problem and then apply a wavelet block thresholding procedure to construct an estimator of the regression function. It is shown that the estimator simultaneously attains the optimal rate of convergence over a wide range of the Besov classes, without prior knowledge of the smoothness of the underlying functions or prior knowledge of the error distribution. The estimator also automatically adapts to the local smoothness of the underlying function, and attains the local adaptive minimax rate for estimating functions at a point. A key technical result in our development is a quantile coupling theorem which gives a tight bound for the quantile coupling between the sample medians and a normal variable. This median coupling inequality may be of independent interest.

Abstract

We consider wavelet block thresholding method for density estimation. A block-thresholded density estimator is proposed and is shown to achieve optimal global rate of convergence over Besov spaces and simultaneously attain the optimal adaptive pointwise convergence rate as well. These results are obtained in part through the determination of an optimal block length.

Abstract

Function estimation over the Besov spaces under pointwise \(\ell^r\ (1\le r<\infty)\) risks is considered. Minimax rates of convergence are derived using a constrained risk inequality and wavelets. Adaptation under pointwise risks is also considered. Sharp lower bounds on the cost of adaptation are obtained and are shown to be attainable by a wavelet estimator. The results demonstrate important differences between the minimax properties under pointwise and global risk measures. The minimax rates and adaptation for estimating derivatives under pointwise risks are also presented. A general \(\ell^r\)-risk oracle inequality is developed for the proofs of the main results.

Abstract

We consider a block thresholding and vaguelet-wavelet approach to certain statistical linear inverse problems. Based on an oracle inequality, an adaptive block thresholding estimator for linear inverse problems is proposed and the asymptotic properties of the estimator are investigated. It is shown that the estimator enjoys a higher degree of adaptivity than the standard term-by-term thresholding methods; it attains the exact optimal rates of convergence over a range of Besov classes. The problem of estimating a derivative is considered in more detail as a test for the general estimation procedure. We show that the derivative estimator is spatially adaptive; it automatically adapts to the local smoothness of the function and attains the local adaptive minimax rate for estimating a derivative at a point.

Abstract

In this article we investigate the asymptotic and numerical properties of a class of block thresholding estimators for wavelet regression. We consider the effect of block size on global and local adaptivity and the choice of thresholding constant. The optimal rate of convergence for block thresholding with a given block size is derived for both the global and local estimation. It is shown that there are conflicting requirements on the block size for achieving the global and local adaptivity. We then consider the choice of thresholding constant for a given block size by treating the block thresholding as a hypothesis testing problem. The combined results lead naturally to an optimal choice of block size and thresholding constant. We conclude with a numerical study which compares the finite-sample performance among block thresholding estimators as well as with other wavelet methods.

Abstract

In standard wavelet methods, the empirical wavelet coefficients are thresholded term by term, on the basis of their individual magnitudes. Information on other coefficients has no influence on the treatment of particular coefficients. We propose and investigate a wavelet shrinkage method that incorporates information on neighboring coefficients into the decision making. The coefficients are considered in overlapping blocks; the treatment of coefficients in the middle of each block depends on the data in the whole block. Both the asymptotic and numerical performances of two particular versions of the estimator are considered. In numerical comparisons with various methods, both versions of the estimator perform excellently; on the theoretical side, we show that one of the versions achieves the exact optimal rates of convergence over a range of Besov classes.

Abstract

Multiresolution wavelet transformation (MWT) and block thresholding used to effectively suppress both background and noise inteference while minimally distorting Raman spectral features. The performance of MWT as a spectral pre-processing algorithm is demonstrated using both synthetic spectra and experimental hyper-spectral Raman images with large background and noise components. The results are quantified by comparing correlation coefficients between synthetic spectra with either the same or different backgrounds. The improved chemical imaging performance obtained using MWT is demonstrated by comparing Principal Component Analysis (PCA) channel images and Spectral Angle Mapping (SAM) classified images before and after MWT pre-processing.

Abstract

We show that for nonparametric regression if the samples have random uniform design, the wavelet method with universal thresholding can be applied directly to the samples as if they were equispaced. The resulting estimator achieves within a logarithmic factor from the minimax rate of convergence over a family of Hölder classes. Simulation result is also discussed.

Abstract

We study wavelet function estimation via the approach of block thresholding and ideal adaptation with oracle. Oracle inequalities are derived and serve as guides for the selection of smoothing parameters. Based on an oracle inequality and motivated by the data compression and localization properties of wavelets, an adaptive wavelet estimator for nonparametric regression is proposed and the optimality of the procedure is investigated. We show that the estimator achieves simultaneously three objectives: adaptivity, spatial adaptivity, and computational efficiency. Specifically, it is proved that the estimator attains the exact optimal rates of convergence over a range of Besov classes and the estimator achieves adaptive local minimax rate for estimating functions at a point. The estimator is easy to implement, at the computational cost of \(O(n)\). Simulation shows that the estimator has excellent numerical performance relative to more traditional wavelet estimators.

Abstract

Standard wavelet shrinkage procedures for nonparametric regression are restricted to equispaced samples. There, data are transformed into empirical wavelet coefficients and threshold rules are applied to the coefficients. The estimators are obtained via the inverse transform of the denoised wavelet coefficients. In many applications, however, the samples are nonequispaced. It can be shown that these procedures would produce suboptimal estimators if they were applied directly to nonequispaced samples. We propose a wavelet shrinkage procedure for nonequispaced samples. We show that the estimate is adaptive and near optimal. For global estimation, the estimate is within a logarithmic factor of the minimax risk over a wide range of piecewise Hölder classes, indeed with a number of discontinuities that grows polynomially fast with the sample size. For estimating a target function at a point, the estimate is optimally adaptive to unknown degree of smoothness within a constant. In addition, the estimate enjoys a smoothness property: if the target function is the zero function, then with probability tending to 1 the estimate is also the zero function.

Abstract

Let \(W_t\) be a Brownian motion and let \(\epsilon_{i,n}\) be iid \(N(0,1)\), \(i=1,\ldots,n\), and independent of \(W_t\). \(\sigma,\tau>0\) are real, unknown parameters. Suppose we observe \(Y_{i,n}=\sigma W_{i/n}+\tau\epsilon_{i,n}\). In this paper we will establish sharp estimators for \(\sigma^2\) and \(\tau^2\) in minimax sense, i.e. they attain asymptotically the minimax constant. These estimators are based on a spectral decomposition of the underlying process \(Y_{i,n}\) and can be computed explicitly in \(O(n\log n)\) operations. A proof for the minimax lower bound is given. Further we show that these estimators are asymptotically normal.

Abstract

Variance function estimation in multivariate nonparametric regression is considered and the minimax rate of convergence is established. Our work uses the approach that generalizes the one used in Munk, Bissantz, Wagner and Freitag (2005, JRSSB) for the constant variance case. As is the case when the number of dimensions \(d=1\), and very much contrary to the common practice, it is often not desirable to base the estimator of the variance function on the residuals from an optimal estimator of the mean. Instead it is desirable to use estimators of the mean with minimal bias. Another important conclusion is that the first order difference-based estimator that achieves minimax rate of convergence in one-dimensional case does not do the same in the high dimensional case. Instead, the optimal order of differences depends on the number of dimensions.

Abstract

We consider a wavelet thresholding approach to adaptive variance function estimation in heteroscedastic nonparametric regression. A data-driven estimator is constructed by applying wavelet thresholding to the squared first-order differences of the observations. We show that the variance function estimator is nearly optimally adaptive to the smoothness of both the mean and variance functions. The estimator is shown to achieve the optimal adaptive rate of convergence under the pointwise squared error simultaneously over a range of smoothness classes. The estimator is also adaptively within a logarithmic factor of the minimax risk under the global mean integrated squared error over a collection of spatially inhomogeneous function classes. Numerical implementation and simulation results are also discussed.

Abstract

Variance function estimation in nonparametric regression is considered and the minimax rate of convergence is derived. We are particularly interested in the effect of the unknown mean on the estimation of the variance function. Our results indicate that, contrary to the common practice, it is often not desirable to base the estimator of the variance function on the residuals from an optimal estimator of the mean. Instead it is desirable to use estimators of the mean with minimal bias. In addition the results also correct the optimal rate claimed in the previous literature.