  We propose an efficient way to sample from a class of structured multivariate
Gaussian distributions which routinely arise as conditional posteriors of model
parameters that are assigned a conditionally Gaussian prior. The proposed
algorithm only requires matrix operations in the form of matrix multiplications
and linear system solutions. We exhibit that the computational complexity of
the proposed algorithm grows linearly with the dimension unlike existing
algorithms relying on Cholesky factorizations with cubic orders of complexity.
The algorithm should be broadly applicable in settings where Gaussian scale
mixture priors are used on high dimensional model parameters. We provide an
illustration through posterior sampling in a high dimensional regression
setting with a horseshoe prior on the vector of regression coefficients.
