Fast and Stable Multiple Smoothing Parameter Selection in Smoothing Spline Analysis of Variance Models With Large Samples

The current parameterization and algorithm used to fit a smoothing spline analysis of variance (SSANOVA) model are computationally expensive, making a generalized additive model (GAM) the preferred method for multivariate smoothing. In this article, we propose an efficient reparameterization of the smoothing parameters in SSANOVA models, and a scalable algorithm for estimating multiple smoothing parameters in SSANOVAs. To validate our approach, we present two simulation studies comparing our reparameterization and algorithm to implementations of SSANOVAs and GAMs that are currently available in R. Our simulation results demonstrate that (a) our scalable SSANOVA algorithm outperforms the currently used SSANOVA algorithm, and (b) SSANOVAs can be a fast and reliable alternative to GAMs. We also provide an example with oceanographic data that demonstrates the practical advantage of our SSANOVA framework. Supplementary materials that are available online can be used to replicate the analyses in this article.     

Sparse Multivariate Regression With Covariance Estimation

We propose a procedure for constructing a sparse estimator of a multivariate regression coefficient matrix that accounts for correlation of the response variables. This method, which we call multivariate regression with covariance estimation (MRCE), involves penalized likelihood with simultaneous estimation of the regression coefficients and the covariance structure. An efficient optimization algorithm and a fast approximation are developed for computing MRCE. Using simulation studies, we show that the proposed method outperforms relevant competitors when the responses are highly correlated. We also apply the new method to a finance example on predicting asset returns. An R-package containing this dataset and code for computing MRCE and its approximation are available online.     

A Fast Algorithm for S-Regression Estimates

Equivariant high-breakdown point regression estimates are computationally expensive, and the corresponding algorithms become unfeasible for moderately large number of regressors. One important advance to improve the computational speed of one such estimator is the fast-LTS algorithm. This article proposes an analogous algorithm for computing S-estimates. The new algorithm, that we call “fast-S”, is also based on a “local improvement” step of the resampling initial candidates. This allows for a substantial reduction of the number of candidates required to obtain a good approximation to the optimal solution. We performed a simulation study which shows that S-estimators computed with the fast-S algorithm compare favorably to the LTS-estimators computed with the fast-LTS algorithm.     

A combined splitting—cross entropy method for rare-event probability estimation of queueing networks

We present a fast algorithm for the efficient estimation of rare-event (buffer overflow) probabilities in queueing networks. Our algorithm presents a combined version of two well known methods: the splitting and the cross-entropy (CE) method. We call the new method SPLITCE. In this method, the optimal change of measure (importance sampling) is determined adaptively by using the CE method. Simulation results for a single queue and queueing networks of the ATM-type are presented. Our numerical results demonstrate higher efficiency of the proposed method as compared to the original splitting and CE methods. In particular, for a single server queue example we demonstrate numerically that both the splitting and the SPLITCE methods can handle our buffer overflow example problems with both light and heavy tails efficiently. Further research must show the full potential of the proposed method.     

Adaptive Splines and Genetic Algorithms

Most existing algorithms for fitting adaptive splines are based on nonlinear optimization and/or stepwise selection. Stepwise knot selection, although computationally fast, is necessarily suboptimal while determining the best model over the space of adaptive knot splines is a very poorly behaved nonlinear optimization problem. A possible alternative is to use a genetic algorithm to perform knot selection. An adaptive modeling technique referred to as adaptive genetic splines (AGS) is introduced which combines the optimization power of a genetic algorithm with the flexibility of polynomial splines. Preliminary simulation results comparing the performance of AGS to those of existing methods such as HAS, SUREshrink and automatic Bayesian curve fitting are discussed. A real data example involving the application of these methods to a fMRI dataset is presented.     

Fast Hamiltonian Monte Carlo Using GPU Computing

