# Derivation of the LMMSE (Linear Minimum Mean Squared Error) Estimate and the MMSE Under Gaussian Prior

I am learning estimation theory through Steven M. Kay - Fundamentals of Statistical Signal Processing, Volume 1: Estimation Theory. In the Chapter 12 (Linear Bayesian Estimator), Theorem 12.1 (Bayesian Gauss-Markov Theorem) gives the LMMSE estimation of the signal based on the linear noisy measurement under the Gaussian prior assumption:

If the data are described by the Bayesian linear model form $$$$\boldsymbol{x}=\boldsymbol{H\theta}+\boldsymbol{w} \tag{12.25}$$$$ where $$\boldsymbol{x}$$ is an $$N \times 1$$ data vector, $$\boldsymbol{H}$$ is a known $$N\times p$$ observation matrix, $$\boldsymbol{\theta}$$ is a $$p \times 1$$ random vector of parameters whose realization is to be estimated and has mean $$E(\boldsymbol{\theta})$$ and covariance matrix $$\boldsymbol{C}_{\theta\theta}$$, and $$\boldsymbol{w}$$ is an $$N \times 1$$ random vector with zero mean and covariance matrix $$\boldsymbol{C}_w$$ and is uncorrelated with $$\boldsymbol{\theta}$$ (the joint PDF $$p(\boldsymbol{w},\boldsymbol{\theta})$$ is otherwise arbitrary), then the LMMSE estimator of $$\boldsymbol{\theta}$$ is \begin{align} \hat{\boldsymbol{\theta}} & = E(\boldsymbol{\theta})+\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T(\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T+\boldsymbol{C}_w)^{-1}(\boldsymbol{x}-\boldsymbol{H}E(\boldsymbol{\theta})) \tag{12.26} \\ & = E(\boldsymbol{\theta})+(\boldsymbol{C}_{\theta\theta}^{-1}+\boldsymbol{H}^T\boldsymbol{C}_w^{-1}\boldsymbol{H})^{-1}\boldsymbol{H}^T\boldsymbol{C}_w^{-1}(\boldsymbol{x}-\boldsymbol{H}E(\boldsymbol{\theta})) \tag{12.27} \end{align} The performance of the estimatior is measured by the error $$\boldsymbol{\epsilon}=\boldsymbol{\theta}-\hat{\boldsymbol{\theta}}$$ whose mean is zero and whose covariance matrix is \begin{align} \boldsymbol{C}_\boldsymbol{\epsilon} &= E_{\boldsymbol{x},\boldsymbol{\theta}}(\boldsymbol{\epsilon}\boldsymbol{\epsilon}^T) \\ & = \boldsymbol{C}_{\theta\theta} - \boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T(\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T+\boldsymbol{C}_w)^{-1}\boldsymbol{H}\boldsymbol{C}_{\theta\theta} \tag{12.28} \\ & = (\boldsymbol{C}_{\theta\theta}^{-1}+\boldsymbol{H}^T\boldsymbol{C}_w^{-1}\boldsymbol{H})^{-1} \tag{12.29} \end{align}

Since the prior of $$\boldsymbol{\theta}$$ is Gaussian, the LMMSE estimate $$\hat{\boldsymbol{\theta}}_{LMMSE}$$ is equivalent to the MMSE estimate $$\hat{\boldsymbol{\theta}}_{MMSE}$$, and $$\hat{\boldsymbol{\theta}}_{MMSE}$$ is equal to the posterior mearn $$E(\boldsymbol{\theta}|\boldsymbol{x})$$. Since the prior and likelihood are both Gaussian, the posterior distribution $$p(\boldsymbol{\theta}|\boldsymbol{x})$$ is also Gaussian.

