First of all, lets simplify things. assume, that we have a zero-mean Gaussian Markov field with an unknown covariance matrix $\Sigma$ $$X\mid \Sigma\sim N\left ( 0,\Sigma \right )$$
Conditional independence of two variables $X_{i}$ and $X_{j}$, given the rest, is garantued iff the $(i,j)$ element of the precision matrix (sometimes refered as concentraition matrix) $Q=\Sigma^{-1}$ is equal to zero.
This precision matrix has to be invertable and positive definite. So, in order to have a posterior samples of $\Sigma$ or $Q$ we have to have a way to randomly draw such structured matrices.
Håvard Rue have advised to look for information about so called G-Wishart distribution and a. Lenkoski work on it. The expression for this distribution is as follows (up to a normalizing constant):
$$f\left ( K|D \right )\propto \left ( det(K)) \right )^\frac{\delta-2}{2}exp\left \{ -tr\left ( KD \right ) \right \}$$
K is a matrix with the structure we need and is positive definite. Now lets put G-Wishart distribution as prior on $Q$, so that posterior will be also G-Wishart distribution:
$$Q|X\sim W_{G}\left ( \delta +n, D+\sum_{k=1}^{n}X_{k}^{T}X_{k} \right )$$
To obtain random samples I have employed the iterated proportional scaling algorithm and Gibbs sampler as described in the paper of Lentoski and Dobra (see reference below). So that each random draw is a positive definetely matrix with zero entries where variables are conditionaly independent.
The sampling algorithm is quite easy to implement.
Lenkoski, A. and Dobra, A. (2010). Computational aspects related to inference in
Gaussian graphical models with the G-Wishart prior. J. Comput. Graph. Statist.
