Importance Link Function Estimation for Markov Chain Monte Carlo Methods

Abstract

This article focuses on improving estimation for Markov chain Monte Carlo simulation. The proposed methodology is based upon the use of importance link functions. With the help of appropriate importance sampling weights, effective estimates of functionals are developed. The method is most easily applied to irreducible Markov chains, where application is typically immediate. An important conceptual point is the applicability of the method to reducible Markov chains through the use of many-to-many importance link functions. Applications discussed include estimation of marginal genotypic probabilities for pedigree data, estimation for models with and without influential observations, and importance sampling for a target distribution with thick tails.     

Explicit error bounds for lazy reversible Markov chain Monte Carlo

We prove explicit, i.e., non-asymptotic, error bounds for Markov Chain Monte Carlo methods, such as the Metropolis algorithm. The problem is to compute the expectation (or integral) of f$f$ with respect to a measure π$\pi$ which can be given by a density ?$?$ with respect to another measure. A straight simulation of the desired distribution by a random number generator is in general not possible. Thus it is reasonable to use Markov chain sampling with a burn-in. We study such an algorithm and extend the analysis of Lovasz and Simonovits [L. Lovász, M. Simonovits, Random walks in a convex body and an improved volume algorithm, Random Structures Algorithms 4 (4) (1993) 359–412] to obtain an explicit error bound.     

Markov Chain Monte Carlo With Mixtures of Mutually Singular Distributions

Markov chain Monte Carlo (MCMC) methods for Bayesian computation are mostly used when the dominating measure is the Lebesgue measure, the counting measure, or a product of these. Many Bayesian problems give rise to distributions that are not dominated by the Lebesgue measure or the counting measure alone. In this article we introduce a simple framework for using MCMC algorithms in Bayesian computation with mixtures of mutually singular distributions. The idea is to find a common dominating measure that allows the use of traditional Metropolis-Hastings algorithms. In particular, using our formulation, the Gibbs sampler can be used whenever the full conditionals are available. We compare our formulation with the reversible jump approach and show that the two are closely related. We give results for three examples, involving testing a normal mean, variable selection in regression, and hypothesis testing for differential gene expression under multiple conditions. This allows us to compare the three methods considered: Metropolis-Hastings with mutually singular distributions, Gibbs sampler with mutually singular distributions, and reversible jump. In our examples, we found the Gibbs sampler to be more precise and to need considerably less computer time than the other methods. In addition, the full conditionals used in the Gibbs sampler can be used to further improve the estimates of the model posterior probabilities via Rao-Blackwellization, at no extra cost.     

Evolutionary Monte Carlo Methods for Clustering

The problem of clustering a group of observations according to some objective function (e.g., K-means clustering, variable selection) or a density (e.g., posterior from a Dirichlet process mixture model prior) can be cast in the framework of Monte Carlo sampling for cluster indicators. We propose a new method called the evolutionary Monte Carlo clustering (EMCC) algorithm, in which three new “crossover moves,” based on swapping and reshuffling sub cluster intersections, are proposed. We apply the EMCC algorithm to several clustering problems including Bernoulli clustering, biological sequence motif clustering, BIC based variable selection, and mixture of normals clustering. We compare EMCC's performance both as a sampler and as a stochastic optimizer with Gibbs sampling, “split-merge” Metropolis–Hastings algorithms, K-means clustering, and the MCLUST algorithm.     

Benchmark Estimation for Markov chain Monte Carlo Samples

While studying various features of the posterior distribution of a vector-valued parameter using an MCMC sample, a subsample is often all that is available for analysis. The goal of benchmark estimation is to use the best available information, that is, the full MCMC sample, to improve future estimates made on the basis of the subsample. We discuss a simple approach to do this and provide a theoretical basis for the method. The methodology and benefits of benchmark estimation are illustrated using a well-known example from the literature. We obtain nearly a 90% reduction in MSE with the technique based on a 1-in-10 subsample and show that greater benefits accrue with the thinner subsamples that are often used in practice.     

Design Ideas for Markov chain Monte Carlo Software

This article discusses design ideas useful in the development of Markov chain Monte Carlo (MCMC) software. Goals of the design are to facilitate analysis of as many statistical models as possible, and to enable users to experiment with different MCMC algorithms as a research tool. These ideas have been used in YADAS, a system written in the Java language, but are also applicable in other object-oriented languages.     

On solving integral equations using Markov chain Monte Carlo methods

In this paper, we propose an original approach to the solution of Fredholm equations of the second kind. We interpret the standard Von Neumann expansion of the solution as an expectation with respect to a probability distribution defined on a union of subspaces of variable dimension. Based on this representation, it is possible to use trans-dimensional Markov chain Monte Carlo (MCMC) methods such as Reversible Jump MCMC to approximate the solution numerically. This can be an attractive alternative to standard Sequential Importance Sampling (SIS) methods routinely used in this context. To motivate our approach, we sketch an application to value function estimation for a Markov decision process. Two computational examples are also provided.     

Bayesian Inference for the One-Factor Copula Model

