Let ${n}$ be a large integer, and let ${M_n}$ be the Gaussian Unitary Ensemble (GUE), i.e. the random Hermitian matrix with probability distribution

$\displaystyle C_n e^{-\hbox{tr}(M_n^2)/2} dM_n$

where ${dM_n}$ is a Haar measure on Hermitian matrices and ${C_n}$ is the normalisation constant required to make the distribution of unit mass. The eigenvalues ${\lambda_1 < \ldots < \lambda_n}$ of this matrix are then a coupled family of ${n}$ real random variables. For any ${1 \leq k \leq n}$, we can define the ${k}$-point correlation function ${\rho_k( x_1,\ldots,x_k )}$ to be the unique symmetric measure on ${{\bf R}^k}$ such that

$\displaystyle \int_{{\bf R}^k} F(x_1,\ldots,x_k) \rho_k(x_1,\ldots,x_k) = {\bf E} \sum_{1 \leq i_1 < \ldots < i_k \leq n} F( \lambda_{i_1}, \ldots, \lambda_{i_k} ).$

A standard computation (given for instance in these lecture notes of mine) gives the Ginebre formula

$\displaystyle \rho_n(x_1,\ldots,x_n) = C'_n (\prod_{1 \leq i < j \leq n} |x_i-x_j|^2) e^{-\sum_{j=1}^n |x_j|^2/2}.$

