In regression analysis of multivariate data, it is tacitly assumed that response and predictor variables in each observed response-predictor pair correspond to the same entity or unit. In this paper, we consider the situation of “permuted data” in which this basic correspondence has been lost. Several recent papers have considered this situation without further assumptions on the underlying permutation. In applications, the latter is often to known to have additional structure that can be leveraged. Specifically, we herein consider the common scenario of “sparsely permuted data” in which only a small fraction of the data is affected by a mismatch between response and predictors. However, an adverse effect already observed for sparsely permuted data is that the least squares estimator as well as other estimators not accounting for such partial mismatch are inconsistent. One approach studied in detail herein is to treat permuted data as outliers which motivates the use of robust regression formulations to estimate the regression parameter. The resulting estimate can subsequently be used to recover the permutation. A notable benefit of the proposed approach is its computational simplicity given the general lack of procedures for the above problem that are both statistically sound and computationally appealing.

In this paper we study the frequentist convergence rate for the Latent Dirichlet Allocation (Blei, Ng and Jordan, 2003) topic models. We show that the maximum likelihood estimator converges to one of the finitely many equivalent parameters in Wasserstein’s distance metric at a rate of $n^{-1/4}$ without assuming separability or non-degeneracy of the underlying topics and/or the existence of more than three words per document, thus generalizing the previous works of Anandkumar et al. (2012, 2014) from an information-theoretical perspective. We also show that the $n^{-1/4}$ convergence rate is optimal in the worst case.

Let us assume that $X$, $Y$ and $U$ are observed and that the conditional mean of $U$ given $X$ and $Y$ can be expressed via an additive dependency of $X$, $\lambda(X)Y$ and $X+Y$ for some unspecified function $\lambda$. This structured regression model can be transferred to a hazard model or a density model when applied on some appropriate grid, and has important forecasting applications via structured marker dependent hazards models or structured density models including age-period-cohort relationships. The structured regression model is also important when the severity of the dependent variable has a complicated dependency on waiting times $X$, $Y$ and the total waiting time $X+Y$. In case the conditional mean of $U$ approximates a density, the regression model can be used to analyse the age-period-cohort model, also when exposure data are not available. In case the conditional mean of $U$ approximates a marker dependent hazard, the regression model introduces new relevant age-period-cohort time scale interdependencies in understanding longevity. A direct use of the regression relationship introduced in this paper is the estimation of the severity of outstanding liabilities in non-life insurance companies. The technical approach taken is to use B-splines to capture the underlying one-dimensional unspecified functions. It is shown via finite sample simulation studies and an application for forecasting future asbestos related deaths in the UK that the B-spline approach works well in practice. Special consideration has been given to ensure identifiability of all models considered.

We propose an adaptive confidence interval procedure (CIP) for the coefficients in the normal linear regression model. This procedure has a frequentist coverage rate that is constant as a function of the model parameters, yet provides smaller intervals than the usual interval procedure, on average across regression coefficients. The proposed procedure is obtained by defining a class of CIPs that all have exact $1-\alpha $ frequentist coverage, and then selecting from this class the procedure that minimizes a prior expected interval width. We describe an adaptive approach for estimating the prior distribution from the data, so that the potential risk of a poorly specified prior is reduced. The resulting adaptive confidence intervals maintain exact non-asymptotic $1-\alpha $ coverage if two conditions are met - that the design matrix is full rank (which will be known) and that the errors are normally distributed (which can be checked empirically). No assumptions on the unknown parameters are necessary to maintain exact coverage. Additionally, in a “$p$ growing with $n$” asymptotic scenario, this adaptive FAB procedure is asymptotically Bayes-optimal among $1-\alpha $ frequentist CIPs.

We study the empirical measure associated to a sample of size $n$ and modified by $N$ iterations of the raking-ratio method. This empirical measure is adjusted to match the true probability of sets in a finite partition which changes each step. We establish asymptotic properties of the raking-ratio empirical process indexed by functions as $n\rightarrow +\infty $, for $N$ fixed. We study nonasymptotic properties by using a Gaussian approximation which yields uniform Berry-Esseen type bounds depending on $n,N$ and provides estimates of the uniform quadratic risk reduction. A closed-form expression of the limiting covariance matrices is derived as $N\rightarrow +\infty $. In the two-way contingency table case the limiting process has a simple explicit formula.

High-dimensional data, where the number of variables exceeds or is comparable to the sample size, is now pervasive in many scientific applications. In recent years, Bayesian shrinkage models have been developed as effective and computationally feasible tools to analyze such data, especially in the context of linear regression. In this paper, we focus on the Normal-Gamma shrinkage model developed by Griffin and Brown [7]. This model subsumes the popular Bayesian lasso model, and a three-block Gibbs sampling algorithm to sample from the resulting intractable posterior distribution has been developed in [7]. We consider an alternative two-block Gibbs sampling algorithm, and rigorously demonstrate its advantage over the three-block sampler by comparing specific spectral properties. In particular, we show that the Markov operator corresponding to the two-block sampler is trace class (and hence Hilbert-Schmidt), whereas the operator corresponding to the three-block sampler is not even Hilbert-Schmidt. The trace class property for the two-block sampler implies geometric convergence for the associated Markov chain, which justifies the use of Markov chain CLT’s to obtain practical error bounds for MCMC based estimates. Additionally, it facilitates theoretical comparisons of the two-block sampler with sandwich algorithms which aim to improve performance by inserting inexpensive extra steps in between the two conditional draws of the two-block sampler.

We consider the problem of detecting a sparse mixture as studied by Ingster (1997) and Donoho and Jin (2004). We consider a wide array of base distributions. In particular, we study the situation when the base distribution has polynomial tails, a situation that has not received much attention in the literature. Perhaps surprisingly, we find that in the context of such a power-law distribution, the higher criticism does not achieve the detection boundary. However, the scan statistic does.

