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?
