Bayesian Variable Selection via Particle Stochastic Search

We focus on Bayesian variable selection in regression models. One challenge is to search the huge model space adequately, while identifying high posterior probability regions. In the past decades, the main focus has been on the use of Markov chain Monte Carlo (MCMC) algorithms for these purposes. In this article, we propose a new computational approach based on sequential Monte Carlo (SMC), which we refer to as particle stochastic search (PSS). We illustrate PSS through applications to linear regression and probit models.     

Markov Chain Monte Carlo Convergence Assessment via Two-Way Analysis of Variance

Abstract

In this article we discuss the problem of assessing the performance of Markov chain Monte Carlo (MCMC) algorithms on the basis of simulation output. In essence, we extend the original ideas of Gelman and Rubin and, more recently, Brooks and Gelman, to problems where we are able to split the variation inherent within the MCMC simulation output into two distinct groups. We show how such a diagnostic may be useful in assessing the performance of MCMC samplers addressing model choice problems, such as the reversible jump MCMC algorithm. In the model choice context, we show how the reversible jump MCMC simulation output for parameters that retain a coherent interpretation throughout the simulation, can be used to assess convergence. By considering various decompositions of the sampling variance of this parameter, we can assess the performance of our MCMC sampler in terms of its mixing properties both within and between models and we illustrate our approach in both the graphical Gaussian models and normal mixtures context. Finally, we provide an example of the application of our diagnostic to the assessment of the influence of different starting values on MCMC simulation output, thereby illustrating the wider utility of our method beyond the Bayesian model choice and reversible jump MCMC context.     

Adaptive approximate Bayesian computation for complex models

We propose a new approximate Bayesian computation (ABC) algorithm that aims at minimizing the number of model runs for reaching a given quality of the posterior approximation. This algorithm automatically determines its sequence of tolerance levels and makes use of an easily interpretable stopping criterion. Moreover, it avoids the problem of particle duplication found when using a MCMC kernel. When applied to a toy example and to a complex social model, our algorithm is 2–8 times faster than the three main sequential ABC algorithms currently available.     

On Dynamic Generalized Linear Models with Applications

In this paper we combine the idea of ‘power steady model’, ‘discount factor’ and ‘power prior’, for a general class of filter model, more specifically within a class of dynamic generalized linear models (DGLM). We show an optimality property for our proposed method and present the particle filter algorithm for DGLM as an alternative to Markov chain Monte Carlo method. We also present two applications; one on dynamic Poisson models for hurricane count data in Atlantic ocean and the another on the dynamic Poisson regression model for longitudinal count data.     

GPU-Powered Shotgun Stochastic Search for Dirichlet Process Mixtures of Gaussian Graphical Models

Gaussian graphical models (GGMs) are popular for modeling high-dimensional multivariate data with sparse conditional dependencies. A mixture of GGMs extends this model to the more realistic scenario where observations come from a heterogenous population composed of a small number of homogeneous subgroups. In this article, we present a novel stochastic search algorithm for finding the posterior mode of high-dimensional Dirichlet process mixtures of decomposable GGMs. Further, we investigate how to harness the massive thread-parallelization capabilities of graphical processing units to accelerate computation. The computational advantages of our algorithms are demonstrated with various simulated data examples in which we compare our stochastic search with a Markov chain Monte Carlo (MCMC) algorithm in moderate dimensional data examples. These experiments show that our stochastic search largely outperforms the MCMC algorithm in terms of computing-times and in terms of the quality of the posterior mode discovered. Finally, we analyze a gene expression dataset in which MCMC algorithms are too slow to be practically useful.     

Variational Approximation for Mixtures of Linear Mixed Models

Mixtures of linear mixed models (MLMMs) are useful for clustering grouped data and can be estimated by likelihood maximization through the Expectation–Maximization algorithm. A suitable number of components is then determined conventionally by comparing different mixture models using penalized log-likelihood criteria such as Bayesian information criterion. We propose fitting MLMMs with variational methods, which can perform parameter estimation and model selection simultaneously. We describe a variational approximation for MLMMs where the variational lower bound is in closed form, allowing for fast evaluation and develop a novel variational greedy algorithm for model selection and learning of the mixture components. This approach handles algorithm initialization and returns a plausible number of mixture components automatically. In cases of weak identifiability of certain model parameters, we use hierarchical centering to reparameterize the model and show empirically that there is a gain in efficiency in variational algorithms similar to that in Markov chain Monte Carlo (MCMC) algorithms. Related to this, we prove that the approximate rate of convergence of variational algorithms by Gaussian approximation is equal to that of the corresponding Gibbs sampler, which suggests that reparameterizations can lead to improved convergence in variational algorithms just as in MCMC algorithms. Supplementary materials for the article are available online.     

