Smallest eigenvalue distributions for two classes of $beta$-Jacobi ensembles

Smallest eigen v alue distributions for t w o classes of β -Jacobi ensem bles Ioana Dumitriu ∗ August 12, 2021 Abstract W e compute the exact and limiting smallest eigen v alue distributions for a class of β -Jacobi ensembles not cov ered by previous studies. In the general β case, these distributions are giv en by multiv ariate h yp ergeometric 2 F 2 /β 1 functions, whose b eha vior can b e analyzed asymptotically for sp ecial v alues of β whic h include β ∈ 2 N + as well as for β = 1. Interest in these ob jects stems from their connections (in the β = 1 , 2 cases) to principal submatrices of Haar-distributed (orthogonal, unitary) matrices app earing in randomized, communication-optimal, fast, and stable algorithms for eigen v alue computations [8], [4]. 1 In tro duction This pap er was b orn from the author’s in terest in the following problem. Motiv ating Problem (MP) . Given an n × n Haar-distributed orthogonal or unitary matrix U n , c ho ose a principal submatrix M r formed b y r rows and r columns, r ≤ n/ 2. Since the Haar (uniform) distribution is in v ariant under p ermutations of ro ws and columns, w e migh t as well assume that M r is the upp er-left corner principal submatrix (in Matlab notation, M r = U n (1 : r , 1 : r )). 1. What is the distribution of the smallest singular v alue of M r ? 2. If r , n → ∞ , how do es the singular v alue asymptotically dep end on n and r ? 3. In particular, what kind of asymptotics do we get when r = Θ( n )? The MP arose in the author’s work in scientific computing, namely on (parallelizable, fast, and stable) randomized algorithms [8, 4], since the quality of the randomized rank-rev ealing decomp osition R UR V presen ted there dep ends on the order of the said smallest singular v alue distribution. The main in terest in that work is for the case when n and r are extremely large (the matrix U n itself does not fit into the processor’s cac he memory); hence Question 2. The reason why Question 3 is of particular imp ortance is b ecause, in that con text, r represen ts the split following a “Divide-and-Conquer” step in an eigen v alue algorithm (and w e must expect that r will be a fraction of the total n umber of eigenv alues n ). The distribution, naturally , dep ends on whether the matrix is orthogonal or unitary . Question 1 relates to w ork b y Collins [6, 7]. Particular instances of Question 2 ha ve b een solv ed by Jiang [15, 16], who, when r = o ( n/ log n ), has shown the m uch more general result that M r M T r con verges entry-b y- en try to a square Wishart distribution (for whic h the extremal eigenv alue asymptotics are known, see [11]). W e set out to examine Question 3, via Question 2; in the process, it transpired that Questions 2 and 3 b oth can b e dealt with simultaneously . W e hav e obtained thus exact answers for all r , n , as ∗ Departmen t of Mathematics, Universit y of W ashington, Seattle, W A. E-mail: dumitriu@math.washington.edu 1 w ell as asymptotics for b oth r fixed, n growing to infinity , and for r , n gro wing to infinity in any way , including when r = Θ( n ); then we noted that all results generalize in a wa y which we explain b elow to yield v arious asymptotics of extremal eigen v alues of sp ecial β -Jacobi ensembles. The β -Jacobi ensem ble, defined b elow for general β > 0, was the k ey to solving the MP . Definition 1.1. The numb ers λ 1 , λ 2 , . . . , λ m ∈ [0 , 1] ar e β -Jac obi distribute d with p ar ameters a, b > − 1 if their joint distribution is given by f β ,a,b,m ( λ 1 , . . . , λ m ) = 1 c β ,a,b,m m Y i =1 λ β 2 ( a +1) − 1 i (1 − λ i ) β 2 ( b +1) − 1 Y i 0 , ( a ) β κ ≡ m Y i =1  a − β 2 ( i − 1)  k i , wher e ( x ) k i is the classic al rising factorial, ( x ) k i = Γ( x + k i ) / Γ( x ) . In addition, given the diagr am of the p artition κ (se e Figur e 1), we define for every squar e s the “arm-length” a κ ( s ) as the numb er of squar es to the right of s , and the “le g-length” of s as the numb er of squar es b elow s . Then j β κ ≡ Y s ∈ κ  l κ ( s ) + 2 β (1 + a κ ( s )   l κ ( s ) + 1 + 2 β a κ ( s )  . s κ κ κ = (3, 1) a (s) = 1 l (s) = 0 κ β β β β β β j = (2/ )(1)(1 + 6/ )(2 + 4/ )(4/ )(1 + 2/ )(2/ )(1) β Figure 1: Arm-length and leg-length for the square s , along with j β [3 , 1] . The Jac k p olynomials hav e v arious normalizations; the one that w e use here has the following tw o prop erties: for any integer k and for κ ` k a partition of k , and denoting b y I m the vector of m 1s, X κ ` k C β κ ( X ) = ( x 1 + . . . x m ) k , C β κ ( I m ) =  2 β  2 k k !  mβ 2  κ j β κ . (2) Finally , we can no w give the definition of the h yp egeometric function of multiple v ariables, following [14]. Definition 2.2. The hyp er ge ometric series p F q is given by p F β q ( a 1 , . . . , a p ; b 1 , . . . , b q ; X ) ≡ ∞ X k =0 X κ ` k 1 k ! ( a 1 ) β κ . . . ( a p ) β κ ( b 1 ) β κ . . . ( b q ) β κ C β κ ( X ) . F or the case p ≤ q + 1, the series conv erges everywhere in the h yp ercub e {| x i | < 1 , 1 ≤ i ≤ n } . W e will mostly focus on the cases p = 2 , q = 1 (sometimes with m = 1, in which case we recov er the classical 2 F 1 h yp ergeometric function of one v ariable), and sometimes p = q = 1 or p = 0 , q = 1. It is w orth men tioning tw o more facts ab out h yp ergeometric series, the first one whic h w e quote from [14] (form ula 13.5 with p = 2 , q = 1), and the second one b eing an immediate consequence of the first one (giv en the homogeneit y of Jack p olynomials): lim a →∞ 2 F β 1 ( a, b ; c ; X/a ) = 1 F β 1 ( b ; c ; X ) , and (3) lim a,b →∞ 2 F β 1 ( a, b ; c ; X/ ( ab )) = 0 F β 1 ( c ; X ) . (4) 4 3 Smallest eigen v alue distributions 3.1 Most general case W e are in terested in finding an expression for the probability density function for the smallest (resp ec- tiv ely , largest) eigen v alues of the Jacobi ensemble of parameters a and b (see Definition 1.1). W e w ould like to obtain the marginal distribution f min for the smallest eigenv alue; assume λ = λ m . T o that extent, we will in tegrate out all other v ariables but λ m , to obtain f min ( λ ) = m c β ,a,b,m λ β 2 ( a +1) − 1 (1 − λ ) β 2 ( b +1) − 1 × Z [ λ, 1] ( m − 1) m − 1 Y i =1 λ β 2 ( a +1) − 1 i (1 − λ i ) β 2 ( b +1) − 1 ( λ i − λ ) β · Y i − 1, m ≥ 1, and β > 0, 2 m + b − 1 + 2 β ≥ m − 1, the series is w ell-defined (if, for some κ , ( − 2 m − b + 1 − 2 β ) κ = 0, then necessarily (1 − m ) κ = 0 as well, and so the p oten tially “offending” term in the series do es not in fact exist). Under these conditions, using Prop osition 13.1.7 from [14], we can change the v ariables in the h yp ergeometric expression from (1 − λ ) k − 1 to λ k − 1 : 2 F 4 /β 1 (1 − m, − m − b + 1; − 2 m − b + 1 − 2 β ; { 1 − λ } k − 1 ) = 2 F 4 /β 1 (1 − m, − m − b + 1; 2 + 2 β ( k − 1); { λ } k − 1 ) 2 F 4 /β 1 (1 − m, − m − b + 1; 2 + 2 β ( k − 1); { 1 } k − 1 ) . According to equation (13.14) of the same [14], 2 F 4 /β 1 (1 − m, − m − b + 1; 2 + 2 β ( k − 1); { 1 } k − 1 ) = k − 1 Y j =1 Γ  2 + 2 β j  Γ  2 m + b + 2 β j  Γ  m + 1 + 2 β j  Γ  m + b + 1 + 2 β j  := A m,b,β ,k . W e now group the left and, resp ectively , righ t terms in the fraction ab ov e, and rewrite the product as follows: A m,b,β ,k = k − 1 Y j =1 m − 1 Y i =1 1 i + 1 + 2 β j ( m + b + i + 2 β j ) = k − 1 Y j =1 m − 1 Y i =1 1 ( i + 1) β 2 + j (( m + b + i ) β 2 + j ) = m − 1 Y i =1 Γ  1 + ( i + 1) β 2  Γ  k + ( i + 1) β 2  Γ  k + β 2 ( m + b + i )  Γ  1 + β 2 ( m + b + i )  . Putting everything together, we obtain Theorem 3.10. F or the c ase when β / 2( a + 1) = k ∈ N + , the smal lest eigenvalue of the β -Jac obi ensemble is given by f m ( λ ) = W m,b,β ,k λ k − 1 (1 − λ ) β 2 m ( b + m ) − 1 × 2 F 4 /β 1 (1 − m, − m − b + 1; 2 + 2 β ( k − 1); { λ } k − 1 ) , (15) with W m,b,β ,k = m c β ,b, 1+2 /β ,m − 1 c β , 2 k /β − 1 ,b,m A m,b,β ,k . After appr opriate c anc el lations, we obtain that W m,b,β ,k = Γ  1 + β 2  Γ( k )Γ  k + β 2  · m Γ  k + β 2 m  Γ  1 + β 2 m  · Γ  k + β 2 ( b + m )  Γ  β 2 ( b + m )  . (16) 10 As in the previous case, we can no w start to analyze the asymptotics in tw o regimes. Regime 1: m fixed, b → ∞ . Using (3) and (12), we can obtain the distribution of the scaled eigenv alue y = ( b + m ) λ , as follows. Theorem 3.11. As b → ∞ while m is kept fixe d, the limiting distribution for the sc ale d smal lest eigenvalue ( ( b + m ) λ ) of the β -Jac obi ensemble with β ( a + 1) / 2 = k ∈ N + is given by f m ( y ) =  β 2  k Γ  1 + β 2  Γ( k )Γ  k + β 2  · m Γ  k + β 2 m  Γ  1 + β 2 m  y k − 1 e − β my / 2 1 F 4 /β 1 (1 − m ; 2 + 2 β ( k − 1); {− y } k − 1 ) . Note that in the case β = 2, k = 1, which is the case 2 of the MP , the hypergeometric abov e simply b ecomes a constan t whic h cancels the more complicated terms of the normalization, and yields that f m ( y ) = me − my / 2 , whic h is consisten t with the results of [15] and [11]. Regime 2: m, ( b + m ) → ∞ . This cov ers the case when b is actually fixed. Again, using (3), (4), and (12), the distribution of the scaled eigenv alue y = m ( b + m ) λ is given b elow. Theorem 3.12. As m, ( b + m ) → ∞ , the limiting distribution for the sc ale d smal lest eigenvalue ( m ( b + m ) λ ) of the β -Jac obi ensemble with β ( a + 1) / 2 = k ∈ N + is given by f ( y ) =  β 2  2 k − 1 Γ  1 + β 2  Γ( k )Γ  k + β 2  · y k − 1 e − β y / 2 0 F 4 /β 1 (2 + 2 β ( k − 1); { y } k − 1 ) . Once again, when β = 2 and k = 1, the hypergeometric expression in the ab ov e is constan t, and the answer is f ( y ) = e − y . Again, this is consisten t with [15] and [11]. The adv antage ov er [15] is that when b = n − 2 r and m = r , as is the case in the MP , there are no conditions o v er the growth to ∞ of r and n − r ; in [15], r m ust be o ( n/ log n ). 4 Numerical Exp erimen ts W e hav e tested numerically all the formulae obtained in the previous sections. All plots were obtained in MA TLAB; to compute the hypergeometric functions inv olv ed w e hav e used Ko ev’s easy-to-use h y- p ergeometric series pac k age [23, 5]. F or the Mon te Carlo tests, we hav e used the con venien t to work with bidiagonal mo dels for β -Jacobi distributions, given b elow: J β ,n,a,b = B β ,n,a,b · B T β ,n,a,b , where B β ,n,a,b ∼         c n − s n c 0 n − 1 c n − 1 s 0 n − 1 − s n − 1 c 0 n − 2 c n − 2 s 0 n − 2 . . . . . . − s 2 c 0 1 c 1 s 0 1         , where all c i ’s and c 0 i ’s are indep endent v ariables distributed as follo ws: c i ∼ s B eta  β 2 ( a + i ) , β 2 ( b + i )  , c 0 i ∼ s B eta  β 2 i, β 2 ( a + b + 1 + i )  , 11 and s i = p 1 − c 2 i , for all 1 ≤ i ≤ n , while s 0 i = q 1 − c 0 i 2 , for al 1 ≤ i ≤ ( n − 1). Here B eta ( s, t ) stands for the well-kno wn distribution with pdf prop ortional to x s − 1 (1 − x ) t − 1 on [0 , 1]. This matrix model, together with a pro of that the matrix J β ,n,a,b has eigenv alue pdf giv en b y the Jacobi ensemble with given parameters, is a b eautiful result of Sutton’s [28]. Figure 2 is an illustration of Theorem 3.2. The solid red line corresp onds to the exact distribution of the smallest eigenv alue, as given b y the theorem, for the case when β = 1 . 75, m = 4, a = 2 . 3, and b = 2 . 5. The normalized histogram (obtained using the histnorm ) represen ts the results of 10 , 000 Mon te Carlo tests using the bidiagonal matrix models. 0 0.1 0.2 0.3 0.4 0 1 2 3 4 5 6 7 Smallest eigenvalue distribution at m=4, β = 1.75, a = 2.3, b = 2. 5 Figure 2: The solid red line represen ts the theoretical distribution, while the normalized histogram represents the results of a Mon te Carlo exp eriment with 10 , 000 trials, using J β ,n,a,b . Figures 3 and 4 represent a comparison b etw een the theoretical results of Case 1: a = 2 β − 2 , resp ectiv ely , Theorems 3.7 and 3.8 (illustrated by the solid red lines), and Monte Carlo experiments with 10 , 000, resp ectively 5 , 000 trials (represented b y the normalized histograms). Note that conv ergence o ccurs fairly quickly; for Figure 3, whic h considers the case b → ∞ , the histogram of β ( b + m ) / 2 λ min is v ery close to the plot of the asymptotical distribution, even though b = 10. F or Figure 4, w e only need to take b = 5 and m = 5 to see ho w close the histogram of β m ( b + m ) / 2 λ min is to the solid line represen ting the asymptotical distribution. F or this latter case, also note the singularity at 0; this is caused by the singularit y of the Bessel function in the asymptotical formula, as given by Theorem 3.8. 0 0.2 0.4 0.6 0.8 1 0 2 4 6 Monte Carlo vs asymptotical distributions for β (b+m)/2 λ min with m = 4, b = 10 (resp., ∞ ), β = 1/3, and a = 2/ β −2 Figure 3: The solid red line represents the asymptotical ( b = ∞ ) distribution, while the normalized histogram represen ts the results of a Monte Carlo exp erimen t for b = 10, with 10 , 000 trials. 12 0 0.2 0.4 0.6 0.8 1 1.2 0 2 4 6 8 Monte Carlo vs asymptotical distributions for β m(b+m)/2 λ min with m = 5 (resp., ∞ ), b = 5 (resp., ∞ ), β = 1, and a = 2/ β − 2 Figure 4: The solid red line represents the asymptotical ( m, b = ∞ ) distribution, while the normalized histogram represen ts the results of a Monte Carlo exp erimen t for m = 5, b = 5, with 5 , 000 trials. Finally , figures 5 and 6 represen t comparisons betw een the theoretical results of Case 2: a = 2k β − 1 , sp ecifically , Theorems 3.11 and 3.12. The exact distributions given by the theorems are represented by the solid red lines, while the histograms are the results of Mon te Carlo tests with 10 , 000 trials. This time, one has to increase β to 50 to obtain significan t results; ho wev er, since the distributions represent asymptotics, that is still rather small. Figure 5 co vers the case β → ∞ , that is, Theorem 3.11, and the histogram is of the distribution on ( b + m ) λ min , while Figure 6 co vers the case m, β → ∞ (Theorem 3.12), and the histogram reflects the empirical distribution of m ( b + m ) λ min . 0 0.5 1 1.5 2 2.5 3 0 0.5 1 1.5 Monte Carlo vs asymptotical distributions for (b+m) λ min with m=6, b = 50 (resp., ∞ ), β = 1.5, k = 2, and a = 2k/ β − 1 Figure 5: The solid red line represents the asymptotical ( b = ∞ ) distribution, while the normalized histogram represents the results of a Monte Carlo exp erimen t for b = 50, with 10 , 000 trials. Ac kno wledgemen ts Ioana would like to thank James Demmel, Peter F orrester, Tiefeng Jiang, and Eric Rains for useful discussions. Ioana’s work is supp orted by NSF CAREER Aw ard DMS-0847661. 13 0 1 2 3 4 5 6 7 0 0.1 0.2 0.3 0.4 0.5 0.6 Monte Carlo vs asymptotical distributions for m(b+m) λ min with m = 15 (resp., ∞ ), b = 50 (resp., ∞ ), β = 3.33, k = 2, and a = 2k/ β − 1 Figure 6: The solid red line represents the asymptotical ( m, b = ∞ ) distribution, while the normalized histogram represents the results of a Monte Carlo exp erimen t for m = 15, b = 50, with 10 , 000 trials. References [1] A.B.J. Kuijlaars and M. V anlessen. Universalit y for eigenv alue correlations from the mo dified Jacobi unitary ensemble. IMRN , 30:1575–1600, 2002. [2] M. Abramowitz and I.A. Stegun, editors. Handb o ok of Mathematic al F unctions . Dov er Publications, New Y ork, 1970. [3] Ric hard Askey and Donald Ric hards. Selb erg’s second beta integral and an in tegral of Meh ta. In T.W. Anderson et al., editor, Pr ob ability, Statistics, and Mathematics: Pap ers in Honor of Samuel Karlin , pages 27–39. Academic Press, New Y ork, 1989. [4] Grey Ballard, James Demmel, and Ioana Dumitriu. Minimizing Communication for Eigenproblems and the Singular V alue Decomp osition. 2010. Preprint. [5] Cy Chan, V esselin Drensky , Alan Edelman, Raymond Kan, and Plamen Ko ev. On Computing Sc hur F unctions and Series Thereof. 2010. Preprin t. [6] Beno ˆ ıt Collins. Int´ egr ales Matriciel les et Pr ob abilit ´ es non-c ommutatives . PhD thesis, Universit ´ e P aris 6, 2003. [7] Beno ˆ ıt Collins. Product of random pro jections, Jacobi ensem bles and universalit y problems arising from free probability. PTRF , 133:315–344, 2005. [8] James Demmel, Ioana Dumitriu, and Olga Holtz. F ast linear algebra is stable. Numer. Math. , 108:59–91, 2007. [9] Ioana Dumitriu, Alan Edelman, and Gene Shuman. MOPS: Multiv ariate Orthogonal Polynomials (sym b olically). J. Symb olic Computation , 42:587–620, 2007. [10] Ioana Dumitriu and Plamen Koev. Distributions of the extreme eigenv alues of Beta-Jacobi random matrices. SIAM J. Matrix Anal. Appl. , 30:1–6, 2008. [11] Alan Edelman. Eigenvalues and Condition Numb ers of R andom Matric es . PhD thesis, Mas- sac husetts Institute of T ec hnology , 1989. [12] Alan Edelman and Brian Sutton. F rom Random Matrices to Sto chastic Op erators. J. Stat. Phys. , 127(6):1121–1165, 2007. [13] P eter J. F orrester. Quantum conductance problems and the Jacobi ensem ble. J. Phys. A , 39(22), 2006. 10.1088/0305-4470/39/22/004. 14 [14] P eter J. F orrester. L o g-Gases and R andom Matric es . London Mathematical So ciet y Monographs. Princeton Universit y Press, Princeton, NJ, 2010. [15] T. Jiang. Approximation of Haar Distributed Matrices and Limiting Distributions of Eigenv alues of Jacobi Ensembles. PTRF , 144:1–38, 2009. [16] T. Jiang. The En tries of Haar-inv arian t Matrices from the Classical Compact Groups. J. The or etic al Pr ob ability , 2009. DOI 10.1007/s10959-009-0241-7. [17] T. Jiang. Limit Theorems for Beta-Jacobi Ensem bles. 2010. [18] Iain Johnstone. Multiv ariate analysis and Jacobi ensembles: Largest eigen v alue, T racyWidom limits and rates of conv ergence. A nn. Stat. , 36(6):2638–2716, 2008. [19] Jos ´ e A. Ram ´ ırez and Brian Rider. Diffusion at the Random Matrix Hard Edge. Commun. Math. Phys. , 288:887–906, 2009. [20] K. Kadell. The Selberg-Jack p olynomials. A dvanc es in Mathematics , 130:33–102, 1997. [21] Joic hi Kaneko. Selb erg in tegrals and hypergeometric functions associated with Jack polynomials. SIAM J. Math. A nal. , 24:1086–1110, 1993. [22] Ro wan Killip and Irina Nenciu. Matrix Mo dels for Circular Ensembles. International Mathematics R ese ar ch Notic es , 50:2665–2701, 2004. [23] Plamen Ko ev and Alan Edelman. The efficient ev aluation of the hypergeometric function of a matrix argument. Math. Comp, , 75:833–846, 2006. [24] R. Lippert. A matrix mo del for the β -jacobi ensem ble. J. Math. Phys. , 44 (10):4807–4816, 2003. [25] Robb J. Muirhead. Asp e cts of Multivariate Statistic al The ory . John Wiley & Sons, New Y ork, 1982. [26] J. Ramirez, B. Rider, and B. Vir´ ag. Beta ensembles, stochastic Airy spectrum, and a diffusion. T echnical rep ort, 2009. arXiv:math-ph/0607331v4. [27] A. Selb erg. Bemerkninger om et multiplet integral. Norsk Matematisk Tiddskrift , 26:71–78, 1944. [28] Brian Sutton. The Sto chastic Op er ator Appr o ach to R andom Matrix The ory . PhD thesis, Mas- sac husetts Insitute of T ec hnology , 2003. [29] Craig A. T racy and Harold Widom. The distribution of the largest eigen v alue in the Gaussian ensem bles. In Calo ger o-Moser-Sutherland Mo dels, CRM Series in Mathematic al Physics , volume 4, pages 461–472. Springer-V erlag, 2000. 15