We consider importance sampling to estimate the probability $\mu$ of a union of $J$ rare events $H_{j}$ defined by a random variable $\boldsymbol{x}$. The sampler we study has been used in spatial statistics, genomics and combinatorics going back at least to Karp and Luby (1983). It works by sampling one event at random, then sampling $\boldsymbol{x}$ conditionally on that event happening and it constructs an unbiased estimate of $\mu$ by multiplying an inverse moment of the number of occuring events by the union bound. We prove some variance bounds for this sampler. For a sample size of $n$, it has a variance no larger than $\mu(\bar{\mu}-\mu)/n$ where $\bar{\mu}$ is the union bound. It also has a coefficient of variation no larger than $\sqrt{(J+J^{-1}-2)/(4n)}$ regardless of the overlap pattern among the $J$ events. Our motivating problem comes from power system reliability, where the phase differences between connected nodes have a joint Gaussian distribution and the $J$ rare events arise from unacceptably large phase differences. In the grid reliability problems even some events defined by $5772$ constraints in $326$ dimensions, with probability below $10^{-22}$, are estimated with a coefficient of variation of about $0.0024$ with only $n=10{,}000$ sample values.

The paper is concerned with the nonparametric statistical estimation of linear spectral functionals for Lévy-driven continuous-time stationary linear models with tapered data. As an estimator for unknown functional we consider the averaged tapered periodogram. We analyze the bias of the estimator and obtain sufficient conditions assuring the proper rate of convergence of the bias to zero, necessary for asymptotic normality of the estimator. We prove a a central limit theorem for a suitable normalized stochastic process generated by a tapered Toeplitz type quadratic functional of the model. As a consequence of these results we obtain the asymptotic normality of our estimator.

This paper presents a method for fast Bayesian variable selection in the normal linear regression model with high dimensional data. A novel approach is adopted in which an explicit posterior probability for including a covariate is obtained. The method is sequential but not order dependent, one deals with each covariate one by one, and a spike and slab prior is only assigned to the coefficient under investigation. We adopt the well-known spike and slab Gaussian priors with a sample size dependent variance, which achieves strong selection consistency for marginal posterior probabilities even when the number of covariates grows almost exponentially with sample size. Numerical illustrations are presented where it is shown that the new approach provides essentially equivalent results to the standard spike and slab priors, i.e. the same marginal posterior probabilities of the coefficients being nonzero, which are estimated via Gibbs sampling. Hence, we obtain the same results via the direct calculation of $p$ probabilities, compared to a stochastic search over a space of $2^{p}$ elements. Our procedure only requires $p$ probabilities to be calculated, which can be done exactly, hence parallel computation when $p$ is large is feasible.

We consider a mixed moving average (MMA) process $X$ driven by a Lévy basis and prove that it is weakly dependent with rates computable in terms of the moving average kernel and the characteristic quadruple of the Lévy basis. Using this property, we show conditions ensuring that sample mean and autocovariances of $X$ have a limiting normal distribution. We extend these results to stochastic volatility models and then investigate a Generalized Method of Moments estimator for the supOU process and the supOU stochastic volatility model after choosing a suitable distribution for the mean reversion parameter. For these estimators, we analyze the asymptotic behavior in detail.

In this paper optimal designs for regression problems with spherical predictors of arbitrary dimension are considered. Our work is motivated by applications in material sciences, where crystallographic textures such as the misorientation distribution or the grain boundary distribution (depending on a four dimensional spherical predictor) are represented by series of hyperspherical harmonics, which are estimated from experimental or simulated data.

For this type of estimation problems we explicitly determine optimal designs with respect to the $\Phi _{p}$-criteria introduced by Kiefer (1974) and a class of orthogonally invariant information criteria recently introduced in the literature. In particular, we show that the uniform distribution on the $m$-dimensional sphere is optimal and construct discrete and implementable designs with the same information matrices as the continuous optimal designs. Finally, we illustrate the advantages of the new designs for series estimation by hyperspherical harmonics, which are symmetric with respect to the first and second crystallographic point group.

We consider an additive partially linear framework for modelling massive heterogeneous data. The major goal is to extract multiple common features simultaneously across all sub-populations while exploring heterogeneity of each sub-population. We propose an aggregation type of estimators for the commonality parameters that possess the asymptotic optimal bounds and the asymptotic distributions as if there were no heterogeneity. This oracle result holds when the number of sub-populations does not grow too fast and the tuning parameters are selected carefully. A plug-in estimator for the heterogeneity parameter is further constructed, and shown to possess the asymptotic distribution as if the commonality information were available. Furthermore, we develop a heterogeneity test for the linear components and a homogeneity test for the non-linear components accordingly. The performance of the proposed methods is evaluated via simulation studies and an application to the Medicare Provider Utilization and Payment data.

The main focus of the analysts who deal with clustered data is usually not on the clustering variables, and hence the group-specific parameters are treated as nuisance. If a fixed effects formulation is preferred and the total number of clusters is large relative to the single-group sizes, classical frequentist techniques relying on the profile likelihood are often misleading. The use of alternative tools, such as modifications to the profile likelihood or integrated likelihoods, for making accurate inference on a parameter of interest can be complicated by the presence of nonstandard modelling and/or sampling assumptions. We show here how to employ Monte Carlo simulation in order to approximate the modified profile likelihood in some of these unconventional frameworks. The proposed solution is widely applicable and is shown to retain the usual properties of the modified profile likelihood. The approach is examined in two instances particularly relevant in applications, i.e. missing-data models and survival models with unspecified censoring distribution. The effectiveness of the proposed solution is validated via simulation studies and two clinical trial applications.