Particle filtering in the Dempster–Shafer theory

This paper derives a particle filter algorithm within the Dempster–Shafer framework. Particle filtering is a well-established Bayesian Monte Carlo technique for estimating the current state of a hidden Markov process using a fixed number of samples. When dealing with incomplete information or qualitative assessments of uncertainty, however, Dempster–Shafer models with their explicit representation of ignorance often turn out to be more appropriate than Bayesian models.The contribution of this paper is twofold. First, the Dempster–Shafer formalism is applied to the problem of maintaining a belief distribution over the state space of a hidden Markov process by deriving the corresponding recursive update equations, which turn out to be a strict generalization of Bayesian filtering. Second, it is shown how the solution of these equations can be efficiently approximated via particle filtering based on importance sampling, which makes the Dempster–Shafer approach tractable even for large state spaces. The performance of the resulting algorithm is compared to exact evidential as well as Bayesian inference.     

A Generalized Markov Sampler

A recent development of the Markov chain Monte Carlo (MCMC) technique is the emergence of MCMC samplers that allow transitions between different models. Such samplers make possible a range of computational tasks involving models, including model selection, model evaluation, model averaging and hypothesis testing. An example of this type of sampler is the reversible jump MCMC sampler, which is a generalization of the Metropolis–Hastings algorithm. Here, we present a new MCMC sampler of this type. The new sampler is a generalization of the Gibbs sampler, but somewhat surprisingly, it also turns out to encompass as particular cases all of the well-known MCMC samplers, including those of Metropolis, Barker, and Hastings. Moreover, the new sampler generalizes the reversible jump MCMC. It therefore appears to be a very general framework for MCMC sampling. This paper describes the new sampler and illustrates its use in three applications in Computational Biology, specifically determination of consensus sequences, phylogenetic inference and delineation of isochores via multiple change-point analysis.     

Covariance Decompositions for Accurate Computation in Bayesian Scale-Usage Models

Analyses of multivariate ordinal probit models typically use data augmentation to link the observed (discrete) data to latent (continuous) data via a censoring mechanism defined by a collection of “cutpoints.” Most standard models, for which effective Markov chain Monte Carlo (MCMC) sampling algorithms have been developed, use a separate (and independent) set of cutpoints for each element of the multivariate response. Motivated by the analysis of ratings data, we describe a particular class of multivariate ordinal probit models where it is desirable to use a common set of cutpoints. While this approach is attractive from a data-analytic perspective, we show that the existing efficient MCMC algorithms can no longer be accurately applied. Moreover, we show that attempts to implement these algorithms by numerically approximating required multivariate normal integrals over high-dimensional rectangular regions can result in severely degraded estimates of the posterior distribution. We propose a new data augmentation that is based on a covariance decomposition and that admits a simple and accurate MCMC algorithm. Our data augmentation requires only that univariate normal integrals be evaluated, which can be done quickly and with high accuracy. We provide theoretical results that suggest optimal decompositions within this class of data augmentations, and, based on the theory, recommend default decompositions that we demonstrate work well in practice. This article has supplementary material online.     

Interweaving Markov Chain Monte Carlo Strategies for Efficient Estimation of Dynamic Linear Models

In dynamic linear models (DLMs) with unknown fixed parameters, a standard Markov chain Monte Carlo (MCMC) sampling strategy is to alternate sampling of latent states conditional on fixed parameters and sampling of fixed parameters conditional on latent states. In some regions of the parameter space, this standard data augmentation (DA) algorithm can be inefficient. To improve efficiency, we apply the interweaving strategies of Yu and Meng to DLMs. For this, we introduce three novel alternative DAs for DLMs: the scaled errors, wrongly scaled errors, and wrongly scaled disturbances. With the latent states and the less well known scaled disturbances, this yields five unique DAs to employ in MCMC algorithms. Each DA implies a unique MCMC sampling strategy and they can be combined into interweaving and alternating strategies that improve MCMC efficiency. We assess these strategies using the local level model and demonstrate that several strategies improve efficiency relative to the standard approach and the most efficient strategy interweaves the scaled errors and scaled disturbances. Supplementary materials are available online for this article.     

