# How to normalize an Inverse Wishart random matrix?

Background:

Let $$d\in \mathbb{N}$$. Define the space of (real symmetric) positive definite matrices of size $$d\times d$$ as follows: \begin{align} \mathcal{S}_{++}^d := \big\{\mathbb{M}\in \mathbb{R}^{d\times d} : \text{\mathbb{M} is symmetric and positive definite}\big\}. \label{eq:def:positive.definite.matrices} \end{align}

According to Wikipedia, given $$\nu > d + 1$$ and $$\mathbb{V}\in \mathcal{S}_{++}^d$$, the density function of the $$\mathrm{Inv\hspace{0.2mm}Wishart}_d(\nu,\mathbb{V})$$ distribution is defined by \label{eq:Wishart.density} \begin{aligned} f_{\nu,\mathbb{V}}(\mathbb{X}) := \frac{|\frac{1}{2} \mathbb{V} \, \mathbb{X}^{-1}|^{\nu/2} \exp\big(-\frac{1}{2}\mathrm{tr}(\mathbb{V} \, \mathbb{X}^{-1})\big)}{|\mathbb{X}|^{(d + 1)/2} \pi^{d(d-1)/4} \prod_{i=1}^d \Gamma(\frac{1}{2} (\nu - (d + i)))}, \quad \mathbb{X}\in \mathcal{S}_{++}^d, \end{aligned} where $$\nu$$ is the number of degrees of freedom, $$\mathbb{V}$$ is the scale matrix, and $$$$\Gamma(a) := \int_0^{\infty} t^{a - 1} e^{-t} d t, \quad a > 0,$$$$ denotes the Euler gamma function.

Question :

If $$\mathbb{X} \sim \mathrm{Inv\hspace{0.2mm}Wishart}_d(\nu,\mathbb{V})$$, how can I normalize $$\mathbb{X}$$ to get a $$\mathrm{Inv\hspace{0.2mm}Wishart}_d(\nu,\mathrm{I}_d)$$ random matrix ? Namely, I think the question reduces to finding $$\mathbb{S}\in \mathcal{S}_{++}^d$$ such that $$\mathbb{Y} = \mathbb{S}^{-1/2} \mathbb{X} \mathbb{S}^{-1/2} \sim \mathrm{Inv\hspace{0.2mm}Wishart}_d(\nu,\mathrm{I}_d)$$.

What I know:

The mean and covariance matrix for the vectorization of $$\mathbb{H}\sim \mathrm{Inv\hspace{0.2mm}Wishart}_d(\nu,\mathbb{V})$$, namely $$$$\label{eq:vectorization} \mathrm{vecp}(\mathbb{H}) := (\mathbb{H}_{11}, \mathbb{H}_{12}, \mathbb{H}_{22}, \dots, \mathbb{H}_{1d}, \mathbb{H}_{2d}, \dots, \mathbb{H}_{dd})^{\top},$$$$ ($$\mathrm{vecp}$$ is the operator that stacks the columns of the upper triangular portion of a symmetric matrix on top of each other) are well known to be: $$$$\mathbb{E}[\mathrm{vecp}(\mathbb{H})] = \frac{\mathrm{vecp}(\mathbb{V})}{\nu - d - 1} \quad \text{(alternatively, \mathbb{E}[\mathbb{H}] = \tfrac{1}{\nu - d - 1} \mathbb{V})}$$$$ and (see Theorem~3.3.16~(ii) of the book Gupta & Nagar (1999) - Matrix Variate Distributions) $$$$\label{eq:covariance.explicit.estimate} \mathbb{V}\mathrm{ar}(\mathrm{vecp}(\mathbb{H})) = \frac{2 \, B_d^{\top} \big[\mathrm{vec}(\mathbb{V}) \mathrm{vec}(\mathbb{V})^{\top} + (\nu - d - 1) (\mathbb{V} \otimes \mathbb{V})\big] B_d}{(\nu - d) (\nu - d - 1)^2 (\nu - d - 3)}, \quad (\text{not 100% sure this is correctly written})$$$$ where $$\mathrm{I}_d$$ is the identity matrix of order $$d$$, $$B_d$$ is a $$d^{\hspace{0.2mm}2} \times \frac{1}{2} d(d + 1)$$ transition matrix (see p.11 of the book Gupta & Nagar (1999) - Matrix Variate Distributions - for the precise definition), and $$\otimes$$ denotes the Kronecker product.

Side question :

Can we rewrite $$\mathrm{vec}(\mathbb{V}) \mathrm{vec}(\mathbb{V})^{\top}$$ or $$\mathrm{vec}(\mathbb{V}) \mathrm{vec}(\mathbb{V})^{\top} + (\nu - d - 1) (\mathbb{V} \otimes \mathbb{V})$$ as $$\mathbb{A} \otimes \mathbb{A}$$ for some positive definite matrix $$\mathbb{A}$$ ?

Given that $$f_{\nu,V}(X) = C_{\nu,V}\times|X|^{-(\nu+d+1)/2} \exp\left(-\tfrac{1}{2}\,\text{tr}\,(VX^{-1})\right),$$ with $$V$$ positive definite having square root $$V^{1/2}$$, it follows that the distribution of $$Y=V^{-1/2}XV^{-1/2}$$ has the desired form, $$f_{\nu,I}(Y) = C_{\nu,I}\times|Y|^{-(\nu+d+1)/2} \exp\left(-\tfrac{1}{2}\,\text{tr}\,(Y^{-1})\right),$$ since $$\text{tr}\, (VX^{-1})=\text{tr}\, Y^{-1}$$ and $$|X|=\text{constant}\times |Y|$$. (The constants $$C_{\nu,V}$$ and $$C_{\nu,I}$$ being fixed by the normalization of the distribution. Also note that the measure $$dY=|V|^{-1}dX$$ differs from $$dX$$ by a constant, so all these constants can be absorbed into $$C$$.)
The transformation of the first moment can be checked easily, since $$E(X)=V$$ we have $$(\nu-p-1)E(Y)=(\nu-p-1)V^{-1/2}E(X)V^{-1/2}=I$$, as it should be.
• It should be possible to find a matrix $M = M(\nu)$ such that $M^{-1/2} (X - \frac{V}{\nu - d - 1}) M^{-1/2}$ tends in distribution (as $\nu\to \infty$) to a standard normal random matrix. Commented Jul 19, 2021 at 13:40