Ranking, also known as learning to rank in machine learning community, is to rank a number of items based on their relevance to a specific query. In literature, most ranking methods use a uniform ranking function to evaluate the relevance, which completely ignores the heterogeneity among queries. To admit different ranking functions for various queries, a general $U$-process formulation for query-dependent ranking is developed. It allows to incorporate neighborhood structure among queries via various forms of smoothing weights to improve the ranking performance. One of its salient features is its capability of producing reasonable rankings for novel queries that are absent in the training set, which is commonly encountered in practice but often neglected in the literature. The proposed method is implemented via an inexact alternating direction method of multipliers (ADMM) for each query parallelly. Its asymptotic risk bound is established, showing that it achieves desirable ranking accuracy at a fast rate for any query including the novel ones. Furthermore, simulated examples and a real application to the Yahoo! challenge dataset also support the advantage of the query-dependent ranking method against existing competitors.

In this paper, we consider the problem of multiple testing where the hypotheses are dependent. In most of the existing literature, either Bayesian or non-Bayesian, the decision rules mainly focus on the validity of the test procedure rather than actually utilizing the dependency to increase efficiency. Moreover, the decisions regarding different hypotheses are marginal in the sense that they do not depend upon each other directly. However, in realistic situations, the hypotheses are usually dependent, and hence it is desirable that the decisions regarding the dependent hypotheses are taken jointly.

In this article, we develop a novel Bayesian multiple testing procedure that coherently takes this requirement into consideration. Our method, which is based on new notions of error and non-error terms, substantially enhances efficiency by judicious exploitation of the dependence structure among the hypotheses. We show that our method minimizes the posterior expected loss associated with an additive “0-1” loss function; we also prove theoretical results on the relevant error probabilities, establishing the coherence and usefulness of our method. The optimal decision configuration is not available in closed form and we propose an efficient simulated annealing algorithm for the purpose of optimization, which is also generically applicable to binary optimization problems.

Extensive simulation studies indicate that in dependent situations, our method performs significantly better than some existing popular conventional multiple testing methods, in terms of accuracy and power control. Moreover, application of our ideas to a real, spatial data set associated with radionuclide concentration in Rongelap islands yielded insightful results.

In the present paper we study three geometrical characteristics for the excursion sets of a two-dimensional stationary isotropic random field. First, we show that these characteristics can be estimated without bias if the considered field satisfies a kinematic formula, this is for instance the case of fields given by a function of smooth Gaussian fields or of some shot noise fields. By using the proposed estimators of these geometric characteristics, we describe some inference procedures for the estimation of the parameters of the field. An extensive simulation study illustrates the performances of each estimator. Then, we use the Euler characteristic estimator to build a test to determine whether a given field is Gaussian or not, when compared to various alternatives. The test is based on a sparse information, i.e., the excursion sets for two different levels of the field to be tested. Finally, the proposed test is adapted to an applied case, synthesized 2D digital mammograms.

This paper develops robust confidence intervals in high-dimensional and left-censored regression. Type-I censored regression models, where a competing event makes the variable of interest unobservable, are extremely common in practice. In this paper, we develop smoothed estimating equations that are adaptive to censoring level and are more robust to the misspecification of the error distribution. We propose a unified class of robust estimators, including one-step Mallow’s, Schweppe’s, and Hill-Ryan’s estimator that are adaptive to the left-censored observations. In the ultra-high-dimensional setting, where the dimensionality can grow exponentially with the sample size, we show that as long as the preliminary estimator converges faster than $n^{-1/4}$, the one-step estimators inherit asymptotic distribution of fully iterated version. Moreover, we show that the size of the residuals of the Bahadur representation matches those of the pure linear models – that is, the effects of censoring disappear asymptotically. Simulation studies demonstrate that our method is adaptive to the censoring level and asymmetry in the error distribution, and does not lose efficiency when the errors are from symmetric distributions.

We consider the problem of isotonic regression, where the underlying signal $x$ is assumed to satisfy a monotonicity constraint, that is, $x$ lies in the cone $\{x\in \mathbb{R}^{n}:x_{1}\leq \dots\leq x_{n}\}$. We study the isotonic projection operator (projection to this cone), and find a necessary and sufficient condition characterizing all norms with respect to which this projection is contractive. This enables a simple and non-asymptotic analysis of the convergence properties of isotonic regression, yielding uniform confidence bands that adapt to the local Lipschitz properties of the signal.

In the present paper, we have studied a Dynamic Stochastic Block Model (DSBM) under the assumptions that the connection probabilities, as functions of time, are smooth and that at most $s$ nodes can switch their class memberships between two consecutive time points. We estimate the edge probability tensor by a kernel-type procedure and extract the group memberships of the nodes by spectral clustering. The procedure is computationally viable, adaptive to the unknown smoothness of the functional connection probabilities, to the rate $s$ of membership switching, and to the unknown number of clusters. In addition, it is accompanied by non-asymptotic guarantees for the precision of estimation and clustering.

This paper studies a two step procedure for monotone increasing additive single index models with Gaussian designs. The proposed procedure is simple, easy to implement with existing software, and consists of consecutively applying LASSO and isotonic regression. Aside from formalizing this procedure, we provide theoretical guarantees regarding its performance: 1) we show that our procedure controls the in-sample squared error; 2) we demonstrate that one can use the procedure for predicting new observations, by showing that the absolute prediction error can be controlled with high-probability. Our bounds show a tradeoff of two rates: the minimax rate for estimating high dimensional quadratic loss, and the minimax nonparametric rate for estimating a monotone increasing function.

The paper gives a unified study of the large sample asymptotic theory of penalized splines including the O-splines using B-splines and an integrated squared derivative penalty [22], the P-splines which use B-splines and a discrete difference penalty [13], and the T-splines which use truncated polynomials and a ridge penalty [24]. Extending existing results for O-splines [7], it is shown that, depending on the number of knots and appropriate smoothing parameters, the $L_{2}$ risk bounds of penalized spline estimators are rate-wise similar to either those of regression splines or to those of smoothing splines and could each attain the optimal minimax rate of convergence [32]. In addition, convergence rate of the $L_{\infty }$ risk bound, and local asymptotic bias and variance are derived for all three types of penalized splines.

