Multivariate t-distribution

One common method of construction of a multivariate t-distribution, for the case of ${\displaystyle p}$ dimensions, is based on the observation that if ${\displaystyle \mathbf {y} }$ and ${\displaystyle u}$ are independent and distributed as ${\displaystyle N({\mathbf {0} },{\boldsymbol {\Sigma }})}$ and ${\displaystyle \chi _{\nu }^{2}}$ (i.e. multivariate normal and chi-squared distributions) respectively, the matrix ${\displaystyle \mathbf {\Sigma } \,}$ is a p × p matrix, and ${\displaystyle {\boldsymbol {\mu }}}$ is a constant vector then the random variable ${\textstyle {\mathbf {x} }={\mathbf {y} }/{\sqrt {u/\nu }}+{\boldsymbol {\mu }}}$ has the density[1]

${\displaystyle {\frac {\Gamma \left[(\nu +p)/2\right]}{\Gamma (\nu /2)\nu ^{p/2}\pi ^{p/2}\left|{\boldsymbol {\Sigma }}\right|^{1/2}}}\left[1+{\frac {1}{\nu }}({\mathbf {x} }-{\boldsymbol {\mu }})^{T}{\boldsymbol {\Sigma }}^{-1}({\mathbf {x} }-{\boldsymbol {\mu }})\right]^{-(\nu +p)/2}}$

and is said to be distributed as a multivariate t-distribution with parameters ${\displaystyle {\boldsymbol {\Sigma }},{\boldsymbol {\mu }},\nu }$. Note that ${\displaystyle \mathbf {\Sigma } }$ is not the covariance matrix since the covariance is given by ${\displaystyle \nu /(\nu -2)\mathbf {\Sigma } }$ (for ${\displaystyle \nu >2}$).

The constructive definition of a multivariate t-distribution simultaneously serves as a sampling algorithm:

1. Generate ${\displaystyle u\sim \chi _{\nu }^{2}}$ and ${\displaystyle \mathbf {y} \sim N(\mathbf {0} ,{\boldsymbol {\Sigma }})}$, independently.
2. Compute ${\displaystyle \mathbf {x} \gets {\sqrt {\nu /u}}\mathbf {y} +{\boldsymbol {\mu }}}$.

This formulation gives rise to the hierarchical representation of a multivariate t-distribution as a scale-mixture of normals: ${\displaystyle u\sim \mathrm {Ga} (\nu /2,\nu /2)}$ where ${\displaystyle \mathrm {Ga} (a,b)}$ indicates a gamma distribution with density proportional to ${\displaystyle x^{a-1}e^{-bx}}$, and ${\displaystyle \mathbf {x} \mid u}$ conditionally follows ${\displaystyle N({\boldsymbol {\mu }},u^{-1}{\boldsymbol {\Sigma }})}$.

In the special case ${\displaystyle \nu =1}$, the distribution is a multivariate Cauchy distribution.

There are in fact many candidates for the multivariate generalization of Student's t-distribution. An extensive survey of the field has been given by Kotz and Nadarajah (2004). The essential issue is to define a probability density function of several variables that is the appropriate generalization of the formula for the univariate case. In one dimension (${\displaystyle p=1}$), with ${\displaystyle t=x-\mu }$ and ${\displaystyle \Sigma =1}$, we have the probability density function

${\displaystyle f(t)={\frac {\Gamma [(\nu +1)/2]}{{\sqrt {\nu \pi \,}}\,\Gamma [\nu /2]}}(1+t^{2}/\nu )^{-(\nu +1)/2}}$

and one approach is to write down a corresponding function of several variables. This is the basic idea of elliptical distribution theory, where one writes down a corresponding function of ${\displaystyle p}$ variables ${\displaystyle t_{i}}$ that replaces ${\displaystyle t^{2}}$ by a quadratic function of all the ${\displaystyle t_{i}}$. It is clear that this only makes sense when all the marginal distributions have the same degrees of freedom ${\displaystyle \nu }$. With ${\displaystyle \mathbf {A} ={\boldsymbol {\Sigma }}^{-1}}$, one has a simple choice of multivariate density function

