本文旨在梳理总结学习到的一些知识。由于笔者水平有限,文中难免存在一些不严谨和错误之处,诚请各位批评指正。
最近看了一篇文章,文章从最小二乘法的角度推导了卡尔曼滤波的公式(链接在文末)。看完后震惊不已,很受启发,于是写了这篇文章。一是为了倒逼输出从而达到知新的效果,二是为了记录一下自己的理解以便日后自己翻阅。
最小二乘法
最小二乘法是解决曲线拟合问题最常用的方法,本质是进行参数估计。实现方式是先使用待定系数法设出拟合函数,然后使用拟合函数和观测数据构造出目标函数来评估拟合效果,并取拟合效果最好时的拟合函数。拟合函数由参数和事先选定的一组线性无关的函数(这里我们姑且叫它基函数)构成,基函数已事先确定,求解拟合函数的问题便转变为参数估计的问题。
为拟合函数,其中
为参数(也就是我们求解的重点)、
为事先选定的一组线性无关的函数:
其中除以
的目的是加权,噪声方差越大的观测数据其权值越小,保证整体结果不会被个别噪声较大的数据所影响。同时也满足了高斯-马尔可夫定理的要求:
在线性回归模型中,如果误差满足零均值、同方差且互不相关,则回归系数的最佳线性无偏估计(BLUE, Best Linear Unbiased Estimator)就是普通最小二乘法估计。
我们的目的便是使
取到最小值,将此时得到的
作为最佳的拟合函数。
如马同学的这张图为例,图中的
至
便是观测数据并且每个都对应
轴上的一个坐标。通过观察我们可以取拟合函数为
,其中
便是我们要估计的参数。
分别对
求偏导,当满足下式时目标函数
取到最小值,可求出参数
:
从而得到拟合函数
。当然如果你对拟合出的函数不满意,可以再取其他的函数作为拟合函数并去估计其参数。
一个简单的状态估计问题
还是这张图,不过这次我们给坐标轴一些物理意义。横坐标代表时间,纵坐标代表位置,我们将其视为一维的匀速运动目标的位置-时间图。你可以想象成一个小车沿着直线匀速前进,你每隔一段时间观测它一次并记录位置和时间数据。
分我们定义状态量
,其中
代表目标所处的位置,
代表目标当前的速度。
不过在被赋予了物理意义之后它现在长这样:
,其中
代表函数与纵轴的截距,也就是
时刻时的位置,
代表目标的运动速度。不过我们对目标 0 时刻的位置不感兴趣,我们对它当前的状态估计值感兴趣,于是有:
其中
代表对当前目标位置的估计值,
代表对当前目标速度的估计值。这样的拟合函数
才是我们想要的,它包含了目标在
时刻的位置与速度。换句话说,我们通过这五个观测数据得到了目标在
时刻的状态估计值
。
于是,如果我们能继续获得更多的观测数据(一直到
),那么有:
其中
为
时刻下的观测矩阵,它将目标
时刻下的状态
转化观测值
。
这是个很怪的拟合函数,非常反直觉。它的观测方式似乎是使用最新的状态估计值,去获得其他时间节点上的观测值;如此一来它既是在对状态进行观测,又是实现了不同时刻间状态的转移。但我们根本不关心拟合函数本身,我们只关心构成它的参数。它的参数包含了目标当前时刻的状态估计值
。如果我们不断地继续获得数据,我们也能不断地对拟合函数的参数进行估计,从而得出目标最新的状态估计值。
铺垫
但是我们发现这个方法 太 慢 了。随着迭代推进,每一次迭代都需要用到历史的所有数据来估计目标当前的状态。也就是说随着运行时间增长,积累的历史数据越多,计算出目标当前状态估计值所需要的时间也就越长。
虽然变成了矩阵形式,但最小二乘法的思想没有变,通过类比一维的情况可以快速的理解矩阵形式。
我们不妨来实实在在的求一下矩阵形式下拟合函数的参数。
在当前我们已经有
个观测数据的情况下,我们对当前目标的状态估计量求导,并令其为零:
上式中的
便是我们目前获得到的观测量数量,可见随着迭代的进行我们的运算量会越来越大。
递推最小二乘法
我们的目的就是避免运算量随着时间而增长,所以必须想方设法将 改为递推形式。换句话说,只使用
时刻的各种数据来推算
时刻的状态,而不是像现在这样将所有的历史信息全部用上。但是忘记历史意味着背叛,所以我们选一个折中的办法——“只送大脑”,我们可以递推地求解目标每个时刻状态的协方差矩阵。这样可以将历史信息蕴含在协方差矩阵中,达到我们的目的。
式
的形式非常有趣,它以
为基准,并通过
的方式来估计
时刻观测量的误差(
表示对
时刻观测量的估计),最后乘以
并补偿到
上。我们一般称
这样的系数为观测增益。这种综合了观测数据与状态数据的计算方式已经有些接近卡尔曼滤波了,但是只有更新过程,几乎没有预测过程。比如在
式中的
和
式中的
完全可以替换成由他们自己经过预测后得出的
时刻的估计量,然后再进行更新,这也是卡尔曼滤波的思想之一。
问题的根源在于我们怪异而又反直觉的观测矩阵。它虽然叫观测矩阵,但它耦合了观测与状态转移两个功能。观测指的是将状态量转化为观测量的过程,状态转移指的是将某时刻的状态转化为另一个时刻的状态。与卡尔曼滤波相比,递推最小二乘法在整个递推过程中缺少了过程噪声。所以我们需要对其进行解耦合,重新定义观测矩阵与状态转移矩阵。
从最小二乘法到卡尔曼滤波
为状态转移矩阵,它可以将
时刻的状态转变为
时刻的。其中
为目标在
时刻的先验状态量,也就是说它是我们使用
时刻的状态预测得来的。
为
时刻的后验状态量,也就是经过上轮递推后得到的
时刻的最优状态估计值。
同时我们用
来表示过程噪声,其中
代表位置噪声,
代表速度噪声,且有:
但要注意的是此时(
时刻)计算的状态协方差矩阵为先验状态协方差矩阵
,也就是我们根据
时刻的状态协方差矩阵推算得出的。
发现式
到
中只出现过
,也就是观测矩阵的下标只有
。观测矩阵在下标为
(也就是最新的时刻)时只具有观测意义,不存在状态转移。因为在定义时观测矩阵
的意义是将最新时刻
目标的状态转变为
时刻的观测量(可以看二、中举的例子),当
时便只有观测意义。于是我们不需要重新定义观测矩阵,当前的式子符合我们的要求。
总结
最小二乘法求解最小目标函数的过程与卡尔曼滤波中求解协方差矩阵迹的最小值的过程非常相似,只不过最小二乘法是使用所有数据来进行优化,而卡尔曼滤波通过对迹求导来获得最优卡尔曼增益。我们可以通过构造特殊拟合函数的方式来让最小二乘法来递推运作,从而只传递状态噪声。然后定义状态转移过程,并加入过程噪声来描述状态转移误差。
相比使用贝叶斯公式与正态分布假设来推导卡尔曼滤波,最小二乘法推导的条件更为宽松。使用贝叶斯公式来推导时需要假设所有的噪声服从正态分布,而最小二乘法推导仅需要满足高斯-马尔可夫定理(噪声零均值、同方差、不相关)即可。当然推导卡尔曼滤波的方法并不只这两种,但了解的越多理解越深刻,对学习也更有帮助。
从递推最小二乘法推导出卡尔曼滤波的过程并不严谨,直接使用先验估计来替换掉递推最小二乘法中的一些项。如果以后有更加严谨的方式我会补上证明。
[1] https://zhuanlan.zhihu.com/p/67250500
[2] https://zhuanlan.zhihu.com/p/339118204
[3] https://www.zhihu.com/question/37031188/answer/411760828