We introduce the notion of DTM-signature, a measure on $\mathbb{R}$ that can be associated to any metric-measure space. This signature is based on the function distance to a measure (DTM) introduced in 2009 by Chazal, Cohen-Steiner and Mérigot. It leads to a pseudo-metric between metric-measure spaces, that is bounded above by the Gromov-Wasserstein distance. This pseudo-metric is used to build a statistical test of isomorphism between two metric-measure spaces, from the observation of two $N$-samples.

The test is based on subsampling methods and comes with theoretical guarantees. It is proven to be of the correct level asymptotically. Also, when the measures are supported on compact subsets of $\mathbb{R}^{d}$, rates of convergence are derived for the $L_{1}$-Wasserstein distance between the distribution of the test statistic and its subsampling approximation. These rates depend on some parameter $\rho >1$. In addition, we prove that the power is bounded above by $\exp (-CN^{1/\rho })$, with $C$ proportional to the square of the aforementioned pseudo-metric between the metric-measure spaces. Under some geometrical assumptions, we also derive lower bounds for this pseudo-metric.

An algorithm is proposed for the implementation of this statistical test, and its performance is compared to the performance of other methods through numerical experiments.

In this paper, we propose an invariant test based on the modified correction to the likelihood ratio test (LRT) of the equality of two high-dimensional covariance matrices. It is well-known that the classical log-LRT is not well defined when the dimension is larger than or equal to one of the sample sizes. Or even the log-LRT is well-defined, it is usually perceived as a bad statistic in the high-dimensional cases because of their low powers under some alternatives. In this paper, we will justify the usefulness of the modified log-LRT, and an invariant test that works well in cases where the dimension is larger than the sample sizes. Besides, the test is established under the weakest conditions on the dimensions and the moments of the samples. The asymptotic distribution of the proposed test statistic is also obtained under the null hypothesis. What is more, we also propose a lite version of the modified LRT in the paper. A simulation study and a real data analysis show that the performances of the two proposed statistics are invariant under affine transformations.

Much effort has been done to control the “false discovery rate” (FDR) when $m$ hypotheses are tested simultaneously. The FDR is the expectation of the “false discovery proportion” $\text{FDP}=V/R$ given by the ratio of the number of false rejections $V$ and all rejections $R$. In this paper, we have a closer look at the FDP for adaptive linear step-up multiple tests. These tests extend the well known Benjamini and Hochberg test by estimating the unknown amount $m_{0}$ of the true null hypotheses. We give exact finite sample formulas for higher moments of the FDP and, in particular, for its variance. Using these allows us a precise discussion about the stability of the FDP, i.e., when the FDP is asymptotically close to its mean. We present sufficient and necessary conditions for this stability. They include the presence of a stable estimator for the proportion $m_{0}/m$. We apply our results to convex combinations of generalized Storey type estimators with various tuning parameters and (possibly) data-driven weights. The corresponding step-up tests allow a flexible adaptation. Moreover, these tests control the FDR at finite sample size. We compare these tests to the classical Benjamini and Hochberg test and discuss the advantages of them.

The estimation of the coefficient matrix in a multivariate response linear regression model is considered in situations where we can observe only strictly increasing transformations of the continuous responses and covariates. It is further assumed that the joint dependence between all the observed variables is characterized by an elliptical copula. Penalized estimators of the coefficient matrix are obtained in a high-dimensional setting by assuming that the coefficient matrix is either element-wise sparse or row-sparse, and by incorporating the precision matrix of the error, which is also assumed to be sparse. Estimation of the copula parameters is achieved by inversion of Kendall’s tau. It is shown that when the true coefficient matrix is row-sparse, the estimator obtained via a group penalty outperforms the one obtained via a simple element-wise penalty. Simulation studies are used to illustrate this fact and the advantage of incorporating the precision matrix of the error when the correlation among the components of the error vector is strong. Moreover, the use of the normal-score rank correlation estimator is revisited in the context of high-dimensional Gaussian copula models. It is shown that this estimator remains as the optimal estimator of the copula correlation matrix in this setting.

This paper studies and critically discusses the construction of nonparametric confidence regions for density level sets. Methodologies based on both vertical variation and horizontal variation are considered. The investigations provide theoretical insight into the behavior of these confidence regions via large sample theory. We also discuss the geometric relationships underlying the construction of horizontal and vertical methods, and how finite sample performance of these confidence regions is influenced by geometric or topological aspects. These discussions are supported by numerical studies.

We investigate the asymptotic behavior of the $L_{p}$-distance between a monotone function on a compact interval and a smooth estimator of this function. Our main result is a central limit theorem for the $L_{p}$-error of smooth isotonic estimators obtained by smoothing a Grenander-type estimator or isotonizing the ordinary kernel estimator. As a preliminary result we establish a similar result for ordinary kernel estimators. Our results are obtained in a general setting, which includes estimation of a monotone density, regression function and hazard rate. We also perform a simulation study for testing monotonicity on the basis of the $L_{2}$-distance between the kernel estimator and the smoothed Grenander-type estimator.

We propose an optimal experimental design for a curvilinear regression model that minimizes the band-width of simultaneous confidence bands. Simultaneous confidence bands for curvilinear regression are constructed by evaluating the volume of a tube about a curve that is defined as a trajectory of a regression basis vector (Naiman, 1986). The proposed criterion is constructed based on the volume of a tube, and the corresponding optimal design that minimizes the volume of tube is referred to as the tube-volume optimal (TV-optimal) design. For Fourier and weighted polynomial regressions, the problem is formalized as one of minimization over the cone of Hankel positive definite matrices, and the criterion to minimize is expressed as an elliptic integral. We show that the Möbius group keeps our problem invariant, and hence, minimization can be conducted over cross-sections of orbits. We demonstrate that for the weighted polynomial regression and the Fourier regression with three bases, the tube-volume optimal design forms an orbit of the Möbius group containing D-optimal designs as representative elements.