${\displaystyle f(\mathbf {t} )={\frac {\Gamma ((\nu +p)/2)\left|\mathbf {A} \right|^{1/2}}{{\sqrt {\nu ^{p}\pi ^{p}\,}}\,\Gamma (\nu /2)}}\left(1+\sum _{i,j=1}^{p,p}A_{ij}t_{i}t_{j}/\nu \right)^{-(\nu +p)/2}}$

which is the standard but not the only choice.

An important special case is the standard bivariate t-distribution, p = 2:

${\displaystyle f(t_{1},t_{2})={\frac {\left|\mathbf {A} \right|^{1/2}}{2\pi }}\left(1+\sum _{i,j=1}^{2,2}A_{ij}t_{i}t_{j}/\nu \right)^{-(\nu +2)/2}}$

Note that ${\displaystyle {\frac {\Gamma \left({\frac {\nu +2}{2}}\right)}{\pi \ \nu \Gamma \left({\frac {\nu }{2}}\right)}}={\frac {1}{2\pi }}}$.

Now, if ${\displaystyle \mathbf {A} }$ is the identity matrix, the density is

${\displaystyle f(t_{1},t_{2})={\frac {1}{2\pi }}\left(1+(t_{1}^{2}+t_{2}^{2})/\nu \right)^{-(\nu +2)/2}.}$

The difficulty with the standard representation is revealed by this formula, which does not factorize into the product of the marginal one-dimensional distributions. When ${\displaystyle \Sigma }$ is diagonal the standard representation can be shown to have zero correlation but the marginal distributions do not agree with statistical independence.

The definition of the cumulative distribution function (cdf) in one dimension can be extended to multiple dimensions by defining the following probability (here ${\displaystyle \mathbf {x} }$ is a real vector):

${\displaystyle F(\mathbf {x} )=\mathbb {P} (\mathbf {X} \leq \mathbf {x} ),\quad {\textrm {where}}\;\;\mathbf {X} \sim t_{\nu }({\boldsymbol {\mu }},{\boldsymbol {\Sigma }}).}$

There is no simple formula for ${\displaystyle F(\mathbf {x} )}$, but it can be approximated numerically via Monte Carlo integration.[2][3]

This was demonstrated by Muirhead [4] though previously derived using the simpler ratio representation above, by Cornish.[5] Let vector ${\displaystyle X}$ follow the multivariate t distribution and partition into two subvectors of ${\displaystyle p_{1},p_{2}}$ elements:

${\displaystyle X_{p}={\begin{bmatrix}X_{1}\\X_{2}\end{bmatrix}}\sim t_{p}\left(\mu _{p},\Sigma _{p\times p},\nu \right)}$

where ${\displaystyle p_{1}+p_{2}=p}$, the known mean vector is ${\displaystyle \mu _{p}={\begin{bmatrix}\mu _{1}\\\mu _{2}\end{bmatrix}}}$ and the scale matrix is ${\displaystyle \Sigma _{p\times p}={\begin{bmatrix}\Sigma _{11}&\Sigma _{12}\\\Sigma _{21}&\Sigma _{22}\end{bmatrix}}}$.

Then

${\displaystyle X_{2}|X_{1}\sim t_{p_{2}}\left(\mu _{2|1},{\frac {\nu +d_{1}}{\nu +p_{1}}}\Sigma _{22|1},\nu +p_{1}\right)}$

where

${\displaystyle \mu _{2|1}=\mu _{2}+\Sigma _{21}\Sigma _{11}^{-1}\left(X_{1}-\mu _{1}\right)}$ is the conditional mean where it exists or median otherwise.

${\displaystyle \Sigma _{22|1}=\Sigma _{22}-\Sigma _{21}\Sigma _{11}^{-1}\Sigma _{12}}$ is the Schur complement of ${\displaystyle \Sigma _{11}{\text{ in }}\Sigma .}$

${\displaystyle d_{1}=(X_{1}-\mu _{1})^{T}\Sigma _{11}^{-1}(X_{1}-\mu _{1})}$ is the squared Mahalanobis distance of ${\displaystyle X_{1}}$ from ${\displaystyle \mu _{1}}$ with scale matrix ${\displaystyle \Sigma _{11}}$

