Auxiliary Variable Method Research Articles

Journal of Global Optimization Best Paper Award for a paper published in 2014

McCormick (Math Prog 10(1):147–175, 1976) provides the framework for convex/concave relaxations of factorable functions, via rules for the product of functions and compositions of the form F ◦ f , where F is a univariate function. Herein, the composition theorem is generalized to allow multivariate outer functions F , and theory for the propagation of subgradients is presented. The generalization interprets the McCormick relaxation approach as a decomposition method for the auxiliary variable method. In addition to extending the framework, the new result provides a tool for the proof of relaxations of specific functions. Moreover, a direct consequence is an improved relaxation for the product of two functions, at least as tight as McCormicks result, and often tighter. The result also allows the direct relaxation of multilinear products of functions. Furthermore, the composition result is applied to obtain improved convex underestimators for the minimum/maximum and the division of two functions for which current relaxations are often weak. These cases can be extended to allow composition of a variety of functions for which relaxations have been proposed. Congratulations to the authors of the awardedpapers! Iwould also like to thankSpringer for sponsoring this award and the committee members for their active involvement. Nominations for the 2015 JOGO Best Paper Award should be emailed to butenko@tamu.edu. All JOGO papers published in 2015 (volumes 61–63) are eligible.

Relative contributions of aleatory and epistemic uncertainty sources in time series prediction

This paper develops a novel computational framework to compute the Sobol indices that quantify the relative contributions of various uncertainty sources towards the system response prediction uncertainty. In the presence of both aleatory and epistemic uncertainty, two challenges are addressed in this paper for the model-based computation of the Sobol indices: due to data uncertainty, input distributions are not precisely known; and due to model uncertainty, the model output is uncertain even for a fixed realization of the input. An auxiliary variable method based on the probability integral transform is introduced to distinguish and represent each uncertainty source explicitly, whether aleatory or epistemic. The auxiliary variables facilitate building a deterministic relationship between the uncertainty sources and the output, which is needed in the Sobol indices computation. The proposed framework is developed for two types of model inputs: random variable input and time series input. A Bayesian autoregressive moving average (ARMA) approach is chosen to model the time series input due to its capability to represent both natural variability and epistemic uncertainty due to limited data. A novel controlled-seed computational technique based on pseudo-random number generation is proposed to efficiently represent the natural variability in the time series input. This controlled-seed method significantly accelerates the Sobol indices computation under time series input, and makes it computationally affordable.

Multivariate McCormick relaxations

McCormick (Math Prog 10(1):147---175, 1976) provides the framework for convex/concave relaxations of factorable functions, via rules for the product of functions and compositions of the form $$F\circ f$$F?f, where $$F$$F is a univariate function. Herein, the composition theorem is generalized to allow multivariate outer functions $$F$$F, and theory for the propagation of subgradients is presented. The generalization interprets the McCormick relaxation approach as a decomposition method for the auxiliary variable method. In addition to extending the framework, the new result provides a tool for the proof of relaxations of specific functions. Moreover, a direct consequence is an improved relaxation for the product of two functions, at least as tight as McCormick's result, and often tighter. The result also allows the direct relaxation of multilinear products of functions. Furthermore, the composition result is applied to obtain improved convex underestimators for the minimum/maximum and the division of two functions for which current relaxations are often weak. These cases can be extended to allow composition of a variety of functions for which relaxations have been proposed.

Adaptive Markov chain Monte Carlo for auxiliary variable method and its application to parallel tempering

Auxiliary variable methods such as the Parallel Tempering and the cluster Monte Carlo methods generate samples that follow a target distribution by using proposal and auxiliary distributions. In sampling from complex distributions, these algorithms are highly more efficient than the standard Markov chain Monte Carlo methods. However, their performance strongly depends on their parameters and determining the parameters is critical. In this paper, we proposed an algorithm for adapting the parameters during drawing samples and proved the convergence theorem of the adaptive algorithm. We applied our algorithm to the Parallel Tempering. That is, we developed an adaptive Parallel Tempering that tunes the parameters on the fly. We confirmed the effectiveness of our algorithm through the validation of the adaptive Parallel Tempering, comparing samples from the target distribution by the adaptive Parallel Tempering and samples by conventional algorithms.

Use of SAMC for Bayesian analysis of statistical models with intractable normalizing constants

