What is the fastest way to obtain the eigenvalues of a Wishart matrix?

Question

I would like a fast method of creating a sample of random numbers which corresponds to the eigenvalues of a Wishart matrix:

For M>N the eigenvalues \lambda_i are given by the jpd

Where $K_N$ is a normalisation constant and $\beta = 1,2$ depends on whether $X$ is real or complex. Does anyone have any suggestions for a fast way of implementing this? My goal is to play with multiple realisations of eigenvalues representing 100x100 matrices.

Many thanks in advance!

EDIT:

Ok, here is a very basic example:

n = 1000;
m = 1001;
μ = 0;
σ = 1/Sqrt[n];
AbsoluteTiming[
a = RandomVariate[NormalDistribution[μ, σ], {n, m}];
wish = a.a\[ConjugateTranspose];
ewish = Eigenvalues@wish;
]
ewishPlot = Histogram[ewish, {0.1}, "PDF"];

ηp = n σ^2 (1 + Sqrt[η])^2;
ηm = n σ^2 (1 - Sqrt[η])^2;
η = m/n;
ρmp[λ_] := 1/(2 π n σ^2 λ)Sqrt[(ηp - λ) (λ - ηm)];
ρmpPlot = Plot[ρmp[λ], {λ, 0, 4}];
Show[ρmpPlot, ewishPlot]

I guess speed won't be a problem but I would still very much like to know the best way of doing this.

Answer 1

Since all the entries of $X$ are i.i.d. normal, it can be converted into a band diagonal matrix with $\chi$-distributed components via elementary reflections. A paper discusses about this is "Matrix Models for Beta Ensembles" and can be found in this arXiv preprint.

The Wishart matrix constructed in this way is tridiagonal and can be represented as SparseArray in Mathematica. The eigenvalues can be obtained efficiently by using Eigenvalues with method "Banded".

Answer 2

One can now use MatrixPropertyDistribution[] along with WishartMatrixDistribution[] to generate the eigenvalues of a random Wishart matrix. For instance,

RandomVariate[MatrixPropertyDistribution[Eigenvalues[X], 
              X \[Distributed] WishartMatrixDistribution[4, IdentityMatrix[4]]], 5]

will generate $5$ sets of eigenvalues of random $4\times4$ Wishart matrices.