Extreme eigenvalues of random matrices from Jacobi ensembles

TL;DR

This paper derives two-term asymptotic formulas for the distribution of the smallest and largest eigenvalues in Jacobi beta-ensembles, revealing explicit expressions and correction terms for large matrices.

Contribution

It provides new two-term asymptotic formulas for eigenvalue distributions in Jacobi beta-ensembles with explicit correction terms and special case formulas involving familiar functions.

Findings

Explicit two-term asymptotic formulas derived

First-order corrections proportional to distribution derivatives

Special cases with explicit formulas involving Bessel functions

Abstract

Two-term asymptotic formulae for the probability distribution functions for the smallest eigenvalue of the Jacobi $ \beta $-Ensembles are derived for matrices of large size in the r\'egime where $ \beta > 0 $ is arbitrary and one of the model parameters $ \alpha_1 $ is an integer. By a straightforward transformation this leads to corresponding results for the distribution of the largest eigenvalue. The explicit expressions are given in terms of multi-variable hypergeometric functions, and it is found that the first-order corrections are proportional to the derivative of the leading order limiting distribution function. In some special cases $ \beta = 2 $ and/or small values of $ \alpha_1 $, explicit formulae involving more familiar functions, such as the modified Bessel function of the first kind, are presented.

Equations

$$\frac{1}{S_N({\alpha }_1+1,{\alpha }_2+1,\beta /2)}\prod _{i=1}^N{x}_i^{{\alpha }_1}(1-x_i{)}^{{\alpha }_2}|\Delta (\mathbf{x}){|}^{\beta },0⩽x_i⩽1,$$

$$S_N(a,b,c)\text{raise}0.34444pt:=\prod _{i=0}^{N-1}\frac{\Gamma (1+(i+1)c)\Gamma (a+ic)\Gamma (b+ic)}{\Gamma (1+c)\Gamma (a+b+(N+i-1)c)},$$

$$\Delta (\mathbf{x})\text{raise}0.34444pt:=\prod _{1⩽i<j⩽N}(x_j-x_i).$$

$${\alpha }_1>-1,{\alpha }_2>-1,\beta >-\frac{1}{2}\min \left\{\frac{1}{N},\frac{{\alpha }_1}{N-1},\frac{{\alpha }_2}{N-1}\right\}.$$

$$\frac{N}{\pi }\frac{1}{\sqrt{x(1-x)}},0<x<1.$$

$$F_{N^2{\varphi }_1}(x)=\mathbb{P}(N^2{\varphi }_1⩽x)=\mathbb{P}\left({\varphi }_1⩽\frac{x}{N^2}\right)=F_{{\varphi }_1}\left(\frac{x}{N^2}\right).$$

$$F_{N^2{\varphi }_1}(x)=1-{\mathrm{e}}^{-x}\det (I_{j-i}(2\sqrt{x}))+\frac{{\alpha }_1+{\alpha }_2}{N}x{\mathrm{e}}^{-x}\det (I_{2+j-i}(2\sqrt{x}))+\mathrm{O}\left(\frac{1}{N^2}\right)$$

$$I_{\nu }(z)\text{raise}0.34444pt:=\frac{z^{\nu }}{2^{\nu }\sqrt{\pi }\Gamma (\nu +\frac{1}{2})}{\int }_{-1}^1{\mathrm{e}}^{-zt}(1-t^2{)}^{\nu -1/2}\mathrm{d}t,\mathrm{Re}\{\nu \}>-\frac{1}{2},$$

$$\underset{N\to \infty }{\lim }{F}_{N^2{\varphi }_1}(x)=1-{\mathrm{e}}^{-\beta x/2}\text{fourIdx}0(\beta /2)1F\left(;\frac{2{\alpha }_1}{\beta };x{1}^{{\alpha }_1}\right),$$

$$\begin{gathered} F_{N^2{\varphi }_1}(x)=1-{\mathrm{e}}^{-\beta x/2}\text{fourIdx}0(\beta /2)1F\left(;\frac{2{\alpha }_1}{\beta };x{1}^{{\alpha }_1}\right) \\ +\frac{x^{1+{\alpha }_1}}{N}\left(({\alpha }_1+{\alpha }_2+1)-\frac{\beta }{2}\right){\left(\frac{\beta }{2}\right)}^{2{\alpha }_1}\frac{\Gamma (1+\beta /2)}{\Gamma (1+{\alpha }_1)\Gamma (1+{\alpha }_1+\beta /2)} \\ \times {\mathrm{e}}^{-\beta x/2}\text{fourIdx}0(\beta /2)1F\left(;\frac{2{\alpha }_1}{\beta }+2;x{1}^{{\alpha }_1}\right)+\mathrm{O}\left(\frac{1}{N^2}\right). \end{gathered}$$

$$\begin{gathered} \mathbb{P}({\varphi }_N⩽1-x/N^2)={\mathrm{e}}^{-\beta x/2}\text{fourIdx}0(\beta /2)1F\left(;\frac{2{\alpha }_2}{\beta };x{1}^{{\alpha }_2}\right) \\ -\frac{x^{1+{\alpha }_2}}{N}\left(({\alpha }_1+{\alpha }_2+1)-\frac{\beta }{2}\right){\left(\frac{\beta }{2}\right)}^{2{\alpha }_2}\frac{\Gamma (1+\beta /2)}{\Gamma (1+{\alpha }_2)\Gamma (1+{\alpha }_2+\beta /2)} \\ \times {\mathrm{e}}^{-\beta x/2}\text{fourIdx}0(\beta /2)1F\left(;\frac{2{\alpha }_2}{\beta }+2;x{1}^{{\alpha }_2}\right)+\mathrm{O}\left(\frac{1}{N^2}\right). \end{gathered}$$

$${\mathrm{D}}_k\text{raise}0.34444pt:=\sum _{i=1}^n{x}_i^k\frac{{\partial }^2}{\partial x_i^2}+\frac{2}{\sigma }\sum _{i\ne j}\frac{x_i^k}{x_i-x_j}\frac{\partial }{\partial x_i}$$

$${\mathrm{E}}_k\text{raise}0.34444pt:=\sum _{i=1}^n{x}^k\frac{\partial }{\partial x_i},$$

$${\mathrm{E}}_1{C}_{\lambda }^{(\sigma )}(\mathbf{x})=|\lambda |C_{\lambda }^{(\sigma )}(\mathbf{x})$$

$${\mathrm{D}}_2{C}_{\lambda }^{(\sigma )}(\mathbf{x})=\left({\rho }_{\lambda }+\frac{2}{\sigma }|\lambda |(n-1))\right)C_{\lambda }^{(\sigma )}(\mathbf{x}),$$

$${\rho }_{\lambda }\text{raise}0.34444pt:=\sum _{i=1}^n{\lambda }_i\left({\lambda }_i-1-\frac{2}{\sigma }(i-1)\right).$$

$$C_{\lambda }^{(\sigma )}(\mathbf{x})=\sum _{\mu }{b}_{\mu \lambda }{m}_{\mu }(\mathbf{x})$$

$$\sum _{|\lambda |=k}{C}_{\lambda }^{(\sigma )}(\mathbf{x})=(x_1+\cdots +x_n{)}^k,k\in \mathbb{N}.$$

$${\mathrm{D}}_2-{\mathrm{D}}_1+(a+b+2){\mathrm{E}}_1-(a+1){\mathrm{E}}_0$$

$$J_{\lambda }^{\sigma ,a,b}(\mathbf{x})=\sum _{\mu \subseteq \lambda }{c}_{\mu ,\lambda }{C}_{\mu }^{(\sigma )}(\mathbf{x}),$$

