递推最小二乘法的推导和理解
本文的框架如下：
- 首先回忆一些最小二乘法的概念，如果很熟悉可以直接跳到递推最小二乘法，评判标准就是可以理解 ( X k T X k ) − 1 X k T Y k (X_k^{T}X_k)^{-1}X_k^{T}Y_k (XkTXk)−1XkTYk这个公式的推导。
- 之后介绍在线实时预测问题，引出递推最小二乘法并进行一些简单的理解。
最小二乘法
快速回顾最小二乘法的推导
最小二乘法解决的是给定一组输入数据和输出数据 ( x i , y i ) (x_i,y_i) (xi,yi)，对其进行参数估计的问题。达到的效果是估计出的参数可以使得这个方程很好拟合当前的数据。其实就是拟合。这里面最重要的思想就是误差平方最小。根据这个思想，推导的思路总体可以概括为：
- 根据误差平方生成损失函数
- 要最小化损失函数
建立误差平方
y
i
y_i
yi是现实中真正测量的值，而我们希望通过
f
(
x
i
)
f(x^i)
f(xi)求出的
y
y
y能够很好的和
y
i
y_i
yi吻合(这里
x
i
x^i
xi可以表示多个值)。现在已知有一组数据{
(
x
1
,
y
1
)
,
(
x
2
,
y
2
)
,
.
.
.
,
(
x
n
,
y
n
)
(x^1,y_1),(x^2,y_2),...,(x^n,y_n)
(x1,y1),(x2,y2),...,(xn,yn)}，假设对于每一组点都希望满足：
y
i
=
f
(
x
i
)
=
(
x
1
i
⋯
x
n
i
)
(
θ
1
⋮
θ
n
)
y_i=f(x^i)= \left( \begin{array}{ccc} x_{1}^i & \cdots & x_{n}^i \end{array} \right)\left( \begin{array}{c} θ_{1} \\ \vdots \\ θ_{n} \end{array} \right)
yi=f(xi)=(x1i⋯xni)⎝⎜⎛θ1⋮θn⎠⎟⎞
但问题是，不可能所有的
f
(
x
i
)
f(x_i)
f(xi)都和测量值相等。那我们自然想到，如果计算每个测量点和理论计算值的误差，让他们总和积累量最小，就可以近似认为满足需求，但是这个误差有正有负，所以我们对其进行平方处理。这里直接写成矩阵的形式：
J
(
θ
)
=
1
2
(
X
θ
−
Y
)
T
(
X
θ
−
Y
)
J(θ) =\frac{1}{2}(Xθ-Y)^T(Xθ-Y)
J(θ)=21(Xθ−Y)T(Xθ−Y)
注意这里
X
,
Y
X,Y
X,Y都是测量出来的数据。未知量只有
θ
θ
θ。
将其最小化
让损失对于
θ
θ
θ最小，自然就是对其进行求导等于零：
∂
∂
θ
J
(
θ
)
=
X
T
(
X
θ
−
Y
)
=
0
\frac{\partial}{\partial θ}J(θ)=X^T(Xθ-Y)=0
∂θ∂J(θ)=XT(Xθ−Y)=0
整理得出：
θ
=
(
X
T
X
)
−
1
X
T
Y
θ=(X^TX)^{-1}X^TY
θ=(XTX)−1XTY
一种对最小二乘法理解的视角
对于最小二乘问题，其实是一个求解方程组的问题：
X
θ
=
Y
Xθ=Y
Xθ=Y
最理想的情况是
X
X
X满秩，那么这样
θ
θ
θ就可以直接解出，为：
θ
=
X
−
1
Y
θ=X^{-1}Y
θ=X−1Y
但是一般情况下，
X
X
X是一个长条的矩阵。这在数学上有个处理，等式左右同时乘
X
T
X^T
XT(不做解释，MIT线性代数中老爷子讲的很好，大家如果想深入理解，可以去看，其实是最小二乘的思想，在Chapter4)，这就变成了最小二乘法：
θ
=
(
X
T
X
)
−
1
X
T
Y
θ=(X^TX)^{-1}X^TY
θ=(XTX)−1XTY
但是其本身想求解的还是：
θ
=
X
−
1
Y
θ=X^{-1}Y
θ=X−1Y
所以之后对递推最小二乘法的结果分析中，我们也可以进行这种简化的理解，即就是要方程组的解。
在进行了简单的回顾之后，下面我们引出递推最小二乘法。
递推最小二乘法
在线实时预测问题
现在我们改变一下需求，假如上述的数据不是一次性给出的，而是隔一段时间给一个数据，且需要根据之前的数据和现在多加进来的一个数据重新进行最小二乘的预测。或者换句话说，如果数据是实时在线给出的，我们需要怎么进行求解？
一种简单的想法就是每一时刻都进行一次最小二乘法的计算，这个当然可以，但是这是相当消耗内存和时间的。还有一种想法是采用迭代：既然我们在上一时刻已经计算出过一组参数，那么下一时刻能否只用上一时刻计算出来的参数，加上这一时刻得到的新的数据，计算出这一时刻最小二乘的结果。这样不仅不用重复计算，而且也不需要记忆数据，大大降低了数据的存储量。主要的思想其实就是找到一种迭代格式,让其满足：
θ
k
+
1
=
θ
k
+
ε
(
θ
k
,
x
k
+
1
)
θ_{k+1} = θ_{k}+ ε(θ_{k},x_{k+1})
θk+1=θk+ε(θk,xk+1)
推导思路与详细过程
首先需要知道，一组数据是一个列向量，在第 k k k时刻，一共有 k k k组数据，这个时候这些列向量共同组成了 X k X_k Xk。这里的下标可以理解成有 k k k组数据。
将k时刻的表达式写成k-1时刻表达式加某一个量
对于最小二乘法的解中的
X
k
T
X
k
X_k^TX_k
XkTXk：
X
k
T
X
k
=
[
x
1
.
.
.
x
k
]
[
x
1
T
⋮
x
k
T
]
=
[
x
1
.
.
.
x
k
−
−
1
]
[
x
1
T
⋮
x
k
−
1
T
]
+
x
k
x
k
T
=
X
k
−
1
T
X
k
−
1
+
x
k
x
k
T
X_k^TX_k=\begin{bmatrix} x_1 & ... & x_k \end{bmatrix}\begin{bmatrix} x_1^T\\\\ \vdots \\\\ x_k^T \end{bmatrix}=\begin{bmatrix} x_1 & ... & x_{k--1} \end{bmatrix}\begin{bmatrix} x_1^T\\\\ \vdots \\\\ x_{k-1}^T \end{bmatrix}+x_kx_k^T=X_{k-1}^TX_{k-1}+x_kx_k^T
XkTXk=[x1...xk]⎣⎢⎢⎢⎢⎢⎡x1T⋮xkT⎦⎥⎥⎥⎥⎥⎤=[x1...xk−−1]⎣⎢⎢⎢⎢⎢⎡x1T⋮xk−1T⎦⎥⎥⎥⎥⎥⎤+xkxkT=Xk−1TXk−1+xkxkT
同理，对于
X
k
T
Y
k
X_k^TY_k
XkTYk做相同的处理，得出：
X
k
T
Y
k
=
X
k
−
1
T
Y
k
−
1
+
x
k
y
k
X_k^TY_k=X_{k-1}^TY_{k-1}+x_ky_k
XkTYk=Xk−1TYk−1+xkyk
写出k-1时刻满足的最小二乘表达式
在
k
−
1
k-1
k−1时刻，也应该满足最小二乘的表达式：
θ
k
−
1
=
(
X
k
−
1
T
X
k
−
1
)
−
1
X
k
−
1
T
Y
k
−
1
θ_{k-1}=(X_{k-1}^{T}X_{k-1})^{-1}X_{k-1}^{T}Y_{k-1}
θk−1=(Xk−1TXk−1)−1Xk−1TYk−1
做一步变换：
(
X
k
−
1
T
X
k
−
1
)
θ
k
−
1
=
X
k
−
1
T
Y
k
−
1
(X_{k-1}^{T}X_{k-1})θ_{k-1}=X_{k-1}^{T}Y_{k-1}
(Xk−1TXk−1)θk−1=Xk−1TYk−1
将前两步的公式带入第k时刻的最小二乘表达式中
θ
k
=
(
X
k
T
X
k
)
−
1
X
k
T
Y
k
=
(
X
k
T
X
k
)
−
1
(
X
k
−
1
T
Y
k
−
1
+
x
k
y
k
)
=
(
X
k
T
X
k
)
−
1
(
(
X
k
−
1
T
X
k
−
1
)
θ
k
−
1
+
x
k
y
k
)
=
(
X
k
T
X
k
)
−
1
(
(
X
k
T
X
k
−
x
k
x
k
T
)
θ
k
−
1
+
x
k
y
k
)
=
θ
k
−
1
−
(
X
k
T
X
k
)
−
1
x
k
x
k
T
θ
k
−
1
+
(
X
k
T
X
k
)
−
1
x
k
y
k
=
θ
k
−
1
+
(
X
k
T
X
k
)
−
1
x
k
(
y
k
−
x
k
T
θ
k
−
1
)
\begin{aligned} θ_{k}&=(X_k^{T}X_k)^{-1}X_k^{T}Y_k\\[2mm] &=(X_k^{T}X_k)^{-1}(X_{k-1}^TY_{k-1}+x_ky_k)\\[2mm] &=(X_k^{T}X_k)^{-1}((X_{k-1}^{T}X_{k-1})θ_{k-1}+x_ky_k)\\[2mm] &=(X_k^{T}X_k)^{-1}((X_{k}^{T}X_{k}-x_kx_k^T)θ_{k-1}+x_ky_k)\\[2mm] &=θ_{k-1}-(X_k^{T}X_k)^{-1}x_kx_k^Tθ_{k-1}+(X_k^{T}X_k)^{-1}x_ky_k\\[2mm] &=θ_{k-1}+(X_k^{T}X_k)^{-1}x_k(y_k-x_k^Tθ_{k-1}) \end{aligned}
θk=(XkTXk)−1XkTYk=(XkTXk)−1(Xk−1TYk−1+xkyk)=(XkTXk)−1((Xk−1TXk−1)θk−1+xkyk)=(XkTXk)−1((XkTXk−xkxkT)θk−1+xkyk)=θk−1−(XkTXk)−1xkxkTθk−1+(XkTXk)−1xkyk=θk−1+(XkTXk)−1xk(yk−xkTθk−1)
至此其实递推最小二乘法的算法已经推到结束了。不过这里其实还有值也可以通过迭代来，就是
(
X
k
T
X
k
)
−
1
(X_k^{T}X_k)^{-1}
(XkTXk)−1这个量。这个我们之前已经推到过了，所以个人感觉最好是写成：
θ
k
=
θ
k
−
1
+
(
X
k
−
1
T
X
k
−
1
+
x
k
x
k
T
)
−
1
x
k
(
y
k
−
x
k
T
θ
k
−
1
)
θ_{k}=θ_{k-1}+(X_{k-1}^TX_{k-1}+x_kx_k^T)^{-1}x_k(y_k-x_k^Tθ_{k-1})
θk=θk−1+(Xk−1TXk−1+xkxkT)−1xk(yk−xkTθk−1)
因为这样，就可以只用上一时刻已经求出来的参数和当前得到的数据计算当前时刻的参数了。
公式的简单理解
我们主要看一下多出来的这一项
(
X
k
T
X
k
)
−
1
x
k
(
y
k
−
x
k
T
θ
k
−
1
)
(X_k^{T}X_k)^{-1}x_k(y_k-x_k^Tθ_{k-1})
(XkTXk)−1xk(yk−xkTθk−1)是什么：
对于
y
k
−
x
k
T
θ
k
−
1
y_k-x_k^Tθ_{k-1}
yk−xkTθk−1这一项，其实可以理解成，使用上一时刻的参数和当前获得的数据计算出来的预测值和当前的测量值之间的误差。即上一时刻参数对当前时刻的误差值。我们写成
ε
ε
ε。此时公式变成：
θ
k
=
θ
k
−
1
+
(
X
k
T
X
k
)
−
1
x
k
ε
θ_{k}=θ_{k-1}+(X_k^{T}X_k)^{-1}x_kε
θk=θk−1+(XkTXk)−1xkε
角度一：回归在线预测问题
我们将
(
X
k
T
X
k
)
−
1
x
k
(X_k^{T}X_k)^{-1}x_k
(XkTXk)−1xk看成加上当前时刻的信息之后产生的价值，这样可以将其简写成
K
k
K_k
Kk。这样就变成了：
θ
k
=
θ
k
−
1
+
K
k
ε
θ_{k}=θ_{k-1}+K_kε
θk=θk−1+Kkε
这个公式其实就回到了当时我们提出在线预测问题那里，我们希望使用上一时刻已经计算出来的参数加上某些值就能算出当前的参数，这里加的值就是使用上一时刻参数造成的误差乘上比例，这个比例可以理解为加上当前数据产生的价值。
角度二：梯度下降视角
之前在推导最小二乘法的时候我们计算过梯度为：
∂
∂
θ
J
(
θ
)
=
X
T
(
X
θ
−
Y
)
=
0
\frac{\partial}{\partial θ}J(θ)=X^T(Xθ-Y)=0
∂θ∂J(θ)=XT(Xθ−Y)=0
如果我们希望更新
θ
θ
θ，按照传统的梯度下降会怎么做：
θ
k
=
θ
k
−
1
+
α
X
T
(
X
θ
−
Y
)
θ_{k}=θ_{k-1}+αX^T(Xθ-Y)
θk=θk−1+αXT(Xθ−Y)
但是这里需要进行一个修改，即将
X
T
(
X
θ
−
Y
)
X^T(Xθ-Y)
XT(Xθ−Y)按照乘法对应关系，变成
x
k
(
y
k
−
x
k
T
θ
k
−
1
)
x_k(y_k-x_k^Tθ_{k-1})
xk(yk−xkTθk−1)。其实就是将
X
T
(
X
θ
−
Y
)
X^T(Xθ-Y)
XT(Xθ−Y)提取出了一部分，只不过这里的
θ
θ
θ不是每一个
k
k
k时刻都一样的(当然我们期望是一样的)。然后步长设置为
(
X
k
T
X
k
)
−
1
(X_k^{T}X_k)^{-1}
(XkTXk)−1,这个时候可能会有疑问，原来梯度下降，梯度是一个列向量，步长是一个标量，为什么这里步长变成了一个矩阵。写到这里我也思考了半天，于是有了第三个视角
角度三：状态方程视角下的 ( X k T X k ) − 1 (X_k^{T}X_k)^{-1} (XkTXk)−1：
状态方程在控制领域里面很常见，比如一个简单的基本运动学方程为：
x
(
k
+
1
)
=
x
(
k
)
+
B
k
u
(
k
)
x(k+1)=x(k)+B_ku(k)
x(k+1)=x(k)+Bku(k)
这是一个动态的变化过程，根据运动学方程，我们就可以根据当前的状态和给出的控制，知道下一时刻的状态。我们将状态变量类比乘递推最小二乘中的参数，控制类比成梯度，那么相对应的矩阵
B
k
B_k
Bk对应的就是
(
X
k
T
X
k
)
−
1
(X_k^{T}X_k)^{-1}
(XkTXk)−1。也就是说，从状态方程的角度来讲，
(
X
k
T
X
k
)
−
1
(X_k^{T}X_k)^{-1}
(XkTXk)−1起到的作用是如何将所谓的梯度，或者控制，作用到参数上面的。这里因为
X
k
X_k
Xk一直在变动，所以通俗的解释为：控制的变动是由数据量的变动引起的。
数据量太大： 矩阵求逆公式
这里主要是对
(
X
k
−
1
T
X
k
−
1
+
x
k
x
k
T
)
−
1
(X_{k-1}^TX_{k-1}+x_kx_k^T)^{-1}
(Xk−1TXk−1+xkxkT)−1进行一点补充，其实如果求这个东西的逆，如果数据量不大，正常直接求就好。但是如果处理的数据量特别大的话，求逆会变成一件很耗时的事情，所以这里有一个数学公式可以进行处理：
[
A
+
B
C
D
]
−
1
=
A
−
1
−
A
−
1
B
[
C
−
1
+
D
A
−
1
B
]
−
1
D
A
−
1
[A+BCD]^{-1}=A^{-1}-A^{-1}B[C^{-1}+DA^{-1}B]^{-1}DA^{-1}
[A+BCD]−1=A−1−A−1B[C−1+DA−1B]−1DA−1
对应这个公式，可以把
(
X
k
−
1
T
X
k
−
1
+
x
k
x
k
T
)
−
1
(X_{k-1}^TX_{k-1}+x_kx_k^T)^{-1}
(Xk−1TXk−1+xkxkT)−1写成：
(
X
k
−
1
T
X
k
−
1
+
x
k
x
k
T
)
−
1
=
X
k
−
1
T
X
k
−
1
−
X
k
−
1
T
X
k
−
1
x
k
x
k
T
X
k
−
1
T
X
k
−
1
1
+
x
k
T
X
k
−
1
T
X
k
−
1
x
k
(X_{k-1}^TX_{k-1}+x_kx_k^T)^{-1}=X_{k-1}^TX_{k-1}-\frac{X_{k-1}^TX_{k-1}x_kx_k^TX_{k-1}^TX_{k-1}}{1+x_k^TX_{k-1}^TX_{k-1}x_k}
(Xk−1TXk−1+xkxkT)−1=Xk−1TXk−1−1+xkTXk−1TXk−1xkXk−1TXk−1xkxkTXk−1TXk−1
9485
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
