A discrepancy bound for deterministic acceptance-rejection samplers beyond $$N^{-1/2}$$ in dimension 1

In this paper we consider an acceptance-rejection (AR) sampler based on deterministic driver sequences. We prove that the discrepancy of an N element sample set generated in this way is bounded by \(\mathcal {O} (N^{-2/3}\log N)\), provided that the target density is twice continuously differentiable with non-vanishing curvature and the AR sampler uses the driver sequence \(\mathcal {K}_M= \{( j \alpha , j \beta ) ~~ mod~~1 \mid j = 1,\ldots ,M\},\) where \(\alpha ,\beta \) are real algebraic numbers such that \(1,\alpha ,\beta \) is a basis of a number field over \(\mathbb {Q}\) of degree 3. For the driver sequence \(\mathcal {F}_k= \{ ({j}/{F_k}, \{{jF_{k-1}}/{F_k}\} ) \mid j=1,\ldots , F_k\},\) where \(F_k\) is the k-th Fibonacci number and \(\{x\}=x-\lfloor x \rfloor \) is the fractional part of a non-negative real number x, we can remove the \(\log \) factor to improve the convergence rate to \(\mathcal {O}(N^{-2/3})\), where again N is the number of samples we accepted. We also introduce a criterion for measuring the goodness of driver sequences. The proposed approach is numerically tested by calculating the star-discrepancy of samples generated for some target densities using \(\mathcal {K}_M\) and \(\mathcal {F}_k\) as driver sequences. These results confirm that achieving a convergence rate beyond \(N^{-1/2}\) is possible in practice using \(\mathcal {K}_M\) and \(\mathcal {F}_k\) as driver sequences in the acceptance-rejection sampler.     

Random projections for Bayesian regression

This article deals with random projections applied as a data reduction technique for Bayesian regression analysis. We show sufficient conditions under which the entire d-dimensional distribution is approximately preserved under random projections by reducing the number of data points from n to \(k\in O({\text {poly}}(d/\varepsilon ))\) in the case \(n\gg d\). Under mild assumptions, we prove that evaluating a Gaussian likelihood function based on the projected data instead of the original data yields a \((1+O(\varepsilon ))\)-approximation in terms of the \(\ell _2\) Wasserstein distance. Our main result shows that the posterior distribution of Bayesian linear regression is approximated up to a small error depending on only an \(\varepsilon \)-fraction of its defining parameters. This holds when using arbitrary Gaussian priors or the degenerate case of uniform distributions over \(\mathbb {R}^d\) for \(\beta \). Our empirical evaluations involve different simulated settings of Bayesian linear regression. Our experiments underline that the proposed method is able to recover the regression model up to small error while considerably reducing the total running time.     

Hypotheses testing about the drift parameter in linear stochastic differential equation driven by stable processes

In this paper, we consider the problem of hypotheses testing about the drift parameter \(\theta \) in the process \(\text {d}Y^{\delta }_{t} = \theta \dot{f}(t)Y^{\delta }_{t}\text {d}t + b(t)\text {d}L^{\delta }_{t}\) driven by symmetric \(\delta \)-stable Lévy process \(L^{\delta }_{t}\) with \(\dot{f}(t)\) being the derivative of a known increasing function f(t) and b(t) being known as well. We consider the hypotheses testing \(H_{0}: \theta \le 0\) and \(K_{0}: \theta =0\) against the alternatives \(H_{1}: \theta >0\) and \(K_{1}: \theta \ne 0\), respectively. For these hypotheses, we propose inverse methods, which are motivated by sequential approach, based on the first hitting time of the observed process (or its absolute value) to a pre-specified boundary or two boundaries until some given time. The applicability of these methods is illustrated. For the case \(Y^{\delta }_{0}=0\), we are able to calculate the values of boundaries and finite observed times more directly. We are able to show the consistencies of proposed tests for \(Y^{\delta }_{0}\ge 0\) with \(\delta \in (1,2]\) and for \(Y^{\delta }_{0}=0\) with \(\delta \in (0,2]\) under quite mild conditions.     

