Sequential Monte Carlo Scheme Research Articles

Accelerating astronomical and cosmological inference with preconditioned Monte Carlo

ABSTRACT We introduce preconditioned Monte Carlo (PMC), a novel Monte Carlo method for Bayesian inference that facilitates efficient sampling of probability distributions with non-trivial geometry. PMC utilizes a Normalizing Flow (NF) in order to decorrelate the parameters of the distribution and then proceeds by sampling from the preconditioned target distribution using an adaptive Sequential Monte Carlo (SMC) scheme. The results produced by PMC include samples from the posterior distribution and an estimate of the model evidence that can be used for parameter inference and model comparison, respectively. The aforementioned framework has been thoroughly tested in a variety of challenging target distributions achieving state-of-the-art sampling performance. In the cases of primordial feature analysis and gravitational wave inference, PMC is approximately 50 and 25 times faster, respectively, than nested sampling (NS). We found that in higher dimensional applications, the acceleration is even greater. Finally, PMC is directly parallelisable, manifesting linear scaling up to thousands of CPUs.

Unscented Dual Quaternion Particle Filter for SE(3) Estimation

We present a novel dual quaternion filter for recursive estimation of rigid body motions. Based on the sequential Monte Carlo scheme, particles are deployed on the manifold of unit dual quaternions. This allows non-parametric modeling of arbitrary distributions underlying on the SE(3) group. The proposal distribution for importance sampling is estimated particle-wise by a novel dual quaternion unscented Kalman filter (DQ-UKF). It is adapted to the manifold geometric structure and drives the prior particles towards high-likelihood regions on the manifold. The resultant unscented dual quaternion particle filter (U-DQPF) incorporates the most recently observed evidence, raising the particle efficiency considerably for nonlinear pose estimation tasks. Compared with ordinary particle filters and other parametric model-based dual quaternion filtering schemes, the proposed U-DQPF shows superior performance in nonlinear SE(3) estimation.

Efficient hbox {SMC}^2 schemes for stochastic kinetic models

Fitting stochastic kinetic models represented by Markov jump processes within the Bayesian paradigm is complicated by the intractability of the observed-data likelihood. There has therefore been considerable attention given to the design of pseudo-marginal Markov chain Monte Carlo algorithms for such models. However, these methods are typically computationally intensive, often require careful tuning and must be restarted from scratch upon receipt of new observations. Sequential Monte Carlo (SMC) methods on the other hand aim to efficiently reuse posterior samples at each time point. Despite their appeal, applying SMC schemes in scenarios with both dynamic states and static parameters is made difficult by the problem of particle degeneracy. A principled approach for overcoming this problem is to move each parameter particle through a Metropolis-Hastings kernel that leaves the target invariant. This rejuvenation step is key to a recently proposed hbox {SMC}^2 algorithm, which can be seen as the pseudo-marginal analogue of an idealised scheme known as iterated batch importance sampling. Computing the parameter weights in hbox {SMC}^2 requires running a particle filter over dynamic states to unbiasedly estimate the intractable observed-data likelihood up to the current time point. In this paper, we propose to use an auxiliary particle filter inside the hbox {SMC}^2 scheme. Our method uses two recently proposed constructs for sampling conditioned jump processes, and we find that the resulting inference schemes typically require fewer state particles than when using a simple bootstrap filter. Using two applications, we compare the performance of the proposed approach with various competing methods, including two global MCMC schemes.

Population annealing simulations of a binary hard-sphere mixture.

Population annealing is a sequential Monte Carlo scheme well suited to simulating equilibrium states of systems with rough free energy landscapes. Here we use population annealing to study a binary mixture of hard spheres. Population annealing is a parallel version of simulated annealing with an extra resampling step that ensures that a population of replicas of the system represents the equilibrium ensemble at every packing fraction in an annealing schedule. The algorithm and its equilibration properties are described, and results are presented for a glass-forming fluid composed of a 50/50 mixture of hard spheres with diameter ratio of 1.4:1. For this system, we obtain precise results for the equation of state in the glassy regime up to packing fractions φ≈0.60 and study deviations from the Boublik-Mansoori-Carnahan-Starling-Leland equation of state. For higher packing fractions, the algorithm falls out of equilibrium and a free volume fit predicts jamming at packing fraction φ≈0.667. We conclude that population annealing is an effective tool for studying equilibrium glassy fluids and the jamming transition.

Fusing Loop and GPS Probe Measurements to Estimate Freeway Density