See [6] for a simple proof of the above conditional distribution.

The use of such distributions[7] is enjoying renewed interest due to applications in mathematical finance, especially through the use of the Student's t copula.

Constructed as an elliptical distribution,[8] take the simplest centralised case with spherical symmetry and no scaling, ${\displaystyle \Sigma =\operatorname {I} \,}$, then the multivariate t-PDF takes the form

${\displaystyle f_{X}(X)=g(X^{T}X)={\frac {\Gamma {\big (}{\frac {1}{2}}(\nu +p)\,{\big )}}{(\nu \pi )^{\,p/2}\Gamma {\big (}{\frac {1}{2}}\nu {\big )}}}{\bigg (}1+\nu ^{-1}X^{T}X{\bigg )}^{-(\nu +p)/2}}$

where ${\displaystyle X=(x_{1},\cdots ,x_{p})^{T}{\text{ is a }}p{\text{-vector}}}$ and ${\displaystyle \nu }$ = degrees of freedom. Muirhead (section 1.5) refers to this as a multivariate Cauchy distribution. The covariance of ${\displaystyle X}$ is

${\displaystyle \int _{-\infty }^{\infty }\cdots \int _{-\infty }^{\infty }f_{X}(x_{1},\dots ,x_{p})XX^{T}\,dx_{1}\dots dx_{p}={\frac {\nu }{\nu -2}}\operatorname {E} (XX^{T})}$

The aim is to convert the Cartesian PDF to a radial one. Kibria and Joarder,[9] in a tutorial-style paper, define radial measure ${\displaystyle r_{2}=R^{2}={\frac {X^{T}X}{p}}}$ such that

${\displaystyle \operatorname {E} [r_{2}]=\int _{-\infty }^{\infty }\cdots \int _{-\infty }^{\infty }f_{X}(x_{1},\dots ,x_{p}){\frac {X^{T}X}{p}}\,dx_{1}\dots dx_{p}}$

which is equivalent to the variance of ${\displaystyle p}$-element vector ${\displaystyle X}$ treated as a univariate zero-mean random sequence. They note that ${\displaystyle r_{2}}$ follows the Fisher-Snedecor or ${\displaystyle F}$ distribution:

${\displaystyle r_{2}\sim F_{F}(p,\nu )=B{\bigg (}{\frac {p}{2}},{\frac {\nu }{2}}{\bigg )}^{-1}{\bigg (}{\frac {p}{\nu }}{\bigg )}^{p/2}r_{2}^{p/2-1}{\bigg (}1+{\frac {p}{\nu }}r_{2}{\bigg )}^{-(p+\nu )/2}}$

having mean value ${\displaystyle \operatorname {E} [r_{2}]={\frac {\nu }{\nu -2}}}$. ${\displaystyle F}$-distributions arise naturally in tests of sums of squares of sampled data after normalization by the sample standard deviation.

By a change of random variable to ${\displaystyle y={\frac {p}{\nu }}r_{2}={\frac {X^{T}X}{\nu }}}$ in the equation above, retaining ${\displaystyle p}$-vector ${\displaystyle X}$, we have ${\displaystyle \operatorname {E} [y]=\int _{-\infty }^{\infty }\cdots \int _{-\infty }^{\infty }f_{X}(X){\frac {X^{T}X}{\nu }}\,dx_{1}\dots dx_{p}={\frac {p}{\nu -2}}}$ and probability distribution

${\displaystyle {\begin{aligned}f_{Y}(y|\,p,\nu )&={\frac {\nu }{p}}B{\bigg (}{\frac {p}{2}},{\frac {\nu }{2}}{\bigg )}^{-1}{\big (}{\frac {p}{\nu }}{\big )}^{\,p/2}{\big (}{\frac {p}{\nu }}{\big )}^{-p/2-1}y^{\,p/2-1}{\big (}1+y{\big )}^{-(p+\nu )/2}\\\\&=B{\bigg (}{\frac {p}{2}},{\frac {\nu }{2}}{\bigg )}^{-1}y^{\,p/2-1}(1+y)^{-(\nu +p)/2}\end{aligned}}}$