We develop efficient Bayesian inference for the one-factor copula model with two significant contributions over existing methodologies. First, our approach leads to straightforward inference on dependence parameters and the latent factor; only inference on the former is available under frequentist alternatives. Second, we develop a reversible jump Markov chain Monte Carlo algorithm that averages over models constructed from different bivariate copula building blocks. Our approach accommodates any combination of discrete and continuous margins. Through extensive simulations, we compare the computational and Monte Carlo efficiency of alternative proposed sampling schemes. The preferred algorithm provides reliable inference on parameters, the latent factor, and model space. The potential of the methodology is highlighted in an empirical study of 10 binary measures of socio-economic deprivation collected for 11,463 East Timorese households. The importance of conducting inference on the latent factor is motivated by constructing a poverty index using estimates of the factor. Compared to a linear Gaussian factor model, our model average improves out-of-sample fit. The relationships between the poverty index and observed variables uncovered by our approach are diverse and allow for a richer and more precise understanding of the dependence between overall deprivation and individual measures of well-being.     

Continuous Contour Monte Carlo for Marginal Density Estimation With an Application to a Spatial Statistical Model

The problem of marginal density estimation for a multivariate density function f(x) can be generally stated as a problem of density function estimation for a random vector λ(x) of dimension lower than that of x. In this article, we propose a technique, the so-called continuous Contour Monte Carlo (CCMC) algorithm, for solving this problem. CCMC can be viewed as a continuous version of the contour Monte Carlo (CMC) algorithm recently proposed in the literature. CCMC abandons the use of sample space partitioning and incorporates the techniques of kernel density estimation into its simulations. CCMC is more general than other marginal density estimation algorithms. First, it works for any density functions, even for those having a rugged or unbalanced energy landscape. Second, it works for any transformation λ(x) regardless of the availability of the analytical form of the inverse transformation. In this article, CCMC is applied to estimate the unknown normalizing constant function for a spatial autologistic model, and the estimate is then used in a Bayesian analysis for the spatial autologistic model in place of the true normalizing constant function. Numerical results on the U.S. cancer mortality data indicate that the Bayesian method can produce much more accurate estimates than the MPLE and MCMLE methods for the parameters of the spatial autologistic model.     

On Markov chain Monte Carlo Algorithms for Computing Conditional Expectations Based on Sufficient Statistics

Much work has focused on developing exact tests for the analysis of discrete data using log linear or logistic regression models. A parametric model is tested for a dataset by conditioning on the value of a sufficient statistic and determining the probability of obtaining another dataset as extreme or more extreme relative to the general model, where extremeness is determined by the value of a test statistic such as the chi-square or the log-likelihood ratio. Exact determination of these probabilities can be infeasible for high dimensional problems, and asymptotic approximations to them are often inaccurate when there are small data entries and/or there are many nuisance parameters. In these cases Monte Carlo methods can be used to estimate exact probabilities by randomly generating datasets (tables) that match the sufficient statistic of the original table. However, naive Monte Carlo methods produce tables that are usually far from matching the sufficient statistic. The Markov chain Monte Carlo method used in this work (the regression/attraction approach) uses attraction to concentrate the distribution around the set of tables that match the sufficient statistic, and uses regression to take advantage of information in tables that “almost” match. It is also more general than others in that it does not require the sufficient statistic to be linear, and it can be adapted to problems involving continuous variables. The method is applied to several high dimensional settings including four-way tables with a model of no four-way interaction, and a table of continuous data based on beta distributions. It is powerful enough to deal with the difficult problem of four-way tables and flexible enough to handle continuous data with a nonlinear sufficient statistic.     

A Markov chain Monte Carlo Algorithm for Approximating Exact Conditional Probabilities

Conditional inference eliminates nuisance parameters by conditioning on their sufficient statistics. For contingency tables conditional inference entails enumerating all tables with the same sufficient statistics as the observed data. For moderately sized tables and/or complex models, the computing time to enumerate these tables is often prohibitive. Monte Carlo approximations offer a viable alternative provided it is possible to obtain samples from the correct conditional distribution. This article presents an MCMC extension of the importance sampling algorithm, using a rounded normal candidate to update randomly chosen cells while leaving the remainder of the table fixed. This local approximation can greatly increase the efficiency of the rounded normal candidate. By choosing the number of cells to be updated at random, a balance is struck between dependency in the Markov chain and accuracy of the candidate.     

Q-矩阵第一特征值估计的Monte Carlo方法

本文对给定的可逆马氏链所对应的 Q-矩阵给出了它的第一非零特征值的 Monte Carlo估计方法 .具体做法是通过增加一个状态构造一个新的可逆马氏链 ,然后利用增加状态的击中时分布去估计第一非零特征值 .     

Robustness and applicability of Markov chain Monte Carlo algorithms for eigenvalue problems

In this paper we analyse applicability and robustness of Markov chain Monte Carlo algorithms for eigenvalue problems. We restrict our consideration to real symmetric matrices.

