Precomputing strategy for Hamiltonian Monte Carlo method based on regularity in parameter space

Markov Chain Monte Carlo (MCMC) algorithms play an important role in statistical inference problems dealing with intractable probability distributions. Recently, many MCMC algorithms such as Hamiltonian Monte Carlo (HMC) and Riemannian Manifold HMC have been proposed to provide distant proposals with high acceptance rate. These algorithms, however, tend to be computationally intensive which could limit their usefulness, especially for big data problems due to repetitive evaluations of functions and statistical quantities that depend on the data. This issue occurs in many statistic computing problems. In this paper, we propose a novel strategy that exploits smoothness (regularity) in parameter space to improve computational efficiency of MCMC algorithms. When evaluation of functions or statistical quantities are needed at a point in parameter space, interpolation from precomputed values or previous computed values is used. More specifically, we focus on HMC algorithms that use geometric information for faster exploration of probability distributions. Our proposed method is based on precomputing the required geometric information on a set of grids before running sampling algorithm and approximating the geometric information for the current location of the sampler using the precomputed information at nearby grids at each iteration of HMC. Sparse grid interpolation method is used for high dimensional problems. Tests on computational examples are shown to illustrate the advantages of our method.     

A Computationally Efficient Projection-Based Approach for Spatial Generalized Linear Mixed Models

Inference for spatial generalized linear mixed models (SGLMMs) for high-dimensional non-Gaussian spatial data is computationally intensive. The computational challenge is due to the high-dimensional random effects and because Markov chain Monte Carlo (MCMC) algorithms for these models tend to be slow mixing. Moreover, spatial confounding inflates the variance of fixed effect (regression coefficient) estimates. Our approach addresses both the computational and confounding issues by replacing the high-dimensional spatial random effects with a reduced-dimensional representation based on random projections. Standard MCMC algorithms mix well and the reduced-dimensional setting speeds up computations per iteration. We show, via simulated examples, that Bayesian inference for this reduced-dimensional approach works well both in terms of inference as well as prediction; our methods also compare favorably to existing “reduced-rank” approaches. We also apply our methods to two real world data examples, one on bird count data and the other classifying rock types. Supplementary material for this article is available online.     

An Artificial Allocations Based Solution to the Label Switching Problem in Bayesian Analysis of Mixtures of Distributions

Label switching is a well-known problem occurring in MCMC outputs in Bayesian mixture modeling. In this article we propose a formal solution to this problem by considering the space of the artificial allocation variables. We show that there exist certain subsets of the allocation space leading to a class of nonsymmetric distributions that have the same support with the symmetric posterior distribution and can reproduce it by simply permuting the labels. Moreover, we select one of these distributions as a solution to the label switching problem using the simple matching distance between the artificial allocation variables. The proposed algorithm can be used in any mixture model and its computational cost depends on the length of the simulated chain but not on the parameter space dimension. Real and simulated data examples are provided in both univariate and multivariate settings. Supplemental material for this article is available online.     

A Marginal Sampler for σ-Stable Poisson–Kingman Mixture Models

We investigate the class of σ-stable Poisson–Kingman random probability measures (RPMs) in the context of Bayesian nonparametric mixture modeling. This is a large class of discrete RPMs, which encompasses most of the popular discrete RPMs used in Bayesian nonparametrics, such as the Dirichlet process, Pitman–Yor process, the normalized inverse Gaussian process, and the normalized generalized Gamma process. We show how certain sampling properties and marginal characterizations of σ-stable Poisson–Kingman RPMs can be usefully exploited for devising a Markov chain Monte Carlo (MCMC) algorithm for performing posterior inference with a Bayesian nonparametric mixture model. Specifically, we introduce a novel and efficient MCMC sampling scheme in an augmented space that has a small number of auxiliary variables per iteration. We apply our sampling scheme to a density estimation and clustering tasks with unidimensional and multidimensional datasets, and compare it against competing MCMC sampling schemes. Supplementary materials for this article are available online.     

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.     

NEWTON-REGULARIZED METHOD BASED ON A BASIC FUNCTION FOR SOLVING THE INVERSE PROBLEM OF CONVECTION-DIFFUSION EQUATIONS

