Sequential Parameter Estimation

Sequential Parameter Estimation#

Often it is desired that an estimate be made as observed data occur rather than waiting until all of the observations are available.

For AWGN the noise samples are independent and an estimate can be obtained in terms of the prior estimate and the current observation.

MAP Estimation#

The estimation procedure based on the current observation and the prior estimates leads to one form of a Kalman filter.

However, an intuitive derivation for a sequential estimate is developed using the estimate given before, repeated here, i.e.,

The Maximum A Posteriori (MAP) estimate is given as:

\[ \vec{\hat{\alpha}}_{MAP} = \left( H^T R_n^{-1} H + A_\alpha^{-1} \right)^{-1} \left( H^T R_n^{-1} \vec{y} + A_\alpha^{-1} \vec{\mu}_\alpha \right) \]

An alternative form of this equation is derived by introducing a zero term using \(H^T R_n^{-1} H \vec{\mu}_\alpha\) within the second parenthesis:

\[ \vec{\hat{\alpha}}_{MAP} = \left( H^T R_n^{-1} H + A_\alpha^{-1} \right)^{-1} \left( H^T R_n^{-1} \vec{y} + A_\alpha^{-1} \vec{\mu}_\alpha + H^T R_n^{-1} H \vec{\mu}_\alpha - H^T R_n^{-1} H \vec{\mu}_\alpha \right) \]

Simplifying:

\[ \vec{\hat{\alpha}}_{MAP} = \left( H^T R_n^{-1} H + A_\alpha^{-1} \right)^{-1} \left[ H^T R_n^{-1} \left( \vec{y} - H \vec{\mu}_\alpha \right) + \left( H^T R_n^{-1} H + A_\alpha^{-1} \right) \vec{\mu}_\alpha \right] \]

Finally:

\[ \vec{\hat{\alpha}}_{MAP} = \vec{\mu}_\alpha + \left( H^T R_n^{-1} H + A_\alpha^{-1} \right)^{-1} H^T R_n^{-1} \left( \vec{y} - H \vec{\mu}_\alpha \right) \]

Since \( \mathbb{E}\{\vec{y}\} = \mathbb{E}\{H \vec{\alpha}\} = H \vec{\mu}_{\alpha} \), it can be seen that the expression \( \vec{y} - H \vec{\mu}_{\alpha} \) is in fact the difference between the observed value and its expected value.

Since the estimate is unbiased, i.e., \( \mathbb{E}\{\vec{\alpha}_e\} = 0 \), the variance of the estimation error is

\[ \text{Var} \left\{ \vec{\alpha}_e \right\} = R_{\alpha_e} = \left( H^T R_n^{-1} H + A_\alpha^{-1} \right)^{-1} \]

Rewriting \(\vec{\hat{\alpha}}_{MAP}\) in compact form as

\[ \boxed{ \vec{\hat{\alpha}}_{MAP} = \vec{\mu}_\alpha + K \left( \vec{y} - H \vec{\mu}_\alpha \right) } \]

where the Kalman gain \(K\) is defined as:

\[ K = R_{\alpha_e} H^T R_n^{-1} \]

Kalman Gain#

The attained expression of \(\vec{\hat{\alpha}}_{MAP}\) above is actually a simple form of the Kalman filter where the variable \(K\) is referred to as the Kalman gain.

In this form, the computations proceed by using the \(k \times L\)-channel matrix \(H\), the \(L \times L\) initial covariance of the estimate \(A_\alpha\), and the \(k \times k\) noise covariance \(R_n\) to calculate the \(L \times L\) covariance \(R_{\alpha_e}\).

Then, the \(L \times k\) Kalman gain is computed as \(K\) above.

Sequential Estimation#

A sequential formulation is desired where the estimate \(\vec{\hat{\alpha}}_{MAP,i}\) at time \(t_i\) is determined from its prior estimate.

Before any observation occurs, the best initial estimate \(\vec{\hat{\alpha}}_{MAP_0}\) is the mean.

\[ \boxed{ \vec{\hat{\alpha}}_{MAP_0} = \vec{\mu}_\alpha } \]

with an initial covariance:

\[ R_{\alpha_{e_0}} = A_\alpha \]

When the first observation arrives, it is represented by:

\[ y_1 = H_1 \vec{\alpha} + n_1 \]

where \(H_1\) is the first row of the channel matrix \(H\) and is given by:

\[ H_1 = [h_{11}, \dots, h_{1k}] \]

The MAP estimate given by solving the differential equation in previous section can be restated as finding the value of \(\vec{\alpha}\) that maximizes \(p(y_1 | \vec{\alpha}) p(\vec{\alpha})\).

From this

\[ \vec{\hat{\alpha}}_{MAP} = \vec{\mu}_\alpha + \left( H^T R_n^{-1} H + A_\alpha^{-1} \right)^{-1} H^T R_n^{-1} \left( \vec{y} - H \vec{\mu}_\alpha \right) \]

we can compute

\[ \boxed{ \vec{\hat{\alpha}}_{MAP_1} = \vec{\hat{\alpha}}_{MAP_0} + K_1 \left( y_1 - H_1 \vec{\hat{\alpha}}_{MAP_0} \right) } \]

with:

\[ K_1 = R_{\alpha_{e_1}} H_1^T R_{n_1}^{-1} \]

and:

\[ R_{\alpha_{e_1}} = \left( H_1^T R_{n_1}^{-1} H_1 + R_{\alpha_{e_0}}^{-1} \right)^{-1} \]