which is a regular Beta-prime distribution ${\displaystyle y\sim \beta \,'{\bigg (}y;{\frac {p}{2}},{\frac {\nu }{2}}{\bigg )}}$ having mean value ${\displaystyle {\frac {{\frac {1}{2}}p}{{\frac {1}{2}}\nu -1}}={\frac {p}{\nu -2}}}$. The cumulative distribution function of ${\displaystyle y}$ is thus ${\displaystyle F_{Y}(y)\sim I\,{\bigg (}{\frac {y}{1+y}};\,{\frac {p}{2}},{\frac {\nu }{2}}{\bigg )}}$

where ${\displaystyle I}$ is the incomplete Beta function.

These results can also be derived via a straightforward coordinate transformation from cartesian to spherical. A constant radius surface at ${\displaystyle R=(X^{T}X)^{1/2}}$ with PDF ${\displaystyle p_{X}(X)\propto {\bigg (}1+\nu ^{-1}R^{2}{\bigg )}^{-(\nu +p)/2}}$ is an iso-density surface. Given this density value, the quantum of probability on a shell of surface area ${\displaystyle A_{R}}$ and thickness ${\displaystyle \delta R}$ at ${\displaystyle R}$ is ${\displaystyle \delta P=p_{X}(R)\,A_{R}\delta R}$.

The enclosed ${\displaystyle p}$-sphere of radius ${\displaystyle R}$ has surface area ${\displaystyle A_{R}={\frac {2\pi ^{p/2}R^{\,p-1}}{\Gamma (p/2)}}}$, and substitution into ${\displaystyle \delta P}$ shows that the shell has element of probability ${\displaystyle \delta P=p_{X}(R){\frac {2\pi ^{p/2}R^{p-1}}{\Gamma (p/2)}}\delta R}$ which is equivalent to radial density function

${\displaystyle f_{R}(R)={\frac {\Gamma {\big (}{\frac {1}{2}}(\nu +p)\,{\big )}}{\nu ^{\,p/2}\pi ^{\,p/2}\Gamma {\big (}{\frac {1}{2}}\nu {\big )}}}{\frac {2\pi ^{p/2}R^{p-1}}{\Gamma (p/2)}}{\bigg (}1+{\frac {R^{2}}{\nu }}{\bigg )}^{-(\nu +p)/2}}$

which further simplifies to ${\displaystyle f_{R}(R)={\frac {2}{\nu ^{1/2}B{\big (}{\frac {1}{2}}p,{\frac {1}{2}}\nu {\big )}}}{\bigg (}{\frac {R^{2}}{\nu }}{\bigg )}^{(p-1)/2}{\bigg (}1+{\frac {R^{2}}{\nu }}{\bigg )}^{-(\nu +p)/2}}$ where ${\displaystyle B(*,*)}$ is the Beta function.

Changing the radial variable to ${\displaystyle y=R^{2}/\nu }$ returns the previous Beta Prime distribution ${\displaystyle f_{Y}(y)={\frac {1}{B{\big (}{\frac {1}{2}}p,{\frac {1}{2}}\nu {\big )}}}y^{\,p/2-1}{\bigg (}1+y{\bigg )}^{-(\nu +p)/2}}$

To scale the radial variables without changing the radial shape function, define scale matrix ${\displaystyle \Sigma =\alpha \operatorname {I} }$ , yielding a 3-parameter Cartesian density function, ie. the probability ${\displaystyle \Delta _{P}}$ in volume element ${\displaystyle dx_{1}\dots dx_{p}}$ is

${\displaystyle \Delta _{P}{\big (}f_{X}(X\,|\alpha ,p,\nu ){\big )}={\frac {\Gamma {\big (}{\frac {1}{2}}(\nu +p)\,{\big )}}{(\nu \pi )^{\,p/2}\alpha ^{\,p/2}\Gamma {\big (}{\frac {1}{2}}\nu {\big )}}}{\bigg (}1+{\frac {X^{T}X}{\alpha \nu }}{\bigg )}^{-(\nu +p)/2}\;dx_{1}\dots dx_{p}}$

or, in terms of scalar radial variable ${\displaystyle R}$,

${\displaystyle f_{R}(R\,|\alpha ,p,\nu )={\frac {2}{\alpha ^{1/2}\;\nu ^{1/2}B{\big (}{\frac {1}{2}}p,{\frac {1}{2}}\nu {\big )}}}{\bigg (}{\frac {R^{2}}{\alpha \,\nu }}{\bigg )}^{(p-1)/2}{\bigg (}1+{\frac {R^{2}}{\alpha \,\nu }}{\bigg )}^{-(\nu +p)/2}}$

The moments of all the radial variables can be derived from the Beta Prime distribution. If ${\displaystyle Z\sim \beta '(a,b)}$ then ${\displaystyle \operatorname {E} (Z^{m})={\frac {B(a+m,b-m)}{B(a,b)}}}$, a known result. Thus, for variable ${\displaystyle y}$, proportional to ${\displaystyle R^{2}}$, we have

${\displaystyle \operatorname {E} (y^{m})={\frac {B({\frac {1}{2}}p+m,{\frac {1}{2}}\nu -m)}{B({\frac {1}{2}}p,{\frac {1}{2}}\nu )}}={\frac {\Gamma {\big (}{\frac {1}{2}}p+m{\big )}\;\Gamma {\big (}{\frac {1}{2}}\nu -m{\big )}}{\Gamma {\big (}{\frac {1}{2}}p{\big )}\;\Gamma {\big (}{\frac {1}{2}}\nu {\big )}}}}$

The moments of ${\displaystyle r_{2}=\nu \,y}$ are

${\displaystyle \operatorname {E} (r_{2}^{m})=\nu ^{m}\operatorname {E} (y^{m})}$

while introducing the scale matrix ${\displaystyle \alpha \operatorname {I} }$ yields

${\displaystyle \operatorname {E} (r_{2}^{m}|\alpha )=\alpha ^{m}\nu ^{m}\operatorname {E} (y^{m})}$

Moments relating to radial variable ${\displaystyle R}$ are found by setting ${\displaystyle R=(\alpha \nu y)^{1/2}}$ and ${\displaystyle M=2m}$ whereupon

${\displaystyle \operatorname {E} (R^{M})=\operatorname {E} {\big (}(\alpha \nu y)^{1/2}{\big )}^{2m}=(\alpha \nu )^{M/2}\operatorname {E} (y^{M/2})=(\alpha \nu )^{M/2}{\frac {B{\big (}{\frac {1}{2}}(p+M),{\frac {1}{2}}(\nu -M){\big )}}{B({\frac {1}{2}}p,{\frac {1}{2}}\nu )}}}$

During affine transformations, the degrees of freedom parameter ${\displaystyle \nu }$ remains invariant throughout and all vectors must ultimately derive from one initial isotropic spherical vector ${\displaystyle Z}$ whose elements remain 'entangled' and are not statistically independent. Adding two sample multivariate t vectors generated with independent Chi-squared samples and different ${\displaystyle \nu }$ values: ${\textstyle {1}/{\sqrt {u_{1}/\nu _{1}}},\;\;{1}/{\sqrt {u_{2}/\nu _{2}}}}$ , as defined in the leading paragraph, will not produce internally consistent distributions, though they will yield a Behrens-Fisher problem.[10]

In univariate statistics, the Student's t-test makes use of Student's t-distribution. Hotelling's T-squared distribution is a distribution that arises in multivariate statistics. The matrix t-distribution is a distribution for random variables arranged in a matrix structure.

• Multivariate normal distribution, which is the limiting case of the multivariate Student's t-distribution when ${\displaystyle \nu \uparrow \infty }$.
• Chi distribution, the pdf of the scaling factor in the construction the Student's t-distribution and also the 2-norm (or Euclidean norm) of a multivariate normally distributed vector (centered at zero).
  • Rayleigh distribution#Student's t, random vector length of multivariate t-distribution
• Mahalanobis distance

• Copula Methods vs Canonical Multivariate Distributions: the multivariate Student T distribution with general degrees of freedom
• Multivariate Student's t distribution