本章解决这样一个问题：对于一个 LTI 离散状态空间方程模型（注意与系统建模中的区别）
{ x k + 1 = A x k + B u k y k = C x k + D u k \left\{ \begin{aligned} x_{k+1}&= Ax_k+Bu_k\\ y_k&= Cx_k+Du_k \end{aligned} \right. { x k + 1 y k = A x k + B u k = C x k + D u k
给定系统若干个输入 u 0 , ⋯ , u s − 1 u_0,\cdots,u_{s-1}{} u 0 , ⋯ , u s − 1 s s s y k y_k y k m m m l l l
系统的阶数 n n n 四个矩阵 A A A B B B C C C D D D # 子空间方程# 推导假定当前时刻 k = i k=i k = i k = 0 ∼ ( i − 1 ) k=0\sim (i-1) k = 0 ∼ ( i − 1 ) i i i
∣ x 1 = A x 0 + B u 0 x 2 = A x 1 + B u 1 = A 2 x 0 + A B u 0 + B u 1 x 3 = A x 2 + B u 2 = A 3 x 0 + A 2 B u 0 + A B u 1 + B u 2 ⋯ x i = A i x 0 + ∑ k = 0 ( i − 1 ) A ( i − 1 ) − k B u k ∣ y 0 = C x 0 + D u 0 y 1 = C x 1 + D u 1 = C A x 0 + C B u 0 + D u 1 y 2 = C x 2 + D u 2 = C A 2 x 0 + C A B u 0 + C B u 1 + D u 2 ⋯ y i − 1 = C A i − 1 x 0 + ∑ k = 0 ( i − 2 ) C A ( i − 2 ) − k B u k + D u i − 1 \left|\ \begin{aligned} x_1&= Ax_0+Bu_0\\ x_2&= Ax_1+Bu_1\\ &= A^2x_0+ABu_0+Bu_1\\ x_3&= Ax_2+Bu_2\\ &= A^3x_0+A^2Bu_0+ABu_1+Bu_2\\ &\cdots\\ x_i&= A^ix_0+\sum\limits_{k=0}^{(i-1)}A^{(i-1)-k}Bu_k \end{aligned} \right. \quad \left|\ \begin{aligned} y_0 &= C x_0 + D u_0\\ y_1 &= Cx_1+Du_1\\ &= C A x_0 + C B u_0 + D u_1\\ y_2 &= Cx_2+Du_2\\ &= C A^2 x_0 + C A B u_0 + C B u_1 + D u_2\\ &\cdots\\ y_{i-1} &= CA^{i-1}x_0+\sum\limits_{k=0}^{(i-2)}CA^{(i-2)-k}Bu_k+Du_{i-1} \end{aligned} \right. ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ x 1 x 2 x 3 x i = A x 0 + B u 0 = A x 1 + B u 1 = A 2 x 0 + A B u 0 + B u 1 = A x 2 + B u 2 = A 3 x 0 + A 2 B u 0 + A B u 1 + B u 2 ⋯ = A i x 0 + k = 0 ∑ ( i − 1 ) A ( i − 1 ) − k B u k ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ ∣ y 0 y 1 y 2 y i − 1 = C x 0 + D u 0 = C x 1 + D u 1 = C A x 0 + C B u 0 + D u 1 = C x 2 + D u 2 = C A 2 x 0 + C A B u 0 + C B u 1 + D u 2 ⋯ = C A i − 1 x 0 + k = 0 ∑ ( i − 2 ) C A ( i − 2 ) − k B u k + D u i − 1
把 y y y
[ y 0 y 1 ⋮ y i − 2 y i − 1 ] = [ C C A ⋮ C A i − 2 C A i − 1 ] x 0 + [ D 0 ⋯ 0 0 C B D ⋯ 0 0 ⋮ ⋮ ⋱ ⋮ ⋮ C A i − 3 B C A i − 4 B ⋯ D 0 C A i − 2 B C A i − 3 B ⋯ C B D ] [ u 0 u 1 ⋮ u i − 2 u i − 1 ] \begin{bmatrix} y_0\\y_1\\\vdots\\y_{i-2}\\y_{i-1} \end{bmatrix} = \begin{bmatrix} C\\ CA\\ \vdots\\ CA^{i-2}\\ CA^{i-1} \end{bmatrix} x_0 + \begin{bmatrix} D & 0 & \cdots & 0 & 0\\ CB & D & \cdots & 0 & 0\\ \vdots & \vdots & \ddots & \vdots & \vdots \\ CA^{i-3}B & CA^{i-4}B & \cdots & D & 0\\ CA^{i-2}B & CA^{i-3}B & \cdots & CB & D \end{bmatrix} \begin{bmatrix} u_0\\u_1\\\vdots\\u_{i-2}\\u_{i-1} \end{bmatrix} ⎣ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎡ y 0 y 1 ⋮ y i − 2 y i − 1 ⎦ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎤ = ⎣ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎡ C C A ⋮ C A i − 2 C A i − 1 ⎦ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎤ x 0 + ⎣ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎡ D C B ⋮ C A i − 3 B C A i − 2 B 0 D ⋮ C A i − 4 B C A i − 3 B ⋯ ⋯ ⋱ ⋯ ⋯ 0 0 ⋮ D C B 0 0 ⋮ 0 D ⎦ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎤ ⎣ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎡ u 0 u 1 ⋮ u i − 2 u i − 1 ⎦ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎤
记 Γ i \Gamma_i Γ i x 0 x_0 x 0 H i H_i H i u u u
[ y 0 y 1 ⋮ y i − 1 ] = Γ i x 0 + H i [ u 0 u 1 ⋮ u i − 1 ] \begin{bmatrix} y_0\\y_1\\ \vdots\\ y_{i-1} \end{bmatrix} = \Gamma_ix_0 + H_i \begin{bmatrix} u_0\\u_1\\\vdots\\u_{i-1} \end{bmatrix} ⎣ ⎢ ⎢ ⎢ ⎢ ⎡ y 0 y 1 ⋮ y i − 1 ⎦ ⎥ ⎥ ⎥ ⎥ ⎤ = Γ i x 0 + H i ⎣ ⎢ ⎢ ⎢ ⎢ ⎡ u 0 u 1 ⋮ u i − 1 ⎦ ⎥ ⎥ ⎥ ⎥ ⎤
刚才我们的 k k k 0 ∼ ( i − 1 ) 0\sim(i-1) 0 ∼ ( i − 1 ) k k k 1 ∼ i 1\sim i 1 ∼ i y y y y 1 ∼ y i y_1\sim y_i y 1 ∼ y i x 1 x_1 x 1 u u u u 1 ∼ u i u_1\sim u_i u 1 ∼ u i Γ i \Gamma_i Γ i H i H_i H i i i i x 1 x_1 x 1 i i i y y y
从 0 0 0 j − 1 j-1 j − 1 j j j
[ y 0 y 1 ⋯ y j − 1 y 1 y 2 ⋯ y j ⋮ ⋮ ⋱ ⋮ y i − 1 y i ⋯ y i + j − 2 ] = Γ i [ x 0 x 1 ⋯ x j − 1 ] + H i [ u 0 u 1 ⋯ u j − 1 u 1 u 2 ⋯ u j ⋮ ⋮ ⋱ ⋮ u i − 1 u i ⋯ u i + j − 2 ] \begin{bmatrix} y_0 & y_1 & \cdots & y_{j-1} \\ y_1 & y_2 & \cdots & y_{j} \\ \vdots & \vdots & \ddots & \vdots \\ y_{i-1} & y_{i} & \cdots & y_{i+j-2} \end{bmatrix} = \Gamma_i \begin{bmatrix} x_0 & x_1 & \cdots & x_{j-1} \end{bmatrix} + H_i \begin{bmatrix} u_0 & u_1 & \cdots & u_{j-1} \\ u_1 & u_2 & \cdots & u_{j} \\ \vdots & \vdots & \ddots & \vdots \\ u_{i-1} & u_{i} & \cdots & u_{i+j-2} \end{bmatrix} ⎣ ⎢ ⎢ ⎢ ⎢ ⎡ y 0 y 1 ⋮ y i − 1 y 1 y 2 ⋮ y i ⋯ ⋯ ⋱ ⋯ y j − 1 y j ⋮ y i + j − 2 ⎦ ⎥ ⎥ ⎥ ⎥ ⎤ = Γ i [ x 0 x 1 ⋯ x j − 1 ] + H i ⎣ ⎢ ⎢ ⎢ ⎢ ⎡ u 0 u 1 ⋮ u i − 1 u 1 u 2 ⋮ u i ⋯ ⋯ ⋱ ⋯ u j − 1 u j ⋮ u i + j − 2 ⎦ ⎥ ⎥ ⎥ ⎥ ⎤
我们再次压缩表达式为
Y p = Γ i X p + H i U p Y_p=\Gamma_iX_p+H_iU_p Y p = Γ i X p + H i U p
这里的 p p p
现在考虑 “未来状态”，特别地，是未来状态相对于当前状态（即 k = i k=i k = i i i i i i i i + j − 1 i+j-1 i + j − 1
[ y i y i + 1 ⋯ y i + j − 1 y i + 1 y i + 2 ⋯ y i + j ⋮ ⋮ ⋱ ⋮ y 2 i − 1 y 2 i ⋯ y 2 i + j − 2 ] = Γ i [ x i x i + 1 ⋯ x i + j − 1 ] + H i [ u i u i + 1 ⋯ u i + j − 1 u i + 1 u i + 2 ⋯ u i + j ⋮ ⋮ ⋱ ⋮ u 2 i − 1 u 2 i ⋯ u 2 i + j − 2 ] \begin{bmatrix} y_i & y_{i+1} & \cdots & y_{i+j-1} \\ y_{i+1} & y_{i+2} & \cdots & y_{i+j} \\ \vdots & \vdots & \ddots & \vdots \\ y_{2i-1} & y_{2i} & \cdots & y_{2i+j-2} \end{bmatrix} = \Gamma_i \begin{bmatrix} x_i & x_{i+1} & \cdots & x_{i+j-1} \end{bmatrix} + H_i \begin{bmatrix} u_i & u_{i+1} & \cdots & u_{i+j-1} \\ u_{i+1} & u_{i+2} & \cdots & u_{i+j} \\ \vdots & \vdots & \ddots & \vdots \\ u_{2i-1} & u_{2i} & \cdots & u_{2i+j-2} \end{bmatrix} ⎣ ⎢ ⎢ ⎢ ⎢ ⎡ y i y i + 1 ⋮ y 2 i − 1 y i + 1 y i + 2 ⋮ y 2 i ⋯ ⋯ ⋱ ⋯ y i + j − 1 y i + j ⋮ y 2 i + j − 2 ⎦ ⎥ ⎥ ⎥ ⎥ ⎤ = Γ i [ x i x i + 1 ⋯ x i + j − 1 ] + H i ⎣ ⎢ ⎢ ⎢ ⎢ ⎡ u i u i + 1 ⋮ u 2 i − 1 u i + 1 u i + 2 ⋮ u 2 i ⋯ ⋯ ⋱ ⋯ u i + j − 1 u i + j ⋮ u 2 i + j − 2 ⎦ ⎥ ⎥ ⎥ ⎥ ⎤
压缩表达式为
Y f = Γ i X f + H i U f Y_f=\Gamma_iX_f+H_iU_f Y f = Γ i X f + H i U f
这里的 f f f
过去和未来的情况都有了，我们使用一开始算出来的各个 x x x i i i 0 0 0 x i = A i x 0 + ∑ k = 0 ( i − 1 ) A ( i − 1 ) − k B u k x_i= A^ix_0+\sum\limits_{k=0}^{(i-1)}A^{(i-1)-k}Bu_k x i = A i x 0 + k = 0 ∑ ( i − 1 ) A ( i − 1 ) − k B u k
x i = A i x 0 + [ A i − 1 B A i − 2 B ⋯ A B B ] [ u 0 u 1 ⋮ u i − 2 u i − 1 ] = A i x 0 + Δ i [ u 0 u 1 ⋮ u i − 1 ] x_i= A^ix_0+\begin{bmatrix}A^{i-1}B&A^{i-2}B&\cdots&AB&B\end{bmatrix}\begin{bmatrix}u_0\\u_1\\\vdots\\u_{i-2}\\u_{i-1}\end{bmatrix}=A^ix_0+\Delta_i\begin{bmatrix}u_0\\u_1\\\vdots\\u_{i-1}\end{bmatrix} x i = A i x 0 + [ A i − 1 B A i − 2 B ⋯ A B B ] ⎣ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎡ u 0 u 1 ⋮ u i − 2 u i − 1 ⎦ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎤ = A i x 0 + Δ i ⎣ ⎢ ⎢ ⎢ ⎢ ⎡ u 0 u 1 ⋮ u i − 1 ⎦ ⎥ ⎥ ⎥ ⎥ ⎤
同样地，依然可以把下标整体移动，而两个系数矩阵不变。然后再把 0 0 0 j − 1 j-1 j − 1
[ x i x i + 1 ⋯ x i + j − 1 ] = A i [ x 0 x 1 ⋯ x j − 1 ] + Δ i [ u 0 u 1 ⋯ u j − 1 u 1 u 2 ⋯ u j ⋮ ⋮ ⋱ ⋮ u i − 1 u i ⋯ u i + j − 2 ] \begin{bmatrix}x_i&x_{i+1}&\cdots&x_{i+j-1}\end{bmatrix}=A^i\begin{bmatrix}x_0&x_1&\cdots&x_{j-1}\end{bmatrix}+\Delta_i\begin{bmatrix}u_0&u_1&\cdots&u_{j-1}\\u_1&u_2&\cdots&u_j\\\vdots&\vdots&\ddots&\vdots\\u_{i-1}&u_i&\cdots&u_{i+j-2}\end{bmatrix} [ x i x i + 1 ⋯ x i + j − 1 ] = A i [ x 0 x 1 ⋯ x j − 1 ] + Δ i ⎣ ⎢ ⎢ ⎢ ⎢ ⎡ u 0 u 1 ⋮ u i − 1 u 1 u 2 ⋮ u i ⋯ ⋯ ⋱ ⋯ u j − 1 u j ⋮ u i + j − 2 ⎦ ⎥ ⎥ ⎥ ⎥ ⎤
自然地，过去与未来 x x x
X f = A i X p + Δ i U p X_f=A^iX_p+\Delta_iU_p X f = A i X p + Δ i U p
最后一步。我们实际上只想知道输入 Y Y Y U U U
第①式移项后代入第③式消 X p X_p X p X f = A i Γ i † ( Y p − H i U p ) + Δ i U p = [ A i Γ i † , Δ i − A i Γ i † H i ] [ Y p U p ] X f = L p W p \begin{aligned}X_f&= A^i\Gamma_i^\dagger(Y_p-H_iU_p)+\Delta_iU_p\\&= \begin{bmatrix}A^i\Gamma_i^\dagger&,&\Delta_i-A^i\Gamma_i^\dagger H_i\end{bmatrix}\begin{bmatrix}Y_p\\U_p\end{bmatrix}\\X_f&= L_pW_p\end{aligned} X f X f = A i Γ i † ( Y p − H i U p ) + Δ i U p = [ A i Γ i † , Δ i − A i Γ i † H i ] [ Y p U p ] = L p W p
相应地第②式改为 Y f = Γ i X f + H i U f = Γ i L p W p + H i U f Y_f=\Gamma_iX_f+H_iU_f=\Gamma_iL_pW_p+H_iU_f Y f = Γ i X f + H i U f = Γ i L p W p + H i U f # 小结对一个 m m m l l l
{ x k + 1 = A x k + B u k y k = C x k + D u k \left\{ \begin{aligned} x_{k+1}&= Ax_k+Bu_k\\ y_k&= Cx_k+Du_k \end{aligned} \right. { x k + 1 y k = A x k + B u k = C x k + D u k
我们有三条子空间辨识方程：
将过去状态的 y y y Y p = Γ i X p + H i U p Y_p=\Gamma_iX_p+H_iU_p Y p = Γ i X p + H i U p 将未来状态的 y y y Y f = Γ i X f + H i U f Y_f=\Gamma_iX_f+H_iU_f Y f = Γ i X f + H i U f 相隔 i i i x x x X f = A i X p + Δ i U p X_f=A^iX_p+\Delta_iU_p X f = A i X p + Δ i U p 其中
X p = [ x 0 x 1 ⋯ x j − 1 ] X f = [ x i x i + 1 ⋯ x i + j − 1 ] \begin{aligned} X_p&= \begin{bmatrix}x_0&x_1&\cdots&x_{j-1}\end{bmatrix}\\ X_f&= \begin{bmatrix}x_i&x_{i+1}&\cdots&x_{i+j-1}\end{bmatrix} \end{aligned} X p X f = [ x 0 x 1 ⋯ x j − 1 ] = [ x i x i + 1 ⋯ x i + j − 1 ]
Y p = [ y 0 y 1 ⋯ y j − 1 y 1 y 2 ⋯ y j ⋮ ⋮ ⋱ ⋮ y i − 1 y i ⋯ y i + j − 2 ] Y f = [ y i y i + 1 ⋯ y i + j − 1 y i + 1 y i + 2 ⋯ y i + j ⋮ ⋮ ⋱ ⋮ y 2 i − 1 y 2 i ⋯ y 2 i + j − 2 ] Y_p=\begin{bmatrix} y_0 & y_1 & \cdots & y_{j-1} \\ y_1 & y_2 & \cdots & y_{j} \\ \vdots & \vdots & \ddots & \vdots \\ y_{i-1} & y_{i} & \cdots & y_{i+j-2} \end{bmatrix}\quad Y_f=\begin{bmatrix} y_i & y_{i+1} & \cdots & y_{i+j-1} \\ y_{i+1} & y_{i+2} & \cdots & y_{i+j} \\ \vdots & \vdots & \ddots & \vdots \\ y_{2i-1} & y_{2i} & \cdots & y_{2i+j-2} \end{bmatrix} Y p = ⎣ ⎢ ⎢ ⎢ ⎢ ⎡ y 0 y 1 ⋮ y i − 1 y 1 y 2 ⋮ y i ⋯ ⋯ ⋱ ⋯ y j − 1 y j ⋮ y i + j − 2 ⎦ ⎥ ⎥ ⎥ ⎥ ⎤ Y f = ⎣ ⎢ ⎢ ⎢ ⎢ ⎡ y i y i + 1 ⋮ y 2 i − 1 y i + 1 y i + 2 ⋮ y 2 i ⋯ ⋯ ⋱ ⋯ y i + j − 1 y i + j ⋮ y 2 i + j − 2 ⎦ ⎥ ⎥ ⎥ ⎥ ⎤
U p = [ u 0 u 1 ⋯ u j − 1 u 1 u 2 ⋯ u j ⋮ ⋮ ⋱ ⋮ u i − 1 u i ⋯ u i + j − 2 ] U f = [ u i u i + 1 ⋯ u i + j − 1 u i + 1 u i + 2 ⋯ u i + j ⋮ ⋮ ⋱ ⋮ u 2 i − 1 u 2 i ⋯ u 2 i + j − 2 ] U_p=\begin{bmatrix}u_0&u_1&\cdots&u_{j-1}\\u_1&u_2&\cdots&u_j\\\vdots&\vdots&\ddots&\vdots\\u_{i-1}&u_i&\cdots&u_{i+j-2}\end{bmatrix}\quad U_f=\begin{bmatrix} u_i & u_{i+1} & \cdots & u_{i+j-1} \\ u_{i+1} & u_{i+2} & \cdots & u_{i+j} \\ \vdots & \vdots & \ddots & \vdots \\ u_{2i-1} & u_{2i} & \cdots & u_{2i+j-2} \end{bmatrix} U p = ⎣ ⎢ ⎢ ⎢ ⎢ ⎡ u 0 u 1 ⋮ u i − 1 u 1 u 2 ⋮ u i ⋯ ⋯ ⋱ ⋯ u j − 1 u j ⋮ u i + j − 2 ⎦ ⎥ ⎥ ⎥ ⎥ ⎤ U f = ⎣ ⎢ ⎢ ⎢ ⎢ ⎡ u i u i + 1 ⋮ u 2 i − 1 u i + 1 u i + 2 ⋮ u 2 i ⋯ ⋯ ⋱ ⋯ u i + j − 1 u i + j ⋮ u 2 i + j − 2 ⎦ ⎥ ⎥ ⎥ ⎥ ⎤
Γ i = [ C C A ⋮ C A i − 2 C A i − 1 ] H i = [ D 0 ⋯ 0 0 C B D ⋯ 0 0 ⋮ ⋮ ⋱ ⋮ ⋮ C A i − 3 B C A i − 4 B ⋯ D 0 C A i − 2 B C A i − 3 B ⋯ C B D ] \Gamma_i=\begin{bmatrix} C\\ CA\\ \vdots\\ CA^{i-2}\\ CA^{i-1} \end{bmatrix}\quad H_i=\begin{bmatrix} D & 0 & \cdots & 0 & 0\\ CB & D & \cdots & 0 & 0\\ \vdots & \vdots & \ddots & \vdots & \vdots \\ CA^{i-3}B & CA^{i-4}B & \cdots & D & 0\\ CA^{i-2}B & CA^{i-3}B & \cdots & CB & D \end{bmatrix} Γ i = ⎣ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎡ C C A ⋮ C A i − 2 C A i − 1 ⎦ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎤ H i = ⎣ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎡ D C B ⋮ C A i − 3 B C A i − 2 B 0 D ⋮ C A i − 4 B C A i − 3 B ⋯ ⋯ ⋱ ⋯ ⋯ 0 0 ⋮ D C B 0 0 ⋮ 0 D ⎦ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎤
Δ i = [ A i − 1 B A i − 2 B ⋯ A B B ] \Delta_i=\begin{bmatrix}A^{i-1}B&A^{i-2}B&\cdots&AB&B\end{bmatrix} Δ i = [ A i − 1 B A i − 2 B ⋯ A B B ]
消去状态变量的形式：
Y f = Γ i L p W p + H i U f Y_f=\Gamma_iL_pW_p+H_iU_f Y f = Γ i L p W p + H i U f
其中
L p = [ A i Γ i † , Δ i − A i Γ i † H i ] W p = [ Y p U p ] L_p=\begin{bmatrix}A^i\Gamma_i^\dagger&,&\Delta_i-A^i\Gamma_i^\dagger H_i\end{bmatrix}\quad\quad W_p=\begin{bmatrix}Y_p\\U_p\end{bmatrix} L p = [ A i Γ i † , Δ i − A i Γ i † H i ] W p = [ Y p U p ]
# 三条假设与两条推论为了让辨识算法有效，需要假定系统足够好，即满足：
r ( X p ) = n {\rm r}(X_p)=n r ( X p ) = n 因为是 n n n r ( [ U p U f ] ) = 2 i m {\rm r}\left(\begin{bmatrix}U_p\\U_f\end{bmatrix}\right)=2im r ( [ U p U f ] ) = 2 i m 即输入样本要满秩，保证输入样本可以探索到整个输入空间. 2 i m 2im 2 i m m m m 2 i 2i 2 i i i i r o w ( U p ) ∩ r o w ( X p ) = ∅ {\rm row}(U_p)\,\cap\,{\rm row}(X_p)=\varnothing r o w ( U p ) ∩ r o w ( X p ) = ∅ 输入空间与状态空间独立，保证系统开环，排除反馈信号的可能 若系统满足上述三个条件，有两个推论
推论 1 ：r ( X f ) = n {\rm r}(X_f)=n r ( X f ) = n
证：X f = A i X p + Δ i U p = [ A i Δ i ] [ X p U p ] X_f=A^iX_p+\Delta_iU_p=\begin{bmatrix}A^i & \Delta_i\end{bmatrix}\begin{bmatrix}X_p \\U_p\end{bmatrix} X f = A i X p + Δ i U p = [ A i Δ i ] [ X p U p ] r ( U p ) = i m {\rm r}(U_p) = im r ( U p ) = i m [ X p U p ] \begin{bmatrix}X_p \\U_p\end{bmatrix} [ X p U p ]
r ( X f ) = r ( [ A i Δ i ] ) ≥ r ( Δ i ) = n {\rm r}(X_f) = {\rm r}(\begin{bmatrix}A^i&\Delta_i\end{bmatrix})\geq {\rm r}(\Delta_i)=n r ( X f ) = r ( [ A i Δ i ] ) ≥ r ( Δ i ) = n
最后一步是因为 Δ i \Delta_i Δ i s p a n ( X f ) {\rm span}(X_f) s p a n ( X f ) r ( X f ) ≤ n {\rm r}(X_f) \leq n r ( X f ) ≤ n r ( X f ) = n {\rm r}(X_f) = n r ( X f ) = n
推论 2 ：s p a n ( X f ) = s p a n ( W p ) ∩ s p a n ( W f ) {\rm span}(X_f)={\rm span}(W_p)\cap {\rm span}(W_f) s p a n ( X f ) = s p a n ( W p ) ∩ s p a n ( W f ) W p = [ Y p U p ] W_p=\begin{bmatrix}Y_p\\U_p\end{bmatrix} W p = [ Y p U p ] W f = [ Y f U f ] W_f=\begin{bmatrix}Y_f\\U_f\end{bmatrix} W f = [ Y f U f ]
这意味着，X f X_f X f
证：(1) 证明 span ( X f ) ⊆ span ( W p ) ∩ span ( W f ) \text{span}(X_f) \subseteq \text{span}(W_p) \cap \text{span}(W_f) span ( X f ) ⊆ span ( W p ) ∩ span ( W f )
由于系统是可达和可观测，系统的状态 X f X_f X f X f X_f X f W p W_p W p W f W_f W f X f ∈ span ( W p ) X_f \in \text{span}(W_p) X f ∈ span ( W p ) X f ∈ span ( W f ) X_f \in \text{span}(W_f) X f ∈ span ( W f )
span ( X f ) ⊆ span ( W p ) ∩ span ( W f ) \text{span}(X_f) \subseteq \text{span}(W_p) \cap \text{span}(W_f) span ( X f ) ⊆ span ( W p ) ∩ span ( W f )
(2) 证明 span ( W p ) ∩ span ( W f ) ⊆ span ( X f ) \text{span}(W_p) \cap \text{span}(W_f) \subseteq \text{span}(X_f) span ( W p ) ∩ span ( W f ) ⊆ span ( X f )
∀ v ∈ span ( W p ) ∩ span ( W f ) \forall\,v \in \text{span}(W_p) \cap \text{span}(W_f) ∀ v ∈ span ( W p ) ∩ span ( W f ) ∃ a , b s.t. \exists\,a,\,b\,\ \text{s.t.}{} ∃ a , b s.t.
v = W p ⋅ a = W f ⋅ b v = W_p\cdot a = W_f\cdot b v = W p ⋅ a = W f ⋅ b
而 W p W_p W p W f W_f W f X f X_f X f v v v X f ⋅ c X_f\cdot c X f ⋅ c v ∈ span ( X f ) v \in \text{span}(X_f) v ∈ span ( X f )
span ( W p ) ∩ span ( W f ) ⊆ span ( X f ) \text{span}(W_p) \cap \text{span}(W_f) \subseteq \text{span}(X_f) span ( W p ) ∩ span ( W f ) ⊆ span ( X f )
# 确定性子空间辨识# N4SIDN umerical algorithm for S ubspace S tate S pace S ystem ID entification 子空间状态空间系统辨识的数值算法
使用消去状态变量的形式 Y f = Γ i L p W p + H i U f Y_f=\Gamma_iL_pW_p+H_iU_f Y f = Γ i L p W p + H i U f
Y f / U f W p = Γ i L p W p / U f W p + H i U f / U f W p = Γ i L p W p = Γ i X f \begin{aligned} Y_f/\!_{U_f}W_p&= \Gamma_iL_pW_p/\!_{U_f}W_p+H_iU_f/\!_{U_f}W_p\\ &= \Gamma_iL_pW_p\\ &= \Gamma_iX_f \end{aligned} Y f / U f W p = Γ i L p W p / U f W p + H i U f / U f W p = Γ i L p W p = Γ i X f
对 O = Y f / U f W p \mathcal{O}=Y_f/\!_{U_f}W_p O = Y f / U f W p
O = [ U 1 U 2 ] [ Σ r 0 0 0 ] [ V 1 T V 2 T ] = U 1 Σ r V 1 T ( = Γ i X f ) \mathcal{O}=\left[ U_1 \quad U_2 \right] \begin{bmatrix} \Sigma_r & 0 \\ 0 & 0 \end{bmatrix} \begin{bmatrix} V_1^\mathrm{T} \\ V_2^\mathrm{T}\end{bmatrix}=U_1\Sigma_rV_1^\mathrm{T}\;(\,=\Gamma_iX_f) O = [ U 1 U 2 ] [ Σ r 0 0 0 ] [ V 1 T V 2 T ] = U 1 Σ r V 1 T ( = Γ i X f )
为了得到 Γ i X f \Gamma_iX_f Γ i X f Σ r \Sigma_r Σ r
Γ i = U 1 Σ r 1 / 2 T \Gamma_i = U_1 \Sigma_r^{1/2}\color{red}T Γ i = U 1 Σ r 1 / 2 T X f = T − 1 Σ r 1 / 2 V 1 T X_f = {\color{red}T^{-1} }\Sigma_r^{1/2} V_1^\mathrm{T}{} X f = T − 1 Σ r 1 / 2 V 1 T 下面就可以进行辨识了。由于 T T T V 1 V_1 V 1 n = r ( X f ) = r ( Σ r ) n=\mathrm{r}(X_f)=\mathrm{r}(\Sigma_r) n = r ( X f ) = r ( Σ r ) O \mathcal O O
然后就是求 A B C D ABCD A B C D X f = [ x i , x i + 1 , ⋯ , x i + j − 1 ] X_f=[x_{i},\,x_{i+1},\,\cdots,\,x_{i+j-1}] X f = [ x i , x i + 1 , ⋯ , x i + j − 1 ] x x x
[ X i + 1 , j − 1 Y i , j − 1 ] = [ A B C D ] [ X i , j − 1 U i , j − 1 ] \begin{bmatrix} X_{i+1,\,j-1} \\ Y_{i,\,j-1} \end{bmatrix} = \begin{bmatrix} A & B \\ C & D \end{bmatrix} \begin{bmatrix} X_{i,\,j-1} \\ U_{i,\,j-1} \end{bmatrix} [ X i + 1 , j − 1 Y i , j − 1 ] = [ A C B D ] [ X i , j − 1 U i , j − 1 ]
其中大写的 X X X Y Y Y U U U
X i + 1 , j − 1 = [ x i + 1 , ⋯ , x i + j − 1 ] Y i , j − 1 = [ y i , y i + 1 , ⋯ , y i + j − 2 ] X i , j − 1 = [ x i , x i + 1 , ⋯ , x i + j − 2 ] U i , j − 1 = [ u i , u i + 1 , ⋯ , u i + j − 2 ] \begin{aligned} X_{i+1,\,j-1} &= [x_{i+1},\,\cdots,\,x_{i+j-1}]\\ Y_{i,\,j-1} &= [y_{i},\,y_{i+1},\,\cdots,\,y_{i+j-2}]\\ X_{i,\,j-1} &= [x_{i},\,x_{i+1},\,\cdots,\,x_{i+j-2}]\\ U_{i,\,j-1} &= [u_{i},\,u_{i+1},\,\cdots,\,u_{i+j-2}]\\ \end{aligned} X i + 1 , j − 1 Y i , j − 1 X i , j − 1 U i , j − 1 = [ x i + 1 , ⋯ , x i + j − 1 ] = [ y i , y i + 1 , ⋯ , y i + j − 2 ] = [ x i , x i + 1 , ⋯ , x i + j − 2 ] = [ u i , u i + 1 , ⋯ , u i + j − 2 ]
既然 X f X_f X f A B C D ABCD A B C D
y = Θ x y= \Theta x y = Θ x
∂ ∥ y − Θ x ∥ 2 ∂ Θ = 2 ( y − Θ x ) x T = 0 Θ ^ = y x T ( x x T ) − 1 \begin{aligned} \dfrac{ {\partial}\|y-\Theta x\|^2}{ {\partial}\Theta}&= 2(y-\Theta x)x^\mathrm{T}=0\\ \hat\Theta&= yx^\mathrm{T}(xx^\mathrm{T})^{-1} \end{aligned} ∂ Θ ∂ ∥ y − Θ x ∥ 2 Θ ^ = 2 ( y − Θ x ) x T = 0 = y x T ( x x T ) − 1
现在剩下一个没弄明白的点：为什么要在 Σ r \Sigma_r Σ r
# MOESPM ultivariable O utput E rror S tate sP ace 多变量输出误差状态空间
做 LQ 分解
[ U p Y p ] = [ L 11 0 L 21 L 22 ] [ Q 1 T Q 2 T ] \begin{bmatrix}U_p\\Y_p\end{bmatrix}= \begin{bmatrix}L_{11}&0\\L_{21}&L_{22}\end{bmatrix}\begin{bmatrix}Q_1^{\mathrm T}\\Q_2^{\mathrm T}\end{bmatrix} [ U p Y p ] = [ L 1 1 L 2 1 0 L 2 2 ] [ Q 1 T Q 2 T ]
⇒ { U p = L 11 Q 1 T Y p = L 21 Q 1 T + L 22 Q 2 T \Rightarrow\,\left\{ \begin{aligned} U_p&= L_{11}Q_1^{\mathrm T}\\ Y_p&= L_{21}Q_1^{\mathrm T}+L_{22}Q_2^{\mathrm T}\\ \end{aligned} \right. ⇒ { U p Y p = L 1 1 Q 1 T = L 2 1 Q 1 T + L 2 2 Q 2 T
由子空间方程之 “过去状态拼一起”：Y p = Γ i X p + H i U p Y_p=\Gamma_iX_p+H_iU_p Y p = Γ i X p + H i U p Y p Y_p Y p U p U_p U p
Γ i X p + H i L 11 Q 1 T = L 21 Q 1 T + L 22 Q 2 T (1) \boxed{\Gamma_iX_p+H_iL_{11}Q_1^{\mathrm T}= L_{21}Q_1^{\mathrm T}+L_{22}Q_2^{\mathrm T} }\tag{1} Γ i X p + H i L 1 1 Q 1 T = L 2 1 Q 1 T + L 2 2 Q 2 T ( 1 )
两边右乘 Q 2 Q_2 Q 2
Γ i X p Q 2 = L 22 \Gamma_iX_pQ_2=L_{22} Γ i X p Q 2 = L 2 2
和 N4SID 类似，为了得到 Γ i X p Q 2 \Gamma_iX_pQ_2 Γ i X p Q 2 L 22 L_{22}{} L 2 2
L 22 = [ U 1 U 2 ] [ Σ r 0 0 0 ] [ V 1 T V 2 T ] = U 1 Σ r V 1 T ( = Γ i X p Q 2 ) L_{22}=\left[ U_1 \quad U_2 \right] \begin{bmatrix} \Sigma_r & 0 \\ 0 & 0 \end{bmatrix} \begin{bmatrix} V_1^\mathrm{T} \\ V_2^\mathrm{T}\end{bmatrix}=U_1\Sigma_rV_1^{\mathrm T}\ (\,=\Gamma_iX_pQ_2) L 2 2 = [ U 1 U 2 ] [ Σ r 0 0 0 ] [ V 1 T V 2 T ] = U 1 Σ r V 1 T ( = Γ i X p Q 2 )
然后从 Σ r \Sigma_r Σ r 同样的，没懂，为啥这样劈开？？ ）
Γ i = U 1 Σ r 1 / 2 \Gamma_i= U_1\Sigma_r^{1/2}{} Γ i = U 1 Σ r 1 / 2 X p Q 2 = Σ r 1 / 2 V 1 T X_pQ_2= \Sigma_r^{1/2}V_1^\mathrm{T}{} X p Q 2 = Σ r 1 / 2 V 1 T 然后就可以开始辨识了：n = r ( X p ) = r ( Σ r ) n=\mathrm{r}(X_p)={\mathrm r}(\Sigma_r) n = r ( X p ) = r ( Σ r ) Q 2 Q_2 Q 2
而 Γ i = [ C C A ⋮ C A i − 1 ] \Gamma_i=\begin{bmatrix}C\\CA\\\vdots\\CA^{i-1}\end{bmatrix} Γ i = ⎣ ⎢ ⎢ ⎢ ⎢ ⎡ C C A ⋮ C A i − 1 ⎦ ⎥ ⎥ ⎥ ⎥ ⎤ A A A C C C
Γ i \Gamma_i Γ i l l l C C C ( Γ i (\Gamma_i ( Γ i l l l ) = ( Γ i )=(\Gamma_i ) = ( Γ i l l l ) × A )\times A ) × A A A A 至于 B B B D D D U 2 T U_2^\mathrm{T}{} U 2 T Γ i \Gamma_i Γ i L 22 L_{22}{} L 2 2 U 1 U_1 U 1 U U U U 2 T U_2^\mathrm{T}{} U 2 T
U 2 T Γ i X p + U 2 T H i L 11 Q 1 T = U 2 T L 21 Q 1 T + U 2 T L 22 Q 2 T \begin{aligned} \bcancel{U_2^\mathrm{T}\Gamma_iX_p}+U_2^\mathrm{T}H_iL_{11}Q_1^{\mathrm T}= U_2^\mathrm{T}L_{21}Q_1^{\mathrm T}+\bcancel{U_2^\mathrm{T}L_{22}Q_2^{\mathrm T} } \end{aligned} U 2 T Γ i X p + U 2 T H i L 1 1 Q 1 T = U 2 T L 2 1 Q 1 T + U 2 T L 2 2 Q 2 T
去掉 Q 1 T Q_1^\mathrm{T}{} Q 1 T L 11 L_{11}{} L 1 1 H i H_i H i
U 2 T [ D 0 ⋯ 0 0 C B D ⋯ 0 0 ⋮ ⋮ ⋱ ⋮ ⋮ C A i − 3 B C A i − 4 B ⋯ D 0 C A i − 2 B C A i − 3 B ⋯ C B D ] = U 2 T L 21 L 11 − 1 U_2^\mathrm{T}\begin{bmatrix} D & 0 & \cdots & 0 & 0\\ CB & D & \cdots & 0 & 0\\ \vdots & \vdots & \ddots & \vdots & \vdots \\ CA^{i-3}B & CA^{i-4}B & \cdots & D & 0\\ CA^{i-2}B & CA^{i-3}B & \cdots & CB & D \end{bmatrix}=U_2^\mathrm{T}L_{21}L_{11}^{-1} U 2 T ⎣ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎡ D C B ⋮ C A i − 3 B C A i − 2 B 0 D ⋮ C A i − 4 B C A i − 3 B ⋯ ⋯ ⋱ ⋯ ⋯ 0 0 ⋮ D C B 0 0 ⋮ 0 D ⎦ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎤ = U 2 T L 2 1 L 1 1 − 1
中间这个大矩阵 H i H_i H i U 2 T U_2^\mathrm{T}{} U 2 T l l l m m m H i H_i H i i i i
U 2 T = [ L 1 , L 2 , ⋯ , L i ] U 2 T L 21 L 11 − 1 = [ M 1 , M 2 , ⋯ , M i ] \begin{aligned} U_2^\mathrm{T}&=[\mathcal{L}_1,\,\mathcal{L}_2,\,\cdots,\,\mathcal{L}_i]\\ U_2^\mathrm{T}L_{21}L_{11}^{-1}&=[\mathcal{M}_1,\,\mathcal{M}_2,\,\cdots,\,\mathcal{M}_i] \end{aligned} U 2 T U 2 T L 2 1 L 1 1 − 1 = [ L 1 , L 2 , ⋯ , L i ] = [ M 1 , M 2 , ⋯ , M i ]
分块之后，实际上是阶梯型，而且 B B B D D D
M 1 = L 1 D + ( L 2 C + ⋯ + L i − 1 C A i − 3 + L i C A i − 2 ) B M 2 = L 2 D + ( L 3 C + ⋯ + L i C A i − 3 ) B ⋮ M i − 1 = L i − 1 D + ( L i C ) B M i = L i D \begin{aligned} \mathcal{M}_1&=\mathcal{L}_1 D + (\mathcal{L}_2 C + \cdots + \mathcal{L}_{i-1} C A^{i-3} + \mathcal{L}_i C A^{i-2})B \\ \mathcal{M}_2&=\mathcal{L}_2 D + (\mathcal{L}_3 C + \cdots + \mathcal{L}_i C A^{i-3}) B \\ &\ \ \vdots \\ \mathcal{M}_{i-1}&= \mathcal{L}_{i-1} D + (\mathcal{L}_i C )B \\ \mathcal{M}_{i}&=\mathcal{L}_i D \end{aligned} M 1 M 2 M i − 1 M i = L 1 D + ( L 2 C + ⋯ + L i − 1 C A i − 3 + L i C A i − 2 ) B = L 2 D + ( L 3 C + ⋯ + L i C A i − 3 ) B ⋮ = L i − 1 D + ( L i C ) B = L i D
注意到括号里面很像 Γ i \Gamma_i Γ i L k ˉ = [ L k , L k + 1 , ⋯ , L i ] \bar{\mathcal{L}_k}=\left[\mathcal{L}_k, \mathcal{L}_{k+1}, \cdots, \mathcal{L}_i\right] L k ˉ = [ L k , L k + 1 , ⋯ , L i ]
[ M 1 M 2 ⋮ M i − 1 M i ] = [ L 1 L 2 ˉ Γ i − 1 L 2 L 3 ˉ Γ i − 2 ⋮ ⋮ L i − 1 L i ˉ Γ 1 L i 0 ] [ D B ] \begin{bmatrix} \mathcal{M}_1 \\ \mathcal{M}_2 \\ \vdots \\ \mathcal{M}_{i-1} \\ \mathcal{M}_i \end{bmatrix} = \begin{bmatrix} \mathcal{L}_1 & \bar{\mathcal{L}_2}\Gamma_{i-1} \\ \mathcal{L}_2 & \bar{\mathcal{L}_3}\Gamma_{i-2} \\ \vdots & \vdots \\ \mathcal{L}_{i-1} & \bar{\mathcal{L}_i}\Gamma_1 \\ \mathcal{L}_i & 0 \end{bmatrix} \begin{bmatrix} D \\ B \end{bmatrix} ⎣ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎡ M 1 M 2 ⋮ M i − 1 M i ⎦ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎤ = ⎣ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎡ L 1 L 2 ⋮ L i − 1 L i L 2 ˉ Γ i − 1 L 3 ˉ Γ i − 2 ⋮ L i ˉ Γ 1 0 ⎦ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎤ [ D B ]
然后只需要使用最小二乘法（传统的就行了）就能解 B B B D D D
BD 这一段的构造，逻辑还是不够丝滑
# 随机子空间辨识# 问题描述在状态空间方程中引入一个随机偏移
{ x k + 1 = A x k + B u k + w k y k = C x k + D u k + v k \left\{ \begin{aligned} x_{k+1}&= Ax_k+Bu_k+w_k\\ y_k&= Cx_k+Du_k+v_k \end{aligned} \right. { x k + 1 y k = A x k + B u k + w k = C x k + D u k + v k
其中 w k w_k w k v k v_k v k 0 ⃗ \vec 0 0
Q = E [ w k w k T ] R = E [ v k v k T ] S = E [ w k v k T ] \begin{aligned} Q&= E[w_kw_k^\mathrm{T}]\\ R&= E[v_kv_k^\mathrm{T}]\\ S&= E[w_kv_k^\mathrm{T}] \end{aligned} Q R S = E [ w k w k T ] = E [ v k v k T ] = E [ w k v k T ]
注意中括号里头，两个角标必须一致，角标不一致时协方差为零（因为白噪的意思就是不同时刻互相独立）。并且进一步假定 E [ x k w k T ] = E [ x k v k T ] = 0 E[x_kw_k^\mathrm{T}]=E[x_kv_k^\mathrm{T}]=0 E [ x k w k T ] = E [ x k v k T ] = 0
考试只考没有 B B B D D D
{ x k + 1 = A x k + w k y k = C x k + v k \left\{ \begin{aligned} x_{k+1}&= Ax_k+w_k\\ y_k&= Cx_k+v_k \end{aligned} \right. { x k + 1 y k = A x k + w k = C x k + v k
# 记号定义：
状态协方差 Σ = E [ x k x k T ] \Sigma=E[x_kx_k^\mathrm{T}] Σ = E [ x k x k T ] 输出协方差 Λ i = E [ y k + i y k T ] \Lambda_i=E[y_{k+i}\,y_k^\mathrm{T}] Λ i = E [ y k + i y k T ] i i i 状态输出协方差 G = E [ x k + 1 y k T ] G=E[x_{k+1}y_k^\mathrm{T}] G = E [ x k + 1 y k T ] 这些协方差应当与 k k k
Σ = E [ x k + 1 x k + 1 T ] = E [ ( A x k + w k ) ( A x k + w k ) T ] = A E [ x k x k T ] A T + E [ w k w k T ] = A Σ A T + Q \begin{aligned} \Sigma&= E[x_{k+1}x_{k+1}^\mathrm{T}]\\ &= E\big[(Ax_k+w_k)(Ax_k+w_k)^\mathrm{T}\big]\\ &= AE[x_kx_k^\mathrm{T}]A^\mathrm{T}+E[w_kw_k^\mathrm{T}]\\ &= A\Sigma A^\mathrm{T}+Q \end{aligned} Σ = E [ x k + 1 x k + 1 T ] = E [ ( A x k + w k ) ( A x k + w k ) T ] = A E [ x k x k T ] A T + E [ w k w k T ] = A Σ A T + Q
G = E [ x k + 1 y k T ] = E [ ( A x k + w k ) ( C x k + v k ) T ] = A E [ x k x k T ] C T + E [ w k v k T ] = A Σ C T + S \begin{aligned} G&= E[x_{k+1}y_{k}^\mathrm{T}]\\ &= E\big[(Ax_k+w_k)(Cx_k+v_k)^\mathrm{T}\big]\\ &= AE[x_kx_k^\mathrm{T}]C^\mathrm{T}+E[w_kv_k^\mathrm{T}]\\ &= A\Sigma C^\mathrm{T}+S \end{aligned} G = E [ x k + 1 y k T ] = E [ ( A x k + w k ) ( C x k + v k ) T ] = A E [ x k x k T ] C T + E [ w k v k T ] = A Σ C T + S
Λ 0 = E [ y k y k T ] = E [ ( C x k + v k ) ( C x k + v k ) T ] = C E [ x k x k T ] C T + E [ v k v k T ] = C Σ C T + R \begin{aligned} \Lambda_0&= E[y_{k}y_{k}^\mathrm{T}]\\ &= E\big[(Cx_k+v_k)(Cx_k+v_k)^\mathrm{T}\big]\\ &= CE[x_kx_k^\mathrm{T}]C^\mathrm{T}+E[v_kv_k^\mathrm{T}]\\ &= C\Sigma C^\mathrm{T}+R \end{aligned} Λ 0 = E [ y k y k T ] = E [ ( C x k + v k ) ( C x k + v k ) T ] = C E [ x k x k T ] C T + E [ v k v k T ] = C Σ C T + R
现在考虑 Λ i \Lambda_i Λ i y y y C x + v Cx+v C x + v x x x A x + w Ax+w A x + w x k + 1 x_{k+1}{} x k + 1
Λ i = E [ y k + i y k T ] = E [ ( C x k + i + v k + i ) y k T ] = C E [ x k + i y k T ] + 0 = C E [ ( A x k + i − 1 + w k + i − 1 ) y k T ] = C A E [ x k + i − 1 y k T ] + 0 = C A E [ ( A x k + i − 2 + w k + i − 2 ) y k T ] = C A 2 E [ x k + i − 2 y k T ] + 0 = ⋯ = C A i − 1 E [ x k + 1 y k T ] + 0 = C A i − 1 G \begin{aligned} \Lambda_i&= E[y_{k+i}y_{k}^\mathrm{T}]\\ &= E\big[(Cx_{k+i}+v_{k+i})y_{k}^\mathrm{T}\big]=CE\big[x_{k+i}y_{k}^\mathrm{T}\big]+0\\ &= CE\big[(Ax_{k+i-1}+w_{k+i-1})y_{k}^\mathrm{T}\big]=CAE\big[x_{k+i-1}y_{k}^\mathrm{T}\big]+0\\ &= CAE\big[(Ax_{k+i-2}+w_{k+i-2})y_{k}^\mathrm{T}\big]=CA^2E\big[x_{k+i-2}y_{k}^\mathrm{T}\big]+0\\ &=\ \ \cdots\ \ =CA^{i-1}E\big[x_{k+1}y_{k}^\mathrm{T}\big]+0\\ &= CA^{i-1}G \end{aligned} Λ i = E [ y k + i y k T ] = E [ ( C x k + i + v k + i ) y k T ] = C E [ x k + i y k T ] + 0 = C E [ ( A x k + i − 1 + w k + i − 1 ) y k T ] = C A E [ x k + i − 1 y k T ] + 0 = C A E [ ( A x k + i − 2 + w k + i − 2 ) y k T ] = C A 2 E [ x k + i − 2 y k T ] + 0 = ⋯ = C A i − 1 E [ x k + 1 y k T ] + 0 = C A i − 1 G
以下仍在施工中，实在是不懂怎么写得丝滑一些，，，这段东西太几把逆天了，课件上一点前因后果都看不出来的那种，莫名其妙。我赌它不考
# 核心原理说穿了其实就是 “用样本均值代替期望”。也即，只要观测的次数 j j j
E [ X ] = 1 j ∑ k = 0 j − 1 x k E[X]=\dfrac1j\sum\limits_{k=0}^{j-1}x_k E [ X ] = j 1 k = 0 ∑ j − 1 x k
用这个思路，有如下观察：
Λ i = E [ y k + i y k T ] = 1 j ∑ k = 0 j − 1 y k + i y k T = 1 j [ y i y i + 1 ⋯ y i + j − 1 ] [ y 0 T y 1 T ⋮ y j − 1 T ] = 1 j Y i , j Y 0 , j T \begin{aligned} \Lambda_i&= E[y_{k+i}y_k^\mathrm{T}]\\ &= \dfrac1j\sum\limits_{k=0}^{j-1}y_{k+i}y_k^\mathrm{T}\\ &= \dfrac1j\begin{bmatrix}y_i&y_{i+1}&\cdots&y_{i+j-1}\end{bmatrix}\begin{bmatrix}y_0^\mathrm{T}\\y_{1}^\mathrm{T}\\\vdots\\y_{j-1}^\mathrm{T}\end{bmatrix}\\ &= \dfrac1jY_{i,\,j}Y_{0,\,j}^\mathrm{T} \end{aligned} Λ i = E [ y k + i y k T ] = j 1 k = 0 ∑ j − 1 y k + i y k T = j 1 [ y i y i + 1 ⋯ y i + j − 1 ] ⎣ ⎢ ⎢ ⎢ ⎢ ⎡ y 0 T y 1 T ⋮ y j − 1 T ⎦ ⎥ ⎥ ⎥ ⎥ ⎤ = j 1 Y i , j Y 0 , j T
其中大写 Y Y Y
此时我定义：
Φ [ A , B ] = 1 j A B T \Phi_{[A,\,B]}=\dfrac1jAB^\mathrm{T} Φ [ A , B ] = j 1 A B T
那么
Λ i = Φ [ Y i , j , Y 0 , j ] \Lambda_i=\Phi_{[Y_{i,\,j},\,Y_{0,\,j}]} Λ i = Φ [ Y i , j , Y 0 , j ]
