2.6: Mínimos cuadrados recursivos (opcional)

¿Y si los datos vienen secuencialmente? ¿Tenemos que recalcular todo cada vez que entra un nuevo punto de datos, o podemos escribir nuestra nueva estimación actualizada en términos de nuestra antigua estimación?

Considera el modelo

\[y_{i}=A_{i} x+e_{i}, \quad i=0,1, \ldots\nonumber\]

donde\(y_{i} \in \mathbf{C}^{m \times 1}, A_{i} \in \mathbf{C}^{m \times n}, x \in \mathbf{C}^{n \times 1}, \text { and } e_{i} \in \mathbf{C}^{m \times 1}\). El vector\(e_{k}\) representa el desajuste entre la medición\(y_{k}\) y el modelo para ello\(A_{k}x\), donde\(A_{k}\) se conoce y\(x\) es el vector de parámetros a estimar. En cada momento\(k\), deseamos encontrar

\[\widehat{x}_{k}=\arg \min _{x}\left(\sum_{i=1}^{k}\left(y_{i}-A_{i} x\right)_{i}^{\prime} S_{i}\left(y_{i}-A_{i} x\right)\right)=\arg \min _{x}\left(\sum_{i=1}^{k} e_{i}^{\prime} S_{i} e_{i}\right)\nonumber\]

donde\(S_{i} \in \mathbf{C}^{m \times 1}\) es una matriz hermitiana definitiva positiva de pesos, para que podamos variar la importancia de los\(e_{i}\)'s y componentes de los\(e_{i}\)'s en la determinación\(\widehat{x}_{k}\).

Para calcular\(\widehat{x}_{k+1}\), deje:

\ [\ bar {y} _ {k+1} =\ left [\ begin {array} {c}
y_ {0}\\
y_ {1}\
\\ cdot\
\ cdot\\
y_ {k+1}
\ end {array}\ derecha];\ quad\ bar {A} _ {k+1} =\ left [\ begin {array} {c}
A_ {0}\
A_ {1}\\
\ cdot\\
\ cdot\\
A_ {k+1}
\ end {array}\ derecha];\ quad\ bar {e} _ {k+1} =\ left [\ begin {array} {c}
e_ {0}\\
e_ {1}\\ cdot
\\ cdot
\\ cdot\\
e_ {k+1}
\ end {array}\ derecho]\ nonumber\]

y

\[\bar{S}_{k+1}=\operatorname{diag}\left(S_{0}, S_{1}, \ldots, S_{k+1}\right)\nonumber\]

donde\(S_{i}\) está la matriz de ponderación para\(e_{i}\).

Nuestro problema es entonces equivalente a

\[\min \left(\bar{e}_{k+1}^{\prime} \bar{S}_{k+1} \bar{e}_{k+1}\right)\nonumber\]

sujeto a:\(\bar{y}_{k+1}=\bar{A}_{k+1} x_{k+1}+\bar{e}_{k+1}\)

Por lo tanto, la solución puede escribirse como

\[\left(\bar{A}_{k+1}^{\prime} \bar{S}_{k+1} \bar{A}_{k+1}\right) \widehat{x}_{k+1}=\bar{A}_{k+1}^{\prime} \bar{S}_{k+1} \bar{y}_{k+1}\nonumber\]

o en forma de suma como

\[\left(\sum_{i=0}^{k+1} A_{i}^{\prime} S_{i} A_{i}\right) \widehat{x}_{k+1}=\sum_{i=0}^{k+1} A_{i}^{\prime} S_{i} y_{i}\nonumber\]

Definiendo

\[Q_{k+1}=\sum_{i=0}^{k+1} A_{i}^{\prime} S_{i} A_{i}\nonumber\]

podemos escribir una recursión para\(Q_{k+1}\) lo siguiente:

\[Q_{k+1}=Q_{k}+A_{k+1}^{\prime} S_{k+1} A_{k+1}\nonumber\]

Reordenando la ecuación de forma de suma para\(\widehat{x}_{k}+1\), obtenemos

\ [\ begin {alineado}
\ sombrero ancho {x} _ {k+1} &=Q_ {k+1} ^ {-1}\ left [\ left (\ sum_ {i=0} ^ {k} A_ {i} ^ {\ prime} S_ {i} A_ {i}\ right)\ sombrero ancho {x} _ _ {k} +A_ {k+1} ^ {prime\ S} _ {k+1} y_ {k+1}\ derecha]\\
&=Q_ {k+1} ^ {-1}\ izquierda [Q_ {k}\ sombrero ancho {x} _ {k} +A_ {k+1} ^ {\ prime} S_ {k+1} y_ {k+1}\ derecha]
\ fin { alineado}\ nonumber\]

