Solved – Understanding this expression of the multivariate t-distribution

gamma distributionmultivariate distributionmultivariate normal distributionpythont-distribution

I found this SO post which expresses the PDF of a multivariate t-distribution in terms of the gamma and normal distribution in python as follows

$$
G = \Gamma (k = \nu /2 ; \theta = 2 / \nu)\\
Z = N (\mu; \Sigma)\\
PDF_t = \mu + Z / \sqrt{G}
$$

where $\mu$ is the mean vector of the distribution, $\nu$ is the degrees of freedom of the t-distribution, $\Gamma$ is the gamma distribution with shape $k$ and scale $\theta$, $N$ is the multivariate normal distribution with mean $\mu$ and covariance $\Sigma$.

My dimensions and notation may be slightly wrong in the equations above, but that is why I am asking the question.

Explicitly, the python code is

d = len(Sigma) # d is length of Sigma, the covariance matrix
# g below generates m samples of the univariate gamma distribution
# then copies (np.tile) these d times and takes the transpose to produce a m*d size matrix
g = np.tile(np.random.gamma(nu/2, 2/nu, m), (d,1)).T # nu is the DOF
Z = np.random.multivariate_normal(np.zeros(d), Sigma, m) # generate samples from multivariate normal
t = mu + Z/np.sqrt(g)
  1. How can samples from the univariate gamma distribution be combined with samples from the multivariate normal distribution in the formula $PDF_t$ above?
  2. Is there some general relation between the t-distribution and the gamma and normal distribution? I found a mathworld explanation of the t-distribution and its relation to the normal and chi-squared distribution which led to a relation between chi-squared and gamma, but I haven't been able to reconcile these to get the relation above. How does this work?

Best Answer

Actually the formula in Python is not great, it should be:

d = len(Sigma) # d is length of Sigma, the covariance matrix
# g below generates m samples of the univariate gamma distribution
# then copies (np.tile) these d times and takes the transpose to produce a m*d size matrix
g = np.random.gamma(nu/2, 2/nu, m) # nu is the DOF
Z = np.random.multivariate_normal(np.zeros(d), Sigma, m) # generate samples from multivariate normal
t = mu + Z/np.sqrt(g[:,None])

Be advised that this doesn't give the PDF of a t-student distribution, but how to generate one new sample out of the distribution (which is something entirely).

  1. This solves this part of the question. The tile is a very bad call here as you only need a broadcasting, as you want to divide each element of your Z vector by the g sample.
  2. There is also the wikipedia page that shows a relation between t-student and gamma.
Related Question