Due to the rapid development of various social networks, the spatial autoregressive (SAR) model is becoming an important tool in social network analysis. However, major bottlenecks remain in analyzing large-scale networks (e.g., Facebook has over 700 million active users), including computational scalability, estimation consistency, and proper network sampling. To address these challenges, we propose a novel least squares estimator (LSE) for analyzing large sparse networks based on the SAR model. Computationally, the LSE is linear in the network size, making it scalable to analysis of huge networks. In theory, the LSE is $\sqrt{n}$-consistent and asymptotically normal under certain regularity conditions. A new LSE-based network sampling technique is further developed, which can automatically adjust autocorrelation between sampled and unsampled units and hence guarantee valid statistical inferences. Moreover, we generalize the LSE approach for the classical SAR model to more complex networks associated with multiple sources of social interaction effect. Numerical results for simulated and real data are presented to illustrate performance of the LSE.

The relative performance of competing point forecasts is usually measured in terms of loss or scoring functions. It is widely accepted that these scoring function should be strictly consistent in the sense that the expected score is minimized by the correctly specified forecast for a certain statistical functional such as the mean, median, or a certain risk measure. Thus, strict consistency opens the way to meaningful forecast comparison, but is also important in regression and M-estimation. Usually strictly consistent scoring functions for an elicitable functional are not unique. To give guidance on the choice of a scoring function, this paper introduces two additional quality criteria. Order-sensitivity opens the possibility to compare two deliberately misspecified forecasts given that the forecasts are ordered in a certain sense. On the other hand, equivariant scoring functions obey similar equivariance properties as the functional at hand – such as translation invariance or positive homogeneity. In our study, we consider scoring functions for popular functionals, putting special emphasis on vector-valued functionals, e.g. the pair (mean, variance) or (Value at Risk, Expected Shortfall).

We consider the problem of variable selection in high-dimensional statistical models where the goal is to report a set of variables, out of many predictors $X_{1},\dotsc ,X_{p}$, that are relevant to a response of interest. For linear high-dimensional model, where the number of parameters exceeds the number of samples $(p>n)$, we propose a procedure for variables selection and prove that it controls the directional false discovery rate (FDR) below a pre-assigned significance level $q\in [0,1]$. We further analyze the statistical power of our framework and show that for designs with subgaussian rows and a common precision matrix $\Omega \in{\mathbb{R}} ^{p\times p}$, if the minimum nonzero parameter $\theta_{\min }$ satisfies \[\sqrt{n}\theta_{\min }-\sigma \sqrt{2(\max_{i\in [p]}\Omega_{ii})\log \left(\frac{2p}{qs_{0}}\right)}\to \infty \,,\] then this procedure achieves asymptotic power one.

Our framework is built upon the debiasing approach and assumes the standard condition $s_{0}=o(\sqrt{n}/(\log p)^{2})$, where $s_{0}$ indicates the number of true positives among the $p$ features. Notably, this framework achieves exact directional FDR control without any assumption on the amplitude of unknown regression parameters, and does not require any knowledge of the distribution of covariates or the noise level. We test our method in synthetic and real data experiments to assess its performance and to corroborate our theoretical results.

Vine copulas, or pair-copula constructions, have become an important tool in high-dimensional dependence modeling. Commonly, it is assumed that the data generating copula can be represented by a simplified vine copula (SVC). In this paper, we study the simplifying assumption and investigate the approximation of multivariate copulas by SVCs. We introduce the partial vine copula (PVC) which is a particular SVC where to any edge a $j$-th order partial copula is assigned. The PVC generalizes the partial correlation matrix and plays a major role in the approximation of copulas by SVCs. We investigate to what extent the PVC describes the dependence structure of the underlying copula. We show that, in general, the PVC does not minimize the Kullback-Leibler divergence from the true copula if the simplifying assumption does not hold. However, under regularity conditions, stepwise estimators of pair-copula constructions converge to the PVC irrespective of whether the simplifying assumption holds or not. Moreover, we elucidate why the PVC is often the best feasible SVC approximation in practice.

We analyze control of the familywise error rate (FWER) in a multiple testing scenario with a great many null hypotheses about the distribution of a high-dimensional random variable among which only a very small fraction are false, or “active”. In order to improve power relative to conservative Bonferroni bounds, we explore a coarse-to-fine procedure adapted to a situation in which tests are partitioned into subsets, or “cells”, and active hypotheses tend to cluster within cells. We develop procedures for a non-parametric case based on generalized permutation testing and a linear Gaussian model, and demonstrate higher power than Bonferroni estimates at the same FWER when the active hypotheses do cluster. The main technical difficulty arises from the correlation between the test statistics at the individual and cell levels, which increases the likelihood of a hypothesis being falsely discovered when the cell that contains it is falsely discovered (survivorship bias). This requires sharp estimates of certain quadrant probabilities when a cell is inactive.

Consider a scalar reflected diffusion $(X_{t}:t\geq 0)$, where the unknown drift function $b$ is modelled nonparametrically. We show that in the low frequency sampling case, when the sample consists of $(X_{0},X_{\Delta },...,X_{n\Delta })$ for some fixed sampling distance $\Delta >0$, the model satisfies the local asymptotic normality (LAN) property, assuming that $b$ satisfies some mild regularity assumptions. This is established by using the connections of diffusion processes to elliptic and parabolic PDEs. The key tools used are regularity estimates for certain parabolic PDEs as well as a detailed analysis of the spectral properties of the elliptic differential operator related to $(X_{t}:t\geq 0)$.

Various problems in manifold estimation make use of a quantity called the reach, denoted by $\tau_{M}$, which is a measure of the regularity of the manifold. This paper is the first investigation into the problem of how to estimate the reach. First, we study the geometry of the reach through an approximation perspective. We derive new geometric results on the reach for submanifolds without boundary. An estimator $\hat{\tau }$ of $\tau_{M}$ is proposed in an oracle framework where tangent spaces are known, and bounds assessing its efficiency are derived. In the case of i.i.d. random point cloud $\mathbb{X}_{n}$, $\hat{\tau }(\mathbb{X}_{n})$ is showed to achieve uniform expected loss bounds over a $\mathcal{C}^{3}$-like model. Finally, we obtain upper and lower bounds on the minimax rate for estimating the reach.

