[proofplan]
We first prove the density for the identity scale matrix $\Sigma=I_p$ by representing the Gaussian sample as a rectangular matrix and applying the Gram polar integration formula, which decomposes Lebesgue measure according to the positive definite Gram matrix and the orthonormal-frame fibre. This change of variables produces the determinant power $(\det W)^{(n-p-1)/2}$; the remaining fibre-volume constant is then determined from the requirement that the density integrate to $1$, using the defining integral for $\Gamma_p(n/2)$. Finally, for general $\Sigma$, we reduce in distribution to the identity scale case through a [Gaussian random vector](/page/Gaussian%20Random%20Vector) linear transformation and compute the Jacobian of the congruence map on $\mathbb{S}^p$.
[/proofplan]
[step:Represent the identity scale Wishart matrix as a Gram matrix]
Assume first that $\Sigma=I_p$. Define the random matrix
\begin{align*}
X: \Omega &\to \mathbb{R}^{n \times p} \\
\omega &\mapsto
\begin{pmatrix}
X_1(\omega)^\top \\
\vdots \\
X_n(\omega)^\top
\end{pmatrix}.
\end{align*}
With respect to Lebesgue measure $\mathcal{L}^{np}$ on $\mathbb{R}^{n\times p}$, the density of $X$ is
\begin{align*}
g: \mathbb{R}^{n\times p} &\to [0,\infty) \\
Y &\mapsto (2\pi)^{-np/2}\exp\left(-\frac{1}{2}\operatorname{tr}(Y^\top Y)\right).
\end{align*}
The Wishart matrix is
\begin{align*}
W=X^\top X.
\end{align*}
Since $n\ge p$, a Gaussian $n\times p$ matrix has rank $p$ almost surely: every $p\times p$ minor is a polynomial in the entries of $X$, at least one such polynomial is not identically zero, and the common zero set of these nonzero polynomial minors has $\mathcal{L}^{np}$-measure zero. Hence $W=X^\top X\in\mathbb{S}_{++}^p$ almost surely.
[guided]
We begin with the identity covariance case because all geometry is already present there and no scale factors obscure the computation. Define
\begin{align*}
X: \Omega &\to \mathbb{R}^{n \times p} \\
\omega &\mapsto
\begin{pmatrix}
X_1(\omega)^\top \\
\vdots \\
X_n(\omega)^\top
\end{pmatrix}.
\end{align*}
Thus the $i$-th row of $X$ is the transpose of the Gaussian vector $X_i$. Because the rows are independent and each has distribution $\mathcal{N}_p(0,I_p)$, the joint density of all $np$ scalar entries is
\begin{align*}
g: \mathbb{R}^{n\times p} &\to [0,\infty) \\
Y &\mapsto (2\pi)^{-np/2}\exp\left(-\frac{1}{2}\operatorname{tr}(Y^\top Y)\right)
\end{align*}
with respect to $\mathcal{L}^{np}$. The trace identity
\begin{align*}
\operatorname{tr}(Y^\top Y)=\sum_{i=1}^n\sum_{j=1}^p Y_{ij}^2
\end{align*}
is exactly the squared Euclidean norm of the matrix entries.
The random matrix in the theorem is
\begin{align*}
W=X^\top X.
\end{align*}
For every nonzero vector $v\in\mathbb{R}^p$,
\begin{align*}
v^\top Wv=v^\top X^\top Xv=|Xv|^2.
\end{align*}
Thus $W$ is positive definite precisely when $X$ has rank $p$. Since $n\ge p$, an $n\times p$ matrix can have full column rank. The set of matrices with rank less than $p$ is the common zero set of all $p\times p$ minors. At least one such minor is a nonzero polynomial, and the zero set of a nonzero polynomial in Euclidean space has Lebesgue measure zero. Since the Gaussian density is absolutely continuous with respect to $\mathcal{L}^{np}$, $X$ has rank $p$ almost surely. Therefore $W\in\mathbb{S}_{++}^p$ almost surely.
[/guided]
[/step]
[step:Apply the Gram polar integration formula]
We invoke the [Gram polar integration formula for rectangular matrices](/page/Gram%20Polar%20Integration%20Formula), stated here explicitly because it is the geometric change of variables behind the Wishart density. Let $\mathbb{V}_{p}(\mathbb{R}^n):=\{Q\in\mathbb{R}^{n\times p}:Q^\top Q=I_p\}$ be the [Stiefel manifold](/page/Stiefel%20Manifold), equipped with its invariant [Hausdorff measure](/page/Hausdorff%20Measure) $\nu_{n,p}$. The formula applies to every non-negative Borel map, so no integrability hypothesis is needed; the iterated integrals are interpreted in the extended non-negative sense. For every non-negative Borel map
\begin{align*}
\Phi:\mathbb{R}^{n\times p}\to[0,\infty),
\end{align*}
one has
\begin{align*}
\int_{\mathbb{R}^{n\times p}}\Phi(Y)\,d\mathcal{L}^{np}(Y)
=
c_{n,p}
\int_{\mathbb{S}_{++}^p}
\int_{\mathbb{V}_p(\mathbb{R}^n)}
\Phi(QA^{1/2})(\det A)^{(n-p-1)/2}
\,d\nu_{n,p}(Q)\,d\mathcal{L}^{p(p+1)/2}(A),
\end{align*}
where $c_{n,p}>0$ is a constant depending only on $n$ and $p$. This is the polar-coordinate formula for the Gram map $Y\mapsto Y^\top Y$; its Jacobian factor is $(\det A)^{(n-p-1)/2}$.
Applying this formula to the density $g$ gives, for every non-negative Borel map
\begin{align*}
\varphi:\mathbb{S}_{++}^p\to[0,\infty),
\end{align*}
the identity
\begin{align*}
\mathbb{E}[\varphi(W)]
&=
\int_{\mathbb{R}^{n\times p}}\varphi(Y^\top Y)g(Y)\,d\mathcal{L}^{np}(Y)\\
&=
c_{n,p}
\int_{\mathbb{S}_{++}^p}
\int_{\mathbb{V}_p(\mathbb{R}^n)}
\varphi(A)(2\pi)^{-np/2}\exp\left(-\frac{1}{2}\operatorname{tr}(A)\right)
(\det A)^{(n-p-1)/2}
\,d\nu_{n,p}(Q)\,d\mathcal{L}^{p(p+1)/2}(A).
\end{align*}
Since the integrand is independent of $Q$, define
\begin{align*}
K_{n,p}:=c_{n,p}\nu_{n,p}(\mathbb{V}_p(\mathbb{R}^n))(2\pi)^{-np/2}.
\end{align*}
Then
\begin{align*}
\mathbb{E}[\varphi(W)]
=
K_{n,p}
\int_{\mathbb{S}_{++}^p}
\varphi(A)\exp\left(-\frac{1}{2}\operatorname{tr}(A)\right)(\det A)^{(n-p-1)/2}
\,d\mathcal{L}^{p(p+1)/2}(A).
\end{align*}
[guided]
The map $Y\mapsto Y^\top Y$ forgets the orthogonal orientation of the columns of $Y$ and remembers only their Gram matrix. The correct analogue of polar coordinates is therefore not a radius and a direction, but a positive definite matrix $A=Y^\top Y$ and an orthonormal $p$-frame $Q$ describing the orientation.
Let
\begin{align*}
\mathbb{V}_{p}(\mathbb{R}^n):=\{Q\in\mathbb{R}^{n\times p}:Q^\top Q=I_p\}
\end{align*}
be the Stiefel manifold, and let $\nu_{n,p}$ denote its invariant Hausdorff measure. The Gram polar integration formula says that for every non-negative Borel map
\begin{align*}
\Phi:\mathbb{R}^{n\times p}\to[0,\infty),
\end{align*}
we have
\begin{align*}
\int_{\mathbb{R}^{n\times p}}\Phi(Y)\,d\mathcal{L}^{np}(Y)
=
c_{n,p}
\int_{\mathbb{S}_{++}^p}
\int_{\mathbb{V}_p(\mathbb{R}^n)}
\Phi(QA^{1/2})(\det A)^{(n-p-1)/2}
\,d\nu_{n,p}(Q)\,d\mathcal{L}^{p(p+1)/2}(A).
\end{align*}
Here $A^{1/2}$ is the positive definite square root of $A$, and $c_{n,p}>0$ depends only on $n$ and $p$. The determinant exponent is the crucial Jacobian contribution: it is the rectangular-matrix analogue of the radial factor $r^{m-1}$ in Euclidean polar coordinates. This is precisely the [Gram polar integration formula for rectangular matrices](/page/Gram%20Polar%20Integration%20Formula). Its hypotheses are satisfied because $\Phi$ is assumed to be non-negative and Borel, and the spaces carry the Borel structures inherited from their Euclidean embeddings.
We now apply this formula to
\begin{align*}
\Phi(Y):=\varphi(Y^\top Y)g(Y),
\end{align*}
where
\begin{align*}
\varphi:\mathbb{S}_{++}^p\to[0,\infty)
\end{align*}
is an arbitrary non-negative Borel map. Since $(QA^{1/2})^\top(QA^{1/2})=A^{1/2}Q^\top QA^{1/2}=A$, and since
\begin{align*}
\operatorname{tr}\left((QA^{1/2})^\top(QA^{1/2})\right)=\operatorname{tr}(A),
\end{align*}
the density factor becomes independent of $Q$. Therefore
\begin{align*}
\mathbb{E}[\varphi(W)]
&=
\int_{\mathbb{R}^{n\times p}}\varphi(Y^\top Y)g(Y)\,d\mathcal{L}^{np}(Y)\\
&=
c_{n,p}
\int_{\mathbb{S}_{++}^p}
\int_{\mathbb{V}_p(\mathbb{R}^n)}
\varphi(A)(2\pi)^{-np/2}\exp\left(-\frac{1}{2}\operatorname{tr}(A)\right)
(\det A)^{(n-p-1)/2}
\,d\nu_{n,p}(Q)\,d\mathcal{L}^{p(p+1)/2}(A).
\end{align*}
Because the integrand no longer depends on $Q$, the inner integral contributes only the finite Stiefel volume. Define
\begin{align*}
K_{n,p}:=c_{n,p}\nu_{n,p}(\mathbb{V}_p(\mathbb{R}^n))(2\pi)^{-np/2}.
\end{align*}
Then
\begin{align*}
\mathbb{E}[\varphi(W)]
=
K_{n,p}
\int_{\mathbb{S}_{++}^p}
\varphi(A)\exp\left(-\frac{1}{2}\operatorname{tr}(A)\right)(\det A)^{(n-p-1)/2}
\,d\mathcal{L}^{p(p+1)/2}(A).
\end{align*}
This identity identifies the density up to the single normalising constant $K_{n,p}$.
[/guided]
[/step]
[step:Evaluate the normalising constant using the multivariate gamma integral]
Set $\varphi\equiv 1$ in the preceding identity. Since $\mathbb{E}[1]=1$,
\begin{align*}
1
=
K_{n,p}
\int_{\mathbb{S}_{++}^p}
\exp\left(-\frac{1}{2}\operatorname{tr}(A)\right)(\det A)^{(n-p-1)/2}
\,d\mathcal{L}^{p(p+1)/2}(A).
\end{align*}
Use the change of variables
\begin{align*}
B:\mathbb{S}_{++}^p&\to\mathbb{S}_{++}^p\\
A&\mapsto \frac{1}{2}A.
\end{align*}
The inverse map is $A=2B$, and because $\mathbb{S}^p$ has dimension $p(p+1)/2$,
\begin{align*}
d\mathcal{L}^{p(p+1)/2}(A)=2^{p(p+1)/2}\,d\mathcal{L}^{p(p+1)/2}(B).
\end{align*}
Also
\begin{align*}
\operatorname{tr}(A)=2\operatorname{tr}(B),
\qquad
\det A=2^p\det B.
\end{align*}
Therefore
\begin{align*}
\int_{\mathbb{S}_{++}^p}
\exp\left(-\frac{1}{2}\operatorname{tr}(A)\right)(\det A)^{(n-p-1)/2}
\,d\mathcal{L}^{p(p+1)/2}(A)
&=
2^{p(n-p-1)/2+p(p+1)/2}
\int_{\mathbb{S}_{++}^p}
e^{-\operatorname{tr}(B)}(\det B)^{(n-p-1)/2}
\,d\mathcal{L}^{p(p+1)/2}(B)\\
&=
2^{np/2}\Gamma_p(n/2),
\end{align*}
where the last equality is the defining integral for the [multivariate gamma function](/page/Multivariate%20Gamma%20Function) at $n/2$:
\begin{align*}
\Gamma_p(n/2)
=
\int_{\mathbb{S}_{++}^p}
e^{-\operatorname{tr}(B)}(\det B)^{(n-p-1)/2}
\,d\mathcal{L}^{p(p+1)/2}(B).
\end{align*}
The hypothesis $n>p-1$ is exactly the condition under which this integral is finite.
\begin{align*}
\end{align*}
Thus
\begin{align*}
K_{n,p}=\frac{1}{2^{np/2}\Gamma_p(n/2)}.
\end{align*}
Hence, for $\Sigma=I_p$, $W$ has density
\begin{align*}
f_{I_p}: \mathbb{S}_{++}^p&\to[0,\infty)\\
A&\mapsto
\frac{(\det A)^{(n-p-1)/2}\exp\left(-\frac{1}{2}\operatorname{tr}(A)\right)}
{2^{np/2}\Gamma_p(n/2)}.
\end{align*}
[guided]
We now determine the one remaining constant. The previous step showed that, for every non-negative Borel map $\varphi:\mathbb{S}_{++}^p\to[0,\infty)$,
\begin{align*}
\mathbb{E}[\varphi(W)]
=
K_{n,p}
\int_{\mathbb{S}_{++}^p}
\varphi(A)\exp\left(-\frac{1}{2}\operatorname{tr}(A)\right)(\det A)^{(n-p-1)/2}
\,d\mathcal{L}^{p(p+1)/2}(A).
\end{align*}
Take $\varphi\equiv 1$. Since the expectation of the constant random variable $1$ is $1$, we obtain
\begin{align*}
1
=
K_{n,p}
\int_{\mathbb{S}_{++}^p}
\exp\left(-\frac{1}{2}\operatorname{tr}(A)\right)(\det A)^{(n-p-1)/2}
\,d\mathcal{L}^{p(p+1)/2}(A).
\end{align*}
To compare this integral with the defining integral for $\Gamma_p(n/2)$, define the linear change of variables
\begin{align*}
B:\mathbb{S}_{++}^p&\to\mathbb{S}_{++}^p\\
A&\mapsto \frac{1}{2}A.
\end{align*}
Its inverse is $A=2B$. Since the real [vector space](/page/Vector%20Space) $\mathbb{S}^p$ has dimension $p(p+1)/2$, scalar multiplication by $2$ multiplies Lebesgue measure on the independent symmetric entries by
\begin{align*}
2^{p(p+1)/2}.
\end{align*}
Thus
\begin{align*}
d\mathcal{L}^{p(p+1)/2}(A)=2^{p(p+1)/2}\,d\mathcal{L}^{p(p+1)/2}(B).
\end{align*}
The trace and determinant transform as
\begin{align*}
\operatorname{tr}(A)=2\operatorname{tr}(B),
\qquad
\det A=2^p\det B.
\end{align*}
Therefore
\begin{align*}
&\int_{\mathbb{S}_{++}^p}
\exp\left(-\frac{1}{2}\operatorname{tr}(A)\right)(\det A)^{(n-p-1)/2}
\,d\mathcal{L}^{p(p+1)/2}(A)\\
&\qquad=
2^{p(n-p-1)/2+p(p+1)/2}
\int_{\mathbb{S}_{++}^p}
e^{-\operatorname{tr}(B)}(\det B)^{(n-p-1)/2}
\,d\mathcal{L}^{p(p+1)/2}(B)\\
&\qquad=
2^{np/2}\Gamma_p(n/2),
\end{align*}
where the final equality is the defining integral for the [multivariate gamma function](/page/Multivariate%20Gamma%20Function). The hypothesis $n>p-1$ is exactly the finiteness condition for this integral. Substituting this value into the normalization identity gives
\begin{align*}
K_{n,p}=\frac{1}{2^{np/2}\Gamma_p(n/2)}.
\end{align*}
Hence, in the identity scale case, $W$ has density
\begin{align*}
f_{I_p}: \mathbb{S}_{++}^p&\to[0,\infty)\\
A&\mapsto
\frac{(\det A)^{(n-p-1)/2}\exp\left(-\frac{1}{2}\operatorname{tr}(A)\right)}
{2^{np/2}\Gamma_p(n/2)}.
\end{align*}
[/guided]
[/step]
[step:Transfer the identity scale density by congruence]
Now let $\Sigma\in\mathbb{S}_{++}^p$. Because the Wishart law is determined by the joint distribution of the independent Gaussian sample vectors, it suffices to compute the law for any representative construction of such a sample. Let $\Sigma^{1/2}\in\mathbb{S}_{++}^p$ denote the positive definite square root of $\Sigma$. Let
\begin{align*}
Z_1,\dots,Z_n:(\Omega,\mathcal{F},\mathbb{P})\to\mathbb{R}^p
\end{align*}
be independent random vectors with common distribution $\mathcal{N}_p(0,I_p)$, and set
\begin{align*}
X_i:=\Sigma^{1/2}Z_i,\qquad 1\le i\le n.
\end{align*}
Then $X_i\sim\mathcal{N}_p(0,\Sigma)$ and
\begin{align*}
W=\sum_{i=1}^n X_iX_i^\top
=\Sigma^{1/2}\left(\sum_{i=1}^n Z_iZ_i^\top\right)\Sigma^{1/2}.
\end{align*}
Define
\begin{align*}
U:=\sum_{i=1}^n Z_iZ_i^\top.
\end{align*}
By the identity scale case, $U$ has density $f_{I_p}$ on $\mathbb{S}_{++}^p$.
Define the congruence map
\begin{align*}
T_\Sigma:\mathbb{S}^p&\to\mathbb{S}^p\\
A&\mapsto \Sigma^{1/2}A\Sigma^{1/2}.
\end{align*}
Then $W=T_\Sigma(U)$. The Jacobian determinant of $T_\Sigma$ on the vector space $\mathbb{S}^p$ is
\begin{align*}
|\det T_\Sigma|=(\det\Sigma)^{(p+1)/2}.
\end{align*}
Indeed, diagonalising $\Sigma^{1/2}=O^\top DO$ with $O$ orthogonal and $D=\operatorname{diag}(d_1,\dots,d_p)$, conjugation by $O$ has Jacobian absolute value $1$ on $\mathbb{S}^p$, while $A\mapsto DAD$ multiplies the diagonal coordinate $A_{ii}$ by $d_i^2$ and the off-diagonal coordinate $A_{ij}$ by $d_id_j$ for $i<j$. Hence the product of coordinate scaling factors is
\begin{align*}
\prod_{i=1}^p d_i^2\prod_{1\le i<j\le p}d_id_j
=
\prod_{i=1}^p d_i^{p+1}
=
(\det \Sigma^{1/2})^{p+1}
=
(\det\Sigma)^{(p+1)/2}.
\end{align*}
Thus the inverse transformation $A\mapsto \Sigma^{-1/2}A\Sigma^{-1/2}$ contributes the factor $(\det\Sigma)^{-(p+1)/2}$ to the density.
[guided]
For the general covariance matrix, the idea is to reduce to the identity case by a linear transformation. Since the theorem asserts a distributional statement, it is enough to build a sample with the same joint law as the given independent $\mathcal{N}_p(0,\Sigma)$ sample and compute the law of its Gram sum. Let $\Sigma^{1/2}\in\mathbb{S}_{++}^p$ be the unique positive definite square root of $\Sigma$. If
\begin{align*}
Z_1,\dots,Z_n:(\Omega,\mathcal{F},\mathbb{P})\to\mathbb{R}^p
\end{align*}
are independent random vectors with distribution $\mathcal{N}_p(0,I_p)$ and we define
\begin{align*}
X_i:=\Sigma^{1/2}Z_i,
\end{align*}
then each $X_i$ is Gaussian with mean $0$ and covariance
\begin{align*}
\mathbb{E}[X_iX_i^\top]
=
\Sigma^{1/2}\mathbb{E}[Z_iZ_i^\top]\Sigma^{1/2}
=
\Sigma^{1/2}I_p\Sigma^{1/2}
=
\Sigma.
\end{align*}
Therefore $X_i\sim\mathcal{N}_p(0,\Sigma)$.
Define
\begin{align*}
U:=\sum_{i=1}^n Z_iZ_i^\top.
\end{align*}
Then the corresponding Wishart matrix satisfies
\begin{align*}
W
=
\sum_{i=1}^n X_iX_i^\top
=
\sum_{i=1}^n \Sigma^{1/2}Z_iZ_i^\top\Sigma^{1/2}
=
\Sigma^{1/2}U\Sigma^{1/2}.
\end{align*}
Thus $W$ is the image of the identity-scale Wishart matrix $U$ under the congruence map
\begin{align*}
T_\Sigma:\mathbb{S}^p&\to\mathbb{S}^p\\
A&\mapsto \Sigma^{1/2}A\Sigma^{1/2}.
\end{align*}
We need the Jacobian of this [linear map](/page/Linear%20Map) on the vector space $\mathbb{S}^p$, not on all $p\times p$ matrices. Diagonalise $\Sigma^{1/2}=O^\top DO$, where $O$ is orthogonal and
\begin{align*}
D=\operatorname{diag}(d_1,\dots,d_p)
\end{align*}
with $d_i>0$. The maps $A\mapsto OAO^\top$ and $A\mapsto O^\top AO$ are orthogonal linear maps on $\mathbb{S}^p$ with respect to the Frobenius inner product, so their Jacobian absolute values are $1$. It remains to compute the Jacobian of
\begin{align*}
A\mapsto DAD.
\end{align*}
In the independent coordinates $A_{ii}$ and $A_{ij}$ for $i<j$, this map sends
\begin{align*}
A_{ii}&\mapsto d_i^2A_{ii},\\
A_{ij}&\mapsto d_id_jA_{ij}.
\end{align*}
Hence the absolute Jacobian determinant is
\begin{align*}
\prod_{i=1}^p d_i^2\prod_{1\le i<j\le p}d_id_j
=
\prod_{i=1}^p d_i^{p+1}
=
(\det \Sigma^{1/2})^{p+1}
=
(\det\Sigma)^{(p+1)/2}.
\end{align*}
Therefore the inverse change of variables from $W$ back to $U$ contributes
\begin{align*}
(\det\Sigma)^{-(p+1)/2}
\end{align*}
to the density.
[/guided]
[/step]
[step:Compute the transformed density]
For $A\in\mathbb{S}_{++}^p$, the density of $W=T_\Sigma(U)$ is
\begin{align*}
f_\Sigma(A)
&=
f_{I_p}(\Sigma^{-1/2}A\Sigma^{-1/2})(\det\Sigma)^{-(p+1)/2}.
\end{align*}
Using multiplicativity of determinant,
\begin{align*}
\det(\Sigma^{-1/2}A\Sigma^{-1/2})
=
(\det\Sigma)^{-1}\det A.
\end{align*}
Using cyclic invariance of trace,
\begin{align*}
\operatorname{tr}(\Sigma^{-1/2}A\Sigma^{-1/2})
=
\operatorname{tr}(\Sigma^{-1}A).
\end{align*}
Therefore
\begin{align*}
f_\Sigma(A)
&=
\frac{\left((\det\Sigma)^{-1}\det A\right)^{(n-p-1)/2}
\exp\left(-\frac{1}{2}\operatorname{tr}(\Sigma^{-1}A)\right)}
{2^{np/2}\Gamma_p(n/2)}
(\det\Sigma)^{-(p+1)/2}\\
&=
\frac{(\det A)^{(n-p-1)/2}
\exp\left(-\frac{1}{2}\operatorname{tr}(\Sigma^{-1}A)\right)}
{2^{np/2}(\det\Sigma)^{n/2}\Gamma_p(n/2)}.
\end{align*}
This is the asserted density with respect to the [Lebesgue measure](/page/Lebesgue%20Measure) $\mathcal{L}^{p(p+1)/2}$ on the independent entries of $\mathbb{S}^p$, and the proof is complete.
[guided]
We now apply the ordinary change-of-variables formula for the invertible linear map
\begin{align*}
T_\Sigma:\mathbb{S}^p&\to\mathbb{S}^p\\
A&\mapsto \Sigma^{1/2}A\Sigma^{1/2}.
\end{align*}
The previous step computed its absolute Jacobian determinant on the vector space $\mathbb{S}^p$ as
\begin{align*}
|\det T_\Sigma|=(\det\Sigma)^{(p+1)/2}.
\end{align*}
Therefore the inverse change of variables $A\mapsto\Sigma^{-1/2}A\Sigma^{-1/2}$ contributes the reciprocal factor $(\det\Sigma)^{-(p+1)/2}$. Since $W=T_\Sigma(U)$ and $U$ has identity-scale density $f_{I_p}$, the density of $W$ is
\begin{align*}
f_\Sigma(A)
&=
f_{I_p}(\Sigma^{-1/2}A\Sigma^{-1/2})(\det\Sigma)^{-(p+1)/2}
\end{align*}
for $A\in\mathbb{S}_{++}^p$.
We now simplify this expression. By multiplicativity of determinant,
\begin{align*}
\det(\Sigma^{-1/2}A\Sigma^{-1/2})
=
(\det\Sigma^{-1/2})(\det A)(\det\Sigma^{-1/2})
=
(\det\Sigma)^{-1}\det A.
\end{align*}
By cyclic invariance of trace,
\begin{align*}
\operatorname{tr}(\Sigma^{-1/2}A\Sigma^{-1/2})
=
\operatorname{tr}(A\Sigma^{-1})
=
\operatorname{tr}(\Sigma^{-1}A).
\end{align*}
Substituting these two identities into the formula for $f_\Sigma$ gives
\begin{align*}
f_\Sigma(A)
&=
\frac{\left((\det\Sigma)^{-1}\det A\right)^{(n-p-1)/2}
\exp\left(-\frac{1}{2}\operatorname{tr}(\Sigma^{-1}A)\right)}
{2^{np/2}\Gamma_p(n/2)}
(\det\Sigma)^{-(p+1)/2}\\
&=
\frac{(\det A)^{(n-p-1)/2}
\exp\left(-\frac{1}{2}\operatorname{tr}(\Sigma^{-1}A)\right)}
{2^{np/2}(\det\Sigma)^{(n-p-1)/2+(p+1)/2}\Gamma_p(n/2)}\\
&=
\frac{(\det A)^{(n-p-1)/2}
\exp\left(-\frac{1}{2}\operatorname{tr}(\Sigma^{-1}A)\right)}
{2^{np/2}(\det\Sigma)^{n/2}\Gamma_p(n/2)}.
\end{align*}
This is the asserted density with respect to the [Lebesgue measure](/page/Lebesgue%20Measure) $\mathcal{L}^{p(p+1)/2}$ on the independent entries of $\mathbb{S}^p$.
[/guided]
[/step]