Moments Of The Inverse-Wishart Distribution

Why Moment Existence Needs Care

The inverse-Wishart distribution is used for random positive-definite matrices, including covariance matrices in Gaussian models. Its mean is often used as a summary, and its second moments describe the variability and dependence of its entries. But positive definiteness alone does not guarantee either quantity is finite. The key is to identify which part of the distribution's tail controls the moments.

We use the parameterization $ \Sigma\sim IW_p(\nu,\Psi) $, where $p$ is the matrix dimension, $\nu>p-1$, and $\Psi$ is positive definite. In this convention the density is proportional to $ |\Sigma|^{-(\nu+p+1)/2}\exp\{-\operatorname{tr}(\Psi\Sigma^{-1})/2\} $. In particular, $\Psi$ is a scale matrix, not the mean.

Definition: If $ \Sigma\sim IW_p(\nu,\Psi) $, then $\Sigma$ is a random symmetric positive-definite $p\times p$ matrix with density $$ f(\Sigma)=\frac{|\Psi|^{\nu/2}}{2^{\nu p/2}\Gamma_p(\nu/2)} |\Sigma|^{-(\nu+p+1)/2} \exp\!\left(-\frac12\operatorname{tr}(\Psi\Sigma^{-1})\right), \qquad \Sigma\succ0. $$ Here $\nu>p-1$, $\Psi\succ0$, and $\Gamma_p$ is the multivariate gamma function. Equivalently, $ \Sigma^{-1}\sim W_p(\nu,\Psi^{-1}) $.

Where the Moment Conditions Come From

A matrix approaching singularity has an eigenvalue tending to zero. In the inverse-Wishart density, that boundary is suppressed by the exponential term: positive definiteness of $\Psi$ gives $\operatorname{tr}(\Psi\Sigma^{-1})\geq c/t$ when an eigenvalue $t$ of $\Sigma$ tends to zero, for some $c>0$ if the other eigenvalues remain in a fixed bounded region. The resulting factor $e^{-c/(2t)}$ decays faster than any power of $t$. This boundary does not cause the moments of entries of $\Sigma$ to diverge.

The relevant tail is instead at large eigenvalues of $\Sigma$. Keep $p-1$ eigenvalues in a bounded region away from zero and let the remaining eigenvalue be $t\to\infty$. The exponential factor tends to a finite positive value rather than suppressing the tail. The determinant power contributes $t^{-(\nu+p+1)/2}$, and the eigenvalue-coordinate Jacobian contributes $t^{p-1}$ through the products of differences between $t$ and the other eigenvalues. Thus the marginal tail in this region is proportional to $t^{(p-3-\nu)/2}$.

An entry of a positive-semidefinite matrix is bounded in magnitude by its largest eigenvalue. Conversely, on a set of eigenvector orientations of positive measure, a specified diagonal entry is a fixed positive fraction of that eigenvalue. Consequently, the tail integral for a moment of order $q$ has the criterion $\int^\infty t^{q+(p-3-\nu)/2}\,dt<\infty$. This integral converges precisely when $ \nu>p-1+2q $. Regions with more than one large eigenvalue do not give a stricter condition: their additional eigenvalue factors and volume are accounted for by the same eigenvalue integral, and the one-large-eigenvalue tail gives the binding threshold.

Theorem: For $ \Sigma\sim IW_p(\nu,\Psi) $, the mean entries are finite if and only if $ \nu>p+1 $. All entrywise second moments $E[\Sigma_{ij}\Sigma_{kl}]$ are finite if and only if $ \nu>p+3 $. More generally, absolute moments of order $q>0$ require $ \nu>p-1+2q $.

The strict inequalities matter. At equality, the corresponding tail integral has exponent $-1$ and diverges logarithmically. Since second moments require the stronger condition, a distribution can have a finite mean but infinite entrywise second moments.

Mean and Second-Moment Formulas

The moment formulas follow from inverse-Wishart integration. One route is to write $Y=\Sigma^{-1}$, whose Wishart density is proportional to $|Y|^{(\nu-p-1)/2}\exp\{-\operatorname{tr}(\Psi Y)/2\}$, and integrate by parts on the positive-definite cone. For clarity, the resulting component identities are stated below with their derivation.