Note that the subscript \(1\) has been added to indicate that the variables are treated using only the knowledge available at time instant \(t_1\). Therefore:

\[ R_{n_1} = \sigma_1^2 \]

To continue with the sequential formulation, the next step requires finding the value of \(\vec{\alpha}\) that maximizes \(p(y_1, y_2 | \vec{\alpha}) p(\vec{\alpha})\).

Since \(y_1\) and \(y_2\) are independent in this section:

\[\begin{split} \begin{align*} p(y_2, y_1 | \vec{\alpha}) p(\vec{\alpha}) &= p(y_2 | y_1, \vec{\alpha}) p(y_1 | \vec{\alpha}) p(\vec{\alpha}) \\ &= p(y_2 | \vec{\alpha}) p(\vec{\alpha} | y_1) p(y_1) \end{align*} \end{split}\]

Since the last term does not involve \(\vec{\alpha}\), it can be ignored. The pdf \(p(\vec{\alpha} | y_1)\) is Gaussian with conditional mean \(\vec{\hat{\alpha}}_1\) and conditional covariance \(R_{\alpha_{e_1}}\) after the observation \(y_1\) occurs. This procedure can be continued in a recursive fashion, where the estimate at step \(l\) is computed from the estimate at step \(l-1\), i.e.,

\[ \boxed{ \vec{\hat{\alpha}}_{MAP_l} = \vec{\hat{\alpha}}_{MAP_{l-1}} + K_l \left( y_l - H_l \vec{\hat{\alpha}}_{MAP_{l-1}} \right) } \]

with:

\[ K_l = R_{\alpha_{e_l}} H_l^T R_{n_l}^{-1} \]

and:

\[ R_{\alpha_{e_l}} = \left( H_l^T R_{n_l}^{-1} H_l + R_{\alpha_{e_{l-1}}}^{-1} \right)^{-1} \]

where:

\[ H_l = [h_{l1} \dots h_{lL}] \]

and:

\[ R_{n_l} = \sigma_l^2 \]

Example#

This example [B2 Ex. 12.6] repeats Example 12.5, where sequential processing is used and \(\rho = 0\).

Since \(\alpha\) is Gaussian with zero mean and variance \(\sigma_\alpha^2\)

\[ \boxed{ \vec{\hat{\alpha}}_{MAP_0} = \vec{\mu}_\alpha = 0 } \]

The interim constant \(a\) defined as:

\[ a = \frac{\sigma_\alpha^2}{\sigma^2 + \sigma_\alpha^2}, \]

We have

\[\begin{split} \begin{align*} R_{\alpha_{e_1}} &= \left( H_1^T R_{n_1}^{-1} H_1 + R_{a_{e_0}}^{-1} \right)^{-1} \\ &= \left( (1)\left(\frac{1}{\sigma^2}\right)(1) + \frac{1}{\sigma_\alpha^2} \right)^{-1} \\ &= \frac{\sigma_\alpha^2 \sigma^2}{\sigma^2 + \sigma_\alpha^2} \\ &= a \sigma^2 \end{align*} \end{split}\]

The Kalman gain is:

\[\begin{split} \begin{align*} K_1 &= R_{\alpha_{e_1}} H^T R_n^{-1} \\ &= a\sigma^2(1) \frac{1}{\sigma^2} \\ &= a \end{align*} \end{split}\]

Knowledge of the observation \(y_1\) and initial MAP estimate \(\hat{\alpha}_{MAP_0} = 0\) allows the first MAP estimate to be computed from Eq. (12.78) as:

\[ \boxed{ \hat{\alpha}_{MAP_1} = K_1 y_1 = a y_1 } \]

In the next iteration, the covariance \(R_{\alpha_{e_2}}\) is determined from Eq. (12.80) as:

\[\begin{split} \begin{align*} R_{\alpha_{e_2}} &= \left( H_2^T R_{n_2}^{-1} H_2 + R_{\alpha_{e_1}}^{-1} \right)^{-1} \\ &= \left( (1)\left(\frac{1}{\sigma^2}\right)(1) + \frac{1}{a \sigma^2} \right)^{-1} \\ &= \frac{a \sigma^2}{1 + a} \end{align*} \end{split}\]

The next Kalman gain is then:

\[\begin{split} \begin{align*} K_2 &= R_{\alpha_{e_2}} H_2^T R_{n_2}^{-1} \\ &= \frac{a \sigma^2}{1 + a} (1) \frac{1}{\sigma^2} \\ &= \frac{a}{1 + a} \end{align*} \end{split}\]

The second MAP estimate obtained from Eq. (12.78) is:

\[\begin{split} \begin{align*} \hat{\alpha}_{MAP_2} &= \hat{\alpha}_{MAP_1} + K_2 \left(y_2 - H_2 \hat{\alpha}_{MAP_1}\right) \\ &= a y_1 + \frac{a}{1 + a} \left( y_2 - a y_1 \right) \\ &= \frac{a}{1 + a} (y_1 + y_2) \\ &= \frac{\sigma_\alpha^2}{\sigma^2 + 2 \sigma_\alpha^2} (y_1 + y_2) \end{align*} \end{split}\]

Thus,

\[ \boxed{ \hat{\alpha}_{MAP_2} = \frac{\sigma_\alpha^2}{\sigma^2 + 2 \sigma_\alpha^2} (y_1 + y_2) } \]

Note that the estimate is the same as that obtained in Example 12.5 with \(\rho = 0\).

Matlab Example#