We propose a two-sample test for large-dimensional covariance matrices in generalized elliptical models. The test statistic is based on a U-statistic estimator of the squared Frobenius norm of the difference between the two population covariance matrices. This statistic was originally introduced by Li and Chen (2012) for the independent component model. As a key theoretical contribution, we establish a new central limit theorem for the U-statistics under elliptical data, valid under both the null and alternative hypotheses. This result enables asymptotic control of the test level and facilitates a power analysis. To the best of our knowledge, the proposed test is the first such method to be supported by theoretical guarantees for elliptical data. Our approach imposes only mild assumptions on the covariance matrices and does not require sparsity nor explicit growth conditions on the dimension-to-sample-size ratio. We illustrate our theoretical findings through applications to both synthetic and real-world data.

We introduce a new instrumental variable (IV) estimator for heterogeneous treatment effects in the presence of endogeneity. Our estimator is based on double/debiased machine learning (DML) and uses efficient machine learning instruments (MLIV) and kernel smoothing. We prove consistency and asymptotic normality of our estimator and also construct confidence sets that are more robust towards weak IV. Along the way, we also provide an accessible discussion of the corresponding estimator for the homogeneous treatment effect with efficient machine learning instruments. The methods are evaluated on synthetic and real datasets and an implementation is made available in the R package IVDML.

The crossed random effects model is widely used, finding applications in various fields such as longitudinal studies, e-commerce, and recommender systems, among others. However, these models encounter scalability challenges, as the computational time for standard algorithms grows superlinearly with the number N of observations in the data set, commonly O(N3∕2) or worse. Recent published works present scalable methods for crossed random effects in linear models and some generalized linear models, but those methods only allow for random intercepts. In this paper, we devise scalable algorithms for models that include random slopes. This addition brings substantial difficulty in estimating the random-effect covariance matrices in a scalable way. We address this issue by using a variational EM algorithm. Our proposed approach accommodates both diagonal covariance matrices and cases where no structure is assumed—a scenario common in fields such as psychology and neuroscience. In simulations, the proposed method is substantially faster than standard methods for large N. It is also more efficient than ordinary least squares which has a problem of greatly underestimating the sampling uncertainty in parameter estimates. We illustrate the new method on a MovieLens data set, as well as a large data set (five million observations) from the online retailer Stitch Fix.

High-dimensional datasets with a network structure are increasingly popular in various fields, such as sociology, finance, and bioinformatics. Existing research focuses on the stationary setting. This paper builds a change point model in which the high-dimensional time series is assumed to follow an autoregressive process with network effect within each time segment. In addition to the heterogeneous model with node-specific coefficients, the homogeneous and the grouped model, which assume different levels of homogeneity among model parameters, are also analyzed. A computationally efficient two-step procedure based on the fused LASSO and a general information criterion is developed to estimate the number and locations of change points simultaneously, as well as model parameters in each time segment. For the grouped model, after the identification of change points, a similar two-step procedure is employed to estimate the group structure for each time segment. The consistency of the proposed method for the number and locations of change points, the model parameters, and the group structure are established. Simulation and real data analysis on U.S. stock returns support the effectiveness of the proposed method in modeling high-dimensional network time series with structural breaks.

We consider a Bayesian approach for the recovery of scalar parameters arising in inverse problems. We consider a general signal-in white noise model where we have access to two independent noisy observations of a function, and of a linear transformation of the function. The linear operator is unknown up to a scalar parameter. We present a Bernstein-von Mises theorem for the marginal posterior of the scalar under regularity assumptions of the operator. We further derive Bernstein-von Mises results for different priors and apply them to two concrete examples: the recovery of the thermal diffusivity in a heat equation problem, and the recovery of a location parameter in a semi-blind deconvolution problem.

We use Stein characterisations to derive new moment-type estimators for the parameters of several truncated multivariate distributions in the i.i.d. case; we also derive the asymptotic properties of these estima-tors. Our examples include the truncated multivariate normal distribution and truncated products of independent univariate distributions. The estimators are explicit, and therefore provide an interesting alternative to the maximum likelihood estimator (MLE). The quality of these estimators is assessed through competitive simulation studies, in which we compare their behavior to the performance of the MLE and the score-matching approach.