Here I am trying to derive $$\hat{\boldsymbol{\theta}}_{MMSE}$$ and $$\boldsymbol{C}_\boldsymbol{\epsilon}$$ from the perspective of PDF multiplication, that is, calculate $$p(\boldsymbol{\theta}|\boldsymbol{x}) \propto p(\boldsymbol{x}|\boldsymbol{\theta})p(\boldsymbol{\theta})=\mathcal{N}(\boldsymbol{x};\boldsymbol{H\theta},\boldsymbol{C}_{w})\mathcal{N}(\boldsymbol{\theta};E(\boldsymbol{\theta}),\boldsymbol{C}_{\theta\theta})$$, and formulate the quadratic and firse-order terms of $$\boldsymbol{\theta}$$ at the exponential to form a Gaussian PDF. The covariance matrix of $$p(\boldsymbol{\theta}|\boldsymbol{x})$$ I got matches 12.29, but the posterior mean is the following form: $$$$E(\boldsymbol{\theta}|\boldsymbol{x}) = \boldsymbol{C}_{\boldsymbol{\epsilon}}(\boldsymbol{H}^T\boldsymbol{C}_w^{-1} \boldsymbol{x}+\boldsymbol{C}_{\theta\theta}^{-1}E(\boldsymbol{\theta})) \tag{q1}$$$$

So my question is, is the posterior mean I got in q1 equal to the $$\hat{\boldsymbol{\theta}}$$ given in 12.26 and 12.27? If so, how can I reach that?

By the way, I can't find the way from 12.26 to 12.27 (12.28 to 12.29 either). So can someone give me a hint?

• Hi: The derivation for this is the same as that for the bayesian ( standard) kalman filter. So, if you have or can get your hands on "bayesian forecasting and dynamic linear models" by west and harrison, that will provide all of the gory details. One needs to use the conditional distribution of the multivariate gaussian, rather than the joint. Aug 14 '20 at 2:11
• Are you after a correct derivation of the MMSE for the Gaussian case?
– Royi
Aug 14 '20 at 9:53
• @Royi I think so, I've check it many times.. Aug 14 '20 at 14:51
• @markleeds That's true a conditional distribution should be used, but $p(\theta|x)$ is a function of $\theta$. I think I can ignore the denominator $p(x)$ when formulating the quadratic form of the numerator. Aug 14 '20 at 14:56
• @McZhang: I don't have time to look at it but I imagine that whatever Royi derived below is correct. At a glance, it looks like he may have left out some steps and I remember West and Harrison providing every single step so it would be good to find their textbook and look at that. I'm sorry that I can't be more helpful. Aug 14 '20 at 16:50

In the past I derived it as following:

It is a little different approach.

If it answers your question I will rewrite it in a proper LaTeX.

Regarding your question about the steps in the derivation you presented, it is using the Woodbury Matrix Identity (Both 12.26 to 12.27 and 12.28 to 12.29).

In the answer above you may see some related derivations.

• Hi Royi, thanks for your reply. I still have the following questions. 1) The results you provide matches 12.26 and 12.28, but I can't see the relationship between 12.26 and 12.27, 12.28 and 12.29 either. 2) In fact I use the similar method as yours. I formulate the quadratic form of the product of $N(x;H\theta,C_w)$ and $N(\theta;E(\theta),C_{\theta \theta})$ (without using forms $C_{xx}, C_{xy}, T_{xy}$ ...), but I get the result as in q1. In fact $T_{xx}^{-1}$ is exactly $C_{\epsilon}$ in q1. The question is whether $T_{xy}\hat{y}$ equals the term within the comma of q1? Aug 14 '20 at 14:49
• The steps there are using the Woodbury Matrix Identity. I added that to my answer.
– Royi
Aug 14 '20 at 15:11
• I see 12.29 to 12.28 is straightforward with Woodbury matrix identity, but I still can not arrive 12.26 after plugging 12.28 into 12.27 .. Can you show more details? I can feel it is already very close Aug 14 '20 at 16:17
• 12.28 and 12.29 are about the error, not the predictor. Why do you tie them to 12.26?
– Royi
Aug 14 '20 at 17:19
• Since 12.27 contains the MSE matrix $(\boldsymbol{C}_{\theta\theta}^{-1}+\boldsymbol{H}^T\boldsymbol{C}_w^{-1}\boldsymbol{H})^{-1}$, which is 12.29. Aug 14 '20 at 17:40

With the help of @Royi and @markleeds, I have found the answer is Yes, that q1 is consistent with 12.26 and 12.27. The key to see this is the Woodbury Matrix Identity.

12.29 to 12.28 is straightforward with the Woodbury Matrix Identity.