Statistical inference for the models with intractable normalizing constants has attracted much attention. During the past two decades, various approximation- or simulation-based methods have been proposed for the problem, such as the Monte Carlo maximum likelihood method and the auxiliary variable Markov chain Monte Carlo methods. The Bayesian stochastic approximation Monte Carlo algorithm specifically addresses this problem: It works by sampling from a sequence of approximate distributions with their average converging to the target posterior distribution, where the approximate distributions can be achieved using the stochastic approximation Monte Carlo algorithm. A strong law of large numbers is established for the Bayesian stochastic approximation Monte Carlo estimator under mild conditions. Compared to the Monte Carlo maximum likelihood method, the Bayesian stochastic approximation Monte Carlo algorithm is more robust to the initial guess of model parameters. Compared to the auxiliary variable MCMC methods, the Bayesian stochastic approximation Monte Carlo algorithm avoids the requirement for perfect samples, and thus can be applied to many models for which perfect sampling is not available or very expensive. The Bayesian stochastic approximation Monte Carlo algorithm also provides a general framework for approximate Bayesian analysis.

On Bayesian Curve Fitting via Auxiliary Variables

In this article we revisit the auxiliary variable method introduced by Smith and Kohn (1996) for the fitting of Pth-order spline regression models with an unknown number of knot points. We introduce modifications which allow the location of knot points to be random, and we further consider an extension of the method to handle models with non-Gaussian errors. We provide a new algorithm for the MCMC sampling of such models. Simulated data examples are used to compare the performance of our method with existing ones. Finally, we make a connection with some change-point problems, and show how they can be reparameterized to the variable selection setting.Supplemental materials including R computing codes used in the examples are available online.

Variational Bayes for estimating the parameters of a hidden Potts model

Hidden Markov random field models provide an appealing representation of images and other spatial problems. The drawback is that inference is not straightforward for these models as the normalisation constant for the likelihood is generally intractable except for very small observation sets. Variational methods are an emerging tool for Bayesian inference and they have already been successfully applied in other contexts. Focusing on the particular case of a hidden Potts model with Gaussian noise, we show how variational Bayesian methods can be applied to hidden Markov random field inference. To tackle the obstacle of the intractable normalising constant for the likelihood, we explore alternative estimation approaches for incorporation into the variational Bayes algorithm. We consider a pseudo-likelihood approach as well as the more recent reduced dependence approximation of the normalisation constant. To illustrate the effectiveness of these approaches we present empirical results from the analysis of simulated datasets. We also analyse a real dataset and compare results with those of previous analyses as well as those obtained from the recently developed auxiliary variable MCMC method and the recursive MCMC method. Our results show that the variational Bayesian analyses can be carried out much faster than the MCMC analyses and produce good estimates of model parameters. We also found that the reduced dependence approximation of the normalisation constant outperformed the pseudo-likelihood approximation in our analysis of real and synthetic datasets.

A fourth-order auxiliary variable projection method for zero-Mach number gas dynamics

A fourth-order numerical method for the zero-Mach-number limit of the equations for compressible flow is presented. The method is formed by discretizing a new auxiliary variable formulation of the conservation equations, which is a variable density analog to the impulse or gauge formulation of the incompressible Euler equations. An auxiliary variable projection method is applied to this formulation, and accuracy is achieved by combining a fourth-order finite-volume spatial discretization with a fourth-order temporal scheme based on spectral deferred corrections. Numerical results are included which demonstrate fourth-order spatial and temporal accuracy for non-trivial flows in simple geometries.

An efficient Markov chain Monte Carlo method for distributions with intractable normalising constants

Maximum likelihood parameter estimation and sampling from Bayesian posterior distributions are problematic when the probability density for the parameter of interest involves an intractable normalising constant which is also a function of that parameter. In this paper, an auxiliary variable method is presented which requires only that independent samples can be drawn from the unnormalised density at any particular parameter value. The proposal distribution is constructed so that the normalising constant cancels from the Metropolis-Hastings ratio. The method is illustrated by producing posterior samples for parameters of the Ising model given a particular lattice realisation.

Bayesian auxiliary variable models for binary and multinomial regression

In this paper we discuss auxiliary variable approaches to Bayesian binary and multinomial regression. These approaches are ideally suited to automated Markov chain Monte Carlo simulation. In the first part we describe a simple technique using joint updating that improves the performance of the conventional probit regression algorithm. In the second part we discuss auxiliary variable methods for inference in Bayesian logistic regression, including covariate set uncertainty. Finally, we show how the logistic method is easily extended to multinomial regression models. All of the algorithms are fully automatic with no user set parameters and no necessary Metropolis-Hastings accept/reject steps. © 2006 International Society for Bayesian Analysis.

Recursive Monte Carlo filters: Algorithms and theoretical analysis