A Sequential Algorithm for Fast Fitting of Dirichlet Process Mixture Models

In this article, we propose an improvement on the sequential updating and greedy search (SUGS) algorithm for fast fitting of Dirichlet process mixture models. The SUGS algorithm provides a means for very fast approximate Bayesian inference for mixture data which is particularly of use when datasets are so large that many standard Markov chain Monte Carlo (MCMC) algorithms cannot be applied efficiently, or take a prohibitively long time to converge. In particular, these ideas are used to initially interrogate the data, and to refine models such that one can potentially apply exact data analysis later on. SUGS relies upon sequentially allocating data to clusters and proceeding with an update of the posterior on the subsequent allocations and parameters which assumes this allocation is correct. Our modification softens this approach, by providing a probability distribution over allocations, with a similar computational cost; this approach has an interpretation as a variational Bayes procedure and hence we term it variational SUGS (VSUGS). It is shown in simulated examples that VSUGS can outperform, in terms of density estimation and classification, a version of the SUGS algorithm in many scenarios. In addition, we present a data analysis for flow cytometry data, and SNP data via a three-class Dirichlet process mixture model, illustrating the apparent improvement over the original SUGS algorithm.     

On the Convergence Rate of Random Permutation Sampler and ECR Algorithm in Missing Data Models

Label switching is a well-known phenomenon that occurs in MCMC outputs targeting the parameters’ posterior distribution of many latent variable models. Although its appearence is necessary for the convergence of the simulated Markov chain, it turns out to be a problem in the estimation procedure. In a recent paper, Papastamoulis and Iliopoulos (J Comput Graph Stat 19:313–331, 2010) introduced the Equivalence Classes Representatives (ECR) algorithm as a solution of this problem in the context of finite mixtures of distributions. In this paper, label switching is considered under a general missing data model framework that includes as special cases finite mixtures, hidden Markov models, and Markov random fields. The use of ECR algorithm is extended to this general framework and is shown that the relabelled sequence which it produces converges to its target distribution at the same rate as the Random Permutation Sampler of Frühwirth-Schnatter (2001) and that both converge at least as fast as the Markov chain generated by the original MCMC output.     

Efficient Data Augmentation for Fitting Stochastic Epidemic Models to Prevalence Data

Stochastic epidemic models describe the dynamics of an epidemic as a disease spreads through a population. Typically, only a fraction of cases are observed at a set of discrete times. The absence of complete information about the time evolution of an epidemic gives rise to a complicated latent variable problem in which the state space size of the epidemic grows large as the population size increases. This makes analytically integrating over the missing data infeasible for populations of even moderate size. We present a data augmentation Markov chain Monte Carlo (MCMC) framework for Bayesian estimation of stochastic epidemic model parameters, in which measurements are augmented with subject-level disease histories. In our MCMC algorithm, we propose each new subject-level path, conditional on the data, using a time-inhomogenous continuous-time Markov process with rates determined by the infection histories of other individuals. The method is general, and may be applied to a broad class of epidemic models with only minimal modifications to the model dynamics and/or emission distribution. We present our algorithm in the context of multiple stochastic epidemic models in which the data are binomially sampled prevalence counts, and apply our method to data from an outbreak of influenza in a British boarding school. Supplementary material for this article is available online.     

Hidden Markov Models Training by a Particle Swarm Optimization Algorithm

In this work we consider the problem of Hidden Markov Models (HMM) training. This problem can be considered as a global optimization problem and we focus our study on the Particle Swarm Optimization (PSO) algorithm. To take advantage of the search strategy adopted by PSO, we need to modify the HMM's search space. Moreover, we introduce a local search technique from the field of HMMs and that is known as the Baum–Welch algorithm. A parameter study is then presented to evaluate the importance of several parameters of PSO on artificial data and natural data extracted from images.     

Efficient sensitivity analysis in hidden markov models

Sensitivity analysis in hidden Markov models (HMMs) is usually performed by means of a perturbation analysis where a small change is applied to the model parameters, upon which the output of interest is re-computed. Recently it was shown that a simple mathematical function describes the relation between HMM parameters and an output probability of interest; this result was established by representing the HMM as a (dynamic) Bayesian network. To determine this sensitivity function, it was suggested to employ existing Bayesian network algorithms. Up till now, however, no special purpose algorithms for establishing sensitivity functions for HMMs existed. In this paper we discuss the drawbacks of computing HMM sensitivity functions, building only upon existing algorithms. We then present a new and efficient algorithm, which is specially tailored for determining sensitivity functions in HMMs.     