In this paper, we consider inference problem about the drift parameter vector in generalized mean reverting processes with multiple and unknown change-points. In particular, we study the case where the parameter may satisfy uncertain restriction. As compared to the results in literature, we generalize some findings in five ways. First, we consider the model which incorporates the uncertain prior knowledge. Second, we derive the unrestricted estimator (UE) and the restricted estimator (RE) and we study their asymptotic properties. Third, we derive a test for testing the hypothesized restriction and we derive its asymptotic local power. We also prove that the proposed test is consistent. Fourth, we construct a class of shrinkage type estimators (SEs) which encloses the UE, the RE and classical SEs. Fifth, we derive the relative risk dominance of the proposed estimators. More precisely, we prove that the SEs dominate the UE. Finally, we present some simulation results which corroborate the established theoretical findings.

This paper develops mixed-normal approximations for probabilities that vectors of multiple Skorohod integrals belong to random convex polytopes when the dimensions of the vectors possibly diverge to infinity. We apply the developed theory to establish the asymptotic mixed normality of the realized covariance matrix of a high-dimensional continuous semimartingale observed at a high-frequency, where the dimension can be much larger than the sample size. We also present an application of this result to testing the residual sparsity of a high-dimensional continuous-time factor model.

We consider estimation of the location and the height of the jump in the $\gamma $-th derivative - a kink of order $\gamma $ - of a regression curve, which is assumed to be Hölder smooth of order $s\geq \gamma +1$ away from the kink. Optimal convergence rates as well as the joint asymptotic normal distribution of estimators based on the zero-crossing-time technique are established. Further, we construct joint as well as marginal asymptotic confidence sets for these parameters which are honest and adaptive with respect to the smoothness parameter $s$ over subsets of the Hölder classes. The finite-sample performance is investigated in a simulation study, and a real data illustration is given to a series of annual global surface temperatures.

Accurate specification of a likelihood function is becoming increasingly difficult in many inference problems in astronomy. As sample sizes resulting from astronomical surveys continue to grow, deficiencies in the likelihood function lead to larger biases in key parameter estimates. These deficiencies result from the oversimplification of the physical processes that generated the data, and from the failure to account for observational limitations. Unfortunately, realistic models often do not yield an analytical form for the likelihood. The estimation of a stellar initial mass function (IMF) is an important example. The stellar IMF is the mass distribution of stars initially formed in a given cluster of stars, a population which is not directly observable due to stellar evolution and other disruptions and observational limitations of the cluster. There are several difficulties with specifying a likelihood in this setting since the physical processes and observational challenges result in measurable masses that cannot legitimately be considered independent draws from an IMF. This work improves inference of the IMF by using an approximate Bayesian computation approach that both accounts for observational and astrophysical effects and incorporates a physically-motivated model for star cluster formation. The methodology is illustrated via a simulation study, demonstrating that the proposed approach can recover the true posterior in realistic situations, and applied to observations from astrophysical simulation data.

Survey data are often collected under informative sampling designs where subject inclusion probabilities are designed to be correlated with the response variable of interest. The data modeler seeks to estimate the parameters of a population model they specify from these data. Sampling weights constructed from marginal inclusion probabilities are typically used to form an exponentiated pseudo likelihood as a plug-in estimator in a partially Bayesian pseudo posterior. We introduce the first fully Bayesian alternative, based on a Bayes rule construction, that simultaneously performs weight smoothing and estimates the population model parameters in a construction that treats the response variable(s) and inclusion probabilities as jointly randomly generated from a population distribution. We formulate conditions on known marginal and pairwise inclusion probabilities that define a class of sampling designs where $L_{1}$ consistency of the joint posterior is guaranteed. We compare performances between the two approaches on synthetic data. We demonstrate that the credibility intervals under our fully Bayesian method achieve nominal coverage. We apply our method to data from the National Health and Nutrition Examination Survey to explore the relationship between caffeine consumption and systolic blood pressure.

This paper looks at the strong consistency of the ordinary least squares (OLS) estimator in linear regression models with adaptive learning. It is a companion to Christopeit & Massmann (2018) which considers the estimator’s convergence in distribution and its weak consistency in the same setting. Under constant gain learning, the model is closely related to stationary, (alternating) unit root or explosive autoregressive processes. Under decreasing gain learning, the regressors in the model are asymptotically collinear. The paper examines, first, the issue of strong convergence of the learning recursion: It is argued that, under constant gain learning, the recursion does not converge in any probabilistic sense, while for decreasing gain learning rates are derived at which the recursion converges almost surely to the rational expectations equilibrium. Secondly, the paper establishes the strong consistency of the OLS estimators, under both constant and decreasing gain learning, as well as rates at which the estimators converge almost surely. In the constant gain model, separate estimators for the intercept and slope parameters are juxtaposed to the joint estimator, drawing on the recent literature on explosive autoregressive models. Thirdly, it is emphasised that strong consistency is obtained in all models although the near-optimal condition for the strong consistency of OLS in linear regression models with stochastic regressors, established by Lai & Wei (1982a), is not always met.

This paper considers a deconvolution regression problem in a multivariate setting with anisotropic structure and constructs an estimator of the function of interest using the hyperbolic wavelet basis. The deconvolution structure assumed is an anisotropic version of the smooth type (either regular-smooth or super-smooth). The function of interest is assumed to belong to a Besov space with anisotropic smoothness. Global performances of the presented hyperbolic wavelet estimators is measured by obtaining upper bounds on convergence rates in the $\mathscr{L}^{p}$-risk with $1\le p\le 2$ and $1\le p<\infty $ in the regular-smooth and super-smooth cases respectively. The results are compared and contrasted with existing convergence results in the literature.

