RLS算法

贴一篇我在知乎的文章，这里面也阐述了最小二乘法与卡尔曼滤波的关系。两个算法一起说了，而且学习的话，比较轻量级。

https://zhuanlan.zhihu.com/p/67250500

一般最小二乘法

一般最小二乘法是给定若干观测值，计算一个最有可能的估计。
$$\left[\begin{array}{c}y(1)\\ y(2)\\ ⋮\\ y(k)\end{array}\right]=\left[\begin{array}{cccc}h(1{)}_1& h(1{)}_2& \cdots & h(1{)}_n\\ h(2{)}_1& h(2{)}_2& \cdots & h(2{)}_n\\ ⋮& ⋮& \ddots & ⋮\\ h(k{)}_1& h(k{)}_2& \cdots & h(k{)}_n\end{array}\right]\left[\begin{array}{c}{x}_1\\ x_2\\ ⋮\\ x_n\end{array}\right]$$
计算上述最有可能的$x$取值，是使得代价函数：
$$\begin{gathered} J(k)=(Y(k)-H(k)x{)}^T(Y(k)-H(k)x) \\ =Y^T(k)Y(k)-Y^T(k)H(k)x-x^T{H}^T(x)Y(k)+x^{T}H(k{)}^{T}H(k)x \end{gathered}$$
最小的 $x$的取值。令其为$\hat{x}(k)$,意思是$k$时刻的估计值
上述中，$H(k)$代表$k$时刻以及$k$时刻以前的观测,即为：
$$H(k)=[h(1{)}^T,h(2{)}^T,\cdots h(k{)}^T{]}^T$$
$Y(k)$类似：
$$Y(k)=[y(1),y(2),\cdots y(k){]}^T$$
当$J(k)$最小时，有：
$$\frac{\partial J(k)}{\partial x}=-2Y^T(k)H(k)+2x^T{H}^T(k)H(k)=0$$
此时，$x=\hat{x}(k)$
所以:$\hat{x}(k)=(H^T(k)H(k){)}^{-}{H}^T(k)Y(k)$
只有$H^T(k)H(k)$ 满秩时，才存在逆矩阵。
所以只有当$k\ge n$时，$H^T(k)H(k)$才可能存在逆矩阵。

带权最小二乘法

有时候，测量有好有坏，带权的最小二乘法就比较有必要了。可以人为赋予每个数据的一个置信度$r(k)$，有时候也会令$$r^2(k)=\frac{1}{D(v(k))}=\frac{1}{{\sigma }_k^2}$$
其中，$v(k)$ 是第$k$测量的误差。区别于估计误差。上面的式子也比较好理解。测量的不确定性会削弱这次测量的置信度。如果利用测量误差的话，则不能存在$0$误差的情况。改写$J(k)$的形式：
$$J(k)=\sum _{i=1}^k{r}^2(k)(y(k)-\hat{y}(k){)}^2$$
令：$$R(k)=\left[\begin{array}{cccc}{r}^2(1)& 0& \cdots & 0\\ 0& r^2(2)& \cdots & 0\\ ⋮& ⋮& \ddots & ⋮\\ 0& 0& \cdots & r^2(n)\end{array}\right]$$
那么：$$\begin{gathered} J(k)=(Y(k)-H(k)\hat{x}{)}^{T}R(Y(k)-H(k)\hat{x}) \\ =Y(k{)}^{T}R(k)Y(k)-Y^T(k)R(k)H(k)\hat{x}-{\hat{x}}^{T}H(k{)}^{T}R(k)Y(k)+{\hat{x}}^T{H}^T(k)R(K)H^T(k)\hat{x} \end{gathered}$$
令：$$\frac{\partial J(k)}{\partial \hat{x}}=-2Y^T(k)R(k)+2x^T{H}^T(k)R(k)H(k)=0$$
则：$$\hat{x}(k)=(H^T(k)R(k)H(k){)}^{-}{H}^T(k)R(k)Y(k)$$
这样设置权重，其实是对噪声的归一化。

$RLS$ 递推算法

有时候，我们不可能一次获得所有的观测结果，有时候随着观测结果的增加。从新修改参数使用最小二乘法计算，会很耗费时间。$RLS$ 就是等价于一般最小二乘法的递推计算方法。

