第四章:非参滤波器
非参滤波器是一种流行的替代高斯技术的方案,不用高斯函数等特定形式表示后验概率,而是通过有限个数值来逼近后验概率,其每一个数值都对应着状态空间中的一个区域。 一些非参贝叶斯滤波器依赖于状态空间的分解,每一个分解都对应着状态空间中一个子域的累积概率cumulative probability。 其它一些方法通过对状态空间的随机采样来近似后验概率。大多数方法中,用于近似后验概率的数值点的数量是可以变化的,这一数量决定了对后验概率近似的质量。 在特定的平滑假设下,随着数量趋向于无穷,非参数技术则一定收敛于真实的后验概率。
本章讨论两种在连续空间下使用有限个数值点来近似后验概率的非参数的方法。第一种方法将状态空间分解为很多个有限的区间,并用直方图(histogram)来表示后验概率, 其中每一个直方都对应着一个区域的累积概率。这种方法被认为是对连续概率密度最好的分段定常近似方法。 第二种技术有限多个样本来描述后验概率,称之为粒子滤波particle filter,在机器人学中十分流行。
直方图和粒子滤波这两种非参技术都不需要对后验概率密度做出很强的参数假设。它们尤其适合表示复杂的多峰分布。因此,人们通常都会选择这种方法来处理全局的不确定性, 解决多个不同假设下的复杂的数据关联问题。但是这些技术的表达能力得益于计算复杂度的增加。 幸运的是,本章中描述的这两种非参数技术都可以自适应的调整数值点的数量。当后验概率的复杂度较低时(比如说,只关注一个具有很小不确定度的状态), 它们就只用很少数量的点。对于复杂的后验概率,就多用一些点。
可以根据后验概率自动调整点数量的技术就被称为自适应的adaptive。 如果它们可以根据计算资源来调整点数量则称之为资源自适应resource-adaptive。 资源自适应技术在机器人学中扮演者重要的角色。它使得机器人能够实时地做决策,而不用关心计算资源是否够用。 粒子滤波器的实现通常就是资源自适应算法,根据可用的计算资源在线调整粒子的数量。
4.1 直方图滤波器
直方图滤波器把状态空间拆分为有限个子区域,并用一个概率值表示每个区域的累积后验概率。对于有限的状态空间,这样的滤波器就称为离散贝叶斯滤波器discrete Bayes filter。 在连续的状态空间中,它们就被称为直方图滤波器histogram filter。本章中,先描述离散贝叶斯滤波器接着讨论它们在连续空间中的形式。
4.1.1 离散贝叶斯滤波器
离散贝叶斯滤波器应用于有限状态空间,在有限状态空间中,随机变量\(X_t\)的可能取值是有限的。在第2.4.2节中,我们已经通过估计门开关状态的例子,讨论了离散贝叶斯滤波器。 本书后续章节中讨论的一些建图问题也会引入离散随机变量。比如说,占用栅格地图算法(grid mappingalgorithm)认为环境中的每个区域只有空和被占用两种状态,这样的随机变量是二值的。 因此,有限状态空间在机器人学中扮演着重要的角色。
算法4.1中描述了离散贝叶斯滤波器的伪代码,它用有限项求和的方式替代了算法2.1中的积分运算。 变量\(x_i\)和\(x_k\)是两个不同的状态,可能的状态的数量是有限的。令\(t\)时刻每个状态\(x_k\)的置信度概率为\(p_{k,t}\)。那么该算法的输入包括一个离散的概率分布\(\left\{p_{k,t}\right\}\), 以及最近一次的控制量\(u_t\)和测量值\(z_t\)。第三行中只根据控制量进行预测,在第4行中,根据测量值更新置信度。离散贝叶斯滤波器在很多单信号领域中很流行, 通常称为隐马尔科夫链模型(hidden Markov Model, HMM)。
算法4.1. 离散贝叶斯滤波器. 其中\(x_i, x_k\)分别表示不同的状态。 |
4.1.2 连续状态
本书中将特别关注离散贝叶斯滤波器,它是一种在连续状态空间中用于近似推理的工具。这样的滤波器称为直方图滤波器histogram filters。 图4.1中说明了直方图滤波器是如何表示随机变量以及非线性映射的。图中显示了一个已经直方图化的高斯分布通过一个非线性函数映射。原始的随机变量有10个直方, 而映射之后的概率分布,也有10个直方,但是其中有两个直方近似为0在图中几乎看不到。为了对比,图4.1中还显示了正确的连续分布。
图 4.1. 直方图表示的连续随机变量。右下角子图中灰色阴影部分显示了连续随机变量\(X\)的概率密度,浅灰色的矩形条是对这一概率密度的直方图近似。 随机变量通过右上角所示的函数映射得到新的变量\(Y\),其概率密度函数以及直方图近似如左上角的图像所示。 |
直方图滤波器将连续的状态空间分割为有限多个区域: $$ \begin{equation}\label{f1} \mathrm{dom}(X_t) = X_{1,t} \bigcup X_{2,t} \bigcup \cdots \bigcup X_{K,t} \end{equation} $$ 其中,\(X_t\)就是我们熟悉的随机变量,描述了\(t\)时刻机器人的状态。函数\(\mathrm{dom}(X_t)\)指的是状态空间,覆盖了\(X_t\)所有的可能取值。 每个\(X_{k,t}\)都描述了一个凸区间,这些区间一起构成了对状态空间的一个分割。对于任何\(i ≠ k\),都有\(X_{i,t} \bigcap X_{k,t} = \emptyset\)并且\(\bigcup_k X_{k, t} = \mathrm{dom}(X_t)\)。
对连续状态空间的直接分割就是多维度网格(multi-dimensional grid),每一个\(X_{k,t}\)就是一个网格元素。通过控制网格的粒度,我们可以在计算效率和精度之间做出权衡。 精细的网格划分具有较小的近似误差,但是计算复杂度更高。
离散贝叶斯滤波器为每一个区域\(X_{k,t}\)赋予一个概率值\(p_{k,t}\)。对于每一个区域,离散滤波器不再具有关于置信度的任何额外信息。 因此后验概率就变成了一个分段逼近的概率密度函数,为同一个区域\(X_{k,t}\)中的所有状态\(x_t\)赋予一个相同的概率值,即: $$ \begin{equation}\label{f2} p(x_t) = \frac{p_{k,t}}{|X_{k,t}|} \end{equation} $$ 其中\(|X_{k,t}|\)是区域\(X_{k,t}\)的体积。
如果状态空间本身就是离散的,条件概率\(p(X_{k,t} | u_t, X_{i,t-1})\)和\(p(z_t | X_{k,t})\)就已经有了定义,可以直接应用上述的算法了。 在连续的状态空间中,这两个条件概率通常表示成概率密度函数\(p(x_t | u_t, x_{t-1})\)和\(p(z_t | x_t)\)的形式,它们都是对于各个状态的定义(不是对状态空间中各个区间的定义)。 当每个区间\(X_{k,t}\)都很小而且具有相同的大小,通常用\(X_{k,t}\)代替整个区间来逼近概率密度函数。比如说,我们可能只是简单的用状态的均值\(X_{K,t}\)来表示: $$ \begin{equation}\label{f3} \hat{x}_{k,t} = |X_{k,t}|^{-1} \int_{X_{k,t}}{x_t d x_t} \end{equation} $$ 并进行如下的简单替换: $$ \begin{equation}\label{f4} p(z_t | X_{k,t}) \approx p(z_t | \hat{x}_{k,t}) \end{equation} $$ $$ \begin{equation}\label{f5} p(X_{k, t} | u_t, X_{i,t-1}) \approx \frac{\eta}{|X_{k,t}|} p(\hat{x}_{k,t} | u_t, \hat{x}_{i,t-1}) \end{equation} $$ 这种近似形式就是公式(\(\ref{f2}\))中离散贝叶斯滤波器的一种解释,与EKF中使用的泰勒展开线性近似形式类似。
4.1.3 直方图近似的数学推导
公式(\(\ref{f4}\))是一种理想的近似方法,\(p(z_t | X_{k,t})\)可以通过下面的积分过程得到: $$ \begin{align}\label{f6} p(z_t | X_{k,t}) & = \frac{p(z_t, X_{k,t})}{p(X_{k,t})} \\ & = \frac{\int_{X_{k,t}} p(z_t, x_t)dx_t }{\int_{X_{k,t}} p(x_t)dx_t } \\ & = \frac{\int_{X_{k,t}} p(z_t | x_t)p(x_t)dx_t }{\int_{X_{k,t}} p(x_t) dx_t } \\ & = \frac{\int_{X_{k,t}} p(z_t | x_t)\frac{p_{k,t}}{|X_{k,t}|}dx_t}{\int_{X_{k,t}} \frac{p_{k,t}}{|X_{k,t}|}dx_t } \\ & = \frac{\int_{X_{k,t}} p(z_t | x_t)dx_t}{\int_{X_{k,t}} 1dx_t} \\ & = |X_{k,t}|^{-1} \int_{X_{k,t}} p(z_t | x_t)dx_t \end{align} $$ 这个推导过程完全是根据公式(\(\ref{f2}\))中的piecewise均匀分布模型(uniform distribution model)展开得到的。对于\(x_t \in X_{k,t}\), 如果用\(p(z_t | \hat{x}_{k,t})\)来近似\(p(z_t | x_t)\),则有: $$ \begin{align}\label{f7} p(z_t | X_{k,t}) & \approx |X_{k,t}|^{-1} \int_{X_{k,t}} p(z_t | \hat{x}_{k,t})dx_t \\ & = |X_{k,t}|^{-1} p(z_t | \hat{x}_{k,t}) \int_{X_{k,t}} 1 dx_t \\ & = |X_{k,t}|^{-1} p(z_t | \hat{x}_{k,t}) |X_{k,t}| \\ & = p(z_t | \hat{x}_{k,t}) \end{align} $$ 这正是式(\(\ref{f4}\))中的近似结果。对公式(\(\ref{f5}\))中关于\(p(X_{k,t} | u_t, X_{i, t-1})\)的近似关系推导要稍微复杂一点,因为在条件符号的左右两侧都有分割区域。我们有: $$ \begin{align}\label{f8} p(X_{k,t} | u_t, X_{i, t-1}) & = \frac{p(X_{k,t}, X_{i,t-1} | u_t)}{p(X_{i,t-1} | u_t)} \\ & = \frac{\int_{X_{k,t}}\int_{X_{i,t-1}} p(x_t, x_{t-1} | u_t) dx_t dx_{t-1}}{\int_{X_{i,t-1}} p(x_{t-1} | u_t)dx_{t-1} } \\ & = \frac{\int_{X_{k,t}}\int_{X_{i,t-1}} p(x_t | x_{t-1}, u_t)p(x_{t-1} | u_t) dx_t dx_{t-1}}{\int_{X_{i,t-1}} p(x_{t-1} | u_t)dx_{t-1} } \end{align} $$ 我们做出马尔可夫假设,认为\(x_{t-1}\)与\(u_t\)之间是相对独立的,即\(p(x_{t-1} | u_t) = p(x_{t-1})\): $$ \begin{align}\label{f9} p(X_{k,t} | u_t, X_{i,t-1}) & = \frac{\int_{X_{k,t}} \int_{X_{i,t-1}} p(x_t | u_t, x_{t-1})p(x_{t-1})dx_t dx_{t-1} }{\int_{X_{i,t-1}} p(x_{t-1}) dx_{t-1} } \\ & = \frac{\int_{X_{k,t}} \int_{X_{i,t-1}} p(x_t | u_t, x_{t-1}) \frac{p_{i,t-1}}{|X_{i,t-1}|} dx_t dx_{t-1} }{\int_{X_{i,t-1}} \frac{p_{i,t-1}}{|X_{i,t-1}|} dx_{t-1} } \\ & = \frac{\int_{X_{k,t}} \int_{X_{i,t-1}} p(x_t | u_t, x_{t-1}) dx_t dx_{t-1} }{\int_{X_{i,t-1}} 1 dx_{t-1} } \\ & = |X_{i,t-1}|^{-1} \int_{X_{k,t}} \int_{X_{i, t-1}} p(x_t | u_t, x_{t-1}) dx_t dx_{t-1} \end{align} $$ 如果我们用\(p(\hat{x}_{k,t} | u_t, \hat{x}_{i,t-1})\)近似\(p(x_t | u_t, x_{t-1})\),就可以得到如下的近似过程。为了保证近似后的结果仍然是一个概率值,需要使用\(\eta\)进行归一化。 $$ \begin{align}\label{f10} p(X_{k,t} | u_t, X_{i,t-1}) & \approx \eta |X_{i,t-1}|^{-1} \int_{X_{k,t}} \int_{X_{i, t-1}} p(\hat{x}_{k,t} | u_t, \hat{x}_{i,t-1}) dx_t dx_{t-1} \\ & = \eta |X_{i,t-1}|^{-1} p(\hat{x}_{k,t} | u_t, \hat{x}_{i,t-1}) \int_{X_{k,t}} \int_{X_{i, t-1}} 1 dx_t dx_{t-1} \\ & = \eta |X_{i,t-1}|^{-1} p(\hat{x}_{k,t} | u_t, \hat{x}_{i,t-1}) |X_{k,t}| |X_{i, t-1}| \\ & = \eta |X_{k,t}| p(\hat{x}_{k,t} | u_t, \hat{x}_{i,t-1}) \end{align} $$ 如果所有的区域都有相同的大小(对于所有的\(k\),\(|X_{k,t}|\)都是相同的),我们就可以直接省略掉因子\(|X_{k,t}|\)。 如此得到的离散贝叶斯滤波器就与算法4.1所描述的算法等价。如果按照这里所描述的实现,辅助参数\(\bar{p}_k\)并不构成一个概率分布, 因为它们没有经过归一化。但是,在算法的第四行中进行了归一化,所以输出的参数实际上就是一个合理的概率密度函数值。
4.1.4 分割方法
在机器人学中,连续状态空间的分割方法大体上分为静态和动态两大类。静态的技术总是使用一个固定的分割,与将要近似的后验概率密度函数的形式无关。 动态的技术根据后验概率密度的形式对分割进行调整。静态的技术通常比较容易实现,但是它们会浪费计算资源。
密度树density trees是一种典型的动态分割方法。密度树根据后验概率密度迭代地分割状态空间。这种分割背后地启发式意义在于,分割等级是后验概率密度的一个函数:区域的概率越小, 对其分割的粒度就越粗The less likely a region, the coarser the decomposition.图4.2中对比了静态格点表示和密度树表示两种分割方法。使用相同数量的直方, 密度树方法具有更好的近似效果。像密度树这样的动态分割技术,总是比静态方法的复杂度低几个数量级,但是实现起来相对比较复杂。
图 4.2. 动态分割V.S.静态分割。。上面左边的图中显示了对随机变量\(Y\)的静态直方图,用10个直方来近似表示\(Y\)的分布, 其中有6个直方的概率基本为0,图中无法体现。上面中间的图中是以密度树的形式对相同随机变量的表示。 |
可以通过选择更新技术selective updating获得与动态分割类似的效果。在更新一个用网格表示的后验概率时,选择更新技术只更新一部分网格。 这种方法地一种常见实现方式是,只更新那些后验概率密度超过了用户指定阈值的网格。选择更新技术可以看作是一种混合的分割方法,用一个较细粒度的网格划分状态空间, 并且有一个包含了所有未被选择更新的区域比较大的集合。它也可以看作是一种动态分割技术,因为具体更新哪些网格是根据后验概率的形式在线决定的。 选择更新技术可以降低计算复杂度,它使得网格分割技术可以应用在三维甚至是更高维空间中。
移动机器人中通常需要在空间的矩阵metric表示形式中识别出各种拓扑逻辑topological。虽然没有明确的形式化定义,拓扑逻辑形式的表达通常被认为是一种类似图的表达形式, 图中的一个节点对应着环境中一个重要的地点(或者特征)。对于室内环境,这样的地点可能对应着十字交叉点、T路口、死胡同等等。这样的分割方案依赖于环境的结构。 我们可能想用规则的空间网格来分割状态空间,因为这样的形式不依赖于环境特征的位置和形状。尽管网格形式通常被看作是一个矩阵,严格地说,是在空间中嵌入了一个矩阵,而不是分割。 在移动机器人中,人们更倾向于使用网格形式。比如说在第7章中的一些例子就是使用10平方厘米甚至更小的网格来分割空间。提高精度就会增加计算的复杂度。
4.2 静态二值贝叶斯滤波器
机器人学中的一些特定问题都可以描述成两状态的状态估计问题,这类问题通常都是由一类称为二值贝叶斯滤波器binary Bayes filter来解决的。 如果机器人需要根据传感器序列估计环境中一个固定二值量,就会产生这类问题。比如说,机器人需要检测门是开着的还是关着的,在这个测量上下文中,门的状态并不会因为传感器的测量而改变。 另一种静态二值贝叶斯滤波器的例子就是占用栅格地图occupancy grid maps,我们将在第9章中详细讨论它。
如果状态是静态的,那么置信度就只是测量值的函数: $$ \begin{equation}\label{f11} bel_t(x) = p(x | z_{1:t}, u_{1:t}) = p(x | z_{1:t}) \end{equation} $$ 其中状态有两种可能的取值,标记为\(x\)和\(\neg x\)。特别的,我们有\(bel_t(\neg x) = 1 - bel_t(x)\)。省略了状态符号\(x\)上的时间索引表示状态是不会改变的。
自然,这类二值估计问题可以用算法4.1中描述的离散贝叶斯滤波器搞定。但置信度通常是以对数差异比log odds ratio的形式实现的。 状态\(x\)的差异比odds被定义为正状态与反状态的概率值之比: $$ \begin{equation}\label{f12} \frac{p(x)}{p(\neg x)} = \frac{p(x)}{1 - p(x)} \end{equation} $$ 对数差异比是对上式的对数运算: $$ \begin{equation}\label{f13} l(x) := \log{ \frac{p(x)}{1 - p(x)} } \end{equation} $$ 对数差异值域为\((-\infty, +\infty)\)。这种使用对数差异比的形式更新置信度的贝叶斯滤波器很适合计算,它避免了概率值接近0或1时的截止误差。
算法4.2中描述了这种算法的基本更新过程,本质上就是一个求和运算。任何因为测量导致随机变量概率值的增加或者减少的计算,都可以用对数差异比的形式计算。 这种二值贝叶斯滤波器用逆测量模型\(p(x|z_t)\)inverse measurement model,替代了我们所熟悉的前向模型\(p(z_t | x)\)forward model。 逆测量模型将二值状态变量的分布描述为一个关于测量值\(z_t\)的函数。
算法4.2. 对数差异比形式的二值贝叶斯滤波器。这是一种逆测量模型,其中\(l_t\)是不随时间改变的二值状态的对数差异比。 |
逆模型通常用于测量值比二值状态更复杂的情形。比如说根据相机图片判定门的开关状态,状态十分简单,但是测量空间却十分庞大。 相比于已知门的状态来估计所有图像出现的概率,我们更容易找到一个函数来根据图像判定门是否开着,在这种情况下逆测量模型要比前向模型更容易实现。
根据公式(\(\ref{f13}\))的定义我们可以很容易的通过对数差异比\(l_t\)得到置信度\(bel_t(x)\): $$ \begin{equation}\label{f14} bel_t(x) = 1 - \frac{1}{1 + \exp \{l_t\}} \end{equation} $$ 为了证明二值贝叶斯滤波器算法的正确性,我们重写贝叶斯滤波器公式如下: $$ \begin{align}\label{f15} p(x | z_{1:t}) & = \frac{p(z_t|x, z_{1:t-1}) p(x | z_{1:t-1})}{p(z_t | z_{1:t-1})} \\ & = \frac{p(z_t | x)p(x|z_{1:t-1})}{p(z_t | z_{1:t-1})} \end{align} $$ 将贝叶斯公式应用到测量模型\(p(z_t | x)\): $$ \begin{equation}\label{f16} p(z_t | x) = \frac{p(x|z_t)p(z_t)}{p(x)} \end{equation} $$ 有: $$ \begin{equation}\label{f17} p(x | z_{1:t}) = \frac{p(x|z_t)p(z_t)p(x|z_{1:t-1})}{p(x)p(z_t | z_{1:t-1})} \end{equation} $$ 类似的,对于\(\neg x\)有: $$ \begin{equation}\label{f18} p(\neg x | z_{1:t}) = \frac{p(\neg x|z_t)p(z_t)p(\neg x|z_{1:t-1})}{p(\neg x)p(z_t | z_{1:t-1})} \end{equation} $$ 用公式(\(\ref{f17}\))除以公式(\(\ref{f18}\))得到: $$ \begin{align}\label{f19} \frac{p(x | z_{1:t})}{p(\neg x | z_{1:t})} & = \frac{p(x|z_t)}{p(\neg x | z_t)} \frac{p(x | z_{1:t-1})}{p(\neg x | z_{1:t-1})} \frac{p(\neg x)}{p(x)} \\ & = \frac{p(x|z_t)}{1 - p(x|z_t)} \frac{p(x|z_{1:t-1})}{1 - p(x|z_{1:t-1})} \frac{1 - p(x)}{p(x)} \end{align} $$ 我们用对数差异比\(l_t(x)\)来表示置信度\(bel_t(x)\)。得到在\(t\)时刻的对数差异置信度: $$ \begin{align}\label{f20} l_t(x) & = \log{\frac{p(x|z_t)}{p(\neg x | z_t)}} + \log{\frac{p(x | z_{1:t-1})}{p(\neg x | z_{1:t-1})}} + \log{\frac{p(\neg x)}{p(x)}} \\ & = \log{\frac{p(x|z_t)}{p(\neg x | z_t)}} - \log{\frac{p(x)}{1-p(x)}} + l_{t-1}(x) \end{align} $$ 其中\(p(x)\)是状态\(x\)的先验概率密度。如式(\(\ref{f20}\))所示,每个测量更新都要用到先加上上一个对数差异比。在处理任何传感器数据之前,定义有初始的对数差异形式的置信度: $$ \begin{equation}\label{f21} l_0(x) = \log{\frac{1- p(x)}{p(x)}} \end{equation} $$
4.3 粒子滤波器
4.3.1 基本算法
粒子滤波器是另一种贝叶斯滤波器的非参数实现形式。与直方图滤波器类似,粒子滤波器用有限的样本来近似后验概率密度。但它们生成样本的方式有很大不同, 而且填充状态空间的方式也不一样。粒子滤波器的核心思想就是用一个随机状态样本集合来表达后验概率\(bel(x_t)\)。 图4.3中以一个高斯函数为例描述了粒子采样的思想。 粒子滤波器用一个对后验分布采样后的样本集合来表示概率密度,而不是用参数的形式写出高斯函数的指数形式。 这种表达形式是一种近似,但它是非参数形式的,因此可以比高斯函数这样的参数方法描述更多的分布空间。此外,这种基于样本的表示方式还具有对非线性函数建模的能力。
图 4.3. 粒子滤波器中的"粒子"。右下角的图中描述了对高斯随机变量\(X\)的采样方式,样本经过右上角的非线性函数映射得到变量\(Y\), \(Y\)的概率分布和采样如左上角的图所示。 |
在粒子滤波器中,后验概率分布的样本称为粒子particles,记为: $$ \begin{equation}\label{f22} \mathcal{X}_t := x_t^{[1]}, x_t^{[2]}, \cdots, x_t^{[M]} \end{equation} $$ 其中\(\mathcal{X}_t\)为粒子集合,每一个粒子\(x_t^{[m]}, 1 ≤ m ≤ M\)都是状态在\(t\)时刻的一个具体实例,也就是在\(t\)时刻真实世界可能状态的一个假设。 其中\(M\)表示粒子的数量。在一般的实现中,粒子数\(M\)通常是一个很大的数字。有些实现中\(M\)是一个关于时间\(t\),或者其它关于置信度\(bel(x_t)\)的参数函数。
粒子滤波器背后的思想是用集合\(\mathcal{X}_t\)来近似置信度\(bel(x_t)\)。 理想情况下粒子集合\(\mathcal{X}_t\)中含有状态假设\(x_t\)的似然概率应当与贝叶斯滤波器的后验概率\(bel(x_t)\)成正比关系: $$ \begin{equation}\label{f23} x_t^{[m]} \sim p(x_t | z_{1:t}, u_{1:t}) \end{equation} $$ 如果状态空间中的一个子区域中具有更多的样本,那么真实的状态就更可能位于该区域中。我们在后面还会讨论, 只有当\(M \to \infty\)公式(\(\ref{f23}\))中所描述的粒子滤波算法的这一特性才能够成立。对于有限的\(M\),粒子所描述的分布与真实分布略有不同。 实际上,只要粒子的数量不是很小这种差异就可以忽略不计。
与我们已经讨论的其它贝叶斯滤波器算法类似,粒子滤波器算法也是用上一个周期的结果\(bel(x_{t-1})\)来迭代地计算当前置信度\(bel(x_t)\)。 因为置信度是用粒子集合来表述的,这意味着粒子滤波器根据\(\mathcal{X}_{t-1}\)来更新\(\mathcal{X}_t\)。 粒子滤波器算法的最基本实现如下算法4.3所示。
算法4.3. 粒子滤波器,一种基于重要性采样的贝叶斯滤波器。 |
算法的输入是粒子集合\(\mathcal{X}_{t-1}\)、最近的控制量\(u_t\)和测量值\(z_t\)。首先构建了一个临时的粒子集合\(\bar{\mathcal{X}}\),作为置信度\(\overline{bel}(x_t)\)的一个估计。 在算法中通过如下的步骤处理输入粒子集合\(\mathcal{X}_{t-1}\)中的每一个粒子\(x_{t-1}^{[m]}\):
- 在第4行中,根据粒子\(x_{t-1}^{[m]}\)和控制量\(u_t\)生成\(t\)时刻的一个假设状态\(x_t^{[m]}\)。 得到的样本由\(m\)索引,标识着它是从集合\(\mathcal{X}_{t-1}\)生成的第\(m\)个样本。 这一步涉及了从下一个状态分布\(p(x_t | u_t, x_{t-1})\)采样的过程。经过\(M\)次循环后得到一个粒子集合,就是粒子滤波器算法所表示的\(\overline{bel}(x_t)\)。
- 第5行中,为每一个粒子\(x_t^{[m]}\)计算一个所谓的重要性因子importance factor,标记为\(w_t^{[m]}\)。 重要性因子用于整合测量值\(z_t\)的信息到粒子集合中。因此,重要性是粒子\(x_t^{[m]}\)观测到测量值\(z_t\)的概率, 即\(w_t^{[m]} = p(z_t | x_t^{[m]})\)。如果我们把\(w_t^{[m]}\)解释为粒子的权重weight, 那么这个带有权重的粒子集合就是对贝叶斯滤波器后验概率\(bel(x_t)\)的一个近似。
- 第8到11行才是粒子滤波器算法真正的"技巧"。这些行实现了一个称为重采样resampling, or importance resampling的操作。 该算法从临时的集合\(\bar{\mathcal{X}}_t\)中抽取出\(M\)个粒子构成新的集合,并根据粒子的重要性因子来描述每个粒子的概率。 如此,重采样把一个有\(M\)个粒子的集合转换为另一个具有相同粒子数量的集合,改变了粒子的分布。 在重采样之前,粒子分布是对\(\overline{bel}(x_t)\)的估计;之后,基本上是后验概率\(bel(x_t) = \eta p(z_t | x_t^{[m]})\overline{bel}(x_t)\)的近似。实际上粒子集合中有很多重复的样本, 因为它们是通过简单的替换得到的。更重要的是,那些没有包含在\(\mathcal{X}_t\)的粒子具有很低的权重。
4.3.2 重要性重采样
在推导粒子滤波器之前,需要更详细的解释一下重采样过程。我们所面临的问题是计算一个关于概率密度函数\(f\)的期望,但是我们只能对另一个概率密度函数\(g\)采样。 比如说计算\(x \in A\)的期望,我们可以将之表示为一个关于\(g\)的期望,如下式。其中\(I\)是一个指示函数,如果输入参数为真则输出1,否则输出0。 $$ \begin{align}\label{f25} E_f\left[I(x \in A)\right] & = \int f(x)I(x\in A) dx \\ & = \int \frac{f(x)}{g(x)} g(x) I(x \in A) dx \\ & = E_g\left[ w(x) I(x \in A) \right] \end{align} $$ 其中\(w(x) = \frac{f(x)}{g(x)}\)是一个权重,描述了\(f\)与\(g\)之间的偏差。为了保证比例因子有意义,需要条件\(f(x) > 0 \rightarrow g(x) > 0\)。
下图4.4中说明了重采样操作背后的思想。其中(a)图中,显示了一个概率分布的密度函数\(f\),称之为目标分布target distribution。 虽然我们想要得到的是对\(f\)的样本,却不能够直接对\(f\)进行采样。退而求其次,我们可以对一个与之相关的概率密度\(g\)进行采样,如(b)所示。 我们称\(g\)所描述的概率分布为建议分布proposal distribution。当\(f(x) > 0\)时,必须有\(g(x) > 0\), 这样通过对\(g\)采样,就不会生成\(f\)为零概率的样本了。(b)中生成的粒子集是对\(g\)的采样而不是\(f\)的。 In particular, for any interval \(A \subseteq \mathrm{dom}(X)\) (or more generally, any Borel set A) the empirical count of particles that fall into \(A\) converges to the integral of \(g\) under \(A\): $$ \begin{equation}\label{f26} \frac{1}{M}\sum_{m=1}^M {I\left(x^{[m]} \in A\right)} \to \int_A{g(x)dx} \end{equation} $$
(a) | (b) | (c) |
为了修正\(f\)与\(g\)之间的偏差,粒子\(x^{[m]}\)被赋予了权重: $$ \begin{equation}\label{f27} w^{[m]} = \frac{f(x^{[m]})}{g(x^{[m]})} \end{equation} $$ 如图4.4(c)所示,图中的蓝色竖条代表重要性权重的值。重要性权重是每个粒子没有归一化的概率密度。特别的,我们有: $$ \begin{equation}\label{f28} \left[\sum_{m=1}^M {w^{[m]}}\right]^{-1} \sum_{m=1}^M {I\left(x^{[m]} \in A\right)} \to \int_A{f(x)dx} \end{equation} $$ 式中的第一项是对重要性因子的归一化。也就是说,虽然我们是从概率密度\(g\)生成的粒子,但是加权的粒子却收敛于概率密度\(f\)。 在一些宽松的条件下,这种近似方式总能够使得任意集合\(A\)收敛于\(E_f[I(x \in A)]\)。在大多数情况下,都是以\(\frac{1}{\sqrt{M}}\)的速率收敛。 其中\(M\)为粒子数量。常数因子取决于\(g\)与\(f\)的相似程度。在粒子滤波器中,概率密度\(f\)对应着目标置信度\(bel(x_t)\)。 假设\(\mathcal{X}_{t-1}\)是关于分布\(bel(x_{t-1})\)的粒子集,那么密度\(g\)就对应着如下的分布的乘积: $$ \begin{equation}\label{f29} p(x_t | u_t, x_{t-1})bel(x_{t-1}) \end{equation} $$ 这个分布称为建议分布proposal distribution。
4.3.3 PF的数学推导
为了方便的从数学上推导粒子滤波器,我们场将粒子看作是状态序列的样本
$$ \begin{equation}\label{f30} x_{0:t}^{[m]} = x_0^{[m]}, x_1^{[m]}, \cdots, x_t^{[m]} \end{equation} $$相应的,我们对算法做一些简单的修改,将粒子\(x_t^{[m]}\)扩展为生成\(x_{0:t-1}^{[m]}\)的状态样本序列。这个粒子滤波器在所有的状态序列上计算后验概率, 来替代置信度\(bel(x_t) = p(x_t | u_{1:t}, z_{1:t})\)。
$$ \begin{equation}\label{f31} bel(x_{0:t}) = p(x_{0:t} | u_{1:t}, z_{1:t}) \end{equation} $$诚然,所有状态序列的空间是很大的,构建能够覆盖整个空间的粒子集合是不现实的。但是这里它并不能给我们带来什么困扰, 因为这个定义只用于推导算法4.3所描述的粒子滤波器。
与2.4.3节中推导\(bel(x_t)\)类似的方式可以得到后验概率\(bel(x_{0:t})\)。具体的,我们有
$$ \begin{align}\label{f32} p(x_{0:t} | z_{1:t}, u_{1:t}) & = \eta p(z_t | x_{0:t}, z_{1:t-1}, u_{1:t}) p(x_{0:t} | z_{1:t-1}, u_{1:t}) \\ & = \eta p(z_t | x_t)p(x_{0:t} | z_{1:t-1}, u_{1:t}) \\ & = \eta p(z_t | x_t)p(x_t | x_{0:t-1}, z_{1:t-1}, u_{1:t})p(x_{0:t-1} | z_{1:t-1}, u_{1:t}) \\ & = \eta p(z_t | x_t)p(x_t | x_{t-1}, u_t)p(x_{0:t-1} | z_{1:t-1}, u_{1:t-1}) \end{align} $$注意这里的推导过程中没有积分符号。这是在后验概率中维护所有状态的结果, 而不是2.4.3节那样计算最新的概率。
现在使用归纳法进行推导。初始条件并不重要,这里假设第一个粒子集合是通过对先验概率\(p(x_0)\)进行采样得到的。我们还假设\(t-1\)时刻的粒子集合按照\(bel(x_{0:t-1})\)分布。 对于第\(m\)个粒子\(x_{0:t-1}^{[m]}\),在算法中的第4步中根据建议分布生成样本\(x_t^{[m]}\)。建议分布形式如下:
$$ \begin{equation}\label{f33} p(x_t | x_{t-1}, u_t)bel(x_{0:t-1}) = p(x_t | x_{t-1}, u_t)p(x_{0:t-1} | z_{1:t-1}, u_{1:t-1}) \end{equation} $$其中
$$ \begin{align}\label{f34} w_t^{[m]} & = \frac{\mathrm{target\ distribution}}{\mathrm{proposal\ distribution}} \\ & = \frac{\eta p(z_t | x_t)p(x_t | x_{t-1}, u_t)p(x_{0:t-1} | z_{1:t-1}, u_{1:t-1})}{p(x_t | x_{t-1}, u_t)p(x_{0:t-1}|z_{0:t-1}, u_{0:t-1})} \\ & = \eta p(z_t | x_t) \end{align} $$常数\(\eta\)没有任何作用,因为进行重采样的机率与重要性权重成正比。以正比于\(w_t^{[m]}\)的概率进行重采样后的粒子,其分布服从建议分布与权重\(w_t^{[m]}\)的乘积:
$$ \begin{equation}\label{f35} \eta w_t^{[m]} p(x_t | x_{t-1}, u_t)p(x_{0:t-1} | z_{0:t-1}, u_{0:t-1}) = bel(x_{0:t}) \end{equation} $$注意这里的常数因子\(\eta\)与式\(\ref{f34}\)不同。算法4.4就满足一个简单的现象,如果\(x_{0:t}^{[m]}\)服从\(bel(x_{0:t})\), 那么状态样本\(x_t^{[m]}\)一般都服从分布\(bel(x_t)\)。
后面我们会讨论,这个这个推导过程只有在\(M \to \infty\)的时候成立,due to a laxness in our consideration of the normalization constants。但即使是对于有限的\(M\), 这一推导过程也解释了粒子滤波器背后的含义。
4.3.4 粒子滤波器的一些特性
提取密度Density Extraction
粒子滤波器的样本集合是对连续置信度的离散近似。但是在很多应用中,都需要估计连续空间,也就是说,除了要能够估计粒子所表示的状态外,还需要估计状态空间中其他点。 这种通过样本点来提取连续概率密度的问题称为密度估计density estimation。
图4.5中列举了几种不同的通过粒子集合抽象密度函数的方法。图4.5(a)中显示了通过图4.3中的非线性函数对一个高斯分布映射之后的粒子和密度函数。 图4.5(b)中虚线显示的高斯函数是一种简单高效的从粒子集中提取密度函数的方法。在这种方法中,提取出来的高斯函数基本上与图中实线所示的真实分布的高斯近似一致。
(a) | (b) | (c) | (d) |
图 4.5. 从粒子中提取密度函数的不同方法。(a)密度和样本集近似, (b)高斯近似(均值和方差), (c)直方图近似, (d)核密度估计。近似方法的选择依赖于应用的需要和计算资源。 |
显然,高斯近似只把握了目标分布的基本属性,而且只适用于单峰分布的情况。多峰样本分布则需要,类似k均值聚类,这样的更复杂的技术,使用混合高斯函数来近似目标分布。 图4.5(c)中显示了另一种方法,使用离散直方图来描述状态空间,图中的每一个直方都是对应区域中样本的权重之和。这种直方图滤波器的一个重要缺点就是, 随着维度的增加空间复杂度呈指数增长。另一方面,直方图可以表示多模态分布,它们可以非常高效的计算,可以及时的得到任何状态的密度,而与粒子的数量无关。
从粒子集中生成如4.1.4节中讨论的密度树,并用之生成直方图可以极大的降低空间复杂度。但是在提取状态空间中任意点的密度使,密度树方法需要付出更大的查找代价, 复杂度为树深度的对数级别。
核密度估计是另一种将粒子集合转换为连续密度的方法。每一个用作中心的粒子都成为核Kernel,总体的密度通过混合核密度来得到。 图4.5(d)中显示了在每一个粒子上放置一个高斯核的混合密度。核密度估计的优势在于它的平滑性核计算的简洁性。 到那时计算任意点的密度是粒子(或者核)数量的线性函数。
那么在实际应该选择哪种密度提取方法呢?这依赖于问题本身。比如说,在很多机器人应用中,计算资源是十分有限的,而粒子的均值提供了足够的控制机器人的信息。 一些应用中,比如说主动定位,还依赖于关于状态空间不确定度这样更为复杂的信息。这样的场景中,直方图或者混合高斯方法是更好的选择。融合多个机器人采集的数据时, 常常需要多个不同样本集合下的密度函数的乘积,密度树或者核密度估计更适合解决这类问题。
样本方差Sampling Variance
粒子滤波器的一种重要的误差源与随机采样的方差有关。从一个概率密度分布中采样有限数量个样本时,从这些样本中提取的统计量与原始分布的统计量略有不同。 比如说,如果我们对一个高斯随机变量采样,那么样本的均值和方差就会与原始随机变量的均值和方差不同。这种由于随机采样而带来的差异性,称之为采样器的方差variance。
设想有两个完全相同的机器人,具有完全相同的高斯置信度,并执行着完全相同的无噪声任务。显然,在执行了任务之后,这两个机器人应当具有相同的置信度。 为了模拟这一情况,我们重复地对一个高斯密度分布采样并通过一个非线性函数映射。图4.6中显示了样本,样本的核密度估计,以及真实的置信度。 第一行图像为对一个高斯函数进行25次采样后的结果,与预期的输出对比,可以看出核密度估计与实际的密度函数有细微的差异,并且不同的核密度估计之间存在着较大的差异。 显然,更多的样本意味着更准确的近似效果,更少的差异。实际上,如果样本足够多,机器人就可以获得一个与真实置信度足够接近的基于样本的近似表示。
(a) | (b) | (c) |
(d) | (e) | (f) |
图 4.6. 随机采样的方差。各图都是对一个高斯函数采样并将之通过一个非线性函数映射。 第一行中的图像为对高斯函数进行25次采样后的样本和核估计的结果。第二行中的图像则是250次采样的结果。 |
重采样Resamplint
样本方差是通过重复的重采样放大的。为了方便理解,我们来考虑一个极端情况,在这种场景下机器人的状态没有改变。 我们知道\(x_t = x_{t-1}\)的情况是存在的。在移动机器人定位问题中,没有移动的机器人就是一个很好的例子。更进一步的, 我们假设机器人没有任何传感器,因此它不能估计状态,也就对周围环境一无所知。因此\(t\)时刻的估计与初始位置估计,甚至是任何时刻的位置估计都是一样的。
很不幸,一般的粒子滤波器都不能得到这样的结果。在初始的时候,我们通过先验概率生成粒子集合,该集合将遍布整个状态空间。 但算法4.3的第8到11行的重采样操作,重新生成状态样本\(x^{[m]}\)的时候,偶尔也会失效。因为状态转移是确定性的,第四行中的前向采样过程不能引入任何新的状态。 随着时间的推移,越来越多的例子被重采样操作移除,但不会创建新的粒子。这样一定会出现一种不希望的状态,最后存活下来的\(M\)个粒子都是一样的,因为多次的重采样导致样本的多样性消失了。 对于一个旁观者而言,他可能认为机器人已经完全确认了世界的状态,而实际上机器人没有任何传感器。
这个例子还说明了粒子滤波器的另一个重要限制,就是重采样过程将丢失粒子群体的多样性,这本身就是近似误差的一种表现。 即使粒子集合本身的多样性variance降低了,用于估计真实置信度的粒子集合的多样性却增加了。粒子滤波器的各种实现的本质就是在控制这个variance,或者说是误差。
降低variance主要有两种策略。首先,我们可以降低重采样的频率。当发现状态稳定了(\(x_t = x_{t-1}\))就不要再进行重采样了。 比如说,对于移动机器人定位问题,机器人停下来的时候,就应该暂停重采样(实际上同时暂停整合测量值也是一个很好的想法)。 即使状态发生了改变,降低重采样的频率通常也是一个很好的想法。如上所述,通过重要性因子的乘法运算总是可以整合多个传感器的数据。更一般的, 它在内存中维护了重要性权重,并通过下式更新它们: $$ \begin{equation}\label{f36} w_t^{[m]} = \begin{cases} 1 & \mathrm{if\ resample\ took\ place} \\ p(z_t | x_t^{[m]})w_{t-1}^{[m]} & \mathrm{if\ no\ resampling\ took\ place} \end{cases}\end{equation} $$ 重采样时机的选择,需要实际的经验: 频繁的重采样会增加丢失多样性的风险。如果采样过于频繁,很多样本将可能落在低概率区间中,这是一种浪费。 一个确定是否进行重采样的标准方法是测量重要性权重的方差,权重的方差关系着样本的表达效率。如果所有的权重都是一样的,那么多样性variance为0,就不应当再进行重采样了。 另一方面,如果只在一小部分的样本上的权重比较大,那么variance就会比较高,应该进行重采样。
第二种策略就是所谓的低方差采样low variance sampling,算法4.4中列出了这种方法的一种实现。 这种方法的基本思想是用一系列的随机过程来替代原来重采样过程中独立的采样方法(参考算法4.3中所描述的基本滤波器)。
算法4.4. 粒子滤波器的低方差采样器。在这一方法中,用一个随机数对带有权重\(\mathcal{W}\)的粒子集合\(\mathcal{X}\)进行采样。 重采样后的粒子的概率仍然正比于其权重。此外这种采样器的时间复杂度为\(O(M)\),其中\(M\)为粒子的数量。 |
与其选择\(M\)个随机数在根据这些随机数选择粒子,这种算法只生成了一个随机数并根据这个随机数来进行采样,但是仍然为每一个粒子设定了一个与概率值成正比关系的权重。 这是通过在区间\([0;M^{-1}]\)中随机选择一个数\(r\)获得的,\(M\)是时刻\(t\)的样本数量。接着算法4.4在每个循环中都在\(r\)上累加一个常数\(M^{-1}\), 然后根据这个值选择粒子。区间\([0:1]\)中的任意数值\(u\)都指着一个粒子,称之为粒子\(i\): $$ \begin{equation}\label{f37} i = \mathrm{argmin}_j {\sum_{m=1}^j {w_t^{[m]}}} ≥ u \end{equation} $$ 算法4.4的整个循环都在做两件事情,计算上式的右边并选择第一个权重和超过\(u\)的粒子索引\(i\)。算法中的第12行做出了这一选择, 下图4.7中示例了这个过程。
低方差采样的优势有三点:首先,它以一种比独立随机采样更系统的方式覆盖了样本空间。很明显,非独立采样器系统的便利了所有的粒子,而不是独立的随机选择这些样本。 其次,如果所有的样本都有相同的重要性权重,那么采样后的集合\(\bar{\mathcal{X}}_t\)与\(\mathcal{X}_t\)等价, 因此,即使我们不把测量值整合进\(\mathcal{X}_t\)中,进行重采样也不会丢失样本。最后,低方差采样器的复杂度为\(O(M)\),而独立采样很难有一样的复杂度。 显而易见,选择了随机数之后,还需要一个\(O(\log M)\)复杂度的搜索算法来生成每个粒子,因此整个重采样过程的复杂度为\(O(M\log M)\)。 计算时间是限制粒子滤波器应用的主要因素,通常一个高效的重采样过程就能带来很大的性能提升。因此,机器人学中粒子滤波器的大部分实现都倾向于选择刚刚讨论的这种机制。
关于采样效率的文献有很多。另一种流行的选择是分层采样stratified sampling,在这种算法中粒子们被分为若干个子集。 采样过程被分为两个阶段,首先根据根据各个子集中粒子的权重和确定从各子集中抽取的样本数量,然后用低方差采样等方法随机的从各个子集中抽取样本。 这样的方法具有更低的采样方差,当机器人需要跟踪多个不同的假设时这种方法具有更好的效果。
采样偏差Sampling Bias
只使用有限多个粒子也会导致后验概率估计的系统偏差。设想一种\(M = 1\)的极端情况,算法4.3中的第3到7行的循环就只会执行一次, 并且\(\bar{\mathcal{X}}_t\)中只有一个从运动模型中采样得到的粒子。关键实在第8到11行中的重采样操作,无论重要性权重\(w_t^{[m]}\)为何值,它都一定会接受这个样本。 因此测量模型\(p(z_t | x_t^{[m]})\)在更新过程中没有任何作用。所以,如果\(M = 1\)粒子滤波器就是从概率\(p(x_t | u_{1:t})\)生成的粒子, 而不是期望的后验概率\(p(x_t | u_{1:t}, z_{1:t})\)。它完全忽略了所有的测量值,这是为什么呢?
这都是重采样过程中隐含的归一化惹的祸。在算法的第9行中,根据重要性权重进行采样的时候,如果\(M = 1\),那么\(w_t^{[m]}\)就是它自己的归一化分母: $$ \begin{equation}\label{f38} p(\mathrm{draw}\ x_t^{[m]}\ \mathrm{in}\ \mathrm{Line 9}) = \frac{w_t^{[m]}}{w_t^{[m]}} = 1 \end{equation} $$ 问题出在,没有归一化的值\(w_t^{[m]}\)是从一个\(M\)维的空间中抽取的,但是归一化之后它们就限定在了一个\(M-1\)维的空间中。 幸运的是,对于很大的\(M\)值,降维的影响就越来越小了。
粒子匮乏partical deprivation
即使粒子的数量很多,正确状态附近都可能没有粒子。这一问题就是所谓的粒子匮乏问题particle deprivation problem。当粒子的数量很少, 不足以覆盖所有高似然度关联的区域时,这一问题就很容易出现。有人会说,这个现象在任何粒子滤波器上都会发生,与粒子集合的大小\(M\)无关。
粒子匮乏时随机采样的方差造成的,一个不那么幸运的随机数序列可能使得所有的粒子都不在真值附近。在每一次重采样过程中,这一现象发生的概率大于0(尽管相比于\(M\)而言这个值通常很小)。 因此,我们只要运行粒子滤波器足够长的时间,就会获得一个完全不准确的估计。
实际上,这个问题只有在粒子的数量\(M\)相比于空间中具有高似然度的状态少很多时,才会大几率的暴露出来。解决粒子匮乏问题的一种方法就是, 不管实际的执行和测量指令是什么,都要在每次重采样之后,再随机生成少量的一些粒子添加到集合中。这种方法可以削弱(但是不能解决)粒子匮乏问题, 所付出的代价就是不准确的后验概率估计。添加随机样本的优点就是简单,软件的改动很小。按照经验来说,添加随机样本应当考虑到最终的诉求, 只有在其它的修正粒子匮乏问题的技术都失效的时候才使用。我们将在第8章中介绍机器人定位的方法中,介绍其他处理粒子匮乏问题的方法。
以上的这些讨论说明,基于样本的表示方法的质量随着样本的增加而提高。这就带来了一个问题,在实际估计过程中究竟应该用多少样本。 不幸的是,这个问题没有一个完美的答案,使用者只能自己瞎猜。根据经验,样本的数量严重依赖于状态空间的维度,以及需要近似的概率分布的不确定性。 比如说,均匀分布比起集中在状态空间中一个小范围的分布而言,就需要更多的样本。我们将在本书后续章节中介绍机器人定位和建图的上下文中,更详细的讨论样本数量的问题。
4.4 总结
本章介绍了直方图滤波器和粒子滤波器两种非参数的贝叶斯滤波器。非参数滤波器用有限数量的数值近似后验概率密度。在对系统模型的温和假设下, 随着数值数量的增加趋向于无穷,这两种方法的近似误差都将收敛到零。
- 直方图滤波器把状态空间分解到有限多个凸区域。每个区域都用一个的数值来代表该区域的累积后验概率。
- 在机器人学中有很多分解方法。分解的粒度可能与环境的结构有关系,也可能没有关系。当有关系时,这种分解方法通常称为是拓扑逻辑的。
- 分解方法可以分为静态和动态的。静态的分解方法与置信度的形状无关。动态的分解方法依赖于机器人置信度具体的形状, 常常致力于正比于后验概率密度的方式提高空间分辨率。动态分解方法往往有更好的效果,但是实现起来更为复杂。
- 另一种非参数的方法称为粒子滤波器。粒子滤波器用从后验概率率中抽取的随机样本来表示后验概率。这些样本称为粒子。粒子滤波器很容易实现, 而且它们也是本书中描述的所有贝叶斯滤波器算法中应用最多的一种。
- 还有一些特殊的策略可以降低粒子滤波器的误差。最流行的技术就是降低由于算法的随机特性导致的估计的方差的技术, 以及根据后验概率的复杂程度选择粒子的数量的技术。
本章以及之前章节中所讨论的滤波器算法是本书中将要讨论的大多数概率机器人学算法的基础。这里列出的材料代表着现如今最流行的算法和概率机器人学中的表达形式。
4.5 参考文献
West和Harrison(1997)深度研究了本章和前面章节中讨论的几种技术。直方图已经在统计学中应用了几十年。 Sturges(1926)提供了一种早期用于选择直方图近似选择方案,Freedman和Diaconis(1981)提供了一种更新的方法。 Scott(1992)做出了一个现代的分析。一旦将状态空间映射为离散直方图,由此产生的时间推理问题就变成了一个离散隐马尔可夫模型的例子, Rabiner和Juang(1986)为这类方法的推广做了很大的贡献。MacDonald和Zucchini(1997)以及Elliott等(1995)是两个比较新的文献。
粒子滤波器可以追溯到Metropolis和Ulam(1949)提出的蒙特卡洛方法,Rubinstein(1981)对这类方法做了相对较新的介绍。 作为粒子滤波器的一部分,重要性重采样技术要追溯到两篇意义深远的文章Rubin(1988)以及Smith和Gelfand(1992)。 分层采样最早是由Neyman(1934)提出的。在过去的几年,在贝叶斯统计学的领域中粒子滤波器得到了广泛的研究(Doucet 1998; Kitagawa 1996; Liu and Chen 1998; Pitt and Shephard 1999)。在AI领域,粒子滤波器被(Kanazawa et al. 1995)以survival of the fittest的名义重新提出。 在计算机视觉中,Isard and Blake (1998)开发的称为condensation的算法应用粒子滤波器解决跟踪问题。 Doucet et al. (2001)对粒子滤波器做了较新的介绍。
4.6 练习
暂略