# Conditional Multivariate Gaussian Identity

I'm trying to verify the form of a multivariate Gaussian provided in a paper I'm reading. It should be pretty elementary.

Let $$Y=X+\varepsilon$$ where $$X\sim N(0,C)$$ and $$\varepsilon\sim N(0,\sigma^2\mathbf{I})$$. The authors then claim that $$X|Y,C,\sigma^2 \sim N(\mu,\Sigma),$$ where $$\mu := C(C+\sigma^2\mathbf I)^{-1}Y\\ \Sigma:=\sigma^2C(C+\sigma^2\mathbf I)^{-1}.$$ My first thought was to consider the joint distribution $$\begin{pmatrix} X\\ Y \end{pmatrix}\sim N\Big(\begin{pmatrix} 0\\ 0 \end{pmatrix},\begin{pmatrix} C & C\\ C^\top & \sigma^2\mathbf I+C \end{pmatrix}\Big)$$ and apply the conditional Gaussian identities. Unfortunately this approach gives me the right $$\mu$$, but I can't see how their form of $$\Sigma$$ comes about. Any thoughts?

Since $$\begin{pmatrix} X\\ \epsilon \end{pmatrix}\sim N\Big(\begin{pmatrix} 0\\ 0 \end{pmatrix},\begin{pmatrix} C & \mathbf O\\ \mathbf O & \sigma^2\mathbf I \end{pmatrix}\Big)$$ and $$\begin{pmatrix} X\\ Y \end{pmatrix} = \begin{pmatrix} \mathbf 1^\text{T} &\mathbf 0^\text{T} \\ \mathbf 1^\text{T} &\mathbf 1^\text{T} \end{pmatrix}\begin{pmatrix} X\\ \epsilon \end{pmatrix}$$ the distribution of $$\begin{pmatrix} X\\ Y \end{pmatrix}$$ is $$\begin{pmatrix} X\\ Y \end{pmatrix}\sim N\Big(\begin{pmatrix} 0\\ 0 \end{pmatrix},\underbrace{\begin{pmatrix} \mathbf 1^\text{T} &\mathbf 0^\text{T} \\ \mathbf 1^\text{T} &\mathbf 1^\text{T} \end{pmatrix}\begin{pmatrix} C & \mathbf O\\ \mathbf O & \sigma^2\mathbf I \end{pmatrix}\begin{pmatrix} \mathbf 1 &\mathbf 1 \\ \mathbf 0 &\mathbf 1 \end{pmatrix}}_{\begin{pmatrix} C & C\\ C & \sigma^2\mathbf I \end{pmatrix} }\Big)$$ indeed. With $$\mathbb E[X|Y] = 0 + C (C+\sigma^2\mathbf I)^{-1} Y$$ and $$\text{var}(X|Y) = C - C (C+\sigma^2I)^{-1} C$$ Applying the Woodbury matrix inversion lemma $$(A+B)^{-1}=A^{-1}-A^{-1}(B^{-1}+A^{-1})^{-1}A^{-1}$$ one gets that \begin{align*} C - C (C+\sigma^2I)^{-1} C &= C - C (C^{-1}- C^{-1}(C^{-1}+\sigma^{-2}\mathbf I)^{-1}C^{-1})C\\ &= C - C +(C^{-1}+\sigma^{-2}\mathbf I)^{-1}\\ &= (C^{-1}\mathbf I+\sigma^{-2}C^{-1}C)^{-1}\\ &= \sigma^2 C (\sigma^2\mathbf I+C)^{-1} \end{align*} The apparent lack of symmetry in the expression may sound suspicious but actually$$C (\sigma^2\mathbf I+C)^{-1} = (\sigma^2\mathbf I+C)^{-1} C$$