这里介绍一个精巧的矩阵反演：
在这只要求 $A,D$可逆的情况下给定分块矩阵：
$$\left[\begin{array}{cc}A& B\\ C& D\end{array}\right]$$
令： $E=D-C{A}^{-}B,F=A-B{D}^{-}C$
计算可得：
$$\begin{gathered} \left[\begin{array}{cc}A& B\\ C& D\end{array}\right]\left[\begin{array}{cc}{A}^{-}+A^{-}B{E}^{-}C{A}^{-}& -A^{-}B{E}^{-}\\ -E^{-}C{A}^{-}& E^{-}\end{array}\right] \\ =\left[\begin{array}{cc}I& O\\ C{A}^{-}+C{A}^{-}B{E}^{-}C{A}^{-}-D{E}^{-}C{A}^{-}& -C{A}^{-}B{E}^{-}+D{E}^{-}\end{array}\right] \\ =\left[\begin{array}{cc}I& O\\ C{A}^{-}+(C{A}^{-}B-D)E^{-}C{A}^{-}& (-C{A}^{-}B+D)E^{-}\end{array}\right]=\left[\begin{array}{cc}I& O\\ O& I\end{array}\right] \end{gathered}$$
另一个方向：$$\left[\begin{array}{cc}A& B\\ C& D\end{array}\right]\left[\begin{array}{cc}{F}^{-}& -A^{-}B{E}^{-}\\ -E^{-}C{A}^{-}& E^{-}\end{array}\right]=\left[\begin{array}{cc}I& O\\ O& I\end{array}\right]$$
这个方向可以自行计算。
这也就是说,当矩阵$A,D,E,F$可逆时，:$(A+B{D}^{-}C{)}^{-}=A^{-}-A^{-}B(D+C{A}^{-}B{)}^{-}C{A}^{-}$

令$P(k)=(H^T(k)H(k){)}^{-}$
则有：
$$\begin{gathered} P(k)=(H^T(k)H(k){)}^{-} \\ =(H^T(k-1)H(k-1)+h^T(k)h(k){)}^{-} \\ =(P^{-}(k-1)+h^T(k)h(k){)}^{-} \\ =P(k-1)-\frac{P(k-1)h^T(k)h(k)P(k-1)}{1+h(k)P(k-1)h^T(k)} \end{gathered}$$
又有：$\hat{x}(k)=P(k)H^T(k)Y(k)$
可得：$P^{-}(k)\hat{x}(k)=H^T(k)Y(k)$
又因为：$P(k)=(P^{-}(k-1)+h^T(k)h(k){)}^{-}$
得：$P^{-}(k-1)=P^{-}(k)-h^T(k)h(k)$
所以：$$\begin{gathered} \hat{x}(k)=P(k)H^T(k)Y(k) \\ =P(k)(P^{-}(k-1)\hat{x}(k-1)+h^T(k)y(k)) \\ =P(k)((P^{-}(k)-h^T(k)h(k))\hat{x}(k-1)+h^T(k)y(k)) \\ =\hat{x}(k-1)+P(k)h^T(k)(y(k)-h(k)\hat{x}(k-1)) \end{gathered}$$
令:$$K(k)=\frac{P(k-1)h^T(k)}{1+h(k)P(k-1)h^T(k)}=P(k)h^T(k)$$
则：
$$\begin{gathered} P(k)=(I-K(k)h(k))P(k-1) \\ \hat{x}(k)=\hat{x}(k-1)+K(k)(y(k)-h(k)\hat{x}(k-1)) \end{gathered}$$
如何确定$P(0)$ 呢？
显然：$$P(k)=(\sum _{i=1}^k{h}^T(i)h(i){)}^{-}\approx (\sum _{i=1}^k{h}^T(i)h(i)+\frac{1}{\infty }I{)}^{-}$$
令$P(0)=\infty I$即可。那么此时$y(0)=h^T(0)=\hat{x}(0)=O_{n\times 1}$
这样定义对代价函数影响忽略不计。

带权RLS 递推计算

