递推最小二乘法推导
递推最小二乘法（Recursive Least Squares, RLS）是一种时间序列分析方法，它用于在线更新线性回归模型的参数，而不需要重新拟合整个数据集。这种方法特别适用于数据流或实时系统。
1.最小二乘法原理
最小二乘法是一种数学优化技术，用于拟合数据点的最佳函数模型。在线性回归的情况下，我们试图找到一组参数 X X X ,使得模型预测 $ \hat{Y} = AX$ 与实际数据 Y Y Y 之间的平方误差之和最小。其中， A A A 是特征矩阵， X X X 是参数向量。
Y = ( y 1 y 2 y 3 ⋮ y m ) , A = ( a 11 ⋯ a 1 n ⋮ ⋮ a m 1 ⋯ a m n ) , X = ( x 1 x 2 x 3 ⋮ x n ) Y=\begin{pmatrix} y_1\\y_2\\y_3\\\vdots\\y_m \end{pmatrix}, A=\begin{pmatrix} a_{11} & \cdots & a_{1n}\\ \vdots & &\vdots \\a_{m1} & \cdots & a_{mn} \end{pmatrix}, X=\begin{pmatrix} x_1\\x_2\\x_3\\\vdots\\x_n \end{pmatrix} Y= y1y2y3⋮ym ,A= a11⋮am1⋯⋯a1n⋮amn ,X= x1x2x3⋮xn
式中， A A A 为输入的数据样本， X X X 为拟合的参数， Y Y Y 为输入的数据样本， Y ^ \hat{Y} Y^ 为输出的估计
为了估计出 n n n 个参数，要求 m ≥ n m\ge n m≥n ,当 m = n m=n m=n 时，存在唯一解。 m > n m>n m>n 时，由于数据中存在模型噪声和测量噪声，一般不可能确定一组参数，定义误差矢量 E = Y − A X ^ , E = [ e 1 , ⋯ , e m ] T E=Y-A\hat{X},E=[e_1,\cdots,e_m]^T E=Y−AX^,E=[e1,⋯,em]T
最小二乘法思想：求测量值与估计值的平方和最小，即：
J = E T E = ∑ i = 1 m e i 2 = ( Y − A X ^ ) T ( Y − A X ^ ) = Y T Y − Y T A X ^ − X ^ T A T Y + X ^ T A T A X ^ J=E^TE=\sum_{i=1}^{m}e_i^2=(Y-A\hat{X})^T(Y-A\hat{X})=Y^TY-Y^TA\hat{X}-\hat{X}^TA^TY+\hat{X}^TA^TA\hat{X} J=ETE=∑i=1mei2=(Y−AX^)T(Y−AX^)=YTY−YTAX^−X^TATY+X^TATAX^
对 J J J 求偏导：
X ^ T A T A X ^ \hat{X}^TA^TA\hat{X} X^TATAX^是一个二次项，其偏导数需要使用矩阵的Frechet导数或利用以下规则：如果 f ( X ) = X T B X f(X)=X^TBX f(X)=XTBX，其中 B B B是一个常数矩阵，则 d f d X = 2 B X \frac{df}{dX}=2BX dXdf=2BX在这里， B = A T A B=A^TA B=ATA。
J J J 是标量，X是向量，求导是对X的每一个元素求偏导， d J / d X dJ/dX dJ/dX 的维度是 n × 1 n\times1 n×1
所以， d J d X ^ = − 2 A T Y + 2 A T A X ^ \frac{dJ}{d\hat{X}}=-2A^TY+2A^TA\hat{X} dX^dJ=−2ATY+2ATAX^
令偏导等于0，即求出最优的拟合参数 X ^ \hat{X} X^
− 2 A T Y + 2 A T A X ^ = 0 ⇒ X ^ = ( A T A ) − 1 A T Y -2A^TY+2A^TA\hat{X}=0 \Rightarrow\hat{X}=(A^TA)^{-1}A^TY −2ATY+2ATAX^=0⇒X^=(ATA)−1ATY
这样，当我们有了足够的样本，可以一次性求解出最优的参数向量X
2.递归最小二乘法
设第 k − 1 k-1 k−1 时刻的样本矩阵 A k − 1 , Y k − 1 A_{k-1},Y_{k-1} Ak−1,Yk−1 分别为：
A k − 1 = ( a 11 ⋯ a 1 n ⋮ ⋮ a ( k − 1 ) 1 ⋯ a ( k − 1 ) n ) , Y k − 1 = ( y 1 y 2 y 3 ⋮ y k − 1 ) A_{k-1}=\begin{pmatrix} a_{11} & \cdots & a_{1n}\\ \vdots & &\vdots \\a_{(k-1)1} & \cdots & a_{(k-1)n} \end{pmatrix},Y_{k-1}=\begin{pmatrix} y_1\\y_2\\y_3\\\vdots\\y_{k-1} \end{pmatrix} Ak−1= a11⋮a(k−1)1⋯⋯a1n⋮a(k−1)n ,Yk−1= y1y2y3⋮yk−1
那么，第k时刻的样本矩阵可以表示为： A k = ( A k − 1 a k ) , Y k = ( Y k − 1 y k ) A_k=\begin{pmatrix} A_{k-1} \\ a_k \end{pmatrix},Y_k=\begin{pmatrix} Y_{k-1} \\ y_k \end{pmatrix} Ak=(Ak−1ak),Yk=(Yk−1yk)
计算中间式：
A k T A k = [ A k − 1 T a k T ] [ A k − 1 a k ] = A k − 1 T A k − 1 + a k T a k A_k^TA_k=\begin{bmatrix}A_{k-1}^T & a_k^T\end{bmatrix}\begin{bmatrix}A_{k-1} \\ a_k\end{bmatrix}=A_{k-1}^TA_{k-1}+a_k^Ta_k AkTAk=[Ak−1TakT][Ak−1ak]=Ak−1TAk−1+akTak
A k T Y k = [ A k − 1 T a k T ] [ Y k − 1 y k ] = A k − 1 T Y k − 1 + a k T y k A_k^TY_k=\begin{bmatrix}A_{k-1}^T & a_k^T\end{bmatrix}\begin{bmatrix}Y_{k-1} \\ y_k\end{bmatrix}=A_{k-1}^TY_{k-1}+a_k^Ty_k AkTYk=[Ak−1TakT][Yk−1yk]=Ak−1TYk−1+akTyk
由第一节推导出的最优拟合结果公式： X ^ k = ( A k T A k ) − 1 A k T Y k \hat{X}_{k}=(A_k^TA_k)^{-1}A_k^TY_k X^k=(AkTAk)−1AkTYk
向前递推得： X ^ k − 1 = ( A k − 1 T A k − 1 ) − 1 A k − 1 T Y k − 1 ⇒ A k − 1 T Y k − 1 = ( A k − 1 T A k − 1 ) X ^ k − 1 \hat{X}_{k-1}=(A_{k-1}^TA_{k-1})^{-1}A_{k-1}^TY_{k-1}\Rightarrow A_{k-1}^TY_{k-1}=(A_{k-1}^TA_{k-1})\hat{X}_{k-1} X^k−1=(Ak−1TAk−1)−1Ak−1TYk−1⇒Ak−1TYk−1=(Ak−1TAk−1)X^k−1
那么： X ^ k = ( A k T A k ) − 1 A k T Y k = ( A k T A k ) − 1 ( A k − 1 T Y k − 1 + a k T y k ) = ( A k T A k ) − 1 ( A k − 1 T A k − 1 X ^ k − 1 + a k T y k ) \hat{X}_{k}=(A_k^TA_k)^{-1}A_k^TY_k=(A_k^TA_k)^{-1}(A_{k-1}^TY_{k-1}+a_k^Ty_k)=(A_k^TA_k)^{-1}(A_{k-1}^TA_{k-1}\hat{X}_{k-1}+a_k^Ty_k) X^k=(AkTAk)−1AkTYk=(AkTAk)−1(Ak−1TYk−1+akTyk)=(AkTAk)−1(Ak−1TAk−1X^k−1+akTyk)
将 A k − 1 T A k − 1 A_{k-1}^TA_{k-1} Ak−1TAk−1 用中间式1进行替换： X ^ k = ( A k T A k ) − 1 ( ( A k T A k − a k T a k ) X ^ k − 1 + a k T y k ) = [ I − ( A k T A k ) − 1 a k T a k ] X ^ k − 1 + ( A k T A k ) − 1 a k T y k = X ^ k − 1 + ( A k T A k ) − 1 a k T ( y k − a k X ^ k − 1 ) \hat{X}_{k}=(A_k^TA_k)^{-1}((A_k^TA_k-a_k^Ta_k)\hat{X}_{k-1}+a_k^Ty_k)=[I-(A_k^TA_k)^{-1}a_k^Ta_k]\hat{X}_{k-1}+(A_k^TA_k)^{-1}a_k^Ty_k=\hat{X}_{k-1}+(A_k^TA_k)^{-1}a_k^T(y_k-a_k\hat{X}_{k-1}) X^k=(AkTAk)−1((AkTAk−akTak)X^k−1+akTyk)=[I−(AkTAk)−1akTak]X^k−1+(AkTAk)−1akTyk=X^k−1+(AkTAk)−1akT(yk−akX^k−1)
即迭代公式初步化简为： X ^ k = X ^ k − 1 + ( A k T A k ) − 1 a k T ( y k − a k X ^ k − 1 ) \hat{X}_{k}=\hat{X}_{k-1}+(A_k^TA_k)^{-1}a_k^T(y_k-a_k\hat{X}_{k-1}) X^k=X^k−1+(AkTAk)−1akT(yk−akX^k−1)
考虑到 ( A k T A k ) − 1 (A_k^TA_k)^{-1} (AkTAk)−1 计算复杂，使用谢尔曼公式对矩阵逆运算进行化简。
谢尔曼公式为： ( A + u v ) − 1 = A − 1 − ( A − 1 u ) ( I + v A − 1 u ) − 1 ( v A − 1 ) (A+uv)^{-1}=A^{-1}-(A^{-1}u)(I+vA^{-1}u)^{-1}(vA^{-1}) (A+uv)−1=A−1−(A−1u)(I+vA−1u)−1(vA−1)
结合中间式1求解 ( A k T A k ) − 1 (A_k^TA_k)^{-1} (AkTAk)−1
令 G k − 1 = ( A k T A k ) − 1 G_k^{-1}=(A_k^TA_k)^{-1} Gk−1=(AkTAk)−1 ,带入中间式1则 G k − 1 = ( G k − 1 + a k T a k ) − 1 G_k^{-1}=(G_{k-1}+a_k^Ta_k)^{-1} Gk−1=(Gk−1+akTak)−1
套用谢尔曼公式： G k − 1 = G k − 1 − 1 − （ G k − 1 − 1 a k T ） ( I + a k G k − 1 − 1 a k T ) − 1 ( a k G k − 1 − 1 ) G_k^{-1}=G_{k-1}^{-1}-（G_{k-1}^{-1}a_k^T）(I+a_kG_{k-1}^{-1}a_k^T)^{-1}(a_kG_{k-1}^{-1}) Gk−1=Gk−1−1−（Gk−1−1akT）(I+akGk−1−1akT)−1(akGk−1−1)
不难看出： ( I + a k G k − 1 − 1 a k T ) − 1 = I 1 + a k G k − 1 − 1 a k T (I+a_kG_{k-1}^{-1}a_k^T)^{-1}=\frac{I}{1+a_kG_{k-1}^{-1}a_k^T} (I+akGk−1−1akT)−1=1+akGk−1−1akTI
并且令 P k = G k − 1 P_k=G_k^{-1} Pk=Gk−1
可以得出另一个递推公式： P k = P k − 1 − P k − 1 a k T a k P k − 1 1 + a k P k − 1 a k T P_k=P_{k-1}-\frac{P_{k-1}a_k^Ta_kP_{k-1}}{1+a_kP_{k-1}a_k^T} Pk=Pk−1−1+akPk−1akTPk−1akTakPk−1
所以最小二乘法递推形式表示为：
{ P k = P k − 1 − P k − 1 a k T a k P k − 1 1 + a k P k − 1 a k T X ^ k = X ^ k − 1 + P k a k T ( y k − a k X ^ k − 1 ) \begin{cases} P_k=P_{k-1}-\frac{P_{k-1}a_k^Ta_kP_{k-1}}{1+a_kP_{k-1}a_k^T} \\ \hat{X}_{k}=\hat{X}_{k-1}+P_ka_k^T(y_k-a_k\hat{X}_{k-1}) \end{cases} {Pk=Pk−1−1+akPk−1akTPk−1akTakPk−1X^k=X^k−1+PkakT(yk−akX^k−1)
扫一扫
5734
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