for the ${n}$-point correlation function, where ${C'_n}$ is another normalisation constant. Using Vandermonde determinants, one can rewrite this expression in determinantal form as

$\displaystyle \rho_n(x_1,\ldots,x_n) = C''_n \det(K_n(x_i,x_j))_{1 \leq i, j \leq n}$

where the kernel ${K_n}$ is given by

$\displaystyle K_n(x,y) := \sum_{k=0}^{n-1} \phi_k(x) \phi_k(y)$

where ${\phi_k(x) := P_k(x) e^{-x^2/4}}$ and ${P_0, P_1, \ldots}$ are the (${L^2}$-normalised) Hermite polynomials (thus the ${\phi_k}$ are an orthonormal family, with each ${P_k}$ being a polynomial of degree ${k}$). Integrating out one or more of the variables, one is led to the Gaudin-Mehta formula

$\displaystyle \rho_k(x_1,\ldots,x_k) = \det(K_n(x_i,x_j))_{1 \leq i, j \leq k}. \ \ \ \ \ (1)$

(In particular, the normalisation constant ${C''_n}$ in the previous formula turns out to simply be equal to ${1}$.) Again, see these lecture notes for details.

The functions ${\phi_k(x)}$ can be viewed as an orthonormal basis of eigenfunctions for the harmonic oscillator operator

$\displaystyle L \phi := (-\frac{d^2}{dx^2} + \frac{x^2}{4})\phi;$

indeed it is a classical fact that

$\displaystyle L \phi_k = (k + \frac{1}{2}) \phi_k.$

As such, the kernel ${K_n}$ can be viewed as the integral kernel of the spectral projection operator ${1_{(-\infty,n+\frac{1}{2}]}(L)}$.

From (1) we see that the fine-scale structure of the eigenvalues of GUE are controlled by the asymptotics of ${K_n}$ as ${n \rightarrow \infty}$. The two main asymptotics of interest are given by the following lemmas:

Lemma 1 (Asymptotics of ${K_n}$ in the bulk) Let ${x_0 \in (-2,2)}$, and let ${\rho_{sc}(x_0) := \frac{1}{2\pi} (4-x_0^2)^{1/2}_+}$ be the semicircular law density at ${x_0}$. Then, we have

$\displaystyle K_n( x_0 \sqrt{n} + \frac{y}{\sqrt{n} \rho_{sc}(x_0)}, x_0 \sqrt{n} + \frac{z}{\sqrt{n} \rho_{sc}(x_0)} )$

$\displaystyle \rightarrow \frac{\sin(\pi(y-z))}{\pi(y-z)} \ \ \ \ \ (2)$

as ${n \rightarrow \infty}$ for any fixed ${y,z \in {\bf R}}$ (removing the singularity at ${y=z}$ in the usual manner).

Lemma 2 (Asymptotics of ${K_n}$ at the edge) We have

$\displaystyle K_n( 2\sqrt{n} + \frac{y}{n^{1/6}}, 2\sqrt{n} + \frac{z}{n^{1/6}} )$

$\displaystyle \rightarrow \frac{Ai(y) Ai'(z) - Ai'(y) Ai(z)}{y-z} \ \ \ \ \ (3)$

as ${n \rightarrow \infty}$ for any fixed ${y,z \in {\bf R}}$, where ${Ai}$ is the Airy function

$\displaystyle Ai(x) := \frac{1}{\pi} \int_0^\infty \cos( \frac{t^3}{3} + tx )\ dt$

and again removing the singularity at ${y=z}$ in the usual manner.

The proof of these asymptotics usually proceeds via computing the asymptotics of Hermite polynomials, together with the Christoffel-Darboux formula; this is for instance the approach taken in the previous notes. However, there is a slightly different approach that is closer in spirit to the methods of semi-classical analysis, which was briefly mentioned in the previous notes but not elaborated upon. For sake of completeness, I am recording some notes on this approach here, although to focus on the main ideas I will not be completely rigorous in the derivation (ignoring issues such as convegence of integrals or of operators, or (removable) singularities in kernels caused by zeroes in the denominator).

— 1. The bulk asymptotics —

We begin with the bulk asymptotics, Lemma 1. Fix ${x_0}$ in the bulk region ${(-2,2)}$. Applying the change of variables

$\displaystyle x = x_0 \sqrt{n} + \frac{y}{\sqrt{n} \rho_{sc}(x_0)}$

we see that the harmonic oscillator ${L}$ becomes

$\displaystyle - n \rho_{sc}(x_0)^2 \frac{d^2}{dy^2} + \frac{1}{4} (x_0 \sqrt{n} + \frac{y}{\sqrt{n} \rho_{sc}(x_0)})^2$

Since ${K_n}$ is the integral kernel of the spectral projection to the region ${L \leq n+\frac{1}{2}}$, we conclude that the left-hand side of (2) (as a function of ${y,z}$) is the integral kernel of the spectral projection to the region

$\displaystyle - n \rho_{sc}(x_0)^2 \frac{d^2}{dy^2} + \frac{1}{4} (x_0 \sqrt{n} + \frac{y}{\sqrt{n} \rho_{sc}(x_0)})^2 \leq n + \frac{1}{2}.$

Isolating out the top order terms in ${n}$, we can rearrange this as

$\displaystyle -\frac{d^2}{dy^2} \leq \pi^2 + o(1).$

Thus, in the limit ${n \rightarrow \infty}$, we expect (heuristically, at least) that the left-hand side of (2) to converge as ${n \rightarrow \infty}$ to the integral kernel of the spectral projection to the region

$\displaystyle -\frac{d^2}{dy^2} \leq \pi^2.$

Introducing the Fourier dual variable ${\xi}$ to ${y}$, as manifested by the Fourier transform

$\displaystyle \hat f(\xi) = \int_{\bf R} e^{-2\pi i \xi y} f(y)\ dy$

and its inverse

$\displaystyle \check F(y) = \int_{\bf R} e^{2\pi i \xi y} F(\xi)\ d\xi,$

then we (heuristically) have ${\frac{d}{dy} = 2\pi i \xi}$, and so we are now projecting to the region

$\displaystyle |\xi|^2 \leq 1/4, \ \ \ \ \ (4)$

i.e. we are restricting the Fourier variable to the interval ${[-1/2,1/2]}$. Back in physical space, the associated projection ${P}$ thus takes the form

$\displaystyle P f(y) = \int_{[-1/2,1/2]} e^{2\pi i \xi y} \hat f(\xi)\ d\xi$

$\displaystyle = \int_{\bf R} \int_{[-1/2,1/2]} e^{2\pi i \xi y} e^{-2\pi i \xi z}\ d\xi f(z)\ dz$

$\displaystyle = \int_{\bf R} \frac{\sin(\pi(y-z))}{y-z} f(z)\ dz$

and the claim follows.

Remark 1 From a semiclassical perspective, the original spectral projection ${L \leq n+\frac{1}{2}}$ can be expressed in phase space (using the dual frequency variable ${\eta}$ to ${x}$) as the ellipse

$\displaystyle 4 \pi^2 \eta^2 + \frac{x^2}{4} \leq n+\frac{1}{2} \ \ \ \ \ (5)$

which after the indicated change of variables becomes the elongated ellipse

$\displaystyle \xi^2 + \frac{1}{2n \rho_{sc}(x_0)(4-x_0^2)} y + \frac{1}{4n^2 \rho_{sc}(x_0)^2 (4-x_0^2)} y^2$

$\displaystyle \leq \frac{1}{4} + \frac{1}{2n (4-x_0^2)}$

which converges (in some suitably weak sense) to the strip (4) as ${n \rightarrow \infty}$.

— 2. The edge asymptotics —

A similar (heuristic) argument gives the edge asymptotics, Lemma 2. Starting with the change of variables

$\displaystyle x = 2 \sqrt{n} + \frac{y}{n^{1/6}}$

the harmonic oscillator ${L}$ now becomes

$\displaystyle - n^{1/3} \frac{d^2}{dy^2} + \frac{1}{4} (2 \sqrt{n} + \frac{y}{n^{1/6}})^2.$

Thus, the left-hand side of (3) becomes the kernel of the spectral projection to the region

$\displaystyle - n^{1/3} \frac{d^2}{dy^2} + \frac{1}{4} (2 \sqrt{n} + \frac{y}{n^{1/6}})^2 \leq n + \frac{1}{2}.$

Expanding out, computing all terms of size ${n^{1/3}}$ or larger, and rearranging, this (heuristically) becomes

$\displaystyle - \frac{d^2}{dy^2} + y \leq o(1)$

and so, heuristically at least, we expect (3) to converge to the kernel of the projection to the region

$\displaystyle - \frac{d^2}{dy^2} + y \leq 0. \ \ \ \ \ (6)$

To compute this, we again pass to the Fourier variable ${\xi}$, converting the above to

$\displaystyle 4 \pi^2 \xi^2 + \frac{1}{2\pi i} \frac{d}{d\xi} \leq 0$

using the usual Fourier-analytic correspondences between multiplication and differentiation. If we then use the integrating factor transformation

$\displaystyle F(\xi) \mapsto e^{8 \pi^3 i \xi^3 / 3} F(\xi)$

we can convert the above region to

$\displaystyle \frac{1}{2\pi i} \frac{d}{d\xi} \leq 0$

which on undoing the Fourier transformation becomes

$\displaystyle y \leq 0,$

and the spectral projection operation for this is simply the spatial multiplier ${1_{(-\infty,0]}}$. Thus, informally at least, we see that the spectral projection ${P}$ to the region (6) is given by the formula

$\displaystyle P = M^{-1} 1_{(-\infty,0]} M$

where the Fourier multiplier ${M}$ is given by the formula

$\displaystyle \widehat{Mf}(\xi) := e^{8 \pi^3 i \xi^3 / 3} \hat f(\xi).$

In other words (ignoring issues about convergence of the integrals),

$\displaystyle Mf(y) = \int_{\bf R} (\int_{\bf R} e^{2\pi i y \xi} e^{8 \pi^3 i \xi^3 / 3} e^{-2\pi i z \xi}\ d\xi) f(z)\ dz$

$\displaystyle = 2 \int_{\bf R} (\int_0^\infty \cos( 2\pi (y-z) \xi + 8 \pi^3 \xi^3 / 3 )\ d\xi) f(z)\ dz$

$\displaystyle = \frac{1}{\pi} \int_{\bf R} (\int_0^\infty \cos( t (y-z) + t^3 / 3 )\ dt) f(z)\ dz$

$\displaystyle = \int_{\bf R} Ai(y-z) f(z)\ dz$

and similarly

$\displaystyle M^{-1} f(z) = \int_{\bf R} Ai(y-z) f(y)\ dy$

(this reflects the unitary nature of ${M}$). We thus see (formally, at least) that

$\displaystyle P f(y) = \int_{\bf R} (\int_{(-\infty,0]} Ai(y-w) Ai(z-w)\ dw) f(z)\ dz.$

To simplify this expression we perform some computations closely related to the ones above. From the Fourier representation

$\displaystyle Ai(y) = \frac{1}{\pi} \int_0^\infty \cos(ty + t^3/3)\ dt$

$\displaystyle = \int_{\bf R} e^{2\pi i y \xi} e^{8 \pi i \xi^3/3}\ d\xi$

we see that

$\displaystyle \widehat{Ai}(\xi) = e^{8 \pi^3 i \xi^3/3}$

which means that

$\displaystyle (4 \pi^2 \xi^2 + \frac{1}{2\pi i} \frac{d}{d\xi}) \widehat{Ai}(\xi) = 0$

and thus

$\displaystyle (- \frac{d^2}{dy^2} + y) Ai(y) = 0,$

thus ${Ai}$ obeys the Airy equation

$\displaystyle Ai''(y) = y Ai(y).$

Using this, one soon computes that

$\displaystyle \frac{d}{dw} \frac{Ai(y-w) Ai'(z-w) - Ai'(y-w) Ai(z-w)}{y-z} = Ai(y-w) Ai(z-w).$

Also, stationary phase asymptotics tell us that ${Ai(y)}$ decays exponentially fast as ${y \rightarrow +\infty}$, and hence ${Ai(y-w)}$ decays exponentially fast as ${w \rightarrow -\infty}$ for fixed ${y}$; similarly for ${Ai'(z-w), Ai'(y-w), Ai(z-w)}$. From the fundamental theorem of calculus, we conclude that

$\displaystyle \int_{(-\infty,0]} Ai(y-w) Ai(z-w)\ dw = \frac{Ai(y) Ai'(z) - Ai'(y) Ai(z)}{y-z},$

(this is a continuous analogue of the Christoffel-Darboux formula), and the claim follows.

Remark 2 As in the bulk case, one can take a semi-classical analysis perspective and track what is going on in phase space. With the scaling we have selected, the ellipse (5) has become

$\displaystyle 4 \pi^2 n^{1/3} \xi^2 + \frac{(2\sqrt{n} + y/n^{1/6})^2}{4} \leq n+\frac{1}{2},$

which we can rearrange as the eccentric ellipse

$\displaystyle 4\pi^2 \xi^2 + y \leq \frac{1}{2n^{1/3}} - \frac{y^2}{4 n^{2/3}}$

which is converging as ${n \rightarrow \infty}$ to the parabolic region

$\displaystyle 4\pi^2 \xi^2 + y \leq 0$

which can then be shifted to the half-plane ${y \leq 0}$ by the parabolic shear transformation ${(y,\xi) \mapsto (y+4\pi^2\xi^2,\xi)}$, which is the canonical relation of the Fourier multiplier ${M}$. (The rapid decay of the kernel ${Ai}$ of ${M}$ at ${+\infty}$ is then reflected in the fact that this transformation only shears to the right and not the left.)

Remark 3 Presumably one should also be able to apply the same heuristics to other invariant ensembles, such as those given by probability distributions of the form

$\displaystyle C_n e^{-\hbox{tr}(P(M_n))} dM_n$

for some potential function ${P}$. Certainly one can soon get to an orthogonal polynomial formulation of the determinantal kernel for such ensembles, but I do not know if the projection operators for such kernels can be viewed as spectral projections to a phase space region as was the case for GUE. But if one could do this, this would provide a heuristic explanation as to the universality phenomenon for such ensembles, as Taylor expansion shows that all (reasonably smooth) regions of phase space converge to universal limits (such as a strip or paraboloid) after rescaling around either a non-critical point or a critical point of the region with the appropriate normalisation.