A curved exponential family model for complex networks

Networks are being increasingly used to represent relational data. As the patterns of relations tends to be complex, many probabilistic models have been proposed to capture the structural properties of the process that generated the networks. Two features of network phenomena not captured by the simplest models is the variation in the number of relations individual entities have and the clustering of their relations. In this paper we present a statistical model within the curved exponential family class that can represent both arbitrary degree distributions and an average clustering coefficient. We present two tunable parameterizations of the model and give their interpretation. We also present a Markov Chain Monte Carlo (MCMC) algorithm that can be used to generate networks from this model.     

Twisting the Alive Particle Filter

This work focuses on sampling from hidden Markov models (Cappe et al. 2005) whose observations have intractable density functions. We develop a new sequential Monte Carlo (e.g. Doucet, 2011) algorithm and a new particle marginal Metropolis-Hastings (Andrieu et al J R Statist Soc Ser B 72:269-342, 2010) algorithm for these purposes. We build from Jasra et al (2013) and Whiteley and Lee (Ann Statist 42:115-141, 2014) to construct the sequential Monte Carlo (SMC) algorithm, which we call the alive twisted particle filter. Like the alive particle filter (Amrein and Künsch, 2011, Jasra et al, 2013), our new SMC algorithm adopts an approximate Bayesian computation (Tavare et al. Genetics 145:505-518, 1997) estimate of the HMM. Our alive twisted particle filter also uses a twisted proposal as in Whiteley and Lee (Ann Statist 42:115-141, 2014) to obtain a low-variance estimate of the HMM normalising constant. We demonstrate via numerical examples that, in some scenarios, this estimate has a much lower variance than that of the estimate obtained via the alive particle filter. The low variance of this normalising constant estimate encourages the implementation of our SMC algorithm within a particle marginal Metropolis-Hastings (PMMH) scheme, and we call the resulting methodology “alive twisted PMMH”. We numerically demonstrate, on a stochastic volatility model, how our alive twisted PMMH can converge faster than the standard alive PMMH of Jasra et al (2013).     

Post-processing of Markov chain Monte Carlo output in Bayesian latent variable models with application to multidimensional scaling

In Bayesian analysis of multidimensional scaling model with MCMC algorithm, we encounter the indeterminacy of rotation, reflection and translation of the parameter matrix of interest. This type of indeterminacy may be seen in other multivariate latent variable models as well. In this paper, we propose to address this indeterminacy problem with a novel, offline post-processing method that is easily implemented using easy-to-use Markov chain Monte Carlo (MCMC) software. Specifically, we propose a post-processing method based on the generalized extended Procrustes analysis to address this problem. The proposed method is compared with four existing methods to deal with indeterminacy thorough analyses of artificial as well as real datasets. The proposed method achieved at least as good a performance as the best existing method. The benefit of the offline processing approach in the era of easy-to-use MCMC software is discussed.     

Using Random Quasi-Monte-Carlo Within Particle Filters,With Application to Financial Time Series

This article presents a new particle filter algorithm which uses random quasi-Monte-Carlo to propagate particles. The filter can be used generally, but here it is shown that for one-dimensional state-space models, if the number of particles is N, then the rate of convergence of this algorithm is N^?1. This compares favorably with the N^?1/2 convergence rate of standard particle filters. The computational complexity of the new filter is quadratic in the number of particles, as opposed to the linear computational complexity of standard methods. I demonstrate the new filter on two important financial time series models, an ARCH model and a stochastic volatility model. Simulation studies show that for fixed CPU time, the new filter can be orders of magnitude more accurate than existing particle filters. The new filter is particularly efficient at estimating smooth functions of the states, where empirical rates of convergence are N^?3/2; and for performing smoothing, where both the new and existing filters have the same computational complexity.     

Generalisation of the Viterbi algorithm

The Viterbi algorithm, derived using dynamic programming techniques,is a maximum a posteriori (MAP) decoding method which was developedin the electrical engineering literature to be used in the analysisof hidden Markov models (HMMs). Given a particular HMM, theoriginal algorithm recovers the MAP state sequence underlyingany observation sequence generated from that model. This paperintroduces a generalization of the algorithm to recover, forarbitrary L, the top L most probable state sequences, with specialreference to its use in the area of automatic speech recognition.