In an age of ever-increasing penetration of GPS-enabled mobile devices, the potential of real-time "probe" location information for estimating the state of transportation networks is receiving increasing attention. Much work has been done on using probe data to estimate the current speed of vehicle traffic (or equivalently, trip travel time). While travel times are useful to individual drivers, the state variable for a large class of traffic models and control algorithms is vehicle density. Our goal is to use probe data to supplement traditional, fixed-location loop detector data for density estimation. To this end, we derive a method based on Rao-Blackwellized particle filters, a sequential Monte Carlo scheme. We present a simulation where we obtain a 30\% reduction in density mean absolute percentage error from fusing loop and probe data, vs. using loop data alone. We also present results using real data from a 19-mile freeway section in Los Angeles, California, where we obtain a 31\% reduction. In addition, our method's estimate when using only the real-world probe data, and no loop data, outperformed the estimate produced when only loop data were used (an 18\% reduction). These results demonstrate that probe data can be used for traffic density estimation.

A modeling and simulation framework for the reliability/availability assessment of a power transmission grid subject to cascading failures under extreme weather conditions

Electrical power transmission networks of many developed countries are undergoing deep transformations aimed at (i) facing the challenge offered by the stochastically fluctuating power contributions due to the continuously growing connections of renewable power generating units and (ii) decreasing their vulnerability to failures or malicious attacks and improving their resilience, in order to provide more reliable services, thus increasing both safety and profits. In this complex context, one of the major concerns is that related to the potentially catastrophic consequences of cascading failures triggered by rare and difficult to predict extreme weather events. In this work, we originally propose to combine an extreme weather stochastic model of literature to a realistic cascading failure simulator based on a direct current (DC) power flow approximation and a proportional re-dispatch strategy. The description of the dynamics of the network is completed by the introduction of a novel restoration model accounting for the operating conditions that a repair crew may encounter during an extreme weather event. The resulting model is solved by a customized sequential Monte Carlo scheme in order to quantify the impact of extreme weather events on the reliability/availability performances of the power grid. The approach is demonstrated with reference to the test case of the IEEE14 power transmission network.

Cooperative parallel particle filters for online model selection and applications to urban mobility

We design a sequential Monte Carlo scheme for the dual purpose of Bayesian inference and model selection. We consider the application context of urban mobility, where several modalities of transport and different measurement devices can be employed. Therefore, we address the joint problem of online tracking and detection of the current modality. For this purpose, we use interacting parallel particle filters, each one addressing a different model. They cooperate for providing a global estimator of the variable of interest and, at the same time, an approximation of the posterior density of each model given the data. The interaction occurs by a parsimonious distribution of the computational effort, with online adaptation for the number of particles of each filter according to the posterior probability of the corresponding model. The resulting scheme is simple and flexible. We have tested the novel technique in different numerical experiments with artificial and real data, which confirm the robustness of the proposed scheme.

Marginalized Continuous Time Bayesian Networks for Network Reconstruction from Incomplete Observations

Continuous Time Bayesian Networks (CTBNs) provide a powerful means to model complex network dynamics. How- ever, their inference is computationally demanding — especially if one considers incomplete and noisy time-series data. The latter gives rise to a joint state- and parameter estimation problem, which can only be solved numerically. Yet, finding the exact parameterization of the CTBN has often only secondary importance in practical scenarios. We therefore focus on the structure learning problem and present a way to analytically marginalize the Markov chain underlying the CTBN model with respect its parameters. Since the resulting stochastic process is parameter-free, its inference reduces to an optimal filtering problem. We solve the latter using an efficient parallel implementation of a sequential Monte Carlo scheme. Our framework enables CTBN inference to be applied to incomplete noisy time-series data frequently found in molecular biology and other disciplines.

Bayesian inference for Markov jump processes with informative observations.

In this paper we consider the problem of parameter inference for Markov jump process (MJP) representations of stochastic kinetic models. Since transition probabilities are intractable for most processes of interest yet forward simulation is straightforward, Bayesian inference typically proceeds through computationally intensive methods such as (particle) MCMC. Such methods ostensibly require the ability to simulate trajectories from the conditioned jump process. When observations are highly informative, use of the forward simulator is likely to be inefficient and may even preclude an exact (simulation based) analysis. We therefore propose three methods for improving the efficiency of simulating conditioned jump processes. A conditioned hazard is derived based on an approximation to the jump process, and used to generate end-point conditioned trajectories for use inside an importance sampling algorithm. We also adapt a recently proposed sequential Monte Carlo scheme to our problem. Essentially, trajectories are reweighted at a set of intermediate time points, with more weight assigned to trajectories that are consistent with the next observation. We consider two implementations of this approach, based on two continuous approximations of the MJP. We compare these constructs for a simple tractable jump process before using them to perform inference for a Lotka-Volterra system. The best performing construct is used to infer the parameters governing a simple model of motility regulation in Bacillus subtilis.

Delayed acceptance particle MCMC for exact inference in stochastic kinetic models