Multilevel particle filters: normalizing constant estimation

In this article, we introduce two new estimates of the normalizing constant (or marginal likelihood) for partially observed diffusion (POD) processes, with discrete observations. One estimate is biased but non-negative and the other is unbiased but not almost surely non-negative. Our method uses the multilevel particle filter of Jasra et al. (Multilevel particle lter, arXiv:1510.04977, 2015). We show that, under assumptions, for Euler discretized PODs and a given \(\varepsilon >0\) in order to obtain a mean square error (MSE) of \({\mathcal {O}}(\varepsilon ^2)\) one requires a work of \({\mathcal {O}}(\varepsilon ^{-2.5})\) for our new estimates versus a standard particle filter that requires a work of \({\mathcal {O}}(\varepsilon ^{-3})\). Our theoretical results are supported by numerical simulations.     

Rate of uniform consistency for a class of mode regression on functional stationary ergodic data

The aim of this paper is to study the asymptotic properties of a class of kernel conditional mode estimates whenever functional stationary ergodic data are considered. To be more precise on the matter, in the ergodic data setting, we consider a random elements (X, Z) taking values in some semi-metric abstract space \(E\times F\). For a real function \(\varphi \) defined on the space F and \(x\in E\), we consider the conditional mode of the real random variable \(\varphi (Z)\) given the event “\(X=x\)”. While estimating the conditional mode function, say \(\theta _\varphi (x)\), using the well-known kernel estimator, we establish the strong consistency with rate of this estimate uniformly over Vapnik–Chervonenkis classes of functions \(\varphi \). Notice that the ergodic setting offers a more general framework than the usual mixing structure. Two applications to energy data are provided to illustrate some examples of the proposed approach in time series forecasting framework. The first one consists in forecasting the daily peak of electricity demand in France (measured in Giga-Watt). Whereas the second one deals with the short-term forecasting of the electrical energy (measured in Giga-Watt per Hour) that may be consumed over some time intervals that cover the peak demand.     

Conditional density estimation using the local Gaussian correlation

Let \(\mathbf {X} = (X_1,\ldots ,X_p)\) be a stochastic vector having joint density function \(f_{\mathbf {X}}(\mathbf {x})\) with partitions \(\mathbf {X}_1 = (X_1,\ldots ,X_k)\) and \(\mathbf {X}_2 = (X_{k+1},\ldots ,X_p)\). A new method for estimating the conditional density function of \(\mathbf {X}_1\) given \(\mathbf {X}_2\) is presented. It is based on locally Gaussian approximations, but simplified in order to tackle the curse of dimensionality in multivariate applications, where both response and explanatory variables can be vectors. We compare our method to some available competitors, and the error of approximation is shown to be small in a series of examples using real and simulated data, and the estimator is shown to be particularly robust against noise caused by independent variables. We also present examples of practical applications of our conditional density estimator in the analysis of time series. Typical values for k in our examples are 1 and 2, and we include simulation experiments with values of p up to 6. Large sample theory is established under a strong mixing condition.     

Objective Bayesian transformation and variable selection using default Bayes factors

In this work, the problem of transformation and simultaneous variable selection is thoroughly treated via objective Bayesian approaches by the use of default Bayes factor variants. Four uniparametric families of transformations (Box–Cox, Modulus, Yeo-Johnson and Dual), denoted by T, are evaluated and compared. The subjective prior elicitation for the transformation parameter \(\lambda _T\), for each T, is not a straightforward task. Additionally, little prior information for \(\lambda _T\) is expected to be available, and therefore, an objective method is required. The intrinsic Bayes factors and the fractional Bayes factors allow us to incorporate default improper priors for \(\lambda _T\). We study the behaviour of each approach using a simulated reference example as well as two real-life examples.     

On comparison of dispersion matrices of estimators under a constrained linear model

