Non-Bayesian particle filters

📝 Abstract

Particle filters for data assimilation in nonlinear problems use “particles” (replicas of the underlying system) to generate a sequence of probability density functions (pdfs) through a Bayesian process. This can be expensive because a significant number of particles has to be used to maintain accuracy. We offer here an alternative, in which the relevant pdfs are sampled directly by an iteration. An example is discussed in detail.

💡 Analysis

Particle filters for data assimilation in nonlinear problems use “particles” (replicas of the underlying system) to generate a sequence of probability density functions (pdfs) through a Bayesian process. This can be expensive because a significant number of particles has to be used to maintain accuracy. We offer here an alternative, in which the relevant pdfs are sampled directly by an iteration. An example is discussed in detail.

📄 Content

arXiv:0905.2181v1 [math.NA] 13 May 2009 Non-Bayesian particle filters Alexandre J. Chorin and Xuemin Tu Department of Mathematics, University of California at Berkeley and Lawrence Berkeley National Laboratory Berkeley, CA, 94720 Abstract Particle filters for data assimilation in nonlinear problems use “par- ticles” (replicas of the underlying system) to generate a sequence of probability density functions (pdfs) through a Bayesian process. This can be expensive because a significant number of particles has to be used to maintain accuracy. We offer here an alternative, in which the relevant pdfs are sampled directly by an iteration. An example is discussed in detail. Keywords particle filter, chainless sampling, normalization factor, iter- ation, non-Bayesian 1 Introduction. There are many problems in science in which the state of a system must be identified from an uncertain equation supplemented by a stream of noisy data (see e.g. [1]). A natural model of this situation consists of a stochastic differential equation (SDE): dx = f(x, t) dt + g(x, t) dw, (1) where x = (x1, x2, . . . , xm) is an m-dimensional vector, dw is m-dimensional Brownian motion, f is an m-dimensional vector function, and g is a scalar (i.e., an m by m diagonal matrix of the form gI, where g is a scalar and I is the identity matrix). The Brownian motion encapsulates all the uncertainty in this equation. The initial state x(0) is assumed given and may be random as well. As the experiment unfolds, it is observed, and the values bn of a mea- surement process are recorded at times tn; for simplicity assume tn = nδ, 1 where δ is a fixed time interval and n is an integer. The measurements are related to the evolving state x(t) by bn = h(xn) + GWn, (2) where h is a k-dimensional, generally nonlinear, vector function with k ≤m, G is a diagonal matrix, xn = x(nδ), and Wn is a vector whose components are independent Gaussian variables of mean 0 and variance 1, independent also of the Brownian motion in equation (1). The task is to estimate x on the basis of equation (1) and the observations (2). If the system (1) is linear and the data are Gaussian, the solution can be found via the Kalman-Bucy filter. In the general case, it is natural to try to estimate x as the mean of its evolving probability density. The initial state x is known and so is its probability density; all one has to do is evaluate sequentially the density Pn+1 of xn+1 given the probability density Pn of xn and the data bn+1. This can be done by following “particles” (replicas of the system) whose empirical distribution approximates Pn. In a Bayesian filter (see e.g [2, 3, 4, 5, 6, 7, 8, 9], one uses the pdf Pn and equation (1) to generate a prior density, and then one uses the new data bn+1 to generate a posterior density Pn+1. In addition, one may have to sample backward to take into account the information each measurement provides about the past and avoid having too many identical particles. Evolving particles is typically expensive, and the backward sampling, usually done by Markov chain Monte Carlo (MCMC), can be expensive as well, because the number of particles needed can grow catastrophically (see e.g. [10]). In this paper we offer an alternative to the standard approach, in which Pn+1 is sampled directly without recourse to Bayes’ theorem and backward sampling, if needed, is done by chainless Monte Carlo [11]. Our direct sam- pling is based on a representation of a variable with density Pn+1 by a col- lection of functions of Gaussian variables parametrized by the support of Pn, with parameters found by iteration. The construction is related to chainless sampling as described in [11]. The idea in chainless sampling is to produce a sample of a large set of variables by sequentially sampling a growing se- quence of nested conditionally independent subsets. As observed in [12, 13], chainless sampling for a SDE reduces to interpolatory sampling, as explained below. Our construction will be explained in the following sections through an example where the position of a ship is deduced from the measurements of an azimuth, already used as a test bed in [6, 14, 15]. 2 2 Sampling by interpolation and iteration. First we explain how to sample via interpolation and iteration in a simple example, related to the example and the construction in [12]. Consider the scalar SDE dx = f(x, t)dt + √σ dw; (3) we want to find sample paths x = x(t), 0 ≤t ≤1, subject to the conditions x(0) = 0, x(1) = X. Let N(a, v) denote a Gaussian variable with mean a and variance v. We first discretize equation (3) on a regular mesh t0, t1, . . . , tN, where tn = nδ, δ = 1/N, 0 ≤n ≤N, with xn = x(tn), and, following [12], use a balanced implicit discretization [16, 17]: xn+1 = xn + f(xn, tn)δ + (xn+1 −xn)f ′(xn)δ + W n+1, where f ′(xn, tn) = ∂f ∂xn(xn, tn) and W n+1 is N(0, σ/N). The joint probability density of the variables x1, . . . , xN−1 is Z−1 exp(−PN 0 V i), where Z is the normalization constant and Vi = ((1 −δf

This content is AI-processed based on ArXiv data.