Blocked Gibbs Sampler Research Articles

Bayesian Graphical Lasso Models and Efficient Posterior Computation

Recently, the graphical lasso procedure has become popular in estimating Gaussian graphical models. In this paper, we introduce a fully Bayesian treatment of graphical lasso models. We first investigate the graphical lasso prior that has been relatively unexplored. Using data augmentation, we develop a simple but highly efficient block Gibbs sampler for simulating covariance matrices. We then generalize the Bayesian graphical lasso to the Bayesian adaptive graphical lasso. Finally, we illustrate and compare the results from our approach to those obtained using the standard graphical lasso procedures for real and simulated data. In terms of both covariance matrix estimation and graphical structure learning, the Bayesian adaptive graphical lasso appears to be the top overall performer among a range of frequentist and Bayesian methods.

Semiparametric Bayesian analysis of nonlinear reproductive dispersion mixed models for longitudinal data

In the development of nonlinear reproductive dispersion mixed models, it is commonly assumed that distribution of random effects is normal. The normality assumption is likely violated in many practical applications. In this paper, we assume that distribution of random effects is specified by a Dirichlet process prior for relaxing this limitation. A semiparametric Bayesian approach combining the stick-breaking prior and the blocked Gibbs sampler as well as the Metropolis–Hastings algorithm is developed for simulating observations from the posterior distributions and producing the joint Bayesian estimates of unknown parameters and random effects. Two goodness-of-fit statistics are presented to assess the plausibility of the posited model, and the procedures for computing the Bayes factor, pseudo-Bayes factor and deviance information criterion for model comparison are given. Also, we propose two Bayesian case deletion influence measures including the ϕ-divergence and Cook’s posterior mean distance. Four simulation studies and a real example are presented to illustrate the newly developed Bayesian methodologies.

Statistical inference in high dimensional spaces of genomics: an RNA structural example

The emergence of genome scale data sets leads to increasingly more precise parameter estimates that are ideally suited for maximum likelihood methods and other highest scoring procedures, when the number of unknowns is modest. However, paradoxically just the opposite is becoming increasingly common in genomics. This paradox has emerged because these technologies have simultaneously opened opportunities to draw inferences on previously unanswerable high dimensional questions. In this regime the curse of dimensionality not only denies frequentist methods including maximum likelihood estimation of all their asymptotic advantages, but also often makes these estimates at best misleading if not downright wrong. However, ensemble based Bayesian inferences do not suffer from these afflictions, as they recognize that drawing inferences is an inherently uncertain process and employ the laws of probability to address this uncertainty. This talk will briefly introduce the ideas probabilistic statistical inference using the following example of RNA secondary structure prediction. RNA secondary structures play a crucial role in the function of many RNAs, and structural features are often essential to their interaction with other cellular components. But as we show the Boltzmann weighted space of RNA secondary structures can be very complex. Here we present a new algorithm, RGibbs, to identify RNA motifs in longer unaligned sequence, and predict consensus secondary structures for using the blocked Gibbs sampler, which has theoretical advantage in convergence time. This algorithm iteratively samples from the conditional probability distributions P(Structure | Alignment) and P(Alignment | Structure). We illustrate how these probabilistically drawn samples can characterize these potentially complex spaces using hierarchical clustering method to characterize the shape of the posterior space, γ-centroid estimator to generate a prediction from sampled structures, and credibility limits to characterize the uncertainty. An analysis of 17 RNA families shows substantially improved structural prediction based on PPV-SEN curves comparisons, compactness of sampled structures around their ensemble centroids, at least eleven families with well separated clusters. The fact that the distances between the references structures and the centroid structures were large compared to the variation among structures within an ensemble raises questions the aptness of the term maximum expected accuracy estimator.

A semiparametric Bayesian approach to generalized partial linear mixed models for longitudinal data