We introduce some new mathematical tools in the analysis of dispersion matrices of the two well-known OLSEs and BLUEs under general linear models with parameter restrictions. We first establish some formulas for calculating the ranks and inertias of the differences of OLSEs’ and BLUEs’ dispersion matrices of parametric functions under the general linear model \({\mathscr {M}}= \{\mathbf{y}, \ \mathbf{X }\pmb {\beta }, \ \pmb {\Sigma }\}\) and the constrained model \({\mathscr {M}}_r = \{\mathbf{y}, \, \mathbf{X }\pmb {\beta }\, | \, \mathbf{A }\pmb {\beta }= \mathbf{b}, \ \pmb {\Sigma }\}\), where \(\mathbf{A }\pmb {\beta }= \mathbf{b}\) is a consistent linear matrix equation for the unknown parameter vector \(\pmb {\beta }\) to satisfy. As applications, we derive necessary and sufficient conditions for many equalities and inequalities of OLSEs’ and BLUEs’ dispersion matrices to hold under \({\mathscr {M}}\) and \({\mathscr {M}}_r\).     

Stochastically optimal bootstrap sample size for shrinkage-type statistics

In nonregular problems where the conventional \(n\) out of \(n\) bootstrap is inconsistent, the \(m\) out of \(n\) bootstrap provides a useful remedy to restore consistency. Conventionally, optimal choice of the bootstrap sample size \(m\) is taken to be the minimiser of a frequentist error measure, estimation of which has posed a major difficulty hindering practical application of the \(m\) out of \(n\) bootstrap method. Relatively little attention has been paid to a stronger, stochastic, version of the optimal bootstrap sample size, defined as the minimiser of an error measure calculated directly from the observed sample. Motivated by this stronger notion of optimality, we develop procedures for calculating the stochastically optimal value of \(m\). Our procedures are shown to work under special forms of Edgeworth-type expansions which are typically satisfied by statistics of the shrinkage type. Theoretical and empirical properties of our methods are illustrated with three examples, namely the James–Stein estimator, the ridge regression estimator and the post-model-selection regression estimator.     

Transformation approaches of linear random-effects models

Assume that a linear random-effects model \(\mathbf{y}= \mathbf{X}\varvec{\beta }+ \varvec{\varepsilon }= \mathbf{X}(\mathbf{A}\varvec{\alpha }+ \varvec{\gamma }) + \varvec{\varepsilon }\) is transformed as \(\mathbf{T}\mathbf{y}= \mathbf{T}\mathbf{X}\varvec{\beta }+ \mathbf{T}\varvec{\varepsilon }= \mathbf{T}\mathbf{X}(\mathbf{A}\varvec{\alpha }+ \varvec{\gamma }) + \mathbf{T}\varvec{\varepsilon }\) by pre-multiplying a given matrix \(\mathbf{T}\) of arbitrary rank. These two models are not necessarily equivalent unless \(\mathbf{T}\) is of full column rank, and we have to work with this derived model in many situations. Because predictors/estimators of the parameter spaces under the two models are not necessarily the same, it is primary work to compare predictors/estimators in the two models and to establish possible links between the inference results obtained from two models. This paper presents a general algebraic approach to the problem of comparing best linear unbiased predictors (BLUPs) of parameter spaces in an original linear random-effects model and its transformations, and provides a group of fundamental and comprehensive results on mathematical and statistical properties of the BLUPs. In particular, we construct many equalities for the BLUPs under an original linear random-effects model and its transformations, and obtain necessary and sufficient conditions for the equalities to hold.     

On periodic time-varying bilinear processes: structure and asymptotic inference

This paper is devoted to the bilinear time series models with periodic-varying coefficients \(\left( { PBL}\right) \). So, firstly conditions ensuring the existence of periodic stationary solutions of the \({ PBL}\) and the existence of higher-order moments of such solutions are given. A distribution free approach to the parameter estimation of \({ PBL}\) is presented. The proposed method relies on minimum distance estimator based on the first and second order empirical moments of the observed process. Consistency and asymptotic normality of the estimator are discussed. Examples and Monte Carlo simulation results illustrate the practical relevancy of our general theoretical results are presented.     