In recent years, the Hamiltonian Monte Carlo (HMC) algorithm has been found to work more efficiently compared to other popular Markov chain Monte Carlo (MCMC) methods (such as random walk Metropolis–Hastings) in generating samples from a high-dimensional probability distribution. HMC has proven more efficient in terms of mixing rates and effective sample size than previous MCMC techniques, but still may not be sufficiently fast for particularly large problems. The use of GPUs promises to push HMC even further greatly increasing the utility of the algorithm. By expressing the computationally intensive portions of HMC (the evaluations of the probability kernel and its gradient) in terms of linear or element-wise operations, HMC can be made highly amenable to the use of graphics processing units (GPUs). A multinomial regression example demonstrates the promise of GPU-based HMC sampling. Using GPU-based memory objects to perform the entire HMC simulation, most of the latency penalties associated with transferring data from main to GPU memory can be avoided. Thus, the proposed computational framework may appear conceptually very simple, but has the potential to be applied to a wide class of hierarchical models relying on HMC sampling. Models whose posterior density and corresponding gradients can be reduced to linear or element-wise operations are amenable to significant speed ups through the use of GPUs. Analyses of datasets that were previously intractable for fully Bayesian approaches due to the prohibitively high computational cost are now feasible using the proposed framework.     

Lower Rank Approximation of Matrices Based on Fast and Robust Alternating Regression

We introduce fast and robust algorithms for lower rank approximation to given matrices based on robust alternating regression. The alternating least squares regression, also called criss-cross regression, was used for lower rank approximation of matrices, but it lacks robustness against outliers in these matrices. We use robust regression estimators and address some of the complications arising from this approach. We find it helpful to use high breakdown estimators in the initial iterations, followed by M estimators with monotone score functions in later iterations towards convergence. In addition to robustness, the computational speed is another important consideration in the development of our proposed algorithm, because alternating robust regression can be computationally intensive for large matrices. Based on a mix of the least trimmed squares (LTS) and Huber's M estimators, we demonstrate that fast and robust lower rank approximations are possible for modestly large matrices.     

A Reactive Local Search-Based Algorithm for the Multiple-Choice Multi-Dimensional Knapsack Problem

In this paper, we approximately solve the multiple-choice multi-dimensional knapsack problem. We propose an algorithm which is based upon reactive local search and where an explicit check for the repetition of configurations is added to the local search. The algorithm starts by an initial solution and improved by using a fast iterative procedure. Later, both deblocking and degrading procedures are introduced in order (i) to escape to local optima and, (ii) to introduce diversification in the search space. Finally, a memory list is applied in order to forbid the repetition of configurations. The performance of the proposed approaches has been evaluated on several problem instances. Encouraging results have been obtained.     

A Geometric Approach to Maximum Likelihood Estimation of the Functional Principal Components From Sparse Longitudinal Data

In this article, we consider the problem of estimating the eigenvalues and eigenfunctions of the covariance kernel (i.e., the functional principal components) from sparse and irregularly observed longitudinal data. We exploit the smoothness of the eigenfunctions to reduce dimensionality by restricting them to a lower dimensional space of smooth functions. We then approach this problem through a restricted maximum likelihood method. The estimation scheme is based on a Newton–Raphson procedure on the Stiefel manifold using the fact that the basis coefficient matrix for representing the eigenfunctions has orthonormal columns. We also address the selection of the number of basis functions, as well as that of the dimension of the covariance kernel by a second-order approximation to the leave-one-curve-out cross-validation score that is computationally very efficient. The effectiveness of our procedure is demonstrated by simulation studies and an application to a CD4+ counts dataset. In the simulation studies, our method performs well on both estimation and model selection. It also outperforms two existing approaches: one based on a local polynomial smoothing, and another using an EM algorithm. Supplementary materials including technical details, the R package fpca, and data analyzed by this article are available online.     

Regression Adjustment for Noncrossing Bayesian Quantile Regression

A two-stage approach is proposed to overcome the problem in quantile regression, where separately fitted curves for several quantiles may cross. The standard Bayesian quantile regression model is applied in the first stage, followed by a Gaussian process regression adjustment, which monotonizes the quantile function while borrowing strength from nearby quantiles. The two-stage approach is computationally efficient, and more general than existing techniques. The method is shown to be competitive with alternative approaches via its performance in simulated examples. Supplementary materials for the article are available online.     

Fast Algorithms for Nonparametric Curve Estimation

Abstract

Naive implementations of local polynomial fits and kernel estimators require almost O(n ²) operations. In this article two fast O(n) algorithms for nonparametric local polynomial fitting are presented. They are based on updating normal equations. Numerical stability is guaranteed by controlling ill-conditioned situations for small bandwidths and data-tuned restarting of the updating procedure. Restarting at every output point results in a moderately fast but highly stable O(n ^7/5) algorithm. Applicability of algorithms is evaluated for estimation of regression curves and their derivatives. The idea is also applied to kernel estimators of regression curves and densities.     