Esto muestra claramente la nueva estimación como una combinación ponderada de la estimación anterior y los nuevos datos, por lo que tenemos la recursión deseada. Otra forma útil de este resultado se obtiene sustituyendo de la recursión por\(Q_{k+1}\) arriba para obtener

\[\widehat{x}_{k+1}=\widehat{x}_{k}-Q_{k+1}^{-1}\left(A_{k+1}^{\prime} S_{k+1} A_{k+1} \widehat{x}_{k}-A_{k+1}^{\prime} S_{k+1} y_{k+1}\right)\nonumber\]

lo que finalmente reduce a

\[\widehat{x}_{k+1}=\widehat{x}_{k}+\underbrace{Q_{k+1}^{-1} A_{k+1}^{\prime} S_{k+1}}_{\text {Kalman Filter Gain }} \underbrace{\left(y_{k+1}-A_{k+1} \widehat{x}_{k}\right)}_{\text {innovations }}\nonumber\]

A la cantidad\(Q_{k+1}^{-1} A_{k+1}^{\prime} S_{k+1}\) se le llama ganancia de Kalman, y\(y_{k+1}-A_{k+1} \widehat{x}_{k}\) se le llama las innovaciones, ya que compara la diferencia entre una actualización de datos y la predicción dada la última estimación.

Desafortunadamente, a medida que uno adquiere más y más datos, es decir, a medida\(k\) que crece, la ganancia de Kalman va a cero. Un punto de datos no puede avanzar mucho contra la masa de datos anteriores que ha `endurecido' la estimación. Si dejamos este estimador tal cual -sin modificación- el estimador `va a dormir' después de un tiempo, y así no se adapta bien a los cambios de parámetros. En la tarea se investiga el concepto de “memoria que se desvanece” para que el estimador no se vaya a dormir.

Un problema de implementación

Otro concepto que es importante en la implementación del algoritmo RLS es el cómputo de\(Q_{k+1}^{-1}\). Si la dimensión de\(Q_{k}\) es muy grande, el cálculo de su inverso puede ser computacionalmente costoso, por lo que a uno le gustaría tener una recursión para\(Q_{k+1}^{-1}\).

Esta recursión es fácil de obtener. Aplicando la práctica identidad matricial

\[(A+B C D)^{-1}=A^{-1}-A^{-1} B\left(D A^{-1} B+C^{-1}\right)^{-1} D A^{-1}\nonumber\]

a la recursión para\(Q_{k+1}\) rendimientos

\[Q_{k+1}^{-1}=Q_{k}^{-1}-Q_{k}^{-1} A_{k+1}^{\prime}\left(A_{k+1} Q_{k}^{-1} A_{k+1}^{\prime}+S_{k+1}^{-1}\right)^{-1} A_{k+1} Q_{k}^{-1}\nonumber\]

Al definir

\[P_{k+1}=Q_{k+1}^{-1}\nonumber\]

esto se convierte

\[P_{k+1}=P_{k}-P_{k} A_{k+1}^{\prime}\left(S_{k+1}^{-1}+A_{k+1} P_{k} A_{k+1}^{\prime}\right)^{-1} A_{k+1} P_{k}\nonumber\]

que se llama la ecuación de Riccati (tiempo discreto).

Interpretación

Tenemos\(\widehat{x}_{k}\) y\({y}_{k+1}\) disponible para calcular nuestra estimación actualizada. Interpretando\(\widehat{x}_{k}\) como medida, vemos que nuestro modelo se convierte en

\ [\ left [\ begin {array} {c}
\ anchohat {x} _ {k}\\
y_ {k+1}
\ end {array}\ right] =\ left [\ begin {array} {c}
I\\
A_ {k+1}
\ end {array}\ right] x+\ left [\ begin {array} {c}
e_ {k}\\
e_ {k+1}
\ end {matriz} \ derecha]\ nonumber\]

El criterio, entonces, por el cual elegimos\(\widehat{x}_{k+1}\) es así

\[\widehat{x}_{k+1}=\operatorname{argmin}\left(e_{k}^{\prime} Q_{k} e_{k}+e_{k+1}^{\prime} S_{k+1} e_{k+1}\right)\nonumber\]

En este contexto, se interpreta\({Q}_{k}\) como el factor de ponderación para la estimación anterior.