In this paper, a new iteration algorithm to solve the coefficient inverse problem is described by using a "basic function" which is specially defined and the idea of regularization. The method is simple and clear.The main advantage of the algorithm is that its computing cost is less than other current algorithms, such as PST and Purlerbalion Methods. Since it has uniform scheme, on the other hand, the method can be easily exleded to other kinds of inverse problems of different leal equations, multidimensional inverse problem and multiparameler inverse problems, etc.     

Probabilistic abductive logic programming using Dirichlet priors

Probabilistic programming is an area of research that aims to develop general inference algorithms for probabilistic models expressed as probabilistic programs whose execution corresponds to inferring the parameters of those models. In this paper, we introduce a probabilistic programming language (PPL) based on abductive logic programming for performing inference in probabilistic models involving categorical distributions with Dirichlet priors. We encode these models as abductive logic programs enriched with probabilistic definitions and queries, and show how to execute and compile them to boolean formulas. Using the latter, we perform generalized inference using one of two proposed Markov Chain Monte Carlo (MCMC) sampling algorithms: an adaptation of uncollapsed Gibbs sampling from related work and a novel collapsed Gibbs sampling (CGS). We show that CGS converges faster than the uncollapsed version on a latent Dirichlet allocation (LDA) task using synthetic data. On similar data, we compare our PPL with LDA-specific algorithms and other PPLs. We find that all methods, except one, perform similarly and that the more expressive the PPL, the slower it is. We illustrate applications of our PPL on real data in two variants of LDA models (Seed and Cluster LDA), and in the repeated insertion model (RIM). In the latter, our PPL yields similar conclusions to inference with EM for Mallows models.     

Population Monte Carlo

Importance sampling methods can be iterated like MCMC algorithms, while being more robust against dependence and starting values. The population Monte Carlo principle consists of iterated generations of importance samples, with importance functions depending on the previously generated importance samples. The advantage over MCMC algorithms is that the scheme is unbiased at any iteration and can thus be stopped at any time, while iterations improve the performances of the importance function, thus leading to an adaptive importance sampling. We illustrate this method on a mixture example with multiscale importance functions. A second example reanalyzes the ion channel model using an importance sampling scheme based on a hidden Markov representation, and compares population Monte Carlo with a corresponding MCMC algorithm.     

On the Bayesian Mixture Model and Identifiability

This article is concerned with Bayesian mixture models and identifiability issues. There are two sources of unidentifiability: the well-known likelihood invariance under label switching and the perhaps less well-known parameter identifiability problem. When using latent allocation variables determined by the mixture model, these sources of unidentifiability create arbitrary labeling that renders estimation of the model very difficult. We endeavor to tackle these problems by proposing a prior distribution on the allocations, which provides an explicit interpretation for the labeling by removing gaps with high probability. We propose a Markov chain Monte Carlo (MCMC) estimation method and present supporting illustrations.     

Randomized Algorithms with Splitting: Why the Classic Randomized Algorithms Do Not Work and How to Make them Work

We show that the original classic randomized algorithms for approximate counting in NP-hard problems, like for counting the number of satisfiability assignments in a SAT problem, counting the number of feasible colorings in a graph and calculating the permanent, typically fail. They either do not converge at all or are heavily biased (converge to a local extremum). Exceptions are convex counting problems, like estimating the volume of a convex polytope. We also show how their performance could be dramatically improved by combining them with the classic splitting method, which is based on simulating simultaneously multiple Markov chains. We present several algorithms of the combined version, which we simple call the splitting algorithms. We show that the most advance splitting version coincides with the cloning algorithm suggested earlier by the author. As compared to the randomized algorithms, the proposed splitting algorithms require very little warm-up time while running the MCMC from iteration to iteration, since the underlying Markov chains are already in steady-state from the beginning. What required is only fine tuning, i.e. keeping the Markov chains in steady-state while moving from iteration to iteration. We present extensive simulation studies with both the splitting and randomized algorithms for different NP-hard counting problems.     

Partitioned treatment of surface-coupled problems with application to the fluid-porous-media interaction