Markov chain Monte Carlo Using an Approximation

This article presents a method for generating samples from an unnormalized posterior distribution f(·) using Markov chain Monte Carlo (MCMC) in which the evaluation of f(·) is very difficult or computationally demanding. Commonly, a less computationally demanding, perhaps local, approximation to f(·) is available, say f^**_x(·). An algorithm is proposed to generate an MCMC that uses such an approximation to calculate acceptance probabilities at each step of a modified Metropolis–Hastings algorithm. Once a proposal is accepted using the approximation, f(·) is calculated with full precision ensuring convergence to the desired distribution. We give sufficient conditions for the algorithm to converge to f(·) and give both theoretical and practical justifications for its usage. Typical applications are in inverse problems using physical data models where computing time is dominated by complex model simulation. We outline Bayesian inference and computing for inverse problems. A stylized example is given of recovering resistor values in a network from electrical measurements made at the boundary. Although this inverse problem has appeared in studies of underground reservoirs, it has primarily been chosen for pedagogical value because model simulation has precisely the same computational structure as a finite element method solution of the complete electrode model used in conductivity imaging, or “electrical impedance tomography.” This example shows a dramatic decrease in CPU time, compared to a standard Metropolis–Hastings algorithm.     

Classical Backfitting for Smooth-Backfitting Additive Models

Smooth backfitting has been shown to have better theoretical properties than classical backfitting for fitting additive models based on local linear regression. In this article, we show that the smooth backfitting procedure in the local linear case can be alternatively performed as a classical backfitting procedure with a different type of smoother matrices. These smoother matrices are symmetric and shrinking and some established results in the literature are readily applicable. The connections allow the smooth backfitting algorithm to be implemented in a much simplified way, give new insights on the differences between the two approaches in the literature, and provide an extension to local polynomial regression. The connections also give rise to a new estimator at data points. Asymptotic properties of general local polynomial smooth backfitting estimates are investigated, allowing for different orders of local polynomials and different bandwidths. Cases of oracle efficiency are discussed. Computer-generated simulations are conducted to demonstrate finite sample behaviors of the methodology and a real data example is given for illustration. Supplementary materials for this article are available online.     

A surrogate and Lagrangian approach to constrained network problems

This paper presents the use of surrogate constraints and Lagrange multipliers to generate advanced starting solutions to constrained network problems. The surrogate constraint approach is used to generate a singly constrained network problem which is solved using the algorithm of Glover, Karney, Klingman and Russell [13]. In addition, we test the use of the Lagrangian function to generate advanced starting solutions. In the Lagrangian approach, the subproblems are capacitated network problems which can be solved using very efficient algorithms.The surrogate constraint approach is implemented using the multiplier update procedure of Held, Wolfe and Crowder [16]. The procedure is modified to include a search in a single direction to prevent periodic regression of the solution. We also introduce a reoptimization procedure which allows the solution from thekth subproblem to be used as the starting point for the next surrogate problem for which it is infeasible once the new surrogate constraint is adjoined.The algorithms are tested under a variety of conditions including: large-scale problems, number and structure of the non-network constraints, and the density of the non-network constraint coefficients.The testing clearly demonstrates that both the surrogate constraint and Langrange multipliers generate advanced starting solutions which greatly improve the computational effort required to generate an optimal solution to the constrained network problem. The testing demonstrates that the extra effort required to solve the singly constrained network subproblems of the surrogate constraints approach yields an improved advanced starting point as compared to the Lagrangian approach. It is further demonstrated that both of the relaxation approaches are much more computationally efficient than solving the problem from the beginning with a linear programming algorithm.     

Robust Parametric Classification and Variable Selection by a Minimum Distance Criterion