From 12.27 to 12.26: \begin{align} \hat{\boldsymbol{\theta}} & = E(\boldsymbol{\theta})+(\boldsymbol{C}_{\theta\theta}^{-1}+\boldsymbol{H}^T\boldsymbol{C}_w^{-1}\boldsymbol{H})^{-1}\boldsymbol{H}^T\boldsymbol{C}_w^{-1}(\boldsymbol{x}-\boldsymbol{H}E(\boldsymbol{\theta})) \tag{12.27} \\ & = E(\boldsymbol{\theta}) + (\boldsymbol{C}_{\theta\theta} - \boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T(\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T+\boldsymbol{C}_w)^{-1}\boldsymbol{H}\boldsymbol{C}_{\theta\theta} ) \boldsymbol{H}^T\boldsymbol{C}_w^{-1}(\boldsymbol{x}-\boldsymbol{H}E(\boldsymbol{\theta})) \\ &= E(\boldsymbol{\theta}) +\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T[\boldsymbol{I}-(\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T+\boldsymbol{C}_w)^{-1}\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T]\boldsymbol{C}_{w}^{-1}(\boldsymbol{x}-\boldsymbol{H}E(\boldsymbol{\theta})) \\ & = E(\boldsymbol{\theta})+\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T(\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T+\boldsymbol{C}_w)^{-1}[\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T+\boldsymbol{C}_w-\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T] \boldsymbol{C}_{w}^{-1}(\boldsymbol{x}-\boldsymbol{H}E(\boldsymbol{\theta})) \\ & = E(\boldsymbol{\theta})+\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T(\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T+\boldsymbol{C}_w)^{-1}(\boldsymbol{x}-\boldsymbol{H}E(\boldsymbol{\theta})) \tag{12.26} \\ \end{align}

From q1 to 12.26: \begin{align} E(\boldsymbol{\theta}|\boldsymbol{x}) &= \boldsymbol{C}_{\boldsymbol{\epsilon}}(\boldsymbol{H}^T\boldsymbol{C}_w^{-1} \boldsymbol{x}+\boldsymbol{C}_{\theta\theta}^{-1}E(\boldsymbol{\theta})) \tag{q1} \\ &= (\boldsymbol{C}_{\theta\theta} - \boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T(\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T+\boldsymbol{C}_w)^{-1}\boldsymbol{H}\boldsymbol{C}_{\theta\theta} )(\boldsymbol{H}^T\boldsymbol{C}_w^{-1} \boldsymbol{x}+\boldsymbol{C}_{\theta\theta}^{-1}E(\boldsymbol{\theta})) \\ &= E(\boldsymbol{\theta}) - \boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T(\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T+\boldsymbol{C}_w)^{-1}\boldsymbol{H}E(\boldsymbol{\theta})\\ & \quad + \boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T\boldsymbol{C}_w^{-1} \boldsymbol{x} - \boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T(\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T+\boldsymbol{C}_w)^{-1}\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T\boldsymbol{C}_w^{-1} \boldsymbol{x} \\ &= E(\boldsymbol{\theta}) - \boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T(\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T+\boldsymbol{C}_w)^{-1}\boldsymbol{H}E(\boldsymbol{\theta})\\ & \quad + \boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T [\boldsymbol{I}-(\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T+\boldsymbol{C}_w)^{-1}\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T]\boldsymbol{C}_w^{-1} \boldsymbol{x} \\ &= E(\boldsymbol{\theta}) - \boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T(\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T+\boldsymbol{C}_w)^{-1}\boldsymbol{H}E(\boldsymbol{\theta})\\ & \quad + \boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T (\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T+\boldsymbol{C}_w)^{-1} [\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T+\boldsymbol{C}_w-\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T]\boldsymbol{C}_w^{-1} \boldsymbol{x} \\ & = E(\boldsymbol{\theta}) - \boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T(\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T+\boldsymbol{C}_w)^{-1}\boldsymbol{H}E(\boldsymbol{\theta}) + \boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T (\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T+\boldsymbol{C}_w)^{-1} \boldsymbol{x} \\ & = E(\boldsymbol{\theta})+\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T(\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T+\boldsymbol{C}_w)^{-1}(\boldsymbol{x}-\boldsymbol{H}E(\boldsymbol{\theta})) \tag{12.26} \end{align}

• I'm glad that whatever I said, helped some. You definitely had a a Kalman Filter in disguise there. Still, if you're into bayesian KF's, the West and Harrison yellow text used to be the bible. Don't know if it still is. All the best. Aug 15 '20 at 14:25