Optimal soft edge scaling variables for the Gaussian and Laguerre even $\beta$ ensembles

TL;DR

This paper investigates the large N asymptotics of eigenvalue densities at the soft edge for Gaussian and Laguerre beta ensembles, revealing a universal correction term of order N^{-2/3} when properly centered.

Contribution

It introduces optimal soft edge scaling variables for Gaussian and Laguerre beta ensembles, demonstrating universal correction behavior in the asymptotic density.

Findings

Leading correction term is O(N^{-2/3}) when scaled variables are centered properly.

The correction behavior is consistent across different parameter regimes in Laguerre ensembles.

Numerical evidence supports the theoretical asymptotic results.

Abstract

The $\beta$ ensembles are a class of eigenvalue probability densities which generalise the invariant ensembles of classical random matrix theory. In the case of the Gaussian and Laguerre weights, the corresponding eigenvalue densities are known in terms of certain $\beta$ dimensional integrals. We study the large $N$ asymptotics of the density with a soft edge scaling. In the Laguerre case, this is done with both the parameter $a$ fixed, and with $a$ proportional to $N$. It is found in all these cases that by appropriately centring the scaled variable, the leading correction term to the limiting density is $O(N^{-2/3})$. A known differential-difference recurrence from the theory of Selberg integrals allows for a numerical demonstration of this effect.

Figures

Equations

$$\prod _{l=1}^{N}w(x_l)\prod _{1\le j<k\le N}\text{absolutevalue}{x}_k-x_j$$