Employing a decoupled solution strategy for the numerical treatment of the set of governing equations describing a surface-coupled phenomenon is a common practice. In this regard, many partitioned solution algorithms have been developed, which usually either belong to the family of Schur-complement methods or to the group of staggered integration schemes. To select a decoupled solution strategy over another is, however, a case-dependent process that should be done with special care. In particular, the performances of the algorithms from the viewpoints of stability and accuracy of the results on the one hand, and the solution speed on the other hand should be investigated. In this contribution, two strategies for a partitioned treatment of the surface-coupled problem of fluid-porous-media interaction (FPMI) are considered. These are one parallel solution algorithm, which is based on the method of localised Lagrange multipliers (LLM), and one sequential solution method, which follows the block-Gauss-Seidel (BGS) integration strategy. In order to investigate the performances of the proposed schemes, an exemplary initial-boundary-value problem is considered and the numerical results obtained by employing the solution algorithms are compared. (© 2016 Wiley-VCH Verlag GmbH & Co. KGaA, Weinheim)     

The Alive Particle Filter and Its Use in Particle Markov Chain Monte Carlo

In the following article, we investigate a particle filter for approximating Feynman–Kac models with indicator potentials and we use this algorithm within Markov chain Monte Carlo (MCMC) to learn static parameters of the model. Examples of such models include approximate Bayesian computation (ABC) posteriors associated with hidden Markov models (HMMs) or rare-event problems. Such models require the use of advanced particle filter or MCMC algorithms to perform estimation. One of the drawbacks of existing particle filters is that they may “collapse,” in that the algorithm may terminate early, due to the indicator potentials. In this article, using a newly developed special case of the locally adaptive particle filter, we use an algorithm that can deal with this latter problem, while introducing a random cost per-time step. In particular, we show how this algorithm can be used within MCMC, using particle MCMC. It is established that, when not taking into account computational time, when the new MCMC algorithm is applied to a simplified model it has a lower asymptotic variance in comparison to a standard particle MCMC algorithm. Numerical examples are presented for ABC approximations of HMMs.     

The Optimal Control of a Train

We consider the problem of determining an optimal driving strategy in a train control problem with a generalised equation of motion. We assume that the journey must be completed within a given time and seek a strategy that minimises fuel consumption. On the one hand we consider the case where continuous control can be used and on the other hand we consider the case where only discrete control is available. We pay particular attention to a unified development of the two cases. For the continuous control problem we use the Pontryagin principle to find necessary conditions on an optimal strategy and show that these conditions yield key equations that determine the optimal switching points. In the discrete control problem, which is the typical situation with diesel-electric locomotives, we show that for each fixed control sequence the cost of fuel can be minimised by finding the optimal switching times. The corresponding strategies are called strategies of optimal type and in this case we use the Kuhn–Tucker equations to find key equations that determine the optimal switching times. We note that the strategies of optimal type can be used to approximate as closely as we please the optimal strategy obtained using continuous control and we present two new derivations of the key equations. We illustrate our general remarks by reference to a typical train control problem.     

An efficient label setting/correcting shortest path algorithm

We design a new label shortest path algorithm by applying the concept of a pseudo permanent label. This approach allows an algorithm to partition the set of nodes into two new sets: pseudo permanently labeled nodes and its complementary set. From this point of view, this new label method can be considered as a label setting method. Moreover, at least one node becomes permanently labeled when some nodes which belong to the set of pseudo permanently labeled nodes are scanned in each iteration of the algorithm. In the case of networks with non-negative length arcs it is easy to prove that this node has the minimum distance label among the non-pseudo permanently labeled nodes. On the other hand, it is not known during the computation which pseudo permanently labeled nodes are permanently labeled. Therefore, all distance labels are temporary and the algorithm becomes a label correcting method. Nevertheless, the proposed algorithm exhibits some nice features, such as: (1) the time bound for the running of the algorithm for a network with n nodes and m arcs is O(nm); (2) the number of node scan operations in the algorithm is less than the number of these operations in the previous label correcting algorithms as is observed in the computational experience; (3) the algorithm incorporates two new rules which allow easy detection of a negative cycle in the network; (4) the algorithm is quite simple and very easy to implement, and does not require sophisticated data structures; (5) the algorithm exhibits flexibility in the order in which the new pseudo permanently labeled nodes are scanned. The above features are possible through the application of the pseudo permanent label concept.     