Almost Optimal Monte Carlo (MAO) algorithms for solving eigenvalue problems are formulated. Results for the structure of both – systematic and probability error are presented. It is shown that the values of both errors can be controlled independently by different algorithmic parameters. The results present how the systematic error depends on the matrix spectrum. The analysis of the probability error is presented. It shows that the close (in some sense) the matrix under consideration is to the stochastic matrix the smaller is this error. Sufficient conditions for constructing robust and interpolation Monte Carlo algorithms are obtained. For stochastic matrices an interpolation Monte Carlo algorithm is constructed.

A number of numerical tests for large symmetric dense matrices are performed in order to study experimentally the dependence of the systematic error from the structure of matrix spectrum. We also study how the probability error depends on the balancing of the matrix.     

A Monte Carlo Estimation of the Entropy for Markov Chains

We introduce an estimate of the entropy of the marginal density p ^t of a (eventually inhomogeneous) Markov chain at time t≥1. This estimate is based on a double Monte Carlo integration over simulated i.i.d. copies of the Markov chain, whose transition density kernel is supposed to be known. The technique is extended to compute the external entropy , where the p ₁ ^t s are the successive marginal densities of another Markov process at time t. We prove, under mild conditions, weak consistency and asymptotic normality of both estimators. The strong consistency is also obtained under stronger assumptions. These estimators can be used to study by simulation the convergence of p ^t to its stationary distribution. Potential applications for this work are presented: (1) a diagnostic by simulation of the stability property of a Markovian dynamical system with respect to various initial conditions; (2) a study of the rate in the Central Limit Theorem for i.i.d. random variables. Simulated examples are provided as illustration.      

Adaptive estimation of the transition density of a Markov chain

In this paper a new estimator for the transition density π of an homogeneous Markov chain is considered. We introduce an original contrast derived from regression framework and we use a model selection method to estimate π under mild conditions. The resulting estimate is adaptive with an optimal rate of convergence over a large range of anisotropic Besov spaces . Some simulations are also presented.     

Non-parametric adaptive estimation of the drift for a jump diffusion process

In this article, we consider a jump diffusion process _(Xt)t≥0${\left(X_t\right)}_{t\ge 0}$ observed at discrete times t=0,Δ,…,nΔ$t=0,\Delta ,\dots ,n\Delta$. The sampling interval Δ$\Delta$ tends to 0 and nΔ$n\Delta$ tends to infinity. We assume that _(Xt)t≥0${\left(X_t\right)}_{t\ge 0}$ is ergodic, strictly stationary and exponentially β$\beta$-mixing. We use a penalised least-square approach to compute two adaptive estimators of the drift function b$b$. We provide bounds for the risks of the two estimators.     

A Markov chain Monte Carlo Expectation Maximization Algorithm for Statistical Analysis of DNA Sequence Evolution with Neighbor-Dependent Substitution Rates

The evolution of DNA sequences can be described by discrete state continuous time Markov processes on a phylogenetic tree. We consider neighbor-dependent evolutionary models where the instantaneous rate of substitution at a site depends on the states of the neighboring sites. Neighbor-dependent substitution models are analytically intractable and must be analyzed using either approximate or simulation-based methods. We describe statistical inference of neighbor-dependent models using a Markov chain Monte Carlo expectation maximization (MCMC-EM) algorithm. In the MCMC-EM algorithm, the high-dimensional integrals required in the EM algorithm are estimated using MCMC sampling. The MCMC sampler requires simulation of sample paths from a continuous time Markov process, conditional on the beginning and ending states and the paths of the neighboring sites. An exact path sampling algorithm is developed for this purpose.     

On the adoption of the Monte Carlo method to solve one-dimensional steady state thermal diffusion problems for non-uniform solids

The present paper is focussed on the investigation of the potential adoption of the Monte Carlo method to solve one-dimensional, steady state, thermal diffusion problems for continuous solids characterised by an isotropic, space-dependent conductivity tensor and subjected to non-uniform heat power deposition.     

Asymptotics of an Efficient Monte Carlo Estimation for the Transition Density of Diffusion Processes

Discretized simulation is widely used to approximate the transition density of discretely observed diffusions. A recently proposed importance sampler, namely modified Brownian bridge, has gained much attention for its high efficiency relative to other samplers. It is unclear for this sampler, however, how to balance the trade-off between the number of imputed values and the number of Monte Carlo simulations under a given computing resource. This paper provides an asymptotically efficient allocation of computing resource to the importance sampling approach with a modified Brownian bridge as importance sampler. The optimal trade-off is established by investigating two types of errors: Euler discretization error and Monte Carlo error. The main results are illustrated with two simulated examples.      

Nonparametric estimation of the stationary density and the transition density of a Markov chain

In this paper, we study first the problem of nonparametric estimation of the stationary density f$f$ of a discrete-time Markov chain (_Xi)$\left(X_i\right)$. We consider a collection of projection estimators on finite dimensional linear spaces. We select an estimator among the collection by minimizing a penalized contrast. The same technique enables us to estimate the density g$g$ of (_Xi,_Xi+1)$(X_i,X_{i+1})$ and so to provide an adaptive estimator of the transition density π=g/f$\pi =g/f$. We give bounds in L²$L^2$ norm for these estimators and we show that they are adaptive in the minimax sense over a large class of Besov spaces. Some examples and simulations are also provided.