$$\text{fourIdx}p(\sigma )qF(a_1,\dots ,a_p;b_1,\dots ,b_q;\mathbf{x})\text{raise}0.34444pt:=\sum _{\lambda }\frac{[a_1{]}_{\lambda }^{(\sigma )}\cdots [a_p{]}_{\lambda }^{(\sigma )}}{[b_1{]}_{\lambda }^{(\sigma )}\cdots [b_q{]}_{\lambda }^{(\sigma )}|\lambda |!}{C}_{\lambda }^{(\sigma )}(\mathbf{x}),$$

$$[a{]}_{\lambda }^{(\sigma )}\text{raise}0.34444pt:=\prod _{i=1}^n{\left(a-\frac{i-1}{\sigma }\right)}_{{\lambda }_i},$$

$$(a{)}_n\text{raise}0.34444pt:=\frac{\Gamma (a+n)}{\Gamma (a)}=\prod _{j=0}^{n-1}(a+j),a\in \mathbb{C},n\in \mathbb{N}.$$

$$\text{fourIdx}2(\sigma )1F(a,b;c;\mathbf{x})=\prod _{i=1}^n(1-x_i{)}^{-a}\text{fourIdx}2(\sigma )1F\left(a,c-b;c;\frac{-x_1}{1-x_1},\dots ,\frac{-x_n}{1-x_n}\right).$$

$$\begin{gathered} {\int }_0^1\cdots {\int }_0^1\prod _{i=1}^n{x}_i^{a-1}(1-x_i{)}^{b-1}\prod _{j=1}^m(x_i-t_j)|\Delta (\mathbf{x}){|}^{2/\sigma }{\mathrm{d}}^n\mathbf{x} \\ =S_n(a+m,b,1/\sigma )\text{fourIdx}2(1/\sigma )1F\left(-n,\sigma (a+b+m-1)+n-1;\sigma (a+m-1);\mathbf{t}\right), \end{gathered}$$

$$F_{{\varphi }_1}(\xi )\text{raise}0.34444pt:=\mathbb{P}({\varphi }_1⩽\xi )=1-\mathbb{P}({\varphi }_1>\xi ).$$