Theorem: Let $ \Sigma\sim IW_p(\nu,\Psi) $, with $\Psi\succ0$. If $ \nu>p+1 $, then $$ E[\Sigma_{ij}]=\frac{\Psi_{ij}}{\nu-p-1}. $$ If $ \nu>p+3 $, then $$ E[\Sigma_{ij}\Sigma_{kl}] = \frac{(\nu-p-2)\Psi_{ij}\Psi_{kl} +\Psi_{ik}\Psi_{jl}+\Psi_{il}\Psi_{jk}} {(\nu-p)(\nu-p-1)(\nu-p-3)}. $$

Here is the integration argument behind the identities. In the Wishart variable $Y$, substitute its density into the integral defining each expectation of an entry of $Y^{-1}$. Differentiate the density and the entries of $Y^{-1}$ with respect to the independent symmetric coordinates of $Y$; the derivative of the inverse is $d(Y^{-1})=-Y^{-1}(dY)Y^{-1}$. Integration by parts gives the first-moment identity with denominator $\nu-p-1$. Applying the same calculation to products of two inverse entries gives the three pairings of indices displayed in the second-moment formula. The determinant exponent supplies the successive factors $\nu-p$, $\nu-p-1$, and $\nu-p-3$; the two cross-index pairings supply $\Psi_{ik}\Psi_{jl}$ and $\Psi_{il}\Psi_{jk}$, while the remaining pairing has coefficient $\nu-p-2$. These integrations have no boundary contribution when their respective moments exist: at a zero eigenvalue the exponential suppresses the integrand, and at infinity the tail conditions above make the boundary terms vanish. The necessity of those same conditions follows from the tail test, so the formulas are used only in their finite-moment ranges. $\square$

The second-moment expression is symmetric in $i,j$, in $k,l$, and under exchange of the two entry pairs, as required for products of entries of a symmetric matrix. In particular, setting $k=i$ and $l=j$ gives a finite raw second moment for $\Sigma_{ij}$, not its variance; the variance is obtained by subtracting $E[\Sigma_{ij}]^2$.

Worked Example: The Scalar Inverse-Wishart Case

For $p=1$, the distribution reduces to $ \Sigma=\psi/U $, where $U\sim\chi^2_\nu$ and $\psi>0$. Take $\nu=7$ and $\psi=2$. Here $d=\nu-p-1=5$, so the mean exists and the second moment exists because $\nu>p+3$.

The mean formula gives $E[\Sigma]=\psi/d=2/5$. For the second moment, the general formula has numerator $(d-1)\psi^2+\psi^2+\psi^2=(d+1)\psi^2=6(4)=24$, and denominator $d(d+1)(d-2)=5\cdot6\cdot3=90$. Therefore $E[\Sigma^2]=24/90=4/15$.

Thus $\operatorname{Var}(\Sigma)=4/15-(2/5)^2=4/15-4/25=8/75$. This agrees with the scalar inverse-chi-square calculation $E[U^{-2}]=1/[(\nu-2)(\nu-4)]=1/15$, multiplied by $\psi^2=4$.

Worked Example: Entrywise Moments for a Two-Dimensional Matrix

Let $ \Sigma\sim IW_2(8,\Psi) $ with $$ \Psi=\begin{pmatrix}2&0\\0&1\end{pmatrix}. $$ Then $d=\nu-p-1=8-2-1=5$, and the second-moment denominator is $d(d+1)(d-2)=5\cdot6\cdot3=90$.

The mean is $E[\Sigma]=\Psi/5$, so $E[\Sigma_{11}]=2/5$, $E[\Sigma_{22}]=1/5$, and $E[\Sigma_{12}]=0$. For the first diagonal entry, $$ E[\Sigma_{11}^2] =\frac{(d-1)\Psi_{11}^2+2\Psi_{11}^2}{90} =\frac{(4+2)4}{90} =\frac{4}{15}. $$ Hence $\operatorname{Var}(\Sigma_{11})=4/15-(2/5)^2=8/75$.