Wavelet regression estimations with strong mixing data

Using a wavelet basis, we establish in this paper upper bounds of wavelet estimation on \( L^{p}({\mathbb {R}}^{d}) \) risk of regression functions with strong mixing data for \( 1\le p<\infty \). In contrast to the independent case, these upper bounds have different analytic formulae for \(p\in [1, 2]\) and \(p\in (2, +\infty )\). For \(p=2\), it turns out that our result reduces to a theorem of Chaubey et al. (J Nonparametr Stat 25:53–71, 2013); and for \(d=1\) and \(p=2\), it becomes the corresponding theorem of Chaubey and Shirazi (Commun Stat Theory Methods 44:885–899, 2015).     

Bayes and maximum likelihood for $$L^1$$-Wasserstein deconvolution of Laplace mixtures

We consider the problem of recovering a distribution function on the real line from observations additively contaminated with errors following the standard Laplace distribution. Assuming that the latent distribution is completely unknown leads to a nonparametric deconvolution problem. We begin by studying the rates of convergence relative to the \(L^2\)-norm and the Hellinger metric for the direct problem of estimating the sampling density, which is a mixture of Laplace densities with a possibly unbounded set of locations: the rate of convergence for the Bayes’ density estimator corresponding to a Dirichlet process prior over the space of all mixing distributions on the real line matches, up to a logarithmic factor, with the \(n^{-3/8}\log ^{1/8}n\) rate for the maximum likelihood estimator. Then, appealing to an inversion inequality translating the \(L^2\)-norm and the Hellinger distance between general kernel mixtures, with a kernel density having polynomially decaying Fourier transform, into any \(L^p\)-Wasserstein distance, \(p\ge 1\), between the corresponding mixing distributions, provided their Laplace transforms are finite in some neighborhood of zero, we derive the rates of convergence in the \(L^1\)-Wasserstein metric for the Bayes’ and maximum likelihood estimators of the mixing distribution. Merging in the \(L^1\)-Wasserstein distance between Bayes and maximum likelihood follows as a by-product, along with an assessment on the stochastic order of the discrepancy between the two estimation procedures.     

Some lower bounds of centered $$L_2$$-discrepancy of $$2^{s-k}$$2s-k designs and their complementary designs

The indicator function is an effective tool in studying factorial designs. This paper presents some lower bounds of centered \(L_2\)-discrepancy through indicator function. Some new lower bounds of centered \(L_2\)-discrepancy for \(2^{s-k}\) designs and their complementary designs are given. Numerical results show that our lower bounds are tight and better than the existing results.     

A coordinate descent algorithm for computing penalized smooth quantile regression

The computation of penalized quantile regression estimates is often computationally intensive in high dimensions. In this paper we propose a coordinate descent algorithm for computing the penalized smooth quantile regression (cdaSQR) with convex and nonconvex penalties. The cdaSQR approach is based on the approximation of the objective check function, which is not differentiable at zero, by a modified check function which is differentiable at zero. Then, using the maximization-minimization trick of the gcdnet algorithm (Yang and Zou in, J Comput Graph Stat 22(2):396–415, 2013), we update each coefficient simply and efficiently. In our implementation, we consider the convex penalties \(\ell _1+\ell _2\) and the nonconvex penalties SCAD (or MCP) \(+ \ell _2\). We establishe the convergence property of the csdSQR with \(\ell _1+\ell _2\) penalty. The numerical results show that our implementation is an order of magnitude faster than its competitors. Using simulations we compare the speed of our algorithm to its competitors. Finally, the performance of our algorithm is illustrated on three real data sets from diabetes, leukemia and Bardet–Bidel syndrome gene expression studies.     

On the decomposition by sources of the Zenga 1984 point and synthetic inequality indexes