当考虑带权递推时 $$\hat{x}(k)=(H^T(k)R(k)H(k){)}^{-}{H}^T(k)R(k)Y(k)$$
一般$RLS$递推的$P(k)$重新定义为：$P(k)=(H^T(k)R(k)H(k){)}^{-}=(H^T(k-1)R(k-1)H(k-1)+h^T(k)r(k)h(k){)}^{-}$
令：$\gamma (k)=\frac{1}{r(k)}$
矩阵反演有：
$$P(k)=P(k-1)-\frac{P(k-1)h^T(k)h(k)P(k-1)}{\gamma (k)-h(k)P(k-1)h^T(k)}$$
由：$\hat{x}(k)=P(k)H^T(k)R(k)Y(k)$
有：$P^{-}(k)\hat{x}(k)=H^T(k)R(k)Y(k)$
由：$P(k)=(P^{-}(k-1)+h^T(k)r(k)h(k))$
有：$P^{-}(k-1)=P^{-}(k)-h^T(k)r(k)h(k)$
则:$$\begin{gathered} \hat{x}(k)=P(k)H^T(k)R(k)Y(k) \\ =P(k)(P^{-}(k-1)\hat{x}(k-1)+h^T(k)r(k)y(k)) \\ =P(k)((P^{-}(k)-h^T(k)r(k)h(k))\hat{x}(k-1)+h^T(k)r(k)y(k)) \\ =\hat{x}(k-1)+P(k)h^T(k)r(k)(y(k)-h(k)x(k)) \end{gathered}$$
令 $$K(k)=P(k)h^T(k)r(k)=\frac{P(k-1)h^T(k)}{\gamma (k)-h(k)P(k-1)h^T(k)}$$
则：$$P(k)=(I-K(k)h(k))P(k-1)$$
$$\hat{x}(k)=\hat{x}(k-1)+K(k)(y(k)-h(k)\hat{x}(k-1))$$
同上，当没有很好的初始设置数值时,可以讲，$P(0)=\infty I$，$\hat{x}(0)=h^T(0)=O_{n\times 1}$

带有遗忘的$RLS$

对于带有遗忘版本的$RLS$算法，将对近期数据更为敏感。
设置遗忘因子$\lambda \in (0,1]$
当$\lambda =1$时，不会遗忘。
每次更新，历史数据的权重会被整体缩小$\lambda$
那么此时：
$J(k)=\lambda J(k-1)+(y(k)-h(k)\hat{x}(k){)}^2$
则：
$J(k)=\lambda ||H(k-1)\hat{x}(k)-Y(k-1)|{|}^2+(y(k)-h(k)\hat{x}(k){)}^2$
从新定义
$H(k)=[\beta H(k-1),h(k)]$
$Y(k)=[\beta Y(k-1),y(k)]$
其中，${\beta }^2=\lambda$
$$\begin{gathered} P(k)=(H^T(k)H(k){)}^{-} \\ =(H^T(k-1)H(k-1)\lambda +h^T(k)h(k){)}^{-} \end{gathered}$$
$$=(P^{-}(k-1)\lambda +h^T(k)h(k){)}^{-}$$
$$=P(k-1)\frac{1}{\lambda }-\frac{\frac{1}{{\lambda }^2}P(k-1)h^T(k)h(k)P(k-1)}{1+\frac{1}{\lambda }h(k)P(k-1)h^T(k)}$$
$$=(I-\frac{P(k-1)h^T(k)h(k)}{\lambda +h(k)P(k-1)h^T(k)})P(k-1)\frac{1}{\lambda }$$
另一边：
$$\begin{gathered} \hat{x}(k)=(H^T(k)H(k){)}^{-}{H}^T(k)Y(k) \\ =P(k)H^T(k)Y(k) \end{gathered}$$
$\hat{x}(k)=P(k)H^T(k)Y(k)$
可得：$P^{-}(k)\hat{x}(k)=H^T(k)Y(k)$
同时：
$P(k)=(P^{-}(k-1)\lambda +h^T(k)h(k){)}^{-}$
$P^{-}(k-1)=\frac{1}{\lambda }(P^{-}(k)-h^T(k)h(k))$
$$\begin{gathered} \hat{x}(k)=P(k)(H^T(k-1)Y(k-1)\lambda +h^T(k)y(k)) \\ =P(k)(P^{-}(k-1)\hat{x}(k-1)\lambda +h^T(k)y(k)) \\ =P(k)\left(\right(P^{-}(k)-h^T(k)h(k))\hat{x}(k-1)+h^T(k)y(k)) \\ =\hat{x}(k-1)+P(k)h^T(k)(y(k)-h(k)\hat{x}(k-1)) \end{gathered}$$
令： $K(k)=P(k)h^T(k)$
那么：
$$\hat{x}(k)=\hat{x}(k-1)+K(k)(y(k)-h(k)\hat{x}(k-1))$$
$$\begin{gathered} K(k)=P(k)h(k)=(I-\frac{P(k-1)h^T(k)h(k)}{\lambda +h(k)P(k-1)h^T(k)})P(k-1)\frac{1}{\lambda }{h}^T(k) \\ =\frac{P(k-1)h^T(k)}{\lambda +h(k)P(k-1)h^T(k)} \end{gathered}$$
$$P(k)=(I-K(k)h(k))P(k-1)\frac{1}{\lambda }$$

总结

RLS 算法一个比较精巧的在线滤波算法。 还可以增加遗忘因子，让其对近期近期数据更为敏感，比如笔记平滑种可能需要这种形式的滤波器。