A generalized partial linear mixed model (GPLMM) is a natural extension of generalized linear mixed models (GLMMs) and partial linear models (PLMs). Almost all existing methods for analyzing GPLMMs are developed on the basis of the assumption that random effects are distributed as a fully parametric distribution such as normal distribution. In this paper, we extend the GPLMMs by specifying a Dirichlet process prior for a general distribution of random effects, and propose a semiparametric Bayesian approach by simultaneously utilizing an approximation truncation Dirichlet process prior of the random effects and a P-spline approximation of the smoothing function. By combining the block Gibbs sampler and the Metropolis–Hastings algorithm, a hybrid algorithm is presented for sampling observations from the posterior distribution. A procedure for selecting the degree of the polynomial components in nonparametric approximation using Bayes factor is given via path sampling. Some goodness-of-fit statistics are proposed to evaluate the plausibility of the posited model. Several simulation studies and a real example are presented to illustrate the proposed methodologies.

Efficient Gaussian graphical model determination under G-Wishart prior distributions

This paper proposes a new algorithm for Bayesian model determination in Gaussian graphical models under G-Wishart prior distributions. We first review recent development in sampling from G-Wishart distributions for given graphs, with a particular interest in the efficiency of the block Gibbs samplers and other competing methods. We generalize the maximum clique block Gibbs samplers to a class of flexible block Gibbs samplers and prove its convergence. This class of block Gibbs samplers substantially outperforms its competitors along a variety of dimensions. We next develop the theory and computational details of a novel Markov chain Monte Carlo sampling scheme for Gaussian graphical model determination. Our method relies on the partial analytic structure of G-Wishart distributions integrated with the exchange algorithm. Unlike existing methods, the new method requires neither proposal tuning nor evaluation of normalizing constants of G-Wishart distributions.

RNAG: a new Gibbs sampler for predicting RNA secondary structure for unaligned sequences

Motivation: RNA secondary structure plays an important role in the function of many RNAs, and structural features are often key to their interaction with other cellular components. Thus, there has been considerable interest in the prediction of secondary structures for RNA families. In this article, we present a new global structural alignment algorithm, RNAG, to predict consensus secondary structures for unaligned sequences. It uses a blocked Gibbs sampling algorithm, which has a theoretical advantage in convergence time. This algorithm iteratively samples from the conditional probability distributions P(Structure | Alignment) and P(Alignment | Structure). Not surprisingly, there is considerable uncertainly in the high-dimensional space of this difficult problem, which has so far received limited attention in this field. We show how the samples drawn from this algorithm can be used to more fully characterize the posterior space and to assess the uncertainty of predictions.Results: Our analysis of three publically available datasets showed a substantial improvement in RNA structure prediction by RNAG over extant prediction methods. Additionally, our analysis of 17 RNA families showed that the RNAG sampled structures were generally compact around their ensemble centroids, and at least 11 families had at least two well-separated clusters of predicted structures. In general, the distance between a reference structure and our predicted structure was large relative to the variation among structures within an ensemble.Availability: The Perl implementation of the RNAG algorithm and the data necessary to reproduce the results described in Sections 3.1 and 3.2 are available at http://ccmbweb.ccv.brown.edu/rnag.htmlContact: charles_lawrence@brown.eduSupplementary information: Supplementary data are available at Bioinformatics online.

Wishart distributions for decomposable covariance graph models

Gaussian covariance graph models encode marginal independence among the components of a multivariate random vector by means of a graph G. These models are distinctly different from the traditional concentration graph models (often also referred to as Gaussian graphical models or covariance selection models) since the zeros in the parameter are now reflected in the covariance matrix Σ, as compared to the concentration matrix Ω = Σ−1. The parameter space of interest for covariance graph models is the cone PG of positive definite matrices with fixed zeros corresponding to the missing edges of G. As in Letac and Massam [Ann. Statist. 35 (2007) 1278–1323], we consider the case where G is decomposable. In this paper, we construct on the cone PG a family of Wishart distributions which serve a similar purpose in the covariance graph setting as those constructed by Letac and Massam [Ann. Statist. 35 (2007) 1278–1323] and Dawid and Lauritzen [Ann. Statist. 21 (1993) 1272–1317] do in the concentration graph setting. We proceed to undertake a rigorous study of these “covariance” Wishart distributions and derive several deep and useful properties of this class. First, they form a rich conjugate family of priors with multiple shape parameters for covariance graph models. Second, we show how to sample from these distributions by using a block Gibbs sampling algorithm and prove convergence of this block Gibbs sampler. Development of this class of distributions enables Bayesian inference, which, in turn, allows for the estimation of Σ, even in the case when the sample size is less than the dimension of the data (i.e., when “n < p”), otherwise not generally possible in the maximum likelihood framework. Third, we prove that when G is a homogeneous graph, our covariance priors correspond to standard conjugate priors for appropriate directed acyclic graph (DAG) models. This correspondence enables closed form expressions for normalizing constants and expected values, and also establishes hyper-Markov properties for our class of priors. We also note that when G is homogeneous, the family IWQG of Letac and Massam [Ann. Statist. 35 (2007) 1278–1323] is a special case of our covariance Wishart distributions. Fourth, and finally, we illustrate the use of our family of conjugate priors on real and simulated data.