The Grenander estimator is a well-studied procedure for univariate nonparametric density estimation. It is usually defined as the Maximum Likelihood Estimator (MLE) over the class of all non-increasing densities on the positive real line. It can also be seen as the MLE over the class of all scale mixtures of uniform densities. Using the latter viewpoint, Pavlides and Wellner [33] proposed a multivariate extension of the Grenander estimator as the nonparametric MLE over the class of all multivariate scale mixtures of uniform densities. We prove that this multivariate estimator achieves the univariate cube root rate of convergence with only a logarithmic multiplicative factor that depends on the dimension. The usual curse of dimensionality is therefore avoided to some extent for this multivariate estimator. This result positively resolves a conjecture of Pavlides and Wellner [33] under an additional lower bound assumption. Our proof proceeds via a general accuracy result for the Hellinger accuracy of MLEs over convex classes of densities. We also provide algorithms for computing the estimator, and illustrate performance on real and simulated datasets.

We introduce a new predictive mechanism that operates in the presence of hidden confounding across distributionally diverse data sources while ensuring consistent estimation of causal parameters—despite their recognized suboptimality for prediction in the literature. Our method is based on a novel estimand that captures the dependence structure between response noise and covariates, incorporating causal parameters into a generative model that adaptively replicates the conditional distribution of the test environment. Identifiability is achieved under a straightforward, empirically verifiable assumption. Our approach ensures probabilistic alignment with test distributions uniformly across arbitrary interventions, enabling valid predictions without requiring worst-case optimization or assumptions about the strength of perturbations at test time. Through extensive simulations, we demonstrate that our method outperforms state-of-the-art invariance-based and domain adaptation approaches. Additionally, we validate its practical applicability and superior target risk performance on a cardiovascular disease dataset.

A core problem in statistical network analysis is to develop network analogues of classical techniques. The problem of bootstrapping network data presents a particular challenge, since one typically observes a single network rather than a sample. Here we propose two methods for obtaining bootstrap samples for networks drawn from latent space models. The first method generates bootstrap replicates of network statistics that can be represented as U-statistics in the latent positions, and avoids actually constructing new bootstrapped networks. Many network quantities can be represented as U-statistics, including average degree and subgraph counts, but other equally popular summaries, such as clustering coefficients, are not expressible as U-statistics. Our second bootstrapping method generates replicates of whole networks, and thus can be used for bootstrapping more general network functions. Under the assumption of a random dot product graph, a type of latent space network model, we show consistency of the proposed bootstrap methods. We give motivating examples throughout and demonstrate the effectiveness of our methods on both synthetic and real data and show that they improve upon prior methods.

The remedian uses a k×b matrix to approximate population and sample medians of n≤bk values in a streaming setting. It recursively replaces buffers of b values with their medians, discarding 200⌈b∕2⌉∕bk% of its most extreme inputs. Prior research (Rousseeuw and Bassett, 1990; Chao and Lin, 1993; Chen and Chen, 2005) explored the remedian’s distribution under the limits k→∞ and k,b→∞. The breakdown point of the remedian, (⌈b∕2⌉∕b)k, vanishes as k→∞, yet converges to 1∕2k as b→∞. This study derives the robust-regime multivariate normal distribution for the standardized vector comprising the sample mean, median, remedian, and remedian rank as b→∞. Our analysis reveals that, under these conditions, the remedian provides an asymptotically precise approximation of the sample median. We compare its statistical efficiency to that of the sample mean and median and extend the study to the simultaneous estimation of multiple quantiles. Specifically, we propose an asymptotic distribution for the random vector generated through parallel remedian estimation applied to components of independently and identically distributed (i.i.d.) random vectors. These findings shed light on robust estimation in streaming settings, highlighting the trade-offs between memory requirements, accuracy, and robustness.