Treatment rules based on individual patient characteristics that are easy to interpret and disseminate are important in clinical practice. Properly planned and conducted randomized clinical trials are used to construct individualized treatment rules. However, it is often a concern that trial participants lack representativeness, so it limits the applicability of the derived rules to a target population. In this work, we use data from a single trial study to propose a two-stage procedure to derive a robust and parsimonious rule to maximize the benefit in the target population. The procedure allows a wide range of possible covariate distributions in the target population, with minimal assumptions on the first two moments of the covariate distribution. The practical utility and favorable performance of the methodology are demonstrated using extensive simulations and a real data application.

Shape constrained estimation in discrete settings has received increasing attention in statistics. Among the most important shape constrained models is multiple monotonicity, including $k$-monotonicity, for a given integer $k\in [1,\infty )$, and complete monotonicity. Multiple monotonicity provides a nice generalization of monotonicity and convexity and has been successfully used in applications related to estimation of species richness. Although fully nonparametric, it is of great interest to determine some of the well-known parametric distributions which belong to this model. Among the most important examples are the family of Poisson distributions and mixtures thereof. In Giguelay (2017) $k$-monotonicity of Poisson distributions was connected to the roots of a certain polynomial, but a typographical error occurred while writing its expression. In this note, we correct that typographical error and give a detailed proof that a Poisson distribution with rate $\lambda \in [0,\infty )$ is $k$-monotone if and only if $\lambda \le \lambda _{k}$, where $\lambda _{k}$ is the smallest zero of the $k$-th degree Laguerre polynomial $L_{k}(x)=\sum _{j=0}^{k}(-1)^{j}\binom{k}{j}x^{j}/j!$, $x\ge 0$. This result yields the sufficient condition that a mixture of Poisson distributions is $k$-monotone if the support of the mixing distribution is included in $[0,\lambda _{k}]$. Furthermore, we show that the only complete monotone Poisson distribution is the Dirac distribution at $0$.

We deal with the problem of nonparametric estimation of a multivariate regression function without any assumption on the compacity of the support of the random design. To tackle the problem, we propose to extend a “warping” device to the multivariate framework. An adaptive warped kernel estimator is first defined in the case of known design distribution and proved to be optimal in the oracle sense. Then, a general procedure is carried out: the marginal distributions of the design are estimated by the empirical cumulative distribution functions, and the dependence structure is built using a kernel estimation of the copula density. The copula density estimator is also studied and proved to be optimal in the oracle and in the minimax sense. The plug-in of this estimator in the regression function estimator provides a fully data-driven procedure. A numerical study illustrates the theoretical results.

The utility of a Markov chain Monte Carlo algorithm is, in large part, determined by the size of the spectral gap of the corresponding Markov operator. However, calculating (and even approximating) the spectral gaps of practical Monte Carlo Markov chains in statistics has proven to be an extremely difficult and often insurmountable task, especially when these chains move on continuous state spaces. In this paper, a method for accurate estimation of the spectral gap is developed for general state space Markov chains whose operators are non-negative and trace-class. The method is based on the fact that the second largest eigenvalue (and hence the spectral gap) of such operators can be bounded above and below by simple functions of the power sums of the eigenvalues. These power sums often have nice integral representations. A classical Monte Carlo method is proposed to estimate these integrals, and a simple sufficient condition for finite variance is provided. This leads to asymptotically valid confidence intervals for the second largest eigenvalue (and the spectral gap) of the Markov operator. In contrast with previously existing techniques, our method is not based on a near-stationary version of the Markov chain, which, paradoxically, cannot be obtained in a principled manner without bounds on the spectral gap. On the other hand, it can be quite expensive from a computational standpoint. The efficiency of the method is studied both theoretically and empirically.

We introduce a novel regression framework which simultaneously models the quantile and the Expected Shortfall (ES) of a response variable given a set of covariates. This regression is based on strictly consistent loss functions for the pair consisting of the quantile and the ES, which allow for M- and Z-estimation of the joint regression parameters. We show consistency and asymptotic normality for both estimators under weak regularity conditions. The underlying loss functions depend on two specification functions, whose choices affect the properties of the resulting estimators. We find that the Z-estimator is numerically unstable and thus, we rely on M-estimation of the model parameters. Extensive simulations verify the asymptotic properties and analyze the small sample behavior of the M-estimator for different specification functions. This joint regression framework allows for various applications including estimating, forecasting and backtesting ES, which is particularly relevant in light of the recent introduction of the ES into the Basel Accords. We illustrate this through two exemplary empirical applications in forecasting and forecast combination of the ES.

We introduce a trimmed version of the Hill estimator for the index of a heavy-tailed distribution, which is robust to perturbations in the extreme order statistics. In the ideal Pareto setting, the estimator is essentially finite-sample efficient among all unbiased estimators with a given strict upper break-down point. For general heavy-tailed models, we establish the asymptotic normality of the estimator under second order regular variation conditions and also show that it is minimax rate-optimal in the Hall class of distributions. We also develop an automatic, data-driven method for the choice of the trimming parameter which yields a new type of robust estimator that can adapt to the unknown level of contamination in the extremes. This adaptive robustness property makes our estimator particularly appealing and superior to other robust estimators in the setting where the extremes of the data are contaminated. As an important application of the data-driven selection of the trimming parameters, we obtain a methodology for the principled identification of extreme outliers in heavy tailed data. Indeed, the method has been shown to correctly identify the number of outliers in the previously explored Condroz data set.

We propose a strategy for computing estimators in some non-standard M-estimation problems, where the data are distributed across different servers and the observations across servers, though independent, can come from heterogeneous sub-populations, thereby violating the identically distributed assumption. Our strategy fixes the super-efficiency phenomenon observed in prior work on distributed computing in (i) the isotonic regression framework, where averaging several isotonic estimates (each computed at a local server) on a central server produces super-efficient estimates that do not replicate the properties of the global isotonic estimator, i.e. the isotonic estimate that would be constructed by transferring all the data to a single server, and (ii) certain types of M-estimation problems involving optimization of discontinuous criterion functions where M-estimates converge at the cube-root rate. The new estimators proposed in this paper work by smoothing the data on each local server, communicating the smoothed summaries to the central server, and then solving a non-linear optimization problem at the central server. They are shown to replicate the asymptotic properties of the corresponding global estimators, and also overcome the super-efficiency phenomenon exhibited by existing estimators.

