Joint posterior distribution of differences I have data $x_1,...,x_n$, $y_1,...,y_m$ and $z_1,...,z_p$ where
$$x_1,...,x_n\sim N(\mu_x,\sigma^2_x)$$
and
$$y_1,...,y_m\sim N(\mu_y,\sigma^2_y)$$
and
$$z_1,...,z_p\sim N(\mu_z,\sigma^2_z)$$
Now let's assume I want to take a Bayesian approach and place the following priors: $p(\mu_x,\sigma^2_x)\propto (\sigma^2_x)^{-1}$, $p(\mu_y,\sigma^2_y)\propto (\sigma^2_y)^{-1}$, and $p(\mu_z,\sigma^2_z)\propto (\sigma^2_z)^{-1}$. Given these priors, I know what the posterior distribution is, but more importantly, I know that the conditional marginal distributions are
$$\mu_x|\sigma^2_x,x_1,...,x_n\sim N(\bar{x},\sigma^2_x/n)$$ and
$$\mu_y|\sigma^2_y, y_1,...,y_m\sim N(\bar{y},\sigma^2_y/m)$$ and
$$\mu_z|\sigma^2_z,z_1,...,z_n\sim N(\bar{z},\sigma^2_z/p)$$
where $\bar x$ is the average of the $x$'s. Similarly, for the case of the $y$'s and $z$'s.
I am interested in deriving the joint distribution of $\mu_x-\mu_y$, $\mu_y-\mu_z$, $\mu_y-\mu_z$ . Does the following approach make sense?
We know that the conditional posterior distribution of $\mu_x$ is
$$p(\mu_x|\sigma^2_x,x_1,...,x_n) = N\left(\bar x, \frac{\sigma^2_x}{n}\right)$$
and similarly for $p(\mu_y|\sigma^2_y,y_1,...,y_m)$ and $p(\mu_z|\sigma^2_z,z_1,...,z_p)$.
Now, define
\begin{align}
    \Delta :=  \begin{bmatrix}
        \delta_{xy}\\\\
        \delta_{xz}\\\\
        \delta_{yz}
    \end{bmatrix} = 
    \begin{bmatrix}
        \mu_x-\mu_y\\\\
        \mu_x - \mu_z\\\\
        \mu_y - \mu_z
    \end{bmatrix} = 
    \begin{bmatrix}
        1 & -1 & 0\\\\
        1 & 0 & -1\\\\
        0 & 1 & -1
    \end{bmatrix}\begin{bmatrix}
        \mu_x\\\\
        \mu_y\\\\
        \mu_z
    \end{bmatrix} =: 
    A\mu,
\end{align}
then we have that $$\Delta|\sigma^2_x,\sigma^2_y,\sigma^2_z, x_1,...,x_n, y_1,...,y_m, z_1,...,z_p \sim N_3\left(A \mathbb{E}(\mu), A\mathbb{C}\text{ov}(\mu)A^T\right)$$
where
\begin{align*}
    \mathbb{E}(\mu) = \begin{bmatrix}
        \bar{x}\\\\
        \bar{y}\\\\
        \bar{z}
    \end{bmatrix} \text{    and     } \mathbb{C}\text{ov}(\mu) = \begin{bmatrix}
        \sigma^2_x/n & 0 & 0 \\\\
         0 & \sigma^2_y/m & 0 \\\\
         0 & 0 & \sigma^2_z/p
    \end{bmatrix}.
\end{align*}
So my ultimate questions are the following:

*

*Would the distribution above for $\Delta$ be the correct joint conditional posterior distribution of the differences?


*And if so, would the appropriate strategy to obtain posterior samples from $\Delta$ be to first samples $\sigma^2_x, \sigma^2_y,$ and $\sigma^2_z$ from their joint posterior distribution?


*And if so for number 2) does the joint posterior distribution of $\sigma^2_x, \sigma^2_y,$ and $\sigma^2_z$ have a closed form?

*

*The conditional analysis is correct (assuming posterior conditional independence between $\mu_x,\mu_y,\mu_z$), as already indicated in the answer to your earlier question!


*Simulating $σ_x$, $σ_y$, and $σ_z$ from their joint (marginal) posterior distribution, then $\mu_x,\mu_y,\mu_z$ from the corresponding conditional posterior is correct. However, if simulating from the (marginal) posterior proves complex, a Gibbs sampling approach may be preferable. But...


*...since
$$\int_{-\infty}^{+\infty} \exp\left\{-\left[n(\bar x_n-\mu_x)^2+\sum_{i=1}^n(x_i-\bar x_n)^2\right]\big/2\sigma^2_x\right\}\,\frac{\text d\mu_x}{\sigma_x^{n+2}}\\\propto\exp\left\{-\sum_{i=1}^n(x_i-\bar x_n)^2\big/2\sigma^2_x\right\}\sigma_x^{-n-1}$$
and the equivalent closed forms for $σ^2_y$, and $σ^2_z$, the marginal posterior on $(\sigma^2_x,\sigma^2_y,\sigma^2_z)$ is available.