For the product of the diagonal entries, $$ E[\Sigma_{11}\Sigma_{22}] =\frac{(d-1)\Psi_{11}\Psi_{22} +\Psi_{12}\Psi_{12}+\Psi_{12}\Psi_{12}}{90} =\frac{4\cdot2}{90} =\frac{4}{45}. $$ The covariance is $4/45-(2/5)(1/5)=4/45-2/25=2/225$.

Finally, for the off-diagonal entry squared, $$ E[\Sigma_{12}^2] =\frac{(d-1)\Psi_{12}^2+\Psi_{11}\Psi_{22}+\Psi_{12}^2}{90} =\frac{2}{90} =\frac{1}{45}. $$ Its mean is zero, so this is also its variance. Each calculation uses the same entrywise formula, but the index pairings determine which elements of $\Psi$ appear.

How Degrees of Freedom Change the Moments

At fixed $p$ and $\Psi$, increasing $\nu$ moves the inverse-Wishart distribution toward smaller matrix values and makes more moments finite. This is visible in the mean denominator $\nu-p-1$, but the change in existence is just as important: the mean and second moments do not turn on at the same threshold.

Worked Example: A Mean Without Finite Second Moments

Set $p=2$ and $\Psi=I_2$. For $\nu=4$, the density is valid because $4>p-1=1$, and the mean exists because $4>p+1=3$. Here $d=\nu-p-1=1$, so $E[\Sigma]=I_2$. But the second-moment condition is $\nu>p+3=5$, which fails. The diagonal second moments are infinite.

Raise the degrees of freedom to $\nu=6$. Now $d=3$, the mean is $I_2/3$, and second moments are finite. The formula gives $E[\Sigma_{11}^2]=1/[d(d-2)]=1/3$, so $\operatorname{Var}(\Sigma_{11})=1/3-1/9=2/9$.

At $\nu=10$, $d=7$, the mean becomes $I_2/7$, and $E[\Sigma_{11}^2]=1/[7(5)]=1/35$. The variance is $1/35-1/49=2/245$. The numerical decrease is consistent with increasing $\nu$, while the $\nu=4$ case illustrates that a finite mean alone does not justify variance calculations.

Common Mistake: Looking Near Singularity for the Wrong Moment

A common error is to argue that inverse-Wishart moments of $\Sigma$ fail because $\Sigma$ can be nearly singular, or that the density's exponential factor controls all tails. The first claim confuses $\Sigma$ with $\Sigma^{-1}$: a small eigenvalue of $\Sigma$ can make its inverse large, but does not itself make the entries of $\Sigma$ large. In addition, the density's exponential factor suppresses that small-eigenvalue boundary.

The correction is to inspect the large-eigenvalue tail of $\Sigma$. There, $\Sigma^{-1}$ has a small eigenvalue and the exponential factor need not decay; the density's polynomial tail determines whether the moments exist. A second frequent error is to use the first-moment threshold for a variance. A variance requires a finite second moment, hence the stricter condition $ \nu>p+3 $.

Key takeaway: For $ \Sigma\sim IW_p(\nu,\Psi) $, the mean is $ \Psi/(\nu-p-1) $ when $ \nu>p+1 $. Entrywise second moments require $ \nu>p+3 $ and use the three-pairing formula above. These thresholds come from the polynomial tail at large eigenvalues of $\Sigma$, not from the near-singular boundary.

Check Your Understanding

Use the parameterization and moment conditions stated above.

For $ \Sigma\sim IW_3(\nu,\Psi) $, what strict condition on $\nu$ guarantees that all mean entries are finite?
For the same dimension, what condition guarantees finite entrywise second moments?
Explain why the exponential term suppresses the boundary where an eigenvalue of $\Sigma$ tends to zero.
If $ \Sigma\sim IW_2(7,I_2) $, calculate $E[\Sigma_{11}]$ and $E[\Sigma_{11}^2]$.
Why does a finite mean not, by itself, justify computing an entrywise variance?