Recursive Monte Carlo filters, also called particle filters, are a powerful tool to perform computations in general state space models. We discuss and compare the accept–reject version with the more common sampling importance resampling version of the algorithm. In particular, we show how auxiliary variable methods and stratification can be used in the accept–reject version, and we compare different resampling techniques. In a second part, we show laws of large numbers and a central limit theorem for these Monte Carlo filters by simple induction arguments that need only weak conditions. We also show that, under stronger conditions, the required sample size is independent of the length of the observed series.

Slice sampling for simulation based fitting of spatial data models

An auxiliary variable method based on a slice sampler is shown to provide an attractive simulation-based model fitting strategy for fitting Bayesian models under proper priors. Though broadly applicable, we illustrate in the context of fitting spatial models for geo-referenced or point source data. Spatial modeling within a Bayesian framework offers inferential advantages and the slice sampler provides an algorithm which is essentially "off the shelf". Further potential advantages over importance sampling approaches and Metropolis approaches are noted and illustrative examples are supplied.

Bayesian Variable Selection and the Swendsen-Wang Algorithm

The need to explore model uncertainty in linear regression models with many predictors has motivated improvements in Markov chain Monte Carlo sampling algorithms for Bayesian variable selection. Currently used sampling algorithms for Bayesian variable selection may perform poorly when there are severe multicollinearities among the predictors. This article describes a new sampling method based on an analogy with the Swendsen-Wang algorithm for the Ising model, and which can give substantial improvements over alternative sampling schemes in the presence of multicollinearity. In linear regression with a given set of potential predictors we can index possible models by a binary parameter vector that indicates which of the predictors are included or excluded. By thinking of the posterior distribution of this parameter as a binary spatial field, we can use auxiliary variable methods inspired by the Swendsen-Wang algorithm for the Ising model to sample from the posterior where dependence among parameters is reduced by conditioning on auxiliary variables. Performance of the method is described for both simulated and real data.

Covariantly quantized spinning particle and its possible connection to noncommutative space-time

Covariant quantization of the Nambu-Goto spinning particle in 2+1 dimensions is studied. The model is relevant in the context of recent activities in noncommutative space-time. From a technical point of view also covariant quantization of the model poses an interesting problem: the set of second class constraints (in the Dirac classification scheme) is reducible. The reducibility problem is analyzed from two contrasting approaches: (i) the auxiliary variable method and (ii) the projection operator method. Finally in the former scheme, a Batalin-Tyutin quantization has been done. This induces a mapping between the non-commutative and the ordinary space-time. The Becchi-Rouet-Stora-Tyutin (BRST) quantization program in the latter scheme has also been discussed.

The Art of Data Augmentation

The term data augmentation refers to methods for constructing iterative optimization or sampling algorithms via the introduction of unobserved data or latent variables. For deterministic algorithms, the method was popularized in the general statistical community by the seminal article by Dempster, Laird, and Rubin on the EM algorithm for maximizing a likelihood function or, more generally, a posterior density. For stochastic algorithms, the method was popularized in the statistical literature by Tanner and Wong's Data Augmentation algorithm for posterior sampling and in the physics literature by Swendsen and Wang's algorithm for sampling from the Ising and Potts models and their generalizations; in the physics literature, the method of data augmentation is referred to as the method of auxiliary variables. Data augmentation schemes were used by Tanner and Wong to make simulation feasible and simple, while auxiliary variables were adopted by Swendsen and Wang to improve the speed of iterative simulation. In general, however, constructing data augmentation schemes that result in both simple and fast algorithms is a matter of art in that successful strategies vary greatly with the (observed-data) models being considered. After an overview of data augmentation/auxiliary variables and some recent developments in methods for constructing such efficient data augmentation schemes, we introduce an effective search strategy that combines the ideas of marginal augmentation and conditional augmentation, together with a deterministic approximation method for selecting good augmentation schemes. We then apply this strategy to three common classes of models (specifically, multivariate t, probit regression, and mixed-effects models) to obtain efficient Markov chain Monte Carlo algorithms for posterior sampling. We provide theoretical and empirical evidence that the resulting algorithms, while requiring similar programming effort, can show dramatic improvement over the Gibbs samplers commonly used for these models in practice. A key feature of all these new algorithms is that they are positive recurrent subchains of nonpositive recurrent Markov chains constructed in larger spaces.

On Bayesian analyses and finite mixtures for proportions