Bayesian estimation of semiparametric nonlinear dynamic factor analysis models using the Dirichlet process prior

Parameters in time series and other dynamic models often show complex range restrictions and their distributions may deviate substantially from multivariate normal or other standard parametric distributions. We use the truncated Dirichlet process (DP) as a non-parametric prior for such dynamic parameters in a novel nonlinear Bayesian dynamic factor analysis model. This is equivalent to specifying the prior distribution to be a mixture distribution composed of an unknown number of discrete point masses (or clusters). The stick-breaking prior and the blocked Gibbs sampler are used to enable efficient simulation of posterior samples. Using a series of empirical and simulation examples, we illustrate the flexibility of the proposed approach in approximating distributions of very diverse shapes.

Model Comparison of Bayesian Semiparametric and Parametric Structural Equation Models

Structural equation models have wide applications. One of the most important issues in analyzing structural equation models is model comparison. This article proposes a Bayesian model comparison statistic, namely the L ν-measure for both semiparametric and parametric structural equation models. For illustration purposes, we consider a Bayesian semiparametric approach for estimation and model comparison in the context of structural equation models with fixed covariates. A finite dimensional Dirichlet process is used to model the crucial latent variables, and a blocked Gibbs sampler is implemented for estimation. Empirical performance of the L ν-measure is evaluated through a simulation study. Results obtained indicate that the L ν-measure, which additionally requires very minor computational effort, gives satisfactory performance. Moreover, the methodologies are demonstrated through an example with a real data set on kidney disease. Finally, the application of the L ν-measure to Bayesian semiparametric nonlinear structural equation models is outlined.

Bayesian nonparametric estimation and consistency of mixed multinomial logit choice models

This paper develops nonparametric estimation for discrete choice models based on the mixed multinomial logit (MMNL) model. It has been shown that MMNL models encompass all discrete choice models derived under the assumption of random utility maximization, subject to the identification of an unknown distribution $G$. Noting the mixture model description of the MMNL, we employ a Bayesian nonparametric approach, using nonparametric priors on the unknown mixing distribution $G$, to estimate choice probabilities. We provide an important theoretical support for the use of the proposed methodology by investigating consistency of the posterior distribution for a general nonparametric prior on the mixing distribution. Consistency is defined according to an $L_1$-type distance on the space of choice probabilities and is achieved by extending to a regression model framework a recent approach to strong consistency based on the summability of square roots of prior probabilities. Moving to estimation, slightly different techniques for non-panel and panel data models are discussed. For practical implementation, we describe efficient and relatively easy-to-use blocked Gibbs sampling procedures. These procedures are based on approximations of the random probability measure by classes of finite stick-breaking processes. A simulation study is also performed to investigate the performance of the proposed methods.

Semiparametric Bayes hierarchical models with mean and variance constraints

In parametric hierarchical models, it is standard practice to place mean and variance constraints on the latent variable distributions for the sake of identifiability and interpretability. Because incorporation of such constraints is challenging in semiparametric models that allow latent variable distributions to be unknown, previous methods either constrain the median or avoid constraints. In this article, we propose a centered stick-breaking process (CSBP), which induces mean and variance constraints on an unknown distribution in a hierarchical model. This is accomplished by viewing an unconstrained stick-breaking process as a parameter-expanded version of a CSBP. An efficient blocked Gibbs sampler is developed for approximate posterior computation. The methods are illustrated through a simulated example and an epidemiologic application.