Recently-proposed particle MCMC methods provide a flexible way of performing Bayesian inference for parameters governing stochastic kinetic models defined as Markov (jump) processes (MJPs). Each iteration of the scheme requires an estimate of the marginal likelihood calculated from the output of a sequential Monte Carlo scheme (also known as a particle filter). Consequently, the method can be extremely computationally intensive. We therefore aim to avoid most instances of the expensive likelihood calculation through use of a fast approximation. We consider two approximations: the chemical Langevin equation diffusion approximation (CLE) and the linear noise approximation (LNA). Either an estimate of the marginal likelihood under the CLE, or the tractable marginal likelihood under the LNA can be used to calculate a first step acceptance probability. Only if a proposal is accepted under the approximation do we then run a sequential Monte Carlo scheme to compute an estimate of the marginal likelihood under the true MJP and construct a second stage acceptance probability that permits exact (simulation based) inference for the MJP. We therefore avoid expensive calculations for proposals that are likely to be rejected. We illustrate the method by considering inference for parameters governing a Lotka-Volterra system, a model of gene expression and a simple epidemic process.

Approximate Bayesian Computation for Smoothing

We consider a method for approximate inference in hidden Markov models (HMMs). The method circumvents the need to evaluate conditional densities of observations given the hidden states. It may be considered an instance of Approximate Bayesian Computation (ABC) and it involves the introduction of auxiliary variables valued in the same space as the observations. The quality of the approximation may be controlled to arbitrary precision through a parameter ε > 0. We provide theoretical results which quantify, in terms of ε, the ABC error in approximation of expectations of additive functionals with respect to the smoothing distributions. Under regularity assumptions, this error is , where n is the number of time steps over which smoothing is performed. For numerical implementation, we adopt the forward-only sequential Monte Carlo (SMC) scheme of [14] and quantify the combined error from the ABC and SMC approximations. This forms some of the first quantitative results for ABC methods which jointly treat the ABC and simulation errors, with a finite number of data and simulated samples.

Conciliating synchronicity with spatial discretization, exclusion, interactions, and detailed balance

The construction of a discrete stochastic system of interacting particles that evolves through a fully synchronous evolution rule while satisfying detailed balance is a highly demanding task. As a consequence, the presence of nontrivial interaction fields can make synchronicity and thermodynamic equilibrium look as two conflicting counterparts. We show that, with the proper prescriptions, the process of migration of particles in a lattice of mutually exclusive nodes can be simulated with a fully synchronous algorithm, which we call parallel Kawasaki dynamics (PKD), that incorporates site exclusion, local interactions, and detailed balance without the need of system partitioning schemes. We show that the underlying pseudo-Hamiltonian (which is derived from the PKD dynamics instead of being assumed a priori as usual in a sequential Monte Carlo scheme) is temperature dependent and causes the resulting equilibrium properties to differ substantially from the conventional hopping model when the system is near critical conditions.

Parameter Estimation Methods for Chaotic Intercellular Networks

We have investigated simulation-based techniques for parameter estimation in chaotic intercellular networks. The proposed methodology combines a synchronization–based framework for parameter estimation in coupled chaotic systems with some state–of–the–art computational inference methods borrowed from the field of computational statistics. The first method is a stochastic optimization algorithm, known as accelerated random search method, and the other two techniques are based on approximate Bayesian computation. The latter is a general methodology for non–parametric inference that can be applied to practically any system of interest. The first method based on approximate Bayesian computation is a Markov Chain Monte Carlo scheme that generates a series of random parameter realizations for which a low synchronization error is guaranteed. We show that accurate parameter estimates can be obtained by averaging over these realizations. The second ABC–based technique is a Sequential Monte Carlo scheme. The algorithm generates a sequence of “populations”, i.e., sets of randomly generated parameter values, where the members of a certain population attain a synchronization error that is lesser than the error attained by members of the previous population. Again, we show that accurate estimates can be obtained by averaging over the parameter values in the last population of the sequence. We have analysed how effective these methods are from a computational perspective. For the numerical simulations we have considered a network that consists of two modified repressilators with identical parameters, coupled by the fast diffusion of the autoinducer across the cell membranes.

OASIS: Online Active Semi-Supervised Learning

We consider a learning setting of importance to large scale machine learning: potentially unlimited data arrives sequentially, but only a small fraction of it is labeled. The learner cannot store the data; it should learn from both labeled and unlabeled data, and it may also request labels for some of the unlabeled items. This setting is frequently encountered in real-world applications and has the characteristics of online, semi-supervised, and active learning. Yet previous learning models fail to consider these characteristics jointly. We present OASIS, a Bayesian model for this learning setting. The main contributions of the model include the novel integration of a semi-supervised likelihood function, a sequential Monte Carlo scheme for efficient online Bayesian updating, and a posterior-reduction criterion for active learning. Encouraging results on both synthetic and real-world optical character recognition data demonstrate the synergy of these characteristics in OASIS.