Parallel Markov chain Monte Carlo Simulation by Pre-Fetching

In recent years, parallel processing has become widely available to researchers. It can be applied in an obvious way in the context of Monte Carlo simulation, but techniques for “parallelizing” Markov chain Monte Carlo (MCMC) algorithms are not so obvious, apart from the natural approach of generating multiple chains in parallel. Although generation of parallel chains is generally the easiest approach, in cases where burn-in is a serious problem, it is often desirable to use parallelization to speed up generation of a single chain. This article briefly discusses some existing methods for parallelization of MCMC algorithms, and proposes a new “pre-fetching” algorithm to parallelize generation of a single chain.     

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.     

Tuned preconditioners for inexact two‐sided inverse and Rayleigh quotient iteration

Convergence results are provided for inexact two‐sided inverse and Rayleigh quotient iteration, which extend the previously established results to the generalized non‐Hermitian eigenproblem and inexact solves with a decreasing solve tolerance. Moreover, the simultaneous solution of the forward and adjoint problem arising in two‐sided methods is considered, and the successful tuning strategy for preconditioners is extended to two‐sided methods, creating a novel way of preconditioning two‐sided algorithms. Furthermore, it is shown that inexact two‐sided Rayleigh quotient iteration and the inexact two‐sided Jacobi‐Davidson method (without subspace expansion) applied to the generalized preconditioned eigenvalue problem are equivalent when a certain number of steps of a Petrov–Galerkin–Krylov method is used and when this specific tuning strategy is applied. Copyright © 2014 John Wiley & Sons, Ltd.     

Efficient planewise‐like preconditioners for solving 3D problems

We deal with the numerical solution of large linear systems resulting from discretizations of three‐dimensional boundary value problems. It has been shown recently that, if the use of presently available planewise pre‐conditionings is as pathological as thought by many people, except for some trivial anisotropic problems, linewise preconditionings could fairly outperform pointwise methods of approximately the same computational complexity. We propose here a zebra (or line red–black) like numbering strategy of the grid points that leads to a rate of convergence comparable to the one predicted for ideal planewise preconditionings. The keys to the success of this strategy are threefold. On the one hand, one gets rid of the, time and memory consuming, task of computing some accurate approximation to the inverse of each pivot plane matrix. On the other hand, at each PCG iteration, there is no longer a need to solve linear systems whose matrices have the same structure as a two‐dimensional boundary value problem matrix. Finally, it is well suited to parallel computations. Copyright © 1999 John Wiley & Sons, Ltd.     

求解阵列天线自适应滤波问题的一种调比随机逼近算法

提出了求解阵列天线自适应滤波问题的一种调比随机逼近算法.每一步迭代中,算法选取调比的带噪负梯度方向作为新的迭代方向.相比已有的其他随机逼近算法,这个算法不需要调整稳定性常数,在一定程度上解决了稳定性常数选取难的问题.数值仿真实验表明,算法优于已有的滤波算法,且比经典Robbins-Monro （RM）算法具有更好的稳定性.     

All roads lead to Rome—New search methods for the optimal triangulation problem

To perform efficient inference in Bayesian networks by means of a Junction Tree method, the network graph needs to be triangulated. The quality of this triangulation largely determines the efficiency of the subsequent inference, but the triangulation problem is unfortunately NP-hard. It is common for existing methods to use the treewidth criterion for optimality of a triangulation. However, this criterion may lead to a somewhat harder inference problem than the total table size criterion. We therefore investigate new methods for depth-first search and best-first search for finding optimal total table size triangulations. The search methods are made faster by efficient dynamic maintenance of the cliques of a graph. This problem was investigated by Stix, and in this paper we derive a new simple method based on the Bron-Kerbosch algorithm that compares favourably to Stix’ approach. The new approach is generic in the sense that it can be used with other algorithms than just Bron-Kerbosch. The algorithms for finding optimal triangulations are mainly supposed to be off-line methods, but they may form the basis for efficient any-time heuristics. Furthermore, the methods make it possible to evaluate the quality of heuristics precisely and allow us to discover parts of the search space that are most important to direct randomized sampling to.