Gibbs sampling for a Bayesian hierarchical general linear model

We consider a Bayesian hierarchical version of the normal theory general linear model which is practically relevant in the sense that it is general enough to have many applications and it is not straightforward to sample directly from the corresponding posterior distribution. Thus we study a block Gibbs sampler that has the posterior as its invariant distribution. In particular, we establish that the Gibbs sampler converges at a geometric rate. This allows us to establish conditions for a central limit theorem for the ergodic averages used to estimate features of the posterior. Geometric ergodicity is also a key requirement for using batch means methods to consistently estimate the variance of the asymptotic normal distribution. Together, our results give practitioners the tools to be as confident in inferences based on the observations from the Gibbs sampler as they would be with inferences based on random samples from the posterior. Our theoretical results are illustrated with an application to data on the cost of health plans issued by health maintenance organizations.

Efficient Markov chain Monte Carlo with incomplete multinomial data

We propose a block Gibbs sampling scheme for incomplete multinomial data. We show that the new approach facilitates maximal blocking, thereby reducing serial dependency and speeding up the convergence of the Gibbs sampler. We compare the efficiency of the new method with the standard, non-block Gibbs sampler via a number of numerical examples.

The local Dirichlet process.

As a generalization of the Dirichlet process (DP) to allow predictor dependence, we propose a local Dirichlet process (lDP). The lDP provides a prior distribution for a collection of random probability measures indexed by predictors. This is accomplished by assigning stick-breaking weights and atoms to random locations in a predictor space. The probability measure at a given predictor value is then formulated using the weights and atoms located in a neighborhood about that predictor value. This construction results in a marginal DP prior for the random measure at any specific predictor value. Dependence is induced through local sharing of random components. Theoretical properties are considered and a blocked Gibbs sampler is proposed for posterior computation in lDP mixture models. The methods are illustrated using simulated examples and an epidemiologic application.

Block Gibbs Sampling for Bayesian Random Effects Models With Improper Priors: Convergence and Regeneration

Bayesian versions of the classical one-way random effects model are widely used to analyze data. If the standard diffuse prior is adopted, there is a simple block Gibbs sampler that can be employed to explore the intractable posterior distribution. In this article, theoretical and methodological results are developed that allow one to use this block Gibbs sampler with the same level of confidence that one would have using classical (iid) Monte Carlo. Indeed, a regenerative simulation method is developed that yields simple, asymptotically valid standard errors for the ergodic averages that are used to estimate intractable posterior expectations. These standard errors can be used to choose an appropriate (Markov chain) Monte Carlo sample size. The regenerative method rests on the assumption that the underlying Markov chain converges to its stationary distribution at a geometric rate. Another contribution of this article is a result showing that, unless the dataset is extremely small and unbalanced, the block Gibbs Markov chain is geometrically ergodic. We illustrate the use of the regenerative method with data from a styrene exposure study. R code for the simulation is posted as an online supplement.

Efficient Markov chain Monte Carlo implementation of Bayesian analysis of additive and dominance genetic variances in noninbred pedigrees.

Accurate and fast computation of quantitative genetic variance parameters is of great importance in both natural and breeding populations. For experimental designs with complex relationship structures it can be important to include both additive and dominance variance components in the statistical model. In this study, we introduce a Bayesian Gibbs sampling approach for estimation of additive and dominance genetic variances in the traditional infinitesimal model. The method can handle general pedigrees without inbreeding. To optimize between computational time and good mixing of the Markov chain Monte Carlo (MCMC) chains, we used a hybrid Gibbs sampler that combines a single site and a blocked Gibbs sampler. The speed of the hybrid sampler and the mixing of the single-site sampler were further improved by the use of pretransformed variables. Two traits (height and trunk diameter) from a previously published diallel progeny test of Scots pine (Pinus sylvestris L.) and two large simulated data sets with different levels of dominance variance were analyzed. We also performed Bayesian model comparison on the basis of the posterior predictive loss approach. Results showed that models with both additive and dominance components had the best fit for both height and diameter and for the simulated data with high dominance. For the simulated data with low dominance, we needed an informative prior to avoid the dominance variance component becoming overestimated. The narrow-sense heritability estimates in the Scots pine data were lower compared to the earlier results, which is not surprising because the level of dominance variance was rather high, especially for diameter. In general, the hybrid sampler was considerably faster than the blocked sampler and displayed better mixing properties than the single-site sampler.

