# Compute the correlation of eigenvalues of a Gaussian Random Matrix

## Background

Let $$M$$ be a $$N\times N$$ random symmetric matrix, which is drawn from a Gaussian distribution of mean $$\mu=0$$ and variance $$\sigma^2=1/N$$.

The first level correlation is defined by: $$\begin{equation} \rho(\lambda)=\left \langle \text{Tr}\delta(\lambda-M)\right \rangle\quad\quad\quad\quad\quad(1) \end{equation}$$ Which is nothing more than the distribution of eigenvalues: following my routine below I can compute this quantity. The brackets represents the average over the gaussian distribution.

RandomMatrix[n_] :=
RandomMatrix[n] =
RandomVariate[NormalDistribution[0, 1/Sqrt[n]], {n, n}];
n := 1500
M := 1/Sqrt*(RandomMatrix[n] + Transpose@RandomMatrix[n])
\[Rho] = Eigenvalues[M];
Histogram[\[Rho]] In the large limit $$N>>1$$ equation $$(1)$$ will describe a semi-circular law.

## Question

The two-level correlation function is defined by: $$\begin{equation} \rho^{(2)}(\lambda, \mu)=\left\langle\frac{1}{N} \operatorname{Tr} \delta(\lambda-M) \frac{1}{N} \operatorname{Tr} \delta(\mu-M)\right\rangle\quad\quad\quad\quad\quad(2) \end{equation}$$ I interpret this quantity is as how correlated two given eigenvalues will be. Following the previous routine, my $$\rho$$ will be a list of the $$N$$ eigenvalues of $$M$$. How can I numerically compute equation $$(2)$$? Isn't it impossible to tell how correlated two elements are in one same list? If I generate several lists, it has to be for different matrices, which does not work anymore.

Given one matrix $$M$$, how can I compute the correlation between two eigenvalues?

Any help or remark is appreciated!

• Your definition of the square random matrix gives the same random matrix each time—is that what you intend? Here you are looking at distributions for a single random matrix and not a sample of many random matrices of a given size. Also, you don’t need to use set-delayed in your definition for n. The set-delaleyed for M is the right idea of you want different random matrices each time, but won’t work with your definition of RandomMatrix. Jun 21, 2021 at 8:27
• Yes it was intended. For two reasons: first because when I want to symmetrize M I need RandomMatrix[n] to stay the same when I call it twice: RandomMatrix[n]+Transpose@RandomMatrix[n]. The second reason was to emphasize the fact that when I speak of computing the correlation function, it is for a single matrix $M$ and not a collection of random matrices. For example when I computed the eigenvalues it was for one matrix, I didn't average it over many matrices.
– Matt
Jun 21, 2021 at 9:45
• However I understand from your comment that I should not use $n$ as the set-delayed, I will use a different parameter next time, thank you.
– Matt
Jun 21, 2021 at 9:47
• One point is hard to understand from your post, why don't you compute according to Eq.(1), but do instead a matrix diagonalization? Jun 21, 2021 at 10:08
• One minor point - the diagonal elements of your RandomMatrix don't have the same distribution as the off-diagonal element (Variance[Diagonal[M]] / Variance[Flatten[M]] is about 2). The off-diagonals are the sum of 2 independent RVs, while the diagonal is the sum of the same RV. Jun 21, 2021 at 10:58

The question is: How do you calculate the two-level correlation function for Gaussian Orthogonal Ensemble (GOE)? One way would be to create a 2d histogram. You create bins of some size $$\Delta\lambda \times \Delta\mu$$ and then for each bin add the number of eigenvalues in the interval $$\lambda + \Delta\lambda$$ to the number of eigenvalues in the interval $$\mu + \Delta\mu$$. The bins should be small enough to see the structure you are interested in, but large enough to capture enough eigenvalues to get decent numerics. I'm not sure it would be easy to get good numerics with just one matrix. For the self-averaging property of the GOE to give a smooth two-point correlation function, the matrix would have to be very big, I guess maybe $$10^8$$ by $$10^8$$, but it would take a long time and a lot of RAM to diagonalize such a large matrix. It would be better to do an ensemble average over smaller matrices.
I don't know what 2d histogram features are available in Mathematica. I looked into calculating a 2d histogram of the two point correlation function of the Gaussian Unitary Ensemble (GUE) using Matlab or Octave which have some built-in 2d histogram functions. But in the end I decided to use the formula for $$\rho^{(2)}(\lambda,\mu)$$ explicitly with gaussian-broadened delta functions. The widths of gaussian functions used to represent the delta functions correspond to the bin sizes in the histogram approach. Here is the resulting publication: https://link.springer.com/article/10.1007/s13538-020-00802-6