$$F_{\phi_{1}}(\xi)=\left\{\begin{array}[]{cc}0,&\xi\leqslant 0,\\ 1,&\xi\geqslant 1,\end{array}\right.$$

$$\mathbb{P}({\varphi }_1>\xi )=(1-\xi {)}^{N(1+{\alpha }_2+(N-1)\beta /2)}\text{fourIdx}2(\beta /2)1F\left(-N,1-N-\frac{2}{\beta }({\alpha }_2+1);\frac{2{\alpha }_1}{\beta };\xi 1^{{\alpha }_1}\right).$$

$$\mathbb{P}({\varphi }_1>\xi )$$

$$=\frac{1}{S_N({\alpha }_1+1,{\alpha }_2+1,\beta /2)}{\int }_{\xi }^1\cdots {\int }_{\xi }^1\prod _{i=1}^N{x}_i^{{\alpha }_1}(1-x_i{)}^{{\alpha }_2}|\Delta (\mathbf{x}){|}^{\beta }{\mathrm{d}}^N\mathbf{x}.$$

$$\begin{gathered} \mathbb{P}(\varphi >\xi )=\frac{(1-\xi {)}^{N(1+{\alpha }_1+{\alpha }_2+(N-1)\beta /2)}}{S_N({\alpha }_1+1,{\alpha }_2+1,\beta /2)} \\ \times {\int }_0^1\cdots {\int }_0^1\prod _{i=1}^N{\left(y_i+\frac{\xi }{1-\xi }\right)}^{{\alpha }_1}(1-y_i{)}^{{\alpha }_2}|\Delta (\mathbf{y}){|}^{\beta }{\mathrm{d}}^N\mathbf{y}. \end{gathered}$$

Full text

In some special cases $\beta=2$ and/or small values of $\alpha_{1}$, explicit formulæ involving more familiar functions, such as the modified Bessel function of the first kind, are presented.

1 Introduction

A random matrix is a matrix whose entries are random variables. As eigenvalues of a matrix are continuous functions of its entries, so the eigenvalues of a random matrix are random variables. A random $N\times N$ matrix has the Jacobi $\beta$-Ensemble (JE) distribution if a joint probability density function of its eigenvalues is

[TABLE]

where

[TABLE]

and the Vandermonde determinant is defined by

[TABLE]

That (1.1) is a properly normalised probability density is a consequence of Selberg’s integral [1].

In many situations $\beta$ is a non-negative integer, and we will mostly be assuming that one of $\alpha_{1}$ or $\alpha_{2}$ is a non-negative integer, but (1.1) makes sense for arbitrary real values of these parameters subject to the constraints

[TABLE]

The naming of this ensemble reflects the presence in (1.1) of the factors $x_{i}^{\alpha_{1}}(1-x_{i})^{\alpha_{2}}$ which are a density with respect to which (a certain version of) the classical Jacobi polynomials form an orthogonal collection.

We label by $\phi_{i}$ the sorted eigenvalues, so that $0\leqslant\phi_{1}\leqslant\phi_{2}\leqslant\cdots\leqslant\phi_{N}\leqslant 1$. This article is concerned with the distribution of the extreme eigenvalues $\phi_{1}$ and $\phi_{N}$. In fact, since the change of variables $x_{i}\mapsto 1-x_{i}$, for $i=1,\ldots,N$, in (1.1) leaves the joint probability density invariant, save for the exchange $\alpha_{1}\leftrightarrow\alpha_{2}$, and reverses the order of the eigenvalues, it will not present a loss of generality to focus on the smallest eigenvalue $\phi_{1}$.

The limiting empirical eigenvalue density for Jacobi random matrices was derived in [2]. For fixed $\alpha_{1},\alpha_{2}$, the large $N$ limiting density is

[TABLE]

This means that for large $N$ the number of eigenvalues in the interval $[0,N^{-2}]$ is approximately $N/\pi\int_{0}^{N^{-2}}(x(1-x))^{-1/2}\,{\mathrm{d}}x={\mathrm{O}}(1)$, and it is natural to expect $N^{2}\phi_{1}$ to have a non-trivial limiting distribution.

Our main objects of interest will be the (cumulative) probabilty distribution function $F_{\phi_{1}}(\xi)\mathbin{\hbox{\raise 0.34444pt\hbox{\rm:}}\!\!=}{\mathbb{P}}(\phi_{1}\leqslant\xi)$, and the rescaled version

[TABLE]

Deferring to below a more comprehensive summary of previous work on this problem, we mention a result [3] of Moreno-Pozas, Morales-Jimenez, McKay in the case $\beta=2$ (the Jacobi Unitary Ensemble, JUE). They proved, for $\alpha_{1}=0,1$ and $\alpha_{2}\in{\mathbb{N}}_{0}$; and for $\alpha_{1}=2,\alpha_{2}\in\{0,1,2\}$, the two-term asymptotic result

[TABLE]

where the determinants appearing in (1.7) are of size $\alpha_{1}\times\alpha_{1}$, and $I_{n}(z)$ is the $I$-Bessel function—the modified Bessel function of the first kind,

[TABLE]

and $I_{-n}(z)=(-1)^{n}I_{n}(z)$ for $n\in{\mathbb{Z}}$.

On the other hand, Borodin and Forrester have derived [4] the leading-order distribution of the smallest eigenvalue of the JE for any $\beta>0$ and $\alpha_{1}\in{\mathbb{N}}_{0}$:

[TABLE]

where $\fourIdx{}{0}{(\sigma)}{1}{F}(;c;{\mathbf{x}})$ is a multivariate hypergeometric function that will be defined precisely in Section 2.3, and ${\mathbf{1}}^{n}\mathbin{\hbox{\raise 0.34444pt\hbox{\rm:}}\!\!=}(1,1,\ldots,1)\in{\mathbb{R}}^{n}$. Our principal result is a version of the two-term asymptotic (1.7) valid for $\beta>0$.

Theorem 1.1**.**

Let $\phi_{1}$ be the smallest eigenvalue of the $N\times N$ Jacobi $\beta$-Ensemble, $\beta>0$, with $\alpha_{1}\in{\mathbb{N}}_{0}$ and $\alpha_{2}>-1$. For $x>0$,

[TABLE]

The error estimate can depend on $\alpha_{1},\alpha_{2},\beta$ but is uniform for $x$ in a compact set.

All our results for the distribution of the smallest eigenvalue can be re-cast to give an analogous result for the largest eigenvalue, as indicated earlier. We will not write down these analogues for every result, allowing just the following Corollary of Theorem 1.1.

Corollary 1.2**.**

Let $\phi_{N}$ be the largest eigenvalue of the $N\times N$ Jacobi $\beta$-Ensemble, $\beta>0$, with $\alpha_{2}\in{\mathbb{N}}_{0}$ and $\alpha_{1}>-1$. For $x>0$,

[TABLE]

There are several known random matrix models that lead to JE eigenvalue distributions. Most famous are perhaps the double-Wishart (or Manova) models from Statistics: set $M_{1}$, $M_{2}$ to be independent $n_{1}\times N$ and $n_{2}\times N$ matrices with independent standard normal real random variable entries, $n_{1},n_{2}\geqslant N$. If $A=M_{1}^{\dagger}M_{1}$, $B=M_{2}^{\dagger}M_{2}$, then the matrix $A(A+B)^{-1}$ has eigenvalues distributed according to the JE with $\beta=1$, $\alpha_{1}=(n_{1}-N-1)/2$ and $\alpha_{2}=(n_{2}-N-1)/2$ [5, 6, 7, 8, 9]. Since our results rely on $\alpha_{1}$ being integer, this requires $n_{1}-N$ to be an odd diference.

If we repeat the above construction, with complex normal random variables, then the eigenvalue distribution of $A(A+B)^{-1}$ is JE with $\beta=2$, $\alpha_{1}=n_{1}-N$ and $\alpha_{2}=n_{2}-N$ [10, Section 8].

Another model leading to the joint probability density function (1.1) is the corners process of random matrices from classical compact groups: if $U$ is a random $m\times m$ unitary or orthogonal matrix chosen with respect to Haar measure, $m\geqslant 2N$, and $M$ is the principal $N\times N$ submatrix of $U$ (the upper-left corner matrix), then, letting $s_{1},\ldots,s_{N}$ denote the $N$ eigenvalues of $M^{\dagger}M$, the points $x_{1}=s_{1}/\beta,\ldots,x_{N}=s_{N}/\beta$ are distributed according to (1.1) with $\alpha_{1}=\beta/2-1$ $\alpha_{2}=(m-2N+1)\beta/2-1$ and $\beta=2$ (unitary case) or $\beta=1$ (orthogonal case) [11, 12, §7.2]. In the latter case $\alpha_{1}=-1/2$ which is not an integer, so Theorem 1.1 does not apply, but Corollary 1.2 does apply for the distribution of the largest eigenvalue if $m$ is an odd number (whence $\alpha_{2}$ is an integer).

Random matrix models that allow full exploration of the parameter space, including to arbitrary $\beta>0$, are also known [13, 14, 15].

The JE exhibits two “hard edges” in the spectrum at $x=0$ and $x=1$ since the eigenvalues are strictly confined between these values, which furthermore coincides with the support of the limiting eigenvalue density (1.5). This is in contrast to some other random matrix models such as the Gaussian ensembles [16] which have compactly supported limiting eigenvalue density—the famous Wigner’s semi-circle law [17, 18]—but without any intrisic obstacle to having individual eigenvalues appearing at any point on the real line. Statistics such as the distribution of smallest eigenvalues are expected to be “universal” in the limit $N\to\infty$, in the sense that they ought not to depend on the precise features of the random matrix model in question. In our present context it means that the limiting distribution (1.9) will be valid for other matrix models with a hard spectral edge. Indeed, the same limiting distribution has been proven for a different set of matrix models exhibiting a hard edge—the Laguerre $\beta$-Ensembles (LE; sometimes called Wishart random matrices) [19], as well as modifications of the JUE that preserve the hard edge [20].

The finite $N$ corrections to the leading order derived in Theorem 1.1 are not expected to be universal—indeed the presence of the parameter $\alpha_{2}$ seems to rule that out—but they do exhibit an interesting feature that had already been conjectured for the Laguerre Unitary Ensemble at the hard edge [21] and proved for that model in [22, 23, 24]: the correction term is proportional to the derivative of the main term. This holds for our two-term asymptotic (1.10), although it may not seem immediately apparent: see (5.7) below. Forrester and Trinh [25] have investigated the eigenvalue density for the LE for $\beta>0$, and found two-term asymptotics at the hard edge of the spectrum, and that the correction term is also proportional to the derivative of the leading term. It seems likely that the methods in the present work could also be adapted to study the hard-edge of the LE too.

JE random matrices have a number of known applications. The outage probability of multiple-input/multiple-output (MIMO) systems subject to interference, such as those used in cellular mobile radio networks, can be modelled in terms of the largest eigenvalue of JUE matrices [26]. The conductance eigenvalues in random matrix models for mesoscopic disordered quantum systems are known to be governed by the JE distribution with $\alpha_{2}=0$ [27, 28]. In this context, expressions for the average spectral density have been derived in terms of multi-variable hypergeometric functions [29], somewhat similar to expressions for the smallest eigenvalue derived in Section 3. Finally, some tests in multivariate Statistics are based on the distributions of extreme eigenvalues of JE (generally the parameters $\beta=1$ and $\beta=2$ corresponding to the real and complex underlying fields are most relevant), see Roy [30]. Some of these statistical applications are reviewed in Section 2 of [31]. Roy’s test has practical applications in signal analysis in the presence of coloured noise, for which the distribution of the largest eigenvalue of the JUE is required [32].

Aside from the references [3, 4] mentioned above, theoretical work on the distribution of extreme eigenvalues for Jacobi ensembles goes back at least to [33] for $\beta=2$ and Constantine [34] for $\beta=1$, motivated by the aforementioned applications in Statistics.

In [35] expressions were derived for distribution functions in terms of multivariate hypergeometric functions in $N$ variables, and corresponding formulae for density functions in $N-1$ variables given in [36] and [37]. Algorithms for a numerical evaluation of the distribution of the smallest eigenvalue in the JUE were given in [38] with methods applicable to arbitrary values of the parameters $\alpha_{1},\alpha_{2}>-1$, and furthermore which extend even to non-integer values of $N$.

Johnstone [31] and Jiang [39] have investigated statistics of extreme eigenvalues, and other quantities, in a setting where the parameter values $\alpha_{1}$ and $\alpha_{2}$ are not fixed, but vary as $N\to\infty$, leading to a soft edge in the spectrum. Scaling limits at the hard and soft-edge were treated together in [40].

Forrester and Li [41] have studied eigenvalue correlations for a broader class of unitary ensembles with a hard edge at the spectrum (which includes the JUE) and found $1/N$-correction terms consistent with [3].

In Section 2 we introduce some of the analytic tools that will be used (multi-variable hypergeometric functions and Jacobi polynomials). In Section 3 we collect some exact formulæ for finite $N$. In Section 4 we prove a two-term asymptotic formula for $\fourIdx{}2{(\sigma)}1{F}$ multivariate hypergeometric functions, that is then used to give the proof of Theorem 1.1 in Section 5. A few special cases are treated in Section 6.

2 Multi-variable hypergeometric functions and

Jacobi polynomials

Multi-variable analogues of classical hypergeometric functions and orthogonal polynomials are a relatively recently-developed area of study that have nevertheless proved very useful in Random Matrix Theory, see, e.g. [42, 43, 44, 45, 46, 47] as well as many other articles cited in the present work. They can be defined as series of Jack polynomials which we define first.

2.1 Jack polynomials

Let ${\mathbf{x}}=(x_{1},\ldots,x_{n})$ be a set of variables, $\lambda=(\lambda_{1},\ldots,\lambda_{n})$ be an integer partition111We can assume that the number of parts of $\lambda$ is equal to the number of variables, since if there are more parts than variables the corresponding Jack polynomial is zero; on the other hand, any partition can be padded with [math]s to increase the number of parts to $n$. of size $|\lambda|=\lambda_{1}+\cdots+\lambda_{n}$, and let $\sigma>0$. The Jack polynomials [48] are certain homogeneous, symmetric polynomials $C_{\lambda}^{(\sigma)}({\mathbf{x}})$ of degree $|\lambda|$.

We define the operators

[TABLE]

and

[TABLE]

for $k\in{\mathbb{N}}_{0}$

Jack polynomials are joint eigenfunctions of $\mathrm{E}_{1}$ and $\mathrm{D}_{2}$ [49]. In fact,

[TABLE]

(a relation satisfied by any homogeneous polynomial of degree $|\lambda|$) and

[TABLE]

where

[TABLE]

The definition of $C_{\lambda}^{(\sigma)}({\mathbf{x}})$ is completed by triangularisation: if

[TABLE]

is the expansion of $C_{\lambda}^{(\sigma)}$ in the basis of monomial symmetric functions, then the coefficient $b_{\mu\lambda}=0$ unless $\mu\leqslant\lambda$ in terms of dominance ordering of partitions [49]; and normalisation:

[TABLE]

(That such a normalisation exists is proved in [49, Prop. 2.3], although a different normalisation for the Jack polynomials is actually used throughout [49]. The normalisation leading to (2.7) is commonly-used for applications in Random Matrix Theory.)

2.2 Multi-variable Jacobi polynomials

As with multi-variable hypergeometric functions defined in the next subsection, multi-variable generalisations of the classical Jacobi polynomials were initially studied for Jack parameter $\sigma=2$ [50] with applications in Statistics in mind. Later these were generalised to other values of $\sigma$, with a variety of conventions for normalisation and support of the orthogonality measure [51, 52, 53, 54]. They sometimes go by the name “Jacobi polynomials associated with the root system $BC_{n}$” [55]. In our definitions, we follow [56] with a difference in the choice of normalisation.

For $a,b\in{\mathbb{R}}$ fixed, $J_{\lambda}^{\sigma,a,b}({\mathbf{x}})$ is a symmetric polynomial eigenfunction of the operator

[TABLE]

of the form

[TABLE]

for constants $c_{\mu,\lambda}$ depending on $a,b$ and $\sigma$, and the notation $\mu\subseteq\lambda$ means $\mu_{i}\leqslant\lambda_{i}$ for $i=1,\ldots,n$. We normalise $J_{\lambda}^{\sigma,a,b}$ by requiring $c_{\lambda,\lambda}=1$ (the “monic” choice).

The multi-variable Jacobi polynomials $\{J_{\lambda}^{2/\beta,\alpha_{1},\alpha_{2}}({\mathbf{x}})\}$ are orthogonal with respect to the joint probability density (1.1) of the JE [56, Théorème 2].

2.3 Multivariate hypergeometric functions

Multivariate hypergeometric functions were introduced for general values of the Jack parameter $\sigma$ by Kaneko [57] and Korányi [58], generalising the definition relevant to the case $\sigma=2$ introduced by Herz [59], the Statistics applications of which being studied in [34, 10, 60]. Efficient numerical implementations of multi-variable hypergeometric functions are available [61].

They are defined as a sum over partitions as

[TABLE]

where $[a]_{\lambda}^{(\sigma)}$ is the generalised Pochhammer symbol defined by

[TABLE]

and the classical Pochhammer symbol $(a)_{n}$ is

[TABLE]

The reader familiar with hypergeometric functions of a single variable (recapitulated in (4.17) below) will recognise the generalisation (2.10).

For general values of the parameters $a_{1},\ldots,a_{p},b_{1},\ldots,b_{q}$ the series in (2.10) converges absolutely for all ${\mathbf{x}}\in{\mathbb{C}}^{n}$ if $p\leqslant q$ and for ${\mathbf{x}}$ in some ball if $p=q+1$ [57]. However, if any of the “upper” parameters, $a_{1}$ say, is equal to a negative integer $-m$, $m>0$, then the series contains only finitely-many terms and defines a multi-variable symmetric polynomial of degree $mn$.

2.4 Some useful identities

Yan undertook one of the first systematic studies of multi-variable hypergeometric functions for arbitrary $\sigma>0$ and proved a number of formulæ and identities, including the Pfaff-like formula [62, eq. (35)]

[TABLE]

(The case $\sigma=2$ was derived in [60, Theorem 7.4.3].)

A number of integral representations are also available. We mention here, and will use below, the formula due to Kaneko [57]:

[TABLE]

valid for $\mathop{\mathfrak{Re}}\{a\}>0,\mathop{\mathfrak{Re}}\{b\}>0$, $\mathop{\mathfrak{Re}}\{1/\sigma\}>-\min\{1/n,\mathop{\mathfrak{Re}}\{a\}/(n-1),\mathop{\mathfrak{Re}}\{b\}/(n-1)\}$, and ${\mathbf{t}}=(t_{1},\ldots,t_{m})$. In (2.14) and below we use ${\mathrm{d}}^{n}{\mathbf{x}}$ as a shorthand for ${\mathrm{d}}x_{1}\cdots{\mathrm{d}}x_{n}$.

3 Calculations for finite-size matrices

In this section we collect some formulæ for the distribution and density of the smallest eigenvalue of a JE matrix of fixed finite size $N\times N$.

3.1 Probability distribution of the smallest eigenvalue

If $\phi_{1}$ is the smallest eigenvalue of the JE then, for any constant $\xi\in{\mathbb{R}}$,

[TABLE]

As all eigenvalues are between [math] and $1$, we obviously have

[TABLE]

so it will be sufficient to find expressions for the probability ${\mathbb{P}}(\phi_{1}>\xi)$ in the range $0<\xi<1$.

Proposition 3.1**.**

Let $\phi_{1}$ be the smallest eigenvalue of the joint distribution (1.1) with $\alpha_{1}\in{\mathbb{N}}$. Then, for $0<\xi<1$,

[TABLE]

Proof.XRecalling that $x_{1},\ldots,x_{n}$ are un-ordered eigenvalues, we integrate the joint probability density (1.1), to get

[TABLE]

If we make the substitution $y_{i}=(x_{i}-\xi)/(1-\xi)$, $1\leqslant i\leqslant N$, this maps each of the integrals to an integral over $[0,1]$, and we have

[TABLE]

For $\alpha_{1}\in{\mathbb{N}}$ this integral can be evaluated by means of Kaneko’s integral (2.14) to get

[TABLE]

Up to this point we have followed Borodin and Forrester’s paper [4] (our (3.6) is equation (3.16) therein). The only novel step in the proof is to simplify the argument of the multivariate hypergeometric function in (3.6) by applying the Pfaff-like identity (2.13) to give (3.3).  $\Box$

Based on (3.6), Borodin and Forrester proved the asymptotic scaling limit (1.9) for the smallest eigenvalue.

We have also a formula for ${\mathbb{P}}(\phi_{1}>\xi)$ in terms of multi-variable Jacobi polynomials.

Corollary 3.2**.**

With $\phi_{1}$ and $\alpha_{1}\in{\mathbb{N}}$ as above, an alternative expression for the probability in Proposition 3.1 is, for $0<\xi<1$,

[TABLE]

where $P_{\lambda}^{\sigma,a,b}({\mathbf{x}})$ is the multi-variable Jacobi polynomial, and an explicit expression for the denominator in (3.7) is

[TABLE]

In these expressions ${\mathbf{0}}^{n}$ is a shorthand for $(0,\ldots,0)\in{\mathbb{R}}^{n}$.

Corollary 3.2 will be proved in Section 5.2. We also note that explicitly computable recursions for coefficients in the series expansion in powers of $\xi$ for $F_{\phi_{1}}(\xi)$ have been derived in [63].

3.2 Probability density of the smallest eigenvalue

Our main interest is in the probability distribution function of the smallest eigenvalue $\phi_{1}$ of the JE. However with little effort we can derive a formula for a probabilty density in terms of a multi-variable hypergeometric function, that will also be used to prove a key differentiation identity (Corollary 3.4 below).

Proposition 3.3**.**

If $\alpha_{1}\in{\mathbb{N}}$, a marginal probability density function for the smallest eigenvalue $\phi_{1}$ of the Jacobi $\beta$-Ensemble (1.1) is given by

[TABLE]

for $0\leqslant\phi_{1}\leqslant 1$, where the normalisation constant is

[TABLE]

Proof.XThe joint probability density function of the ordered eigenvalues of the JE is

[TABLE]

This has the same functional form as (1.1), except for the factor $N!$ in the numerator to account for the ordering of the variables. To derive the marginal density function for $\phi_{1}$ we integrate out all the other variables

[TABLE]

un-ordering the $N-1$ integrations. With the change of variables $y_{i}=(\phi_{i}-\phi_{1})/(1-\phi_{1})$ for $i=2,\ldots,N$, this multiple integral becomes

[TABLE]

where, in a slightly unusual notation ${\mathbf{y}}=(y_{2},\ldots,y_{N})\in{\mathbb{R}}^{N-1}$. The $(N-1)$-fold multiple integral may be evaluated by means of Kaneko’s integral (2.14) to give

[TABLE]

By an application of the Pfaff-like identity (2.13), this may be re-written as (3.9).  $\Box$

The formula (3.9) for the probability density function was first derived by Dumitriu [36], with a different method of proof. Slightly different, but equivalent, multivariable hypergeometric function representations for the probability density function have been given in [37].

Corollary 3.4**.**

For $\alpha_{1}\in{\mathbb{N}}$ we have the derivative identity

[TABLE]

for all $\xi\in{\mathbb{C}}$ except possibly $\xi=1$.

Proof.XFrom (3.1) and (3.3) above the probability distribution function of $\phi_{1}$, the smallest eigenvalue, is

[TABLE]

for $0\leqslant\xi\leqslant 1$, and a probability density function is given by (3.9). The result (3.15) follows because the density function agrees with the derivative of the distribution function at points of continuity. By analytic continuation the identity persists outside of the interval $0<\xi<1$.  $\Box$

We remark that the result (3.15) does not seem easy to prove in a direct way starting from the definition (2.10) of the multivariate hypergeometric functions. A similar observation was made by Forrester [19] who found the analogous identity at the level of $\fourIdx{}1{(\sigma)}1{F}$ multivariate hypergeometric functions.

Later, we will want to take the limit $N\to\infty$, so we record here the asymptotic behaviour of $Z_{N}$ in this limit.

Lemma 3.5**.**

As $N\to\infty$ we have

[TABLE]

Proof.XUsing the value (1.2) for the Selberg integrals, and cancelling common factors we get

[TABLE]

Re-writing the factors appearing in the product in (3.18) as

[TABLE]

we realise that many factors cancel in the product over $k$ and we are left with

[TABLE]

Reuniting the product with the prefactors in (3.18) and further cancellation results in

[TABLE]

The asymptotic (3.17) follows by applying the asymptotic formula

[TABLE]

to the $N$-dependent factors. $\Box$

4 Two-term asymptotic formula

Our main analytic tool is going to be a two-term asymptotic formula for the $\fourIdx{}2{(\sigma)}1F$ multi-variable hypergeometric function, stated below, and proved in the following subsections.

Theorem 4.1**.**

Let $a,b,c\in{\mathbb{C}}$ and $\sigma>0$ be fixed, such that $c-(i-1)/\sigma$ is not a negative integer for $1\leqslant i\leqslant n$. Then with $p_{1}({\mathbf{x}})\mathbin{\hbox{\raise 0.34444pt\hbox{\rm:}}\!\!=}x_{1}+\cdots+x_{n}$,

[TABLE]

where the error estimate is uniform for ${\mathbf{x}}$ in compact subsets of ${\mathbb{C}}^{n}$, but may depend on $a,b,c,\sigma,n$. The operator $\mathrm{E}_{1}$ in (4.1) is defined in (2.2).

Our strategy will be to split the sum defining the $\fourIdx{}2{(\sigma)}1F$ multi-variable hypergeometric function as

[TABLE]

recalling that the Jack polynomials $C_{\lambda}^{(\sigma)}$ are homogeneous of order $|\lambda|$. It will turn out that the tail terms (the second sum) do not contribute significantly to the $N\to\infty$ limit.

4.1 Preliminary results

It is a consequence of Stirling’s formula that

[TABLE]

as $|z|\to\infty$. We will need a form of this result with control on how the error depends on the parameters $\alpha,\beta$.

Lemma 4.2**.**

Suppose $\alpha\in{\mathbb{C}}$ is a quantity such that $|\alpha|^{2}/z={\mathrm{o}}(1)$ as $z\to\infty$. Then

[TABLE]

as $z\to\infty$ with $|\arg\{z\}|<\pi$.

Proof.XBy the classical Stirling formula [64, 6.1.37]

[TABLE]

as $z\to\infty$ with $|\arg\{z\}|<\pi$. Therefore,

[TABLE]

using

[TABLE]

and

[TABLE]

provided $|\alpha/z|\leqslant 1/2$.

Now,

[TABLE]

provided $|\alpha|^{2}/|z|\leqslant 1/2$. Putting it together with (4.6) we get (4.4).  $\Box$

Corollary 4.3**.**

Suppose $\alpha$ and $\beta$ satisfy $|\alpha|^{2}/z={\mathrm{o}}(1)$ and $|\beta|^{2}/z={\mathrm{o}}(1)$ as $z\to\infty$. Then

[TABLE]

Proof.XApplying Lemma 4.2 and cancelling common factors,

[TABLE]

by the binomial theorem. Re-writing the second term gives (4.10).  $\Box$

We will also require some bounds on Pochhammer symbols (2.12).

Lemma 4.4**.**

Let $a\in{\mathbb{C}}$ be fixed, and $n,N\in{\mathbb{N}}$.

If $n\leqslant N$ then $|(a-N)_{n}|\leqslant(|a|+N-n+1)_{n}\leqslant(N+|a|)^{n}$; 2. 2.

If $n>N$ then $|(a-N)_{n}|\leqslant(|a|)_{n}$.

Proof.XWe have that

[TABLE]

If $n\leqslant N$ then using the second line of (4.12) and the triangle inequality

[TABLE]

If $n>N$ then re-ordering the product,

[TABLE]

However,

[TABLE]

so from (4.14) we end up with

[TABLE]

$\Box$

We shall test certain sums for convergence by comparison with the classical hypergeometric series

[TABLE]

For generic choice of parameters $a_{1},\ldots,a_{p},b_{1},\ldots,b_{q}\in{\mathbb{C}}$ the power series (4.17) is known to have radius of convergence $r=1$ if $p=q+1$ and infinite radius of convergence if $p\leqslant q$ [65, §2.2]. The exceptions to these rules are when the series has only a finite number of terms and reduces to a polynomial in $z$. This can happen when one of the parameters $a_{1},\ldots,a_{p}$ is a negative integer.

4.2 The main contribution

We start by analysing the first sum on the right-hand side of (4.2).

Proposition 4.5**.**

Let $a,b$ be fixed quantities. Then

[TABLE]

where the implied constant may depend on $n,a,b,c$ and $\sigma$ but is uniform for ${\mathbf{x}}$ in compact subsets of ${\mathbb{C}}^{n}$.

Making use of the representation

[TABLE]

we can write, using Corollary 4.3 for the ratios of gamma functions,

[TABLE]

where $\sigma^{\prime}=\max\{1,\sigma^{-1}\}$. This leads to

[TABLE]

We adopt here a convenient shorthand $\|\lambda\|_{4}^{4}=\lambda_{1}^{4}+\cdots\lambda_{n}^{4}$. Recalling the actions (2.3), (2.4) of the operators $\mathrm{E}_{1}$, $\mathrm{D}_{2}$ on Jack polynomials,

[TABLE]

Substituting (4.22) for the numerator in (4.18) leads via cancellation of the factors $N^{2|\lambda|}$ to the sum on the right-hand side of (4.18). To prove the error estimate it will be sufficient to demonstrate that

[TABLE]

uniformly for ${\mathbf{x}}$ in compact subsets of ${\mathbb{C}}^{n}$. We do this following a method elaborated by Kaneko [57]. Namely: there exists a constant $C$ depending only on $n$ such that

[TABLE]

where $\|{\mathbf{x}}\|_{\infty}=\max\{|x_{1}|,\ldots,|x_{n}|\}$ and $\sigma^{\prime}=\max\{1,\sigma^{-1}\}$ [57, Lemma 1]. Thus there exists a constant $R>0$ depending only on $n,\sigma$ such that $\|\lambda\|_{4}^{4}|C_{\lambda}^{(\sigma)}({\mathbf{x}})|\leqslant CR^{|\lambda|}\|{\mathbf{x}}\|_{\infty}^{|\lambda|}$. We also may observe that

[TABLE]

So

[TABLE]

using the definition (2.11) of the generalised Pochhammer symbol together with (4.25). It can be seen that (4.26) is bounded by comparing each factor to a convergent $\fourIdx{}{0}{}{1}{F}$ hypergeometric series.  $\Box$

4.3 Bounding the tail terms

Using Kaneko’s bound for Jack polynomials from the end of subsubsection 4.2 we get

[TABLE]

for some constant $C$ depending only on $n$ and $R$ depending only on $n$ and $\sigma$. If $\lambda$ is a partition with $|\lambda|\geqslant N^{1/3}$ then, as $\lambda_{1}$ is the largest part, we must have $\lambda_{1}\geqslant N^{1/3}/n$. Factorising the generalised Pochhammer symbols according to (2.11) we achieve the inequality

[TABLE]

Proposition 4.6**.**

For every $\rho>0$ and compact set $K\subseteq{\mathbb{C}}^{n}$, there exists a constant $C_{\rho,K}$ (which may additionally depend on $a,b,c,\sigma,n$) such that for every ${\mathbf{x}}\in K$, and all $N$ sufficiently large, we have

[TABLE]

Proof.XFor brevity, let us define $a^{\prime}=a-(i-1)/\sigma$, $b^{\prime}=b-(i-1)/\sigma$, $c^{\prime}=c-(i-1)/\sigma$. The following three estimates give us the bound we need:

Using part 2. of Lemma 4.4,

[TABLE]

provided we additionally choose $N$ large enough so that $R\|{\mathbf{x}}\|_{\infty}/N<1$, and the absolute convergence of a $\fourIdx{}2{}1{F}$ hypergeometric series in the unit disc. 2. 2.

Using part 1. of Lemma 4.4,

[TABLE]

by comparison with a $\fourIdx{}0{}1{F}$ hypergeometric series. 3. 3.

For the sum over $\lambda_{1}$ between $N^{1/3}/n$ and $N$ we again use part 1. of Lemma 4.4, as in the previous step leading to

[TABLE]

where the series has been compared to a $\fourIdx{}1{}1{F}$ hypergeometric series.

We use 1. and 2. to bound each factor for $i\geqslant 2$ in (4.28) by a constant. We then use 1. and 3. to deduce the rapid decay in $N$ of the remaining sum over $\lambda_{1}$.  $\Box$

With essentially the same method and calculations we can bound similarly the tail of two further series.

Proposition 4.7**.**

For every $\rho>0$ and compact set $K\subseteq{\mathbb{C}}^{n}$, there exists a constant $C_{\rho,K}$ (which may additionally depend on $a,b,c,\sigma,n$) such that for every ${\mathbf{x}}\in K$, and all $N$ sufficiently large, we have

[TABLE]

and

[TABLE]

Proof.XUsing the fact that $C_{\lambda}^{(\sigma)}({\mathbf{x}})$ is an eigenfunction of $\mathrm{E}_{1}$ and $\mathrm{D}_{2}$ with eigenvalues that depend only polynomially on the parts of $\lambda$, we may follow the proof of the preceding Proposition 4.6 making only trivial changes.  $\Box$

4.4 Asymptotic formula

Proposition 4.8**.**

For fixed $a,b$, we have

[TABLE]

where the implied constant may depend on $n,a,b,\sigma$ but is uniform for ${\mathbf{x}}$ in compact subsets of ${\mathbb{C}}^{n}$.

Proof.XStarting from (4.2) and as a consequence of the bound (4.29) of Proposition 4.6, we have

[TABLE]

for any $\rho>0$. By Proposition 4.5 this is

[TABLE]

By Proposition 4.7 we can complete the sums in (4.37) without affecting the error estimate to get

[TABLE]

This is (4.35), recognising

[TABLE]

$\Box$

Given that $\fourIdx{}0{(\sigma)}1{F}(;c;{\mathbf{x}})$ is continuous, this already proves

[TABLE]

locally uniformly in ${\mathbf{x}}\in{\mathbb{C}}^{n}$—the leading-order of Theorem 4.1.

Our final task in this section will be to put the “$1/N$” term of (4.35) into a nicer form.

4.5 Partial Differential Equation satisfied by $\fourIdx{}0{(\sigma)}1{F}$

It is a Theorem of Yan [62, Theorem 2.1] and Kaneko [57, Theorem 2] that if $c-(i-1)/\sigma$ is not a negative integer for any $1\leqslant i\leqslant n$ then the unique solution to system of equations

[TABLE]

$1\leqslant i\leqslant n$, subject to $\varPsi({\mathbf{x}})$ being symmetric in its variables and analytic at ${\mathbf{x}}={\mathbf{0}}$, is

[TABLE]

This result for $\sigma=2$ was first proved by Muirhead [66], having been conjectured, apparently, by A. G. Constantine. Muirhead also shows how the system (4.41) can be degenerated to give the holonomic system of equations for $\fourIdx{}0{(2)}{1}{F}$ multivariate hypergeometric functions, which can easily be generalised for arbitrary $\sigma$ as follows.

Proposition 4.9**.**

Provided that $c-(i-1)/\sigma$ is not a negative integer for any $1\leqslant i\leqslant n$, the multivariate hypergeometric function $\fourIdx{}0{(\sigma)}1{F}(c;{\mathbf{x}})$ is the unique solution of the $n$ differential equations

[TABLE]

$1\leqslant i\leqslant n$, subject to the constraints that $\varPsi({\mathbf{x}})$ is symmetric in its variables, is analytic at ${\mathbf{x}}={\mathbf{0}}$ and satisfies $\varPsi({\mathbf{0}})=1$.

Proof.XSince we now know (4.40) that

[TABLE]

we set $a=b=-N$ and make the change of variables $x_{i}\mapsto x_{i}/N^{2}$ in (4.41) to get

[TABLE]

Dividing through by $N^{2}$ and letting $N\to\infty$ we recover (4.43).  $\Box$

Corollary 4.10**.**

Under the same condition on $c$ as in Proposition 4.9, the function $\varPsi({\mathbf{x}})=\fourIdx{}0{(\sigma)}1{F}(;c;{\mathbf{x}})$ is a solution to the partial differential equation

[TABLE]

where $p_{1}({\mathbf{x}})=x_{1}+\cdots+x_{n}$.

Proof.XWe multiply through the $i$th equation (4.43), satisfied by $\varPsi({\mathbf{x}})=\fourIdx{}0{(\sigma)}1{F}(;c;{\mathbf{x}})$, by $x_{i}$:

[TABLE]

Since

[TABLE]

(4.47) is equivalent to

[TABLE]

Summing over (4.49) for $i=1,\ldots,n$, we arrive at

[TABLE]

Re-arranged, this is (4.46).  $\Box$

Proof of Theorem 4.1.X We use Corollary 4.10 to replace the term $\mathrm{D}_{2}\left\{\fourIdx{}{0}{(\sigma)}1{F}(;c;{\mathbf{x}})\right\}$ in (4.35) by

[TABLE]

The resultant cancellation of terms involving $\sigma$ leads to (4.1).  $\Box$

5 Main Results

5.1 Proof of Theorem 1.1

Looking-back to (3.3) we had

[TABLE]

Standard asymptotic arguments give

[TABLE]

uniformly for ${x}$ in compact sets.

In the result of Theorem 4.1, taking ${\mathbf{x}}$ to be a constant multiple of ${\mathbf{1}}^{n}$ we have

[TABLE]

which may be applied in (5.1) with $a=0$, $b=1-2(\alpha_{2}+1)/\beta$, $c=2\alpha_{1}/\beta$, $\sigma=\beta/2$ and $n=\alpha_{1}$ (so that $c-(i-1)/\sigma=2(\alpha_{1}+1-i)/\beta$ is not a negative integer for $i=1,\ldots,\alpha_{1}$, justifying the use of Proposition 4.9), to give

[TABLE]

Combining with (5.2) we get

[TABLE]

Setting $\xi=x/N^{2}$ in (3.15) of Corollary 3.4

[TABLE]

where $Z_{N}$ was defined in (3.10). A further application of (the leading order of) Theorem 4.1 and equation (5.2), and Lemma 3.5 for the asymptotics of $Z_{N}$, brings (5.6) to

[TABLE]

We use (5.7) to remove the derivative term from (5.5), giving

[TABLE]

This completes the proof.  $\Box$

That the first-order correction term in (5.8) is proportional to the derivative of the leading term implies that the finite-$N$ behaviour can be interpreted, up to an error of order ${\mathrm{O}}(N^{-2})$ as a correction to the width: if we let

[TABLE]

then we may interpret Theorem 1.1 as saying

[TABLE]

By Taylor’s theorem, this is equivalent to the re-centring

[TABLE]

5.2 Connection with Jacobi polynomials

We now prove the formula (3.7) from Corollary 3.2 giving a formula for the distribution of the smallest JE eigenvalue in terms of multi-variable Jacobi polynomials.

Proof of Corollary 3.2.X We have for ${\mathbf{x}}\in{\mathbb{C}}^{n}$

[TABLE]

This is essentially Théorème 5 of [56], incorporating our different choice of normalisation of the Jacobi polynomials. It may be proved by observing that both sides of (5.12) are multivariate symmetric polynomials that satisfy the same partial differential equation (see [62] or [57] for the PDE satisfied by $\fourIdx{}2{(\sigma)}1{F}$) and that the leading term on both sides is proportional to $C_{(N^{n})}^{(\sigma)}({\mathbf{x}})$ with identical constant.

Comparing to (3.6) we find the appropriate parameters for the multivariate Jacobi polynomial are

[TABLE]

From (3.6) and the identity (5.12) with parameters as above we have

[TABLE]

Taking the limit $\xi\to 0^{+}$ we need to have ${\mathbb{P}}(\phi_{1}>0)=1$ and so

[TABLE]

Some of the quantities in (5.15) simplify: we have, starting from (2.11),

[TABLE]

and

[TABLE]

so that

[TABLE]

In (5.15) the ratio $[2\alpha_{1}/\beta]_{(N^{\alpha_{1}})}^{(\beta/2)}/[-N]_{(N^{\alpha_{1}})}^{(\beta/2)}$ is of the form (5.18) and the factor $[N-1+2(\alpha_{1}+\alpha_{2}+1)/\beta]_{(N^{n})}^{(\beta/2)}$ is of the form (5.17). Thus we get

[TABLE]

$\Box$

5.3 Numerical simulations

In order to visualise our results better, we present in this subsection the results of some numerical simulations, and compare them against our theoretical predictions.

In Figures 1 and 2 we present numerical results for two different values of $\beta$: $\beta=1$ (Figure 1) and $\beta=3$ (Figure 2). In both cases we plot the empirical distribution function for 1 000 samples of the scaled smallest eigenvalue of JE random matrices for $N=30$ and $N=20$ respectively, together with the same calculation at $N=1000$. We observe that for the smaller value of $N$, the empirical distribution fits well to the two-term asymptotic prediction (1.10) of Theorem 1.1. At $N=1000$, the data follows the leading order term of (1.10) (which is the formula (1.9) derived by Borodin and Forrester [4]), the order $1/N$ corrections being negligble at such matrix size.

The implementation details differ somewhat between Figures 1 and 2. To construct random matrices from the Jacobi Orthogonal Ensemble ($\beta=1$) to generate the data for Figure 1, we used the double-Wishart matrix construction starting from independent $(N+1)\times N$ and $(N+7)\times N$ matrices as described in the Introduction, which results in samples from the $\alpha_{1}=0,\alpha_{2}=3$ JOE. For the theoretical prediction with $\alpha_{1}=0$ we were able to use the simpler formula (6.54) from Section 6.3 below, rather than (1.10).

To generate the $\beta=3$ data for Figure 2 there is no double-Wishart construction available, so we instead implemented the algorithm of Killip and Nenciu [14] which can generate samples from the JE for arbitrary $\beta>0$. For calculating the values of the multivariate hypergeometric functions needed for the theoretical prediction (1.10) we used the numerical routines of Koev and Edelman [61].

6 Explicit formulæ

In certain situations we are able to derive expressions for $F_{\phi_{1}}$ and asymptotics for $F_{N^{2}\phi_{1}}$ that are more explicit, and these are expounded in the present Section.

To begin with we focus on the case $\beta=2$ corresponding to the Jacobi Unitary Ensemble. In this case we benefit from the fact that the multi-variable Jacobi polynomials enjoy a determinantal structure. In Section 6.3 we record some formulæ (for arbitrary $\beta>0$) for the special cases $\alpha_{1}=0$ and $\alpha_{1}=1$.

6.1 Determinantal identities

In order to state the determinantal identities let us denote by $p_{m}^{a,b}(x)$ the $m$th monic Jacobi polynomial orthogonal with respect to the measure $x^{a}(1-x)^{b}$ on the interval $[0,1]$. In terms of the definition $P_{m}^{(a,b)}$ of Jacobi polynomials given by Szegő [67, Ch. IV], our definition satisfies

[TABLE]

As a hypergeometric function,

[TABLE]

We also need the hook-length $h_{\lambda}$ for a partition, defined by

[TABLE]

with $\ell(\lambda)$ the number of non-zero parts.

If $\sigma=1$ the following Lemma gives alternative expresions for the multivariate Jacobi polynomials.

Lemma 6.1**.**

Let ${\mathbf{x}}\in{\mathbb{C}}^{n}$. Then

[TABLE]

where $h_{\lambda}$ is the hook-length of the partition $\lambda$ and $\Delta({\mathbf{x}})$ is the Vandermonde determinant (1.3). If, furthermore ${\mathbf{x}}=x{\mathbf{1}}^{n}$, $x\in{\mathbb{C}}$, then we have

[TABLE]

Proof.XThat the multivariate Jacobi polynomial at $\sigma=1$ has a determinant evaluation in terms of univariate Jacobi polynomials is known since [56, Théorème 10]:

[TABLE]

Lasalle uses a different version of the Jacobi polynomials to us which changes the numerical value of the constant in (6.6). We can fix the constant by observing that since our Jacobi polynomials are monic,

[TABLE]

where $\mathfrak{s}_{\lambda}$ is a Schur polynomial. The Jack polynomials at $\sigma=1$ are proportional to Schur polynomials [48], and in fact

[TABLE]

The implication is

[TABLE]

In our applications ${\mathbf{x}}$ is a scalar multiple of ${\mathbf{1}}^{n}$. To take the confluent limit ${\mathbf{x}}\to x{\mathbf{1}}^{n}$, we use the formula

[TABLE]

proved in [68, Lemma A.1], where $\mathcal{W}$ denotes the Wronskian

[TABLE]

and the functions $\varphi_{1},\ldots,\varphi_{n}$ must be regular at $x$. (A version of (6.10) valid for polynomials $\varphi_{1},\ldots,\varphi_{n}$ was proved in [69, Theorem 1], which would suffice to handle (6.4), but later we will have reason to apply (6.10) to non-polynomial functions.) In the application we presently have in mind, $\varphi_{i}$ is the Jacobi polynomial $p_{\lambda_{i}+n-i}^{a,b}$ and using the fact that

[TABLE]

which extends to

[TABLE]

we have

[TABLE]

Combining (6.4), (6.10) and (6.14) we get (6.5).  $\Box$

Corollary 6.2**.**

For ${\mathbf{x}}\in{\mathbb{C}}^{n}$, and fixed $c$ such that $c-i$ is not a negative integer for $i=0,\ldots,n-1$, we have

[TABLE]

where $I_{\nu}(z)$ is the $I$-Bessel function. If ${\mathbf{x}}=x{\mathbf{1}}^{n}$, $x\in{\mathbb{C}}$ we have a further formula:

[TABLE]

A special case (with $c=2n$) of (6.16) was given in [70], where it was proved by using the relationships with Painlevé functions [71]. We further remark that (6.15) could be proved in a less direct way through the use of tau functions of hypergeometric type [72].

Proof.XWe know that if $c-(i-1)/\sigma$ is not a negative integer for any $i=1,\ldots,n$ then

[TABLE]

uniformly for ${\mathbf{x}}$ in compact subsets of $\mathbb{C}^{n}$. From (5.12)

[TABLE]

Setting $\sigma=1$ and applying (6.4),

[TABLE]

It is problematic to take the limit $N\to\infty$ here directly. The determinant of Jacobi polynomials tends to [math] as $N\to\infty$ but to find the rate and leading-term some further manipulations are necessary. These involve repeated use of the contiguous identity

[TABLE]

to add successively to each row a multiple of the row below, in a recursive fashion, to get that

[TABLE]

and we have

[TABLE]

We return to the hypergeometric function representation of Jacobi polynomials (6.2),

[TABLE]

and observe that

[TABLE]

so

[TABLE]

Using (4.3) and the fact that

[TABLE]

we have the asymptotic behaviour

[TABLE]

as $N\to\infty$. Putting this into the determinant from (6.22),

[TABLE]

as $N\to\infty$.

We already know (equation (5.16)) that

[TABLE]

and we have

[TABLE]

so

[TABLE]

Putting (6.29) and (6.31) together into the right-hand side of (6.19) we find that all the factors of $N$ cancel, and we recover the $N\to\infty$ limit, which yields

[TABLE]

We can derive an expression involving the more familar Bessel functions by means of the identity [73, §7.8, eq. (1)]

[TABLE]

Inserting this into (6.32), we have

[TABLE]

which is equivalent to (6.15).

To take the confluent limit ${\mathbf{x}}\to x{\mathbf{1}}^{n}$, we prefer to work with (6.32). Using, for a second time, identity (6.10),

[TABLE]

Since

[TABLE]

we have

[TABLE]

and

[TABLE]

Putting this into (6.35),

[TABLE]

reversing the order of the rows in the determinant. Using (6.33) this becomes

[TABLE]

As a final step we use the identity $\det(\theta^{i-j}a_{ij})=\det(a_{ij})$ to reduce this to (6.16).  $\Box$

6.2 Smallest eigenvalue of the Jacobi Unitary Ensemble

Proposition 6.3**.**

Let $F_{\phi_{1}}$ be the probability distribution function of the smallest eigenvalue of the $N\times N$ Jacobi $\beta$-Ensemble with $\beta=2$, $\alpha_{1}\in{\mathbb{N}}_{0}$ and $\alpha_{2}>-1$. For $0<\xi<1$,

[TABLE]

As $N\to\infty$ we have, for $x>0$,

[TABLE]

where the determinants in (6.42) are of size $\alpha_{1}\times\alpha_{1}$.

Proof.XIn view of (3.7) we need to specialise (6.5) to the rectangular partition $\lambda=(N^{n})$, to get

[TABLE]

We have already seen (equation (6.30)) that

[TABLE]

and, after cancellation, (6.43) becomes

[TABLE]

In the determinant in (6.45) there is a factor $(N+n-i)!$ multiplying the $i$th row. If we extract these factors, they cancel the factorials in the denominator. Finally, we reverse the order of the rows producing a factor $(-1)^{\lfloor n/2\rfloor}$ and the end result

[TABLE]

So, if $\beta=2$,

[TABLE]

where, from (5.19),

[TABLE]

Note that since $\alpha_{1}!(1+\alpha_{1})_{N-1}=(N+\alpha_{1}-1)!$ and $(N-1)!(N)_{\alpha_{1}}=(N+\alpha_{1}-1)!$ too we get a simplified formula

[TABLE]

Equations (6.47) and (6.49) yield (6.41).

We now turn to the two-term asymptotic formula (6.42). With $\beta=2$ in (1.10),

[TABLE]

Using the representation (6.16) for the multi-variable hypergeometric functions this becomes

[TABLE]

and, upon cancellation of the gamma function factors,

[TABLE]

$\Box$

The formula (6.52) proves a conjecture made in [3]. (The result is also implicit in the recent work [41].) It seems likely that explicit formulæ along the lines of (6.52) will also be available in the other privileged cases $\beta=1$ and $\beta=4$. The details will appear elsewhere.

6.3 Small values of $\alpha_{1}$

If $\alpha_{1}=0$, then most of the analysis above is quite unnecessary and we already recover from (3.5)

[TABLE]

recognising the value of the Selberg integral (or observing that both sides must be unity as $\xi\to 0^{+}$). So, for $\alpha_{1}=0$ and $x>0$,

[TABLE]

referring-back to (5.2). This generalises, to arbitrary $\beta>0$, Corollary 1 of [3].

If $\alpha_{1}=1$ then the multivariate hypergeometric funtions of $\alpha_{1}$ arguments become ordinary one-variable hypergeometric functions. In this case, for $0\leqslant\xi\leqslant 1$,

[TABLE]

by (3.3). Using Corollary 3.2, and the fact that multi-variable Jacobi polynomials of a single variable coincide with classical single-variable ones, this may be expressed further as

[TABLE]

We apply Theorem 4.1 with $n=1$ to (6.55) to find that

[TABLE]

using (6.36) for the derivative of the hypergeometric function. We may further use (6.33) to replace hypergeometric functions with Bessel functions, yielding

[TABLE]

We note that

[TABLE]

an application of the Bessel function identity [73, §7.11, eq. (23)]

[TABLE]

Putting (6.59) into (6.58) gives

[TABLE]

This result is consistent with (6.52) when $\beta=2$.

Acknowledgements

The author wishes to acknowledge helpful conversations about this work with J. P. Keating and D. Savin.