# How to perform joint inference on multivariate normal variables?

Suppose I have the following model:

\begin{aligned} \text C &\sim \mathcal N \left(\mu, \delta^2\right) \\ \forall i: \text L_i | \text C = c &\sim \mathcal N \left(c, \lambda_i^2 \right) \\ \forall i: \text P_i | \text L_i = l_i &\sim \mathcal N \left(l_i, \sigma_i^2 \right) \end{aligned}

In the above, $$\text C$$ and the $$\text L_i$$ are the hidden variables I'm reasoning about, and the $$\text P_i$$ are the observed variables that serve as evidence. Then:

\begin{aligned} p \left(\text C, \vec {\text L} \middle | \vec {\text P} = \vec p \right) &\propto p \left(\vec {\text P} = \vec p \middle | \text C, \vec {\text L} \right) p \left(\vec {\text L} \middle | \text C \right) p \left(\text C \right) \\ &= p \left(\text C \right) \prod\limits_{i=1}^n p\left(\text P_i = p_i \middle| \text L_i\right) p\left(\text L_i \middle | \text C\right) \\ &= \mathcal N \left(C \middle |\mu, \delta^2\right) \prod\limits_{i=1}^n \mathcal N \left(p_i \middle |\text L_i, \sigma_i^2 \right) \mathcal N \left(\text L_i \middle| \text C, \lambda_i^2 \right) \\ &\propto \mathcal N \left(C \middle |\mu, \delta^2\right) \prod\limits_{i=1}^n \mathcal N \left(\text L_i \middle| \frac {\lambda_i^{-2} \text C + \sigma_i^{-2} p_i}{\lambda_i^{-2} + \sigma_i^{-2}}, \left(\lambda_i^{-2} + \sigma_i^{-2} \right)^{-1} \right) \\ \end{aligned}

Presumably the above is a multivariate normal distribution, with some posterior mean $$\vec \mu$$ and some posterior covariance matrix $$\mathbf \Sigma$$. What are they?

A possible way to figure this out would be to describe the original problem vectorially, but I have no idea how to do that either. I think, since in the multivariate normal definition the mean vector seems to be defined as the vector of the means of each variable unconditionally on all others, we'd have:

\begin{aligned} \text C, \vec {\text L} &\sim \mathcal N \left( \begin{bmatrix} \mu \\ \mu \\ \vdots \\ \mu\end{bmatrix}, \begin{bmatrix} \delta^2 & \delta^2 & \cdots & \delta^2 \\ \delta^2 & \delta^2 + \lambda_1^2 & & \delta^2 \\ \vdots &&\ddots & \\ \delta^2 & \delta^2 & & \delta^2 + \lambda_n^2 \end{bmatrix} \right) \end{aligned}

Although I'm not sure the first row and column have the right values.

If the above is correct, then:

\begin{aligned} \vec {\text P} | \vec {\text L} = \vec l &\sim \mathcal N \left(\vec l, \begin{bmatrix} \sigma_1^2 && 0\\ & \ddots &\\ 0 && \sigma_n^2 \end{bmatrix} \right) \\ \vec {\text L} | \text C = c &\sim \mathcal N \left( \begin{bmatrix}c \\ \vdots \\ c \end{bmatrix}, \begin{bmatrix} \lambda_1^2 && 0 \\ & \ddots & \\ 0 && \lambda_n^2 \end{bmatrix} \right) \end{aligned}

So for the inference step:

\begin{aligned} p \left(\text C, \vec {\text L} \middle| \vec {\text P} = \vec p \right) &\propto \mathcal N \left(\text C \middle |\mu, \delta^2\right) \mathcal N \left(\vec p \middle | \vec {\text L}, \begin{bmatrix} \sigma_1^2 && 0\\ & \ddots &\\ 0 && \sigma_n^2 \end{bmatrix} \right) \mathcal N \left(\vec {\text L} \middle| \begin{bmatrix}\text C \\ \vdots \\ \text C \end{bmatrix}, \begin{bmatrix} \lambda_1^2 && 0 \\ & \ddots & \\ 0 && \lambda_n^2 \end{bmatrix} \right) \\ &= \mathcal N \left(\text C \middle |\mu, \delta^2\right) \mathcal N \left(\vec p \middle | \vec {\text L}, \mathbf \Sigma \right) \mathcal N \left(\vec {\text L} \middle |\vec {\text C}, \mathbf \Lambda \right) \\ &\propto \mathcal N \left(\text C \middle |\mu, \delta^2\right) \mathcal N \left(\vec {\text L} \middle | \left( \mathbf \Sigma ^ {-1} + \mathbf \Lambda ^ {-1}\right) ^ {-1} \left( \mathbf \Sigma ^ {-1} \vec p + \mathbf \Lambda ^ {-1} \vec {\text C} \right), \left( \mathbf \Sigma ^ {-1} + \mathbf \Lambda ^ {-1} \right) ^ {-1} \right) \end{aligned}

Except the above is just... the same thing, but written vectorially. So I'm stumped. What's next, how do I find the posterior mean and covariance matrix for the whole thing?