Semiparametric Bayesian analysis of structural equation models with fixed covariates

Latent variables play the most important role in structural equation modeling. In almost all existing structural equation models (SEMs), it is assumed that the distribution of the latent variables is normal. As this assumption is likely to be violated in many biomedical researches, a semiparametric Bayesian approach for relaxing it is developed in this paper. In the context of SEMs with covariates, we provide a general Bayesian framework in which a semiparametric hierarchical modeling with an approximate truncation Dirichlet process prior distribution is specified for the latent variables. The stick-breaking prior and the blocked Gibbs sampler are used for efficient simulation in the posterior analysis. The developed methodology is applied to a study of kidney disease in diabetes patients. A simulation study is conducted to reveal the empirical performance of the proposed approach. Supplementary electronic material for this paper is available in Wiley InterScience at http://www.mrw.interscience.wiley.com/suppmat/1097-0258/suppmat/.

Sufficient burn-in for Gibbs samplers for a hierarchical random effects model

We consider Gibbs and block Gibbs samplers for a Bayesian hierarchical version of the one-way random effects model. Drift and minorization conditions are established for the underlying Markov chains. The drift and minorization are used in conjunction with results from J. S. Rosenthal [J. Amer. Statist. Assoc. 90 (1995) 558-566] and G. O. Roberts and R. L. Tweedie [Stochastic Process. Appl. 80 (1999) 211-229] to construct analytical upper bounds on the distance to stationarity. These lead to upper bounds on the amount of burn-in that is required to get the chain within a prespecified (total variation) distance of the stationary distribution. The results are illustrated with a numerical example.

A blocking Gibbs sampling method to detect major genes with phenotypic data from a diallel mating.

Diallel mating is a frequently used design for estimating the additive and dominance genetic (polygenic) effects involved in quantitative traits observed in the half- and full-sib progenies generated in plant breeding programmes. Gibbs sampling has been used for making statistical inferences for a mixed-inheritance model (MIM) that includes both major genes and polygenes. However, using this approach it has not been possible to incorporate the genetic properties of major genes with the additive and dominance polygenic effects in a diallel mating population. A parent block Gibbs sampling method was developed in this study to make statistical inferences about the major gene and polygenic effects on quantitative traits for progenies derived from a half-diallel mating design. Using simulated data sets with different major and polygenic effects, the proposed method accurately estimated the major and polygenic effects of quantitative traits, and possible genotypes of parents and progenies. The impact of specifying different prior distributions was examined and was found to have little effect on inference on the posterior distribution. This approach was applied to an experimental data set of Loblolly pine (Pinus taeda L.) derived from a 6-parent half-diallel mating. The result indicated that there might be a recessive major gene affecting height growth in this diallel population.

Computational Methods for Multiplicative Intensity Models Using Weighted Gamma Processes

We develop computational procedures for a class of Bayesian nonparametric and semiparametric multiplicative intensity models incorporating kernel mixtures of spatial weighted gamma measures. A key feature of our approach is that explicit expressions for posterior distributions of these models share many common structural features with the posterior distributions of Bayesian hierarchical models using the Dirichlet process. Using this fact, along with an approximation for the weighted gamma process, we show that with some care, one can adapt efficient algorithms used for the Dirichlet process to this setting. We discuss blocked Gibbs sampling procedures and Pólya urn Gibbs samplers. We illustrate our methods with applications to proportional hazard models, Poisson spatial regression models, recurrent events, and panel count data.