In this work we provide a decomposition by sources of the inequality index \(\zeta \) defined by Zenga (Giornale degli Economisti e Annali di economia 43(5–6):301–326, 1984). The source contributions are obtained with the method proposed in Zenga et al. (Stat Appl X(1):3–31, 2012) and Zenga (Stat Appl XI(2):133–161, 2013), that allows to compare different inequality measures. This method is based on the decomposition of inequality curves. To apply this decomposition to the index \(\zeta \) and its inequality curve, we adapt the method to the “cograduation” table. Moreover, we consider the case of linear transformation of sources and analyse the corresponding results.     

Finie mixures of canonical fundamenal skew $$$$-disribuions

This paper introduces a finite mixture of canonical fundamental skew \(t\) (CFUST) distributions for a model-based approach to clustering where the clusters are asymmetric and possibly long-tailed (in: Lee and McLachlan, arXiv:1401.8182 [statME], 2014b). The family of CFUST distributions includes the restricted multivariate skew \(t\) and unrestricted multivariate skew \(t\) distributions as special cases. In recent years, a few versions of the multivariate skew \(t\) (MST) mixture model have been put forward, together with various EM-type algorithms for parameter estimation. These formulations adopted either a restricted or unrestricted characterization for their MST densities. In this paper, we examine a natural generalization of these developments, employing the CFUST distribution as the parametric family for the component distributions, and point out that the restricted and unrestricted characterizations can be unified under this general formulation. We show that an exact implementation of the EM algorithm can be achieved for the CFUST distribution and mixtures of this distribution, and present some new analytical results for a conditional expectation involved in the E-step.     

Variable selection for survival data with a class of adaptive elastic net techniques

The accelerated failure time (AFT) models have proved useful in many contexts, though heavy censoring (as for example in cancer survival) and high dimensionality (as for example in microarray data) cause difficulties for model fitting and model selection. We propose new approaches to variable selection for censored data, based on AFT models optimized using regularized weighted least squares. The regularized technique uses a mixture of \(\ell _1\) and \(\ell _2\) norm penalties under two proposed elastic net type approaches. One is the adaptive elastic net and the other is weighted elastic net. The approaches extend the original approaches proposed by Ghosh (Adaptive elastic net: an improvement of elastic net to achieve oracle properties, Technical Reports 2007) and Hong and Zhang (Math Model Nat Phenom 5(3):115–133 2010), respectively. We also extend the two proposed approaches by adding censoring observations as constraints into their model optimization frameworks. The approaches are evaluated on microarray and by simulation. We compare the performance of these approaches with six other variable selection techniques-three are generally used for censored data and the other three are correlation-based greedy methods used for high-dimensional data.     

Optimal regular graph designs

A typical problem in optimal design theory is finding an experimental design that is optimal with respect to some criteria in a class of designs. The most popular criteria include the A- and D-criteria. Regular graph designs occur in many optimality results, and if the number of blocks is large enough, an A-optimal (or D-optimal) design is among them (if any exist). To explore the landscape of designs with a large number of blocks, we introduce extensions of regular graph designs. These are constructed by adding the blocks of a balanced incomplete block design repeatedly to the original design. We present the results of an exact computer search for the best regular graph designs and the best extended regular graph designs with up to 20 treatments v, block size \(k \le 10\) and replication r \(\le 10\) and \(r(k-1)-(v-1)\lfloor r(k-1)/(v-1)\rfloor \le 9\).     

Fast parallel $$\alpha $$-stable distribution function evaluation and parameter estimation using OpenCL in GPGPUs

\(\alpha \)-Stable distributions are a family of probability distributions found to be suitable to model many complex processes and phenomena in several research fields, such as medicine, physics, finance and networking, among others. However, the lack of closed expressions makes their evaluation analytically intractable, and alternative approaches are computationally expensive. Existing numerical programs are not fast enough for certain applications and do not make use of the parallel power of general purpose graphic processing units. In this paper, we develop novel parallel algorithms for the probability density function and cumulative distribution function—including a parallel Gauss–Kronrod quadrature—, quantile function, random number generator and maximum likelihood estimation of \(\alpha \)-stable distributions using OpenCL, achieving significant speedups and precision in all cases. Thanks to the use of OpenCL, we also evaluate the results of our library with different GPU architectures.