We investigate a robust penalized logistic regression algorithm based on a minimum distance criterion. Influential outliers are often associated with the explosion of parameter vector estimates, but in the context of standard logistic regression, the bias due to outliers always causes the parameter vector to implode, that is, shrink toward the zero vector. Thus, using LASSO-like penalties to perform variable selection in the presence of outliers can result in missed detections of relevant covariates. We show that by choosing a minimum distance criterion together with an elastic net penalty, we can simultaneously find a parsimonious model and avoid estimation implosion even in the presence of many outliers in the important small n large p situation. Minimizing the penalized minimum distance criterion is a challenging problem due to its nonconvexity. To meet the challenge, we develop a simple and efficient MM (majorization–minimization) algorithm that can be adapted gracefully to the small n large p context. Performance of our algorithm is evaluated on simulated and real datasets. This article has supplementary materials available online.     

On the numerical relaxation of single-slip plasticity in finite strains

The modeling of the elastoplastic behaviour of single crystals with infinite latent hardening leads to a nonconvex energy density, whose minimization produces fine structures. The computation of the quasiconvex envelope of the energy density is faced in this case with huge numerical difficulties caused by the clusters of local minima. By exploiting the structure of the problem, we present a fast and efficient numerical relaxation algorithm as alternative to global optimization techniques usually adopted in literature which are computationally more expensive. (© 2006 WILEY-VCH Verlag GmbH & Co. KGaA, Weinheim)     

Optimal estimation of sensor biases for asynchronous multi-sensor data fusion

An important step in a multi-sensor surveillance system is to estimate sensor biases from their noisy asynchronous measurements. This estimation problem is computationally challenging due to the highly nonlinear transformation between the global and local coordinate systems as well as the measurement asynchrony from different sensors. In this paper, we propose a novel nonlinear least squares formulation for the problem by assuming the existence of a reference target moving with an (unknown) constant velocity. We also propose an efficient block coordinate decent (BCD) optimization algorithm, with a judicious initialization, to solve the problem. The proposed BCD algorithm alternately updates the range and azimuth bias estimates by solving linear least squares problems and semidefinite programs. In the absence of measurement noise, the proposed algorithm is guaranteed to find the global solution of the problem and the true biases. Simulation results show that the proposed algorithm significantly outperforms the existing approaches in terms of the root mean square error.     

Semismooth Newton Coordinate Descent Algorithm for Elastic-Net Penalized Huber Loss Regression and Quantile Regression

We propose an algorithm, semismooth Newton coordinate descent (SNCD), for the elastic-net penalized Huber loss regression and quantile regression in high dimensional settings. Unlike existing coordinate descent type algorithms, the SNCD updates a regression coefficient and its corresponding subgradient simultaneously in each iteration. It combines the strengths of the coordinate descent and the semismooth Newton algorithm, and effectively solves the computational challenges posed by dimensionality and nonsmoothness. We establish the convergence properties of the algorithm. In addition, we present an adaptive version of the “strong rule” for screening predictors to gain extra efficiency. Through numerical experiments, we demonstrate that the proposed algorithm is very efficient and scalable to ultrahigh dimensions. We illustrate the application via a real data example. Supplementary materials for this article are available online.     

The prize-collecting generalized minimum spanning tree problem

We introduce the prize-collecting generalized minimum spanning tree problem. In this problem a network of node clusters needs to be connected via a tree architecture using exactly one node per cluster. Nodes in each cluster compete by offering a payment for selection. This problem is NP-hard, and we describe several heuristic strategies, including local search and a genetic algorithm. Further, we present a simple and computationally efficient branch-and-cut algorithm. Our computational study indicates that our branch-and-cut algorithm finds optimal solutions for networks with up to 200 nodes within two hours of CPU time, while the heuristic search procedures rapidly find near-optimal solutions for all of the test instances.     

A Computationally Efficient Oracle Estimator for Additive Nonparametric Regression with Bootstrap Confidence Intervals

Abstract

This article makes three contributions. First, we introduce a computationally efficient estimator for the component functions in additive nonparametric regression exploiting a different motivation from the marginal integration estimator of Linton and Nielsen. Our method provides a reduction in computation of order n which is highly significant in practice. Second, we define an efficient estimator of the additive components, by inserting the preliminary estimator into a backfitting˙ algorithm but taking one step only, and establish that it is equivalent, in various senses, to the oracle estimator based on knowing the other components. Our two-step estimator is minimax superior to that considered in Opsomer and Ruppert, due to its better bias. Third, we define a bootstrap algorithm for computing pointwise confidence intervals and show that it achieves the correct coverage.