• Since you have two latent variables, why do you want to write the posterior as single multivariate normal distribution? Aug 23, 2019 at 11:21
• Two reasons: one, I want the mean vector more than anything, I wanna know what the expected value of each variable there is; two, this is one slice of a recurrent Bayesian network and it's gonna be way easier to update on the next step if the previous step comes out in a standard format. Aug 23, 2019 at 13:29
• Just to clarify, does this mean you want the posterior of a vector $Z = [C, \vec{L}]$ containing the concatenation of the latent variables? Aug 23, 2019 at 14:04
• Yep, exactly, both the posterior mean and posterior covariance matrix, and also ideally how to even....... get them from what I have, because so far it seems very Mysterious to me how to perform update on data of lower dimensionality than the prior. Aug 23, 2019 at 14:12

Okay, so let me try to continue the vectorial thoughts. If I define $$\vec {\text {CL}}$$ as the joint vector of my latent variables, then there exists a matrix $$\mathbf X$$ such that $$\vec L = \mathbf X \vec{\text {CL}}$$:
$$\mathbf X = \begin{bmatrix} 0 & 1 & & 0 \\ \vdots & & \ddots & \\ 0 & 0 & & 1 \end{bmatrix}$$
\begin{aligned} p \left(\vec {\text {CL}} \middle| \vec {\text P} = \vec p \right) &= \frac {p \left(\vec {\text P} = \vec p \middle| \vec {\text {CL}} \right)p \left(\vec {\text {CL}} \right)}{p \left(\vec {\text P} = \vec p \right)} \\ &\propto \mathcal N \left(\vec p \middle | \mathbf X \vec{\text {CL}}, \mathbf \Sigma\right) \mathcal N \left(\vec{\text {CL}} \middle | \begin{bmatrix} \mu \\ \mu \\ \vdots \\ \mu\end{bmatrix}, \begin{bmatrix} \delta^2 & \delta^2 & \cdots & \delta^2 \\ \delta^2 & \delta^2 + \lambda_1^2 & & \delta^2 \\ \vdots &&\ddots & \\ \delta^2 & \delta^2 & & \delta^2 + \lambda_n^2 \end{bmatrix} \right) \\ &= \mathcal N \left(\vec p \middle | \mathbf X \vec{\text {CL}}, \mathbf \Sigma\right) \mathcal N \left(\vec{\text {CL}} \middle | \vec \mu, \mathbf \Delta \right) \\ &\propto \exp {\left(\left( \vec p - \mathbf X \vec{\text {CL}} \right)^\top \mathbf \Sigma^{-1} \left( \vec p - \mathbf X \vec{\text {CL}} \right) + \left(\vec{\text {CL}} - \vec\mu \right)^\top \mathbf \Delta^{-1} \left(\vec{\text {CL}} - \vec\mu \right) \right)}^{-\frac 1 2} \\ &= \exp {\left( \begin{aligned} &\left(\vec p ^\top - \vec{\text {CL}}^\top \mathbf X ^\top \right) \left(\mathbf \Sigma^{-1} \vec p - \mathbf \Sigma^{-1} \mathbf X \vec{\text {CL}} \right) + \\ &\left(\vec{\text {CL}}^\top - \vec\mu^\top \right) \left(\mathbf \Delta^{-1} \vec{\text {CL}} - \mathbf \Delta^{-1} \vec\mu \right) \end{aligned} \right)}^{- \frac 1 2} \\ &= \exp {\left( \begin{aligned} &\vec p ^\top \mathbf \Sigma^{-1} \vec p - \vec p ^\top \mathbf \Sigma^{-1} \mathbf X \vec {\text {CL}} - \vec {\text {CL}} ^ \top \mathbf X ^\top \mathbf \Sigma^{-1}\vec p + \vec {\text {CL}} ^ \top \mathbf X ^\top \mathbf \Sigma^{-1} \mathbf X \vec{\text {CL}} + \\ & \vec{ \text {CL}} ^\top \mathbf \Delta ^{-1} \vec{\text {CL}} - \vec{\text {CL}} ^\top \mathbf \Delta ^{-1} \vec\mu - \vec \mu ^\top \mathbf \Delta^{-1} \vec{\text {CL}} + \vec \mu ^\top \mathbf \Delta^{-1} \vec \mu \end{aligned} \right)}^{- \frac 1 2} \\ &= \exp {\left( \begin{aligned} &\vec{\text {CL}}^\top \left(\mathbf X^\top \mathbf \Sigma^{-1} \mathbf X + \mathbf \Delta^{-1} \right) \vec {\text {CL}} \\ -&\vec{\text {CL}}^\top \left( \mathbf X^\top \mathbf \Sigma^{-1} \vec p + \mathbf \Delta^{-1} \vec \mu \right) - \left(\vec p ^\top \mathbf \Sigma^{-1} \mathbf X + \vec \mu ^\top\mathbf \Delta^{-1} \right) \vec{\text {CL}} \\ + & \vec p^\top \mathbf \Sigma^{-1} \vec p + \vec \mu ^\top \mathbf \Delta^{-1} \vec \mu \end{aligned} \right)}^{- \frac 1 2} \\ &\propto \mathcal N \left( \vec {\text {CL}} \middle | \left(\mathbf X^\top \mathbf \Sigma^{-1} \mathbf X + \mathbf \Delta^{-1} \right) ^{-1} \left( \mathbf X^\top \mathbf \Sigma^{-1} \vec p + \mathbf \Delta^{-1} \vec \mu \right), \left(\mathbf X^\top \mathbf \Sigma^{-1} \mathbf X + \mathbf \Delta^{-1} \right) ^{-1} \right) \end{aligned}
So, assuming I got my definition of $$\mathbf \Delta$$ and $$\vec \mu$$ right, this should be the answer I'm looking for.