We propose a computationally efficient estimator, formulated as a convex program, for a broad class of nonlinear regression problems that involve difference of convex(DC) nonlinearities. The proposed method can be viewed as a significant extension of the “anchored regression” method formulated and analyzed in [10] for regression with convex nonlinearities. Our main assumption, in addition to other mild statistical and computational assumptions, is availability of a certain approximation oracle for the average of the gradients of the observation functions at a ground truth. Under this assumption and using a PAC-Bayesian analysis we show that the proposed estimator produces an accurate estimate with high probability. As a concrete example, we study the proposed framework in the bilinear regression problem with Gaussian factors and quantify a sufficient sample complexity for exact recovery. Furthermore, we describe a computationally tractable scheme that provably produces the required approximation oracle in the considered bilinear regression problem.

Pretest-posttest studies are commonly used for assessing the effect of a treatment or an intervention. We propose an empirical likelihood based approach to both testing and estimation of the treatment effect in non-randomized pretest-posttest studies where the posttest outcomes are subject to missingness. The proposed empirical likelihood ratio test and the estimation procedure are multiply robust in the sense that multiple working models are allowed for the propensity score of treatment assignment, the missingness probability and the outcome regression, and the validity of the test and the estimation requires only a certain combination of those multiple working models to be correctly specified. An empirical likelihood ratio confidence interval can be constructed for the treatment effect and has better coverage probabilities than confidence intervals based on the Wald statistic. Simulations are conducted to demonstrate the finite-sample performances of the proposed methods.

This article proposes a generalized threshold latent variable model for flexible threshold modeling of time series. The proposed model encompasses several existing models, and allows a discrete valued threshold variable. Sufficient conditions for stationarity and ergodicity are investigated. The minimum description length principle is applied to formulate a criterion function for parameter estimation and model selection. A computationally efficient procedure for optimizing the criterion function is developed based on a genetic algorithm. Consistency and weak convergence of the parameter estimates are established. Moreover, simulation studies and an application for initial public offering data are presented to illustrate the proposed methodology.

We introduce a way of modelling temporal dependence in random functions $X(t)$ in the framework of linear regression. Based on discretised curves $\{X_{i}(t_{0}),X_{i}(t_{1}),\ldots ,X_{i}(t_{T})\}$, the final point $X_{i}(t_{T})$ is predicted from $\{X_{i}(t_{0}),X_{i}(t_{1}),\ldots ,X_{i}(t_{T-1})\}$. The proposed model flexibly reflects the relative importance of predictors by partitioning the regression parameters into a smooth and a rough regime. Specifically, unconstrained (rough) regression parameters are used for observations located close to $X_{i}(t_{T})$, while the set of regression coefficients for the predictors positioned far from $X_{i}(t_{T})$ are assumed to be sampled from a smooth function. This both regularises the prediction problem and reflects the ‘decaying memory’ structure of the time series. The point at which the change in smoothness occurs is estimated from the data via a technique akin to change-point detection. The joint estimation procedure for the smoothness change-point and the regression parameters is presented, and the asymptotic behaviour of the estimated change-point is analysed. The usefulness of the new model is demonstrated through simulations and four real data examples, involving country fertility data, pollution data, stock volatility series and sunspot number data. Our methodology is implemented in the R package srp, available from CRAN.

The estimation of marks for a finite population of points scattered onto a study region is considered when a sample of these points is selected by a probabilistic sampling scheme. At each point, the mark is estimated by means of an inverse distance weighting interpolator. The design-based asymptotic properties of the resulting maps are derived when the study area remains fixed, a sequence of nested populations with increasing size is considered and samples of increasing size are selected. Conditions ensuring design-based asymptotic unbiasedness and consistency are given. They essentially require that marks are the values of a pointwise or uniformly continuous deterministic function, the enlargement of the populations is rather regular and the sequence of sampling designs ensures an asymptotic spatial balance. A computationally simple mean squared error estimator is proposed. A simulation study is performed to assess the theoretical results on artificial populations. Finally, an application for mapping the values of the height of trees in a forest stand located in North Italy is reported.

The likelihood ratio statistic, with its asymptotic $\chi ^{2}$ distribution at regular model points, is often used for hypothesis testing. However, the asymptotic distribution can differ at model singularities and boundaries, suggesting the use of a $\chi ^{2}$ might be problematic nearby. Indeed, its poor behavior for testing near singularities and boundaries is apparent in simulations, and can lead to conservative or anti-conservative tests. Here we develop a new distribution designed for use in hypothesis testing near singularities and boundaries, which asymptotically agrees with that of the likelihood ratio statistic. For two example trinomial models, arising in the context of inference of evolutionary trees, we show the new distributions outperform a $\chi ^{2}$.

In this paper, we propose to construct confidence bands by bootstrapping the debiased kernel density estimator (for density estimation) and the debiased local polynomial regression estimator (for regression analysis). The idea of using a debiased estimator was recently employed by Calonico et al. (2018b) to construct a confidence interval of the density function (and regression function) at a given point by explicitly estimating stochastic variations. We extend their ideas of using the debiased estimator and further propose a bootstrap approach for constructing simultaneous confidence bands. This modified method has an advantage that we can easily choose the smoothing bandwidth from conventional bandwidth selectors and the confidence band will be asymptotically valid. We prove the validity of the bootstrap confidence band and generalize it to density level sets and inverse regression problems. Simulation studies confirm the validity of the proposed confidence bands/sets. We apply our approach to an Astronomy dataset to show its applicability.