STOCHASTIC DESIGN AND CONTROL IN RANDOM HETEROGENEOUS MATERIALS

The present paper discusses a sampling framework that enables the optimization of complex systems characterized by high-dimensional uncertainties and design variables. We are especially concerned with problems relating to random heterogeneous materials where uncertainties arise from the stochastic variability of their properties. In particular, we reformulate topology optimization problems to account for the design of truly random composites. In addition, we address the optimal perscription of input loads/excitations in order to achieve a target response by the random material system. The methodological advances proposed in this paper consist of an adaptive Sequential Monte Carlo scheme that economizes the number of runs of the forward solver and allows the analyst to identify several local maxima that provide important information with regards to the robustness of the design. We further propose a principled manner of introducing information from approximate models that can ultimately lead to further reductions in computational cost.

A two-stage procedure for determining the number of trials in the application of a Monte Carlo method for uncertainty evaluation

Adaptive Monte Carlo schemes can be used to determine the number of Monte Carlo trials (the number of evaluations of the measurement model) necessary for the evaluation of uncertainty according to Supplement 1 to the GUM (GUM S1). The goal is to reach a prescribed numerical accuracy of the Monte Carlo results (the estimate, associated standard uncertainty and coverage interval endpoints) for a chosen confidence level. It is shown that simple sequential adaptive Monte Carlo schemes may not perform well in this regard and an alternative method based on a two-stage procedure due to Stein is proposed. The implementation of this two-stage scheme for GUM S1 is described, and its performance and robustness are demonstrated in terms of simulation results.

Sequential Monte Carlo pricing of American-style options under stochastic volatility models

We introduce a new method to price American-style options on underlying investments governed by stochastic volatility (SV) models. The method does not require the volatility process to be observed. Instead, it exploits the fact that the optimal decision functions in the corresponding dynamic programming problem can be expressed as functions of conditional distributions of volatility, given observed data. By constructing statistics summarizing information about these conditional distributions, one can obtain high quality approximate solutions. Although the required conditional distributions are in general intractable, they can be arbitrarily precisely approximated using sequential Monte Carlo schemes. The drawback, as with many Monte Carlo schemes, is potentially heavy computational demand. We present two variants of the algorithm, one closely related to the well-known least-squares Monte Carlo algorithm of Longstaff and Schwartz [The Review of Financial Studies 14 (2001) 113-147], and the other solving the same problem using a "brute force" gridding approach. We estimate an illustrative SV model using Markov chain Monte Carlo (MCMC) methods for three equities. We also demonstrate the use of our algorithm by estimating the posterior distribution of the market price of volatility risk for each of the three equities.

Sequential Monte Carlo methods for complexity-constrained MAP equalization of dispersive MIMO channels

The ability to perform nearly optimal equalization of multiple input multiple output (MIMO) wireless channels using sequential Monte Carlo (SMC) techniques has recently been demonstrated. SMC methods allow to recursively approximate the a posteriori probabilities of the transmitted symbols, as observations are sequentially collected, using samples from adequate probability distributions. Hence, they are a class of online (adaptive) algorithms, suitable to handle the time-varying channels typical of high speed mobile communication applications. The main drawback of the SMC-based MIMO-channel equalizers so far proposed is that their computational complexity grows exponentially with the number of input data streams and the length of the channel impulse response, rendering these methods impractical. In this paper, we introduce novel SMC schemes that overcome this limitation by the adequate design of proposal probability distribution functions that can be sampled with a lesser computational burden, yet provide a close-to-optimal performance in terms of the resulting equalizer bit error rate and channel estimation error. We show that the complexity of the resulting receivers grows polynomially with the number of input data streams and the length of the channel response, and present computer simulation results that illustrate their performance in some typical scenarios.

Adaptation and Change Detection With a Sequential Monte Carlo Scheme

Given the sequential data from an unknown target system with changing parameters, the first part of this paper discusses online algorithms that adapt to smooth as well as abrupt changes. This paper examines four different parameter/ hyperparameter dynamics for online learning and compares their performance within an online Bayesian learning framework. Using the dynamics that performed best in the first part, the second part of this paper attempts to perform change detection in unknown systems in terms of the time dependence of the marginal likelihood. Because of the sequential nature of the algorithms, a sequential Monte Carlo scheme (particle filter) is a natural means for implementation.

Sequential Monte Carlo scheme for Bayesian estimation in the presence of data outliers

Bayesian inference has been used widely in physics, biology, and engineering for a variety of experiment- or observation-based estimation problems. Sequential Monte Carlo simulations are effective for realizing Bayesian estimations when the system and observational processes are nonlinear. In realistic applications, large disturbances in the observation, or outliers, may be present. We develop a theory and practical strategy to suppress the effect of outliers in the experimental observation and provide numerical support.