particle filter (2)

前一篇介绍贝叶斯滤波的数学原理和蒙特卡洛定位的算法。本篇将介绍基于序贯重要性采样。粒子滤波的思想就是采用一个加权粒子分布去近似后验概率分布p(x)

蒙特卡洛积分

定义一连续随机变量X, 其概率密度分布函数为 p(X); 定义Y=f(X)，则随机变量Y的数学期望：

$E(Y) = \int f(x)p(x) dx \approx \frac{1}{N} \sum_{i=1}^{N} f(x_i) p(x_i)$

实际中，概率密度分布p(X)未知，如何保障所采样的点服从p(X)

直接采样

通过对均匀分布采样，实现对任意分布的采样。

任何未知概率密度分布的累积概率函数cdf都映射在[0-1]区间，通过在[0-1]区间的均匀采样，再函数z = cdf(y)求逆，即是符合真实 y的概率密度分布的采样点。

但如果cdf()函数未知或无法求逆，直接采样不可行。

接受-拒绝采样

用一个已知概率分布函数q(X)去采样，然后按照一定的方法拒绝某些样本，达到近似p(X)分布:

p(x_i) <= k p(x_i)

该采样的限制是确定参数k。

重要性采样

在一定的抽样数量基础上，增加准确度。未知p(x), 在已知概率密度分布的q(x)上采样{x_1, x_2, … x_n}后估计f的期望：

$E(Y) = \int f(x)p(x) dx \approx \frac{1}{N} \sum_{i=1}^{N} f(x_i)$

$E(Y) = \int \frac{p(x)}{q(x)} q(x) f(x) dx$

定义新的随机变量：
$Z = \frac{p(x)}{q(x)} f(x)$

关于原随机变量Y在未知概率分布p(x)下的期望，转化为新的随机变量Z在已知概率分布q(x)下的期望。已知概率分布，即知道如何采样。这里 p(x)/q(x) 就是权值。

so the posterior expectations can be computed as:

$Z(x^i) = \frac{p(x^i | z_{1:t})}{q(x^i | z_{1:t})} f(x^i) \\ E[Y] = \frac{1}{N} \sum_{i=1}^{N} Z(x^i)$

as the importance weights can be defined as:

$w^i = \frac{p(x^i | z_{1:t})}{q(x^i | z_{1:t})}$

the problem is we can’t get p(x|z) , but a loosed (unnormalized) importance weights as:

$\tilde{w^i} = \frac{p(x^i) p(z_{1:t} | x^i)}{q(x^i | z_{1:t})}$

then do normalized from it:

$w^i = \frac{\tilde{w^i}}{\sum_{j} \tilde{w^j} }$

so the posterior expectation is approximated as:

$E[Y | z_{1:t}] \approx \sum_{i=1}^{N} w^i f(x^i)$

sequential importance sampling(SIS)

consider the full posterior distribution of states X_{0:k} given measurements y_{1:k} :

$p(x_{0:k} | z_{1:k}) \approx p(z_k | x_{0:k}, z_{1:k-1}) p(x_{0:k} | z_{1:k-1}) \\ = p(z_k | x_k) p(x_k | x_{0:k-1}, z_{1:k-1}) p(x_{0:k-1}| z_{1:k-1}) \\ = p(z_k | x_k) p(x_k | x_{k-1}) p(x_{0:k-1} | z_{1:k-1})$

consider the sequential of q(x):

$q(x_{0:k}| z_{1:k}) = q(x_{0:k-1} | z_{1:k-1}) q(x_k | x_{0:k-1}, z_{1:k})$

then the unnormalized importance weights can be as:

$w^i \approx \frac{p(z_k|x_k)p(x_k|x_{k-1}) p(x_{0:k-1} | z_{1:k-1})}{ q(x_k | x_{0:k-1}, z_{1:k}) q(x_{0:k-1} | z_{1:k-1}) }$

namely:

$w_k^i \approx \frac{p(z_k|x_k)p(x_k|x_{k-1})}{ q(x_k | x_{0:k-1}, z_{1:k}) } w_{k-1}^i$

the problem in SIS is the algorithm is degenerate, that variance of the weights increases at every step, which means the algorithm will converget to single none-zero (w=1) weight and the rest being zero.

so resampling.

sequential importance resampling(SER)

resampling process:

1) interpret each weight as the probability of obtaining the sample index i in the set x^i

2) draw N samples from that discrete distribtuion and replace the old sample set with the new one.

SER process:

1) draw point x_i from the q distribution.
2) calculate the weights in iteration SIS and normalized the weights to sum to unity
3) if the effective number of particles is too low, go to resampling process