$$w(x)=\{\begin{array}{ll}{e}^{-x^2/2},& \text{Gaussian}\\ x^{a/2}{e}^{-x/2}{\chi }_{x>0},& \text{Laguerre}\\ x^a(1-x{)}^b{\chi }_{0<x<1},& \text{Jacobi}\end{array}$$

$$\prod _{l=1}^N{w}_{\beta }(x_l)\prod _{1\le j<k\le N}\text{absolutevalue}{x_k-x_j}^{\beta }$$

$$w_{\beta }(x)=\{\begin{array}{ll}{e}^{-\beta x^2/2},& \text{Gaussian}\\ x^{\beta a/2}{e}^{-\beta x/2}{\chi }_{x>0},& \text{Laguerre}.\end{array}$$

$${\rho }_{N,\beta }^{\#}(x)=w_{\beta }(x)\frac{N}{Z_{\beta ,N}^{\#}}{\int }_{{\mathbb{R}}^{N-1}}\prod _{l=2}^N{w}_{\beta }(x_l)\text{absolutevalue}{x-x_l}^{\beta }\prod _{2\le j<k\le N}\text{absolutevalue}{x_k-x_j}^{\beta }d{x}_2\cdots d{x}_N,$$

$$Z_{\beta ,N}^{\#}={\int }_{{\mathbb{R}}^N}\prod _{l=1}^N{w}_{\beta }(x_l)\prod _{1\le j<k\le N}\text{absolutevalue}{x_k-x_j}^{\beta }d{x}_2\cdots d{x}_N$$

$${\rho }_{N,\beta }^{\#}(x)=w_{\beta }(x)\frac{N{Z}_{\beta ,N-1}^{\#}}{Z_{\beta ,N}^{\#}}{⟨\prod _{l=2}^N\text{absolutevalue}{x-x_l}^{\beta }⟩}_{N-1}^{\#},$$

$$K_{n,\beta }(x)=-\frac{1}{(2\pi i{)}^n}{\int }_{-i\infty }^{i\infty }d{v}_1\cdots {\int }_{-i\infty }^{i\infty }d{v}_{\beta }\prod _{j=1}^{\beta }{e}^{v_j^3/3-x{v}_j}\prod _{1\le k<l\le n}\text{absolutevalue}{v_k-v_j}^{4/\beta },$$

$$\mathrm{Ai}(x)=\frac{1}{2\pi i}{\int }_{-i\infty }^{i\infty }{e}^{v^3/3-xv}dv.$$

$$\underset{N\to \infty }{\lim }\frac{1}{\sqrt{2}{N}^{1/6}}{\rho }_{N,\beta }^G\left(\sqrt{2N}+\frac{x}{\sqrt{2}{N}^{1/6}}\right)$$

$$=\frac{1}{2\pi }{\left(\frac{4\pi }{\beta }\right)}^{\beta /2}\Gamma (1+\beta /2)\prod _{j=1}^{\beta }\frac{\Gamma (1+2/\beta )}{\Gamma (1+2j/\beta )}{K}_{\beta ,\beta }(x)$$

$$=:{\rho }_{\infty ,\beta ,0}(x).$$

$$\frac{1}{\sqrt{2}{N}^{1/6}}{\rho }_{N,\beta }^G\left(\sqrt{2N}+\frac{1}{\sqrt{2}{N}^{1/6}}\left(x+\left(\frac{1}{2}-\frac{1}{\beta }\right)\frac{1}{N^{1/3}}\right)\right)$$

$$2(2N{)}^{1/3}{\rho }_{N,\beta }^{(L,a)}\left(4N+2(2N{)}^{1/3}\left(x+\frac{a}{(2N{)}^{1/3}}\right)\right)$$

$$b=\left(\frac{1}{\sqrt{1+\alpha }}+1\right){\left(\frac{\sqrt{1+\alpha }+1}{2}\right)}^3.$$

$$2(bN{)}^{1/3}{\rho }_{N,\beta }^{(L,\alpha N)}\left(N(\sqrt{1+\alpha }+1{)}^2+2(bN{)}^{1/3}\left(x+\frac{\alpha }{\sqrt{1+\alpha }}\left(\frac{1}{2}-\frac{1}{\beta }\right)\frac{1}{2(bN{)}^{1/3}}\right)\right)$$

$${\rho }_{(k)}(x_1,...,x_k)=\det [K_N(x_j,x_l){]}_{j,l=1}^k$$

$$s_x=\{\begin{array}{ll}\sqrt{2N}+x/\sqrt{2}{N}^{1/6},& \text{Gaussian}\\ 4N+2a+2(2N{)}^{1/3}x,& \text{Laguerre}\\ N(\sqrt{1+\alpha }+1{)}^2+2(bN{)}^{1/3}x,& \text{Laguerre }a=\alpha N,\end{array}$$

$${\rho }_{\infty ,\beta =2,1}(x)$$

$${\rho }_{\infty ,\beta =2,1}^G(x)$$

$${\rho }_{\infty ,\beta =2,1}^{(L,a)}(x)$$

$${\rho }_{\infty ,\beta =2,1}^{(L,\alpha N)}(x)$$

$$\displaystyle\phantom{b^{-2/3}\bigg{[}\left(\frac{3}{5}\frac{b}{\sqrt{1+\alpha}}-\frac{3}{20}\frac{\alpha^{2}}{1+\alpha}\right)x^{2}\mathrm{Ai}(x)^{2}}+\left(-\frac{2}{5}\frac{b}{\sqrt{1+\alpha}}+\frac{1}{10}\frac{\alpha^{2}}{1+\alpha}\right)x\mathrm{Ai}^{\prime}(x)^{2}\bigg{]}.$$

$$\frac{\partial s_x}{\partial x}{K}_N^{\#}(s_x,s_y)=\frac{\mathrm{Ai}(x){\mathrm{Ai}}^{\prime }(y)-{\mathrm{Ai}}^{\prime }(x)\mathrm{Ai}(y)}{x-y}+O(N^{-2/3})$$

$$\frac{\partial s_x}{\partial x}{\rho }_{N,\beta =2}^{\#}(s_x)={\rho }_{\infty ,\beta =2,0}(x)+\frac{1}{N^{2/3}}{\rho }_{\infty ,\beta =2,1}^{\#}(x)+\cdots .$$

$${⟨\prod _{l=1}^N(x-x_l{)}^n⟩}_{{\mathrm{ME}}_{\beta ,N}(e^{-\beta x^2/2})}={⟨\prod _{l=1}^n(x-i{x}_l{)}^N⟩}_{{\mathrm{ME}}_{4/\beta ,n}(e^{-x^2})}.$$

$${\rho }_{N,\beta }^G(x)=e^{-\beta x^2/2}\frac{N{Z}_{\beta ,N-1}^G}{Z_{\beta ,N}^G}{⟨\prod _{l=1}^{\beta }(x-i{x}_l{)}^{N-1}⟩}_{{\mathrm{ME}}_{4/\beta ,\beta }(e^{-x^2})}.$$

$${\rho }_{N,\beta }^G(\sqrt{2N}x)=\frac{1}{2}\frac{G_{\beta ,N-1}}{G_{\beta ,N}{\mathcal{N}}_{\beta }}(2N{)}^{\beta N/2+\beta }{e}^{-\beta N{x}^2}{R}_{N,\beta }(x),$$

$$G_{\beta ,N}$$

$${\mathcal{N}}_{\beta }$$

$$R_{N,\beta }(x)={\int }_{\mathcal{C}}d{u}_1{e}^{Nf(u_1,x)}\cdots {\int }_{\mathcal{C}}d{u}_{\beta }{e}^{Nf(u_{\beta },x)}\prod _{1\le j<k\le \beta }\text{absolutevalue}{u_k-u_j}^{4/\beta }$$

$$f(u,x)=-2u^2+\log (iu+x)-\frac{1}{N}\log (iu+x).$$

$$x\mapsto 1+\frac{x}{2N^{2/3}},$$

$$\begin{gathered} \frac{1}{2}\frac{G_{\beta ,N-1}}{G_{\beta ,N}{\mathcal{N}}_{\beta }}(2N{)}^{\beta N/2+\beta }={\beta }^{-\beta /2}{2}^{-5/2+2\beta +N\beta }{\pi }^{-1-\beta /2}{N}^{-1/2+\beta }{e}^{N\beta /2}\Gamma (1+\beta /2)\prod _{j=2}^{\beta }\frac{\Gamma (1+2/\beta )}{\Gamma (1+2j/\beta )} \\ \times (1+O(N^{-1})). \end{gathered}$$

Full text

1. Introduction

Universal results, whereby scaled large $N$ properties are independent of many local properties, abound in random matrix theory. Due to this, it is possible to talk of the statistical properties of large sized, bulk scaled complex Hermitian matrices, for example, without having to be more specific about the details of the ensemble. Thus the random Hermitian matrices may be formed entry-wise by sampling from a zero mean, finite standard deviation distribution – this is the Wigner class – or they may be defined by imposing a unitary invariant measure i.e. depending on the eigenvalues but not the eigenvectors; see e.g. [33].

Bulk scaling refers to choosing the origin in the support of the spectrum but away from the edges, then choosing for the length scale units of the mean spacing in the neighbourhood of this point. This is to be contrasted to the soft edge scaling, for which the origin is chosen to be in the neighborhood of the largest (or possibly smallest) eigenvalue, and the length scale is chosen to be of the order of the spacing between the first and the second eigenvalues.

Universality results apply in both the bulk and soft edge scaling regimes [12]. These two cases show themselves in some celebrated applications of random matrices. In the study of the statistical properties of the large Riemann zeros, the Montgomery-Odlyzko law asserts that upon scaling to have mean spacing unity, these coincide with the statistical properties of the bulk scaled eigenvalues of Hermitian matrices; see e.g. [25]. And for growth models on curved interfaces in the KPZ class, fluctuations have been predicted theoretically [32] and shown experimentally [34], to be identical to the fluctuations of the soft edge scaled largest eigenvalue of a random Hermitian matrix. Since many statistical quantities can be computed exactly for bulk and soft edge scaled Hermitian random matrices (see e.g. [14]), there are explicit functional forms available against which to compare data relating to the applications.

Data on the Montgomery-Odlyzko law comes by way of various data sets of Odlyzko – originally the high precision evaluation of the $10^{20}$-th Riemann zero and over $70\times 10^{6}$ neighbours [29], and more recently high precision evaluation of the ($10^{23}+985,531,550$) zero and over $10^{9}$ neighbours as announced in [30]. The latter sample size is sufficiently large that the difference between the empirical statistical quantity of interest (e.g. the nearest neighbour spacing) as determined from the data, and the random matrix prediction implied by the Montgomery-Odlyzko law, exhibits a distinct functional form [30]. In keeping with the proposal of Keating and Snaith [26] that the eigenvalues of random unitary matrices chosen with Haar measure provides the correct statistical description of finite size effects, one can use knowledge of the leading correction to the $N\to\infty$ bulk scaled distribution in this random matrix ensemble to predict the forms seen in the data [2, 17, 4].

Soft edge finite size effects can also be analysed [21, 5, 15, 11]. While different ensembles exhibit different functional forms for the leading correction, there is growing evidence that the optimal order – as controlled by the precise choice of the centring variable – has the leading soft edge correction for complex Hermitian matrices being proportional to $N^{-2/3}$ [20]. Although not the subject of the present work, there is also a literature along these lines with respect to hard edge scaling; see [10, 3, 31].

In the case of real symmetric matrices, the work of Johnstone [23] and collaborators [24, 28] has shown that the optimal choice of the soft edge centring and scaling variables for the classical ensembles again leads to a correction proportional to $N^{-2/3}$. Here the term classical ensembles refers to the particular class of orthogonally invariant real symmetric random matrices with eigenvalue probability density (PDF) proportional to

[TABLE]

chosen with weight function one of

[TABLE]

(the notation $\chi_{A}$ denotes the function taking the value $1$ for $A$ true and [math] otherwise).

The purpose of the present paper is to study the optimal choice of the soft edge scaling for the $\beta$-generalisation of the Gaussian and Laguerre orthogonal ensembles as specified by (1.1) and (1.2). Now the eigenvalue PDF is proportional to

[TABLE]

with weight function

[TABLE]

The matrix ensemble implied by (1.3) in the Gaussian case was first introduced by Dyson [9] as an interpolation of the classical Gaussian orthogonal ($\beta=1$), unitary ($\beta=2$) and sympletic ($\beta=4$) ensembles. Later it was shown that in both the Gaussian and Laguerre cases, (1.3) for general $\beta>0$ can be realised as the eigenvalue PDF for particular real tridiagonal matrices [8].

The quantity to be considered is the eigenvalue density

[TABLE]

where $\#=G,(L,a)$ depending on the weight (1.4) and

[TABLE]

is the normalisation. We note that (1.5) can be written as

[TABLE]

where $\langle\,\cdot\,\rangle_{N-1}^{\#}$ denotes the average with respect to (1.3) in the case of $N-1$ eigenvalues $\{x_{j}\}_{j=2}^{N}$. For $\beta$ even, the average is a polynomial in $\beta$, and it is this feature which gives opportunity to analyse the soft edge scaling limit.

Earlier work of Desrosiers and Forrester [7] has, for $\beta$ even, computed the soft edge scaling limit of $\rho_{N,\beta}^{\#}(x)$ for $\#=G$ and $(L,a)$. This involves the multidimensional integral

[TABLE]

which extends the integral representation of the Airy function,

[TABLE]

It was shown in [7] that

[TABLE]

Here we extend these results to bound the leading correction, by identifying the optimal choice of the centring variable.

Theorem 1**.**

We have that

[TABLE]

and

[TABLE]

are, for large $N$, equal to the RHS of (1.10) with a correction term $O(N^{-2/3})$.

An analogous result can also be established for a soft edge scaling in the Laguerre case with $a=\alpha N$.

Theorem 2**.**

Consider the Laguerre case of (1.3) with $a=\alpha N$. Let

[TABLE]

We have that

[TABLE]

is, for large $N$, equal to the RHS of (1.10) with a correction term $O(N^{-2/3})$.

Moreover, our working below gives the explicit form of the $O(N^{-2/3})$ leading order correction term in each case in terms of generalisations of (1.8). In the Gaussian case this is given by (3.20) while in the Laguerre case the relevant expressions are (4.11) with $a=\alpha N$ and (4.21) with a fixed.

In Section 2 we revise known results for the leading correction to not only the density, but also the general $k$-point correlation function for the Gaussian and Laguerre unitary ensembles ($\beta=2$ case of (1.3) and (1.4)). Specialised to the density, these results provide a check on our subsequent working. We proceed in Section 3 to establish the Gaussian case of Theorem 1. Section 4 establishes the Laguerre case of Theorem 1, as well as Theorem 2. A known differential-difference recurrence for computing the polynomial in (1.7) allows for a numerical and graphical demonstration, given in Section 5, that the identified scaling variables exhibit convergence at the claimed rate.

2. The case $\beta=2$

A distinguishing feature of the case $\beta=2$ is that the general $k$-point correlation is determinantal,

[TABLE]

for some correlation kernel $K_{N}$ independent of $k$. This provides a strategy to establish the leading correction term to the limiting form of $\rho_{(k)}$ under optimal soft edge scaling of this correlation function: thus one should analyse the corresponding asymptotic properties of $K_{N}$. This was done in the Gaussian and Laguerre cases (the latter with $a$ fixed) by Choup [5] and later extended to the Laguerre unitary ensemble with $a=\alpha N$ in [20], giving an $O(N^{-2/3})$ correction in each case.

To present these results, define

[TABLE]

where in the final formula $b$ is given by (1.11). Also set

[TABLE]

The results of [5, 20] tell us that for large $N$

[TABLE]

for $\#=G,(L,a),(L,\alpha N)$ and furthermore tell us the precise functional form of the term $O(N^{-2/3})$. To present the latter, we specialise to the case $x=y$, when the LHS of (2.4) corresponds to the soft edge scaled density. Using the notation (2.3) we have

[TABLE]

3. The Gaussian $\beta$ ensemble

Key to analysing (1.7) in the Gaussian case for $\beta$ even is the duality formula [1]

[TABLE]

Here we have used the notation $\mathrm{ME}_{\beta,N}[w(x)]$ to refer to the PDF corresponding to (1.3). Note the exchange of roles $(N,n,\beta)\to(n,N,4/\beta)$ between the two sides of (3.1). We see that the average in (1.7) is of the form of the LHS of (3.1) for $\beta$ even, telling us that then

[TABLE]

The normalisation $Z_{\beta,N}^{G}$ and that implied by the average in (3.2) are known explicitly, see e.g. [1]. Introducing the notation of this latter reference allows (3.2) to be written, after minor manipulation,

[TABLE]

where

[TABLE]

and

[TABLE]

with

[TABLE]

The contour $\mathcal{C}$ in (3.4) is equal to $\mathbb{R}$ as read off from (3.2), however by first ordering $u_{1}<u_{2}<\cdots<u_{\beta}$ the product of differences in the integrand can be interpreted as an analytic function and deformed into the complex plane by Cauchy’s theorem (see [7]).

Our interest is the large $N$ asymptotics of (3.3) after the replacement

[TABLE]

where the new scaled variable $x$ is $O(1)$. The choice (3.6) is consistent with the first case of (2.2). We begin by noting that upon simplification, then use of Stirling’s formula, one can check

[TABLE]

Also we have

[TABLE]

In relation to the exponent in (3.4), with $g(u)=-2u^{2}+\log(iu+1)$ we see from (3.5) that

[TABLE]

With (3.9) substituted in (3.4), the next step is to perform a saddle point analysis of (3.4). This has been done to leading order in [7] – here we extend the working therein to the next order but appealing to the latter reference for the formal justification of the methodology.

According to (3.9) the leading order term in the exponent of (3.4) is $g(u)$, which we can check has a double saddle point at $u=u_{0}:=i/2$. Expanding about this point by introducing the variable $t=u-u_{0}$, then changing variables $v=2iN^{1/3}t$ – this is in keeping with the direction of steepest descent; see [7] – gives

[TABLE]

Consequently

[TABLE]

where

[TABLE]

As an extension of (1.8), define

[TABLE]

Substituting (3.10) in (3.4), then substituting the result together with the expansion (3.7) in (3.3) with $x$ replaced by (3.6), we deduce that

[TABLE]

Our task now is to obtain a simplification of $\mathbb{K}_{\beta,\beta}[c_{1}^{G}(v;x,\beta)](x)$.

Proposition 1**.**

We have

[TABLE]

Proof.

In (3.12), let is regard $f$ as an operator acting on all terms in the integrand to the right.

It follows from the fundamental theorem of calculus that

[TABLE]

On the other hand, by direct differentiation

[TABLE]

where the second equality follows from the fact that the integrand implied by the final term in the line above is anti-symmetric upon the interchange of any two co-ordinates and thus vanishes. Hence

[TABLE]

It similarly follows from the fundamental theorem of calculus that

[TABLE]

while direct differentiation gives

[TABLE]

Taking the arithmetic mean of the final term with itself after interchanging $v_{j}\leftrightarrow v_{k}$ shows

[TABLE]

Hence

[TABLE]

Recalling the definition of $c_{1}^{G}$ in (3.11), it follows from (3.15) and (3.16) that

[TABLE]

The result (3.14) now follows from the fact that

[TABLE]

∎

The results (3.13), (1.10) and (3.14) tell us that for $\beta$ even

[TABLE]

where

[TABLE]

The first statement in Theorem 1 now follows from (3.18) by recentring $x$ according to $x\mapsto x+\left({1\over 2}-{1\over\beta}\right)$ then expanding the first term according to Taylor’s theorem up to second order, and the second term up to first order; this process also gives the explicit $O(N^{-2/3})$ correction term as

[TABLE]

Remark 1*.*

(A.) In the case $\beta=2$ we have from (3.20), (3.19) and (3.12) that the correction term (3.20) is equal to

[TABLE]

where $c_{2}^{G}$ is given by (3.11). Expanding $|v_{2}-v_{1}|^{2}$ and making use of (1.9) gives agreement with $\rho_{\infty,\beta=2,1}(x)$ as specified in (2.3).

(B.) The factor $({1\over\beta}-{1\over 2})$ in (3.14) is familiar as a factor in the $O(1/N)$ corrections to the bulk smoothed densities in the Gaussian and Laguerre $\beta$ ensembles; see [22] and the recent works [35], [19].

4. The Laguerre $\beta$ ensemble

Introduce the normalisations

[TABLE]

Use the notation CE${}_{\beta,N}^{(c,d)}$ to denote the ensemble of eigen-angles specified by the PDF proportional to

[TABLE]

Combining [14, Prop. 13.2.5 & Exercises 13.1 (4.4)] gives the duality formula

[TABLE]

According to [7] the multi-dimensional integral implied by the RHS of (4.1) can be written in the form (3.4), but with

[TABLE]

and $\mathcal{C}$ a counterclockwise closed contour about $u=1$. Furthermore, as in the Gaussian case, by suitably ordering the contours the absolute value signs about $|u_{k}-u_{j}|$ can be removed and the contours deformed into the complex plane. Given this, after recalling (1.7), and using the notation $R_{N,\beta}$ as is specified in terms of $f(u,x)$ and $\mathcal{C}$ by (3.4), we have in keeping with [7, Eq. (22)] that

[TABLE]

The case $a=\alpha N$

In keeping with the third case of (2.2), we want to expand for large $N$ to the first two orders (4.3) with $a=\alpha N$ and

[TABLE]

We begin by noting

[TABLE]

and

[TABLE]

In (4.6) use was made use of Stirling’s formula and the multiplication formula for the gamma function

[TABLE]

Next, we turn to a saddle point analysis of the integral in (4.3). After making the substitution (4.4) in (4.2) we see that up to terms of order $1/N$, the latter quantity is equal to

[TABLE]

This possesses a double saddle point at $u_{0}=1/(\sqrt{1+\alpha}+1)$. We expand about this point in the direction of steepest descent by the changing of variables $v=-2(bN)^{1/3}(u-u_{0})$. This gives

[TABLE]

Combining the results (4.5), (4.6), (4.7) with (4.3) shows

[TABLE]

where

[TABLE]

Analogous to the result of Proposition 1, the $O(N^{-1/3})$ term in (4.8) can be identified as being proportional to the derivative of the leading term.

Proposition 2**.**

We have

[TABLE]

Proof.

This follows directly from (4.9), (3.15), (3.16) and (3.17). ∎

Combining (4.8) and (4.10) establishes Theorem 2. In addition, with the notation (3.19), we see that the explicit form of the $O(N^{-2/3})$ correction term is

[TABLE]

In the case $\beta=2$, proceeding as in the Gaussian case of Remark 1, we can check that this is consistent with $\rho_{\infty,\beta=2,1}^{(L,\alpha N)}$ as specified in (2.3).

The case $a$ fixed

We fix $a$ and make the replacement

[TABLE]

A simple calculation shows

[TABLE]

Further, by use of Stirling’s formula and the multiplication formula for the gamma function we deduce

[TABLE]

With $a$ fixed we read of that the leading large $N$ term of (4.2) is

[TABLE]

This has a double saddle point at $u_{0}=1/2$. Expanding around this point in the direction of the steepest descent by the change of variables $v=-2(2N)^{1/3}(u-u_{0})$ then shows

[TABLE]

By combining (4.13), (4.14), (4.15) with (4.3), for the Laguerre $\beta$-ensemble with $a$ fixed, the large $N$ expansion of the soft edge scaled density as defined in the variable (4.12) is given by

[TABLE]

where

[TABLE]

Immediate from (4.17) and (3.17) is that

[TABLE]

Applying the translation $x\mapsto x+a/(2N)^{1/3}$ to (4.16) then shows that the scaled soft edge density in the variable listed second in (2.2), up to order $O(N^{-2/3})$, is given by

[TABLE]

where

[TABLE]

This corroborates with Theorem 1 in the Laguerre case, and furthermore gives the leading correction as

[TABLE]

where $r_{\infty,\beta}^{(L,a)}$ is specified according to (3.19) with the substitution (4.20). Note that this correction is independent of $a$. In the case $\beta=2$ this latter property is already evident from the appropriate case of (2.3). The procedure of Remark 1 can be used to reduce (4.21) with $\beta=2$ to this functional form.

5. Numerics

In this section we provide a method to determine the graphical form of the optimal $O(N^{-2/3})$ correction term for the soft edge density of the Gaussian and Laguerre even $\beta$-ensembles. Multiple integral forms have been given in (3.20), (4.11) and (4.21). However these are not well suited to accurate numerical evaluation. Instead we turn to known [13, 16, 18] recursive properties relating to the average in (1.7), recalling that for $\beta$ even the latter is a polynomial of degree $\beta(N-1)$. As such we are making use of the broader theory of the Selberg integral; see [14, Ch. 4].

For fixed parameters $(\lambda_{1},\lambda_{2},\lambda,\alpha)$, introduce the auxiliary function $I_{p}[w(t)](x)=I_{p}[w(t)](x;\alpha,\lambda)$ as the multiple integral

[TABLE]

where

[TABLE]

denotes the elementary symmetric polynomials. In the case $w(t)=t^{\lambda_{1}}(1-t)^{\lambda_{2}}\chi_{0<t<1}$ we know from [13] that these integrals satisfy the differential-difference equation

[TABLE]

where

[TABLE]

The immediate relevance of (5.1) to the average (1.7) comes from the fact that

[TABLE]

where the superscript $(J)$ indicates use of the Jacobi weight in (1.7). In the case $\alpha=1$, $p=0$, (5.1) is independent of $x$, which provides an initial condition to iterate (5.2) up to $p=N$. The crucial property

[TABLE]

then gives the value of the polynomial in the case $\alpha=2$, $p=0$, and the iteration can be continued.

Introduce the notation

[TABLE]

for the Gaussian and Laguerre cases of (5.1). A simple change of variables shows

[TABLE]

Applying this limiting procedure to (5.2) gives

[TABLE]

and

[TABLE]

where $\lambda=\beta/2$ and $\lambda_{1}=a\beta/2$. The constant (i.e. $x$ independent) values for $\alpha=1$, $p=0$ can be read off from the normalisations appearing in (3.3) and (4.1), which being limiting cases of the Selberg integral can be expressed in terms of products of gamma functions. Thus G_{0}\big{\rvert}_{\alpha=1}=G_{2\lambda,N} and L_{0}\big{\rvert}_{\alpha=1}=W_{a,2\lambda,N}.

Figures 1, 2, 3 and 4 were generated by substituting the appropriate soft edge variables

[TABLE]

into (5.5) and (5.6). With $s_{x,\beta}$, $s_{x}^{\prime}$ defined by (5.7), (2.2) (with the exception for the Laguerre case with $a$ fixed, where we take instead $s_{x}^{\prime}=4N+2(2N)^{1/3}x$), Figures 1, 2 and 3 plot

[TABLE]

For some large $N_{1},N_{2}\in\mathbb{N}$ values, the successive differences given by (5.8) admit, as seen graphically, the same functional form and are $O(1)$. (5.8) plots the corrections (3.20), (4.21), (4.11) up to some error, to leading order, proportional to $\left({1\over N_{2}^{2/3}}-{1\over N_{1}^{2/3}}\right)^{-1}\left({1\over N_{2}}-{1\over N_{1}}\right).$ Figure 4 plots

[TABLE]

where

[TABLE]

and shows the approximate functional form of the correctional term of order $N^{-1/3}$ due to the scaling (2.2) which bears resemblance to the derivative of (1.7). This further emphasises the results found in Proposition 1, 2 and (4.18). The small displacement relative to the derivative in the displayed examples – which we take as a finite N artifact –- is not a feature of the previous figures (right captions) due to the fundamentally different nature of what is being plotted.

6. More general invariant ensembles

In our previous study [20], the question of the optimal soft edge scaling for complex Hermitian Wigner matrices was briefly addressed. We recall that an Hermitian random matrix is termed a Wigner ensemble if all its diagonal entries (which must be real) are independently chosen from the same zero mean, finite variance distribution, and all its upper triangular entries (which may be complex) chosen similarly. The Gaussian unitary ensemble, corresponding to the Gaussian case of (1.1) with $\beta=2$, is both a unitary invariant ensemble, and a complex Wigner matrix: matrices from this ensemble can be generated by choosing the diagonal entries from the normal distribution N$[0,1/\sqrt{2}]$ and the upper triangular entries from N$[0,1/2]+i{\rm N}[0,1/2]$. In [20] numerical simulations were performed on the particular complex Wigner matrices $Y={1\over 2}(X+X^{\dagger})$, where all entries of $X$ are chosen independently and uniformly from the set of four values ${1\over\sqrt{2}}(\pm 1\pm i)$. The off diagonal entries of $Y$ then have mean zero and variance one half, as for the GUE, guaranteeing (see e.g. [33]) that to leading order the largest eigenvalue occurs at $\sqrt{2N}$. By the use of simulation, evidence was presented that a scaled and centred variable can be identified such that the leading correction to the limiting distribution of the soft edge scaled largest eigenvalue is $O(N^{-2/3})$. Of course the natural question is to ask if this effect holds more generally.

What then for the case of general invariant ensembles beyond the classical cases? In the case $\beta=2$, and for the family of weights $w_{2}(x)=e^{-Nx^{2m}}$, as noted in [27, Eq. (1.9)], results of [6] give the asymptotic formula for the soft edge scaled correlation kernel (2.1)

[TABLE]

where $b,\gamma$ depend on $m$. Thus once again there is a choice of scaling variables such that the correction to the correlations at the soft edge is $O(N^{-2/3})$.

Acknowledgements

This work is part of a research program supported by the Australian Research Council (ARC) through the ARC Centre of Excellence for Mathematical and Statistical frontiers (ACEMS). PJF also acknowledges partial support from ARC grant DP170102028, and AKT acknowledges the support of a Melbourne postgraduate award.