When the results of biological experiments are tested for a possible difference between treatment and control groups, the inference is only valid if based upon a model that fits the experimental results satisfactorily. In dominant-lethal testing, foetal death has previously been assumed to follow a variety of models, including a Poisson, Binomial, Beta-binomial and various mixture models. However, discriminating between models has always been a particularly difficult problem. In this paper, we consider the data from 6 separate dominant-lethal assay experiments and discriminate between the competing models which could be used to describe them. We adopt a Bayesian approach and illustrate how a variety of different models may be considered, using Markov chain Monte Carlo (MCMC) simulation techniques and comparing the results with the corresponding maximum likelihood analyses. We present an auxiliary variable method for determining the probability that any particular data cell is assigned to a given component in a mixture and we illustrate the value of this approach. Finally, we show how the Bayesian approach provides a natural and unique perspective on the model selection problem via reversible jump MCMC and illustrate how probabilities associated with each of the different models may be calculated for each data set. In terms of estimation we show how, by averaging over the different models, we obtain reliable and robust inference for any statistic of interest.

Seeking efficient data augmentation schemes via conditional and marginal augmentation

Data augmentation, sometimes known as the method of auxiliary variables, is a powerful tool for constructing optimisation and simulation algorithms. In the context of optimisation, Meng & van Dyk (1997, 1998) reported several successes of the 'working parameter' approach for constructing efficient data-augmentation schemes for fast and simple EM-type algorithms. This paper investigates the use of working parameters in the context of Markov chain Monte Carlo, in particular in the context of Tanner & Wong's (1987) data augmentation algorithm, via a theoretical study of two working-parameter approaches, the conditional augmentation approach and the marginal augmentation approach. Posterior sampling under the univariate t model is used as a running example, which particularly illustrates how the marginal augmentation approach obtains a fast-mixing positive recurrent Markov chain by first constructing a nonpositive recurrent Markov chain in a larger space.

Auxiliary Variable Methods for Markov Chain Monte Carlo with Applications

Abstract Suppose that one wishes to sample from the density π(x) using Markov chain Monte Carlo (MCMC). An auxiliary variable u and its conditional distribution π(u|x) can be defined, giving the joint distribution π(x, u) = π(x)π(u|x). A MCMC scheme that samples over this joint distribution can lead to substantial gains in efficiency compared to standard approaches. The revolutionary algorithm of Swendsen and Wang is one such example. Besides reviewing the Swendsen-Wang algorithm and its generalizations, this article introduces a new auxiliary variable method called partial decoupling. Two applications in Bayesian image analysis are considered: a binary classification problem in which partial decoupling out performs Swendsen-Wang and single-site Metropolis methods, and a positron emission tomography (PET) reconstruction that uses the gray level prior of Geman and McClure. A generalized Swendsen–Wang algorithm is developed for this problem, which reduces the computing time to the point where MCMC is a viable method of posterior exploration.

Auxiliary Variable Methods for Markov Chain Monte Carlo with Applications

Suppose that one wishes to sample from the density π(x) using Markov chain Monte Carlo (MCMC). An auxiliary variable u and its conditional distribution π(u|x) can be defined, giving the joint distribution π(x, u) = π(x)π(u|x). A MCMC scheme that samples over this joint distribution can lead to substantial gains in efficiency compared to standard approaches. The revolutionary algorithm of Swendsen and Wang is one such example. Besides reviewing the Swendsen-Wang algorithm and its generalizations, this article introduces a new auxiliary variable method called partial decoupling. Two applications in Bayesian image analysis are considered: a binary classification problem in which partial decoupling out performs Swendsen-Wang and single-site Metropolis methods, and a positron emission tomography (PET) reconstruction that uses the gray level prior of Geman and McClure. A generalized Swendsen–Wang algorithm is developed for this problem, which reduces the computing time to the point where MCMC is a viable method of posterior exploration.

Analytical solution for an infinite elastic power-law hardening plate containing an elastic circular inhomogeneity and subjected to equi-biaxial tension

This paper presents an clasto-plastic analytical solution for the plane stress inclusion problem of an elastic power-law hardening plate containing an elastic circular inhomogeneity and subjected to equi-biaxial far-field tension. Hencky's deformation theory (for compressible materials) and von Mises' yield criterion are applied, and infinitesimal deformations are assumed. The solution is derived by using a stress formulation and with the help of a modified Nadai's auxiliary-variable method and the extended Michell theorem. All expressions for the stress, strain and displacement components are derived in explicit forms in terms of an auxiliary variable and four constant parameters, which are determined from the given boundary conditions. Three specific solutions of practical interest are presented as limiting cases, one of which is the closed-form solution for the plate containing a traction-free circular hole. Numerical results are also provided to demonstrate quantitatively applications of the solution in the opening and reinforcement design of spherical pressure vessels.