English: A reconstruction of Figure 1 of Moments of Wishart-Laguerre and Jacobi ensembles of random matrices (Livan 2011)
https://arxiv.org/pdf/1103.2638.pdf
```python
import numpy as np
import matplotlib.pyplot as plt
- 1 for LOE, 2 for LUE, 4 for LSE
betas = 1, 2, 4
NMs = [(8, 15)]
- Choose number of samples
Nmatr = 100000
Es = {}
for n, m in NMs:
for beta in betas:
if beta == 1: # Wishart Orthogonal Ensemble
X = np.random.randn(Nmatr, n, m)
M = np.einsum('ijk,ilk->ijl', X, X)
E = np.linalg.eigvals(M.reshape(Nmatr, n, n)).flatten()
elif beta == 2: # Wishart Unitary Ensemble
X_real = np.random.randn(Nmatr, n, m)
X_imag = np.random.randn(Nmatr, n, m)
X = X_real + 1j * X_imag
M = np.einsum('ijk,ilk->ijl', X, X.conjugate())
E = np.linalg.eigvals(M.reshape(Nmatr, n, n)).flatten()
elif beta == 4: # Wishart Symplectic Ensemble
A = np.random.randn(Nmatr, n,m) + 1j * np.random.randn(Nmatr, n,m)
B = np.random.randn(Nmatr, n,m) + 1j * np.random.randn(Nmatr, n,m)
X = np.block([[A, B],[-np.conj(B), np.conj(A)]])
M = np.einsum('ijk,ilk->ijl', X, X.conjugate())
E = np.linalg.eigvals(M.reshape(Nmatr, 2 * n, 2 * n)).flatten()
Es[(n, m, beta)] = E
for n, m in NMs:
plt.figure(figsize=(16, 8))
legends = {1: "LOE", 2:"LUE", 4:"LSE"}
colors={1:"blue", 2:"red", 4:"green"}
for beta in betas:
color=colors[beta]
E = Es[(n, m, beta)]
xs = np.real(E) / np.sqrt(beta)
bin_heights, bin_borders, _ = plt.hist(xs, bins=500, density=True, color=color, alpha=0.1)
bin_centers = bin_borders[:-1] + np.diff(bin_borders) / 2
# Compute sliding window average
window_size = 5
window = np.ones(window_size) / window_size
smoothed_heights = np.convolve(bin_heights, window, mode='same')
# Plot sliding window average
plt.plot(bin_centers, smoothed_heights, label=legends[beta], color=color)
# Add plot labels and title
plt.xlabel('x', fontsize=14)
plt.ylabel('ρ(x)', fontsize=14)
plt.title(r'Eigenvalues $/\sqrtTemplate:\beta$, with (N, M) = {}'.format((n, m)), fontsize=18)
plt.grid(True)
plt.legend()
plt.show()
```