I would like to be able to sample the standard deviation of a multidimensional Gaussian distribution of dimension $n$; that is, given some $\phi$, I would like to sample

$P(\sigma | \phi) \propto \frac{1}{\sigma^n} e^{-\frac{\phi}{2\sigma^2}}$

For high $n$, this distribution is sharply peaked; for my purposes $n$ will be on the order of 1000-3000. What is an ``efficient'' method of obtaining a single sample from this distribution? (A single sample as $\phi$ will change between each sample.)


Your problem, in my opinion, can be solved if you are willing to sample from $\sigma^2$ and not $\sigma$.

Look at the inverse Gamma distribution, with parameter $\alpha$ and $\beta$. $\sigma^2 \sim IG(\alpha, \beta)$ if its density is : $$ \frac{\beta^\alpha (\sigma^2)^{-(\alpha+1)} \exp(-\frac{\beta}{\sigma^2}) }{\Gamma(\alpha)} $$ If you set $\alpha = \frac{n}{2}$ and $\beta = \frac{\phi}{2}$, then your distribution is proportional to $IG(\alpha, \beta)$

  • 3
    $\begingroup$ +1 You could also sample from an appropriate Gamma distribution and transform that result (take the reciprocal square root). Since Gamma RNGs are easier to find than inverse-Gamma RNGs this option widens the possible implementation choices. $\endgroup$ – whuber Aug 4 '14 at 12:53

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.