On the distribution of maximum value of the characteristic polynomial of GUE random matrices
aa r X i v : . [ m a t h - ph ] J un ON THE DISTRIBUTION OF THE MAXIMUM VALUE OF THECHARACTERISTIC POLYNOMIAL OF GUE RANDOM MATRICES
Y. V. FYODOROV AND N. J. SIMM
Abstract.
Motivated by recently discovered relations between logarithmically correlatedGaussian processes and characteristic polynomials of large random N × N matrices H fromthe Gaussian Unitary Ensemble (GUE), we consider the problem of characterising the distri-bution of the global maximum of D N ( x ) := − log | det( xI − H ) | as N → ∞ and x ∈ ( − , freezing transition scenario for logarithmically correlatedprocesses. Although the general idea behind the method is the same as for the earlier con-sidered case of the Circular Unitary Ensemble, the present GUE case poses new challenges.In particular we show how the conjectured self-duality in the freezing scenario plays thecrucial role in our selection of the form of the maximum distribution. Finally, we demon-strate a good agreement of the found probability density with the results of direct numericalsimulations of the maxima of D N ( x ). Introduction.
The space of all N × N Hermitian matrices H with probability density function(1.1) P ( H ) ∝ exp( − N Tr( H ))is known as the Gaussian Unitary Ensemble (or GUE)[1, 38, 43]. Here and henceforth thevariance is chosen to ensure that asymptotically for N → ∞ , the limiting mean density ofthe GUE eigenvalues is given by the Wigner semicircle law ρ ( x ) = (2 /π ) √ − x supportedin the interval x ∈ [ − , p N ( x ) = det( xI − H ) of the matrix H constitutes one of the most basic quantities of interest, encoding all eigenvalues of H through the roots of p N ( x ). As one varies the argument x over an interval containing manyeigenvalues for a given realization of the ensemble, the value of the polynomial p N ( x ) showshuge variations by the orders of magnitude for large N , see Figure 1.1 for N = 50 and Figure1.2 for N = 3000. x Figure 1.1.
A plot of a single realization of | p N ( x ) | e − E log | p N ( x ) | for N = 50.The global maximum is marked with a red circle. −1 −0.8 −0.6 −0.4 −0.2 0 0.2 0.4 0.6 0.8 1−30−20−100102030 x Figure 1.2.
A plot of a single realization of 2 log (cid:0) | p N ( x ) | e − E log | p N ( x ) | (cid:1) with N = 3000. The maximum value is marked with a red circle. he purpose of this article is to describe the statistical properties of the highest peak displayed by the modulus of the GUE polynomial | p N ( x ) | , namely the probability densityfor the maximum value attained by | p N ( x ) | over the interval [ − ,
1] on the real line as N → ∞ . Our main result is the following Conjecture 1.1.
Consider the random variable (1.2) M ∗ N := max x ∈ [ − , (cid:26) | p N ( x ) | − E (log | p N ( x ) | ) (cid:27) Then in the limit N → ∞ we have (1.3) M ∗ N = 2 log( N ) −
32 log(log( N )) − (1 + o (1)) y + o (1) where y is a continuous random variable characterized by the two-sided Laplace transform ofits probability density: (1.4) E ( e ys ) = 1 C (2 π ) s Γ( s + 1)Γ( s + 3) G ( s + 7 / G ( s + 6) G ( s + 1) where Γ( z ) and G ( z ) stand for the Euler gamma-function and the Barnes digamma-function,correspondingly. The normalization C can be evaluated explicitly as (1.5) C = e / π / / A where A is the Glaisher-Kinkelin constant A = e / − ζ ′ ( − = 1 . ... . Remark 1.2.
The product form of the Laplace transform (1.4) offers an interesting inter-pretation of the above results. Noting that Γ(1 + s ) is the moment generating function of astandard Gumbel random variable G , we can write(1.6) y = G + y ′ where y ′ is an independent random variable with two-sided Laplace transform(1.7) E ( e y ′ s ) = 1 C (2 π ) s Γ( s + 3) G ( s + 7 / G ( s + 6) G ( s + 1) . In the end of the paper we provide convincing numerical evidence that this Laplace transformdoes indeed define a unique random variable y ′ . This immediately implies that the probabilitydensity of y is the convolution of a Gumbel random variable with y ′ . Such a convolutionstructure is expected to appear universally when studying the extreme value statistics oflogarithmically correlated Gaussian fields, see the discussion around and after Eq. (1.11).In recent years, much interest has accumulated regarding the statistical behaviour ofcharacteristic polynomials of various random matrices as a function of the spectral variable x . To a large extent this interest was stimulated by the established paradigm that manystatistical properties of the Riemann zeta function along the critical line, that is ζ (1 / it ),can be understood by comparison with analogous properties of the characteristic polynomialsof random matrices [32, 30, 12, 14, 29, 2].For invariant ensembles [38, 43] of self-adjoint matrices with real eigenvalues, statisticalcharacteristics of p N ( x ) depend very essentially on the choice of scale spanned by the realvariable x . From that end it is conventional to say that x spans the local (or microscopic ) cale if one considers intervals containing in the limit N → ∞ typically only a finite numberof eigenvalues (the corresponding scale for GUE in (1.1) is of the order of 1 /N ). At suchscales, standard objects of interest are correlation functions containing products and ratios ofcharacteristic polynomials, which show determinantal/Pfaffian structures [6, 28, 45, 5, 8, 33]for Hermitian/real symmetric matrices and tend to universal limits at the local scale. Similarstructures arise for properly defined characteristic polynomials p N ( θ ) = det (cid:0) I − U e − iθ (cid:1) ofcircular ensembles (like CUE, COE, and CSE)[38] of unitary random matrices U uniformlydistributed with respect to the Haar measure on U ( N ) (and other classical groups) [12, 9, 13],whose properties on the local scale are indistinguishable from their Hermitian counterparts.Next, when x spans an interval containing in the limit N → ∞ typically of order of N eigenvalues one speaks of the global (or macroscopic ) scale behaviour. At such a scaleproperties of p N ( x ) display both universal and non-universal features, the latter depending onthe ensemble chosen. The study of characteristic polynomials at such a scale was initiated in[30] where it was shown that the function V N ( θ ) = − | det(1 − U e − iθ ) | , with U belongingto the CUE, converges (in an appropriate sense) to a random Gaussian Fourier series of theform(1.8) V ( θ ) = ∞ X n =1 √ n (cid:0) v n e inθ + v n e − inθ (cid:1) , where the coefficients v n , v n are independent standard complex Gaussian random variables,i.e. E { v n } = 0, E { v n } = 0 and E { v n v n } = 1. The covariance structure associated withsuch a process is given by E { V ( θ ) V ( θ ) } = − | e iθ − e iθ | as long as θ = θ . Such a(generalized) random function V ( θ ) is a representative of random processes known in theliterature under the name of , see [27, 23] for background discussion and furtherreferences.Recently the study of the global scale behaviour was extended to the GUE polynomial p N ( x ) in [24] by using earlier insights from [31] and [35]. That work revealed again a structureanalogous to that of (1.8), though different in detail. Namely, it was shown that the naturallimit of ˜ D N ( x ) := − log | p N ( x ) | + E { log | p N ( x ) |} is given by the random Chebyshev-Fourierseries(1.9) F ( x ) = ∞ X n =1 √ n a n T n ( x ) , x ∈ ( − , , with T n ( x ) = cos( n arccos( x )) being Chebyshev polynomials and real a n being independentstandard Gaussians. A quick computation shows that the covariance structure associatedwith the generalized process F ( x ) is given by an integral operator with kernel(1.10) E { F ( x ) F ( y ) } = ∞ X n =1 n T n ( x ) T n ( y ) = −
12 log(2 | x − y | ) , as long as x = y . Such a limiting process F ( x ) is an example of an aperiodic 1 /f -noise.Finally, one can consider an intermediate, or mesoscopic spectral scales, with intervalstypically containing in the limit N → ∞ the number of eigenvalues growing with N , butrepresenting still a vanishingly small fraction of the total number N of all eigenvalues. Theproperties of the characteristic polynomials at such scales were again addressed in [24] whereit was shown that for the GUE, that object gives rise to a particular (singular) instance of he so-called fractional Brownian motion (fBm) [37, 16] with the Hurst index H = 0, againcharacterized by correlations logarithmic in the spectral parameter.The discussion above serves, in particular, the purpose of pointing to an intimate con-nection between Gaussian random processes with logarithmic correlations and the modulusof characteristic polynomials at global and mesoscopic scales. The relation is importantas logarithmically correlated Gaussian (LCG) random processes and fields attract growingattention in Mathematical Physics and Probability and play an important role in problemsof Quantum Gravity, Turbulence, and Financial Mathematics, see e.g. [18]. In particular,the periodic 1 /f noise (1.8) emerged in constructions of conformally invariant planar ran-dom curves [4]. Among other things, the statistics of the global maximum of LCG fieldsattracted considerable attention, see [15] and references therein. Particularly relevant in thepresent context are the results of Ding, Roy and Zeitouni [15] on the maxima of regular-ized lattice versions of LCG fields which we discuss informally below. Let V N = Z dN be the d − dimensional box of side length N with the left bottom corner located at the origin. Asuitably normalized version of the logarithmically correlated Gaussian field is a collectionof Gaussian variables φ N,v : v ∈ V N with variance E { φ N,v } = 2 log N + f ( v ) and covariancestructure(1.11) E { φ N,v , φ
N,u } = 2 log + N | u − v | + g ( u, v ) , for u = v ∈ V N where ln + ( w ) = max (ln w,
0) and both f ( v ) and g ( u, v ) are continuous bounded functionsfar enough from the boundary of V N . Now set M N = max v ∈ V N φ N,v and m N = √ d log N − d log log N . The limiting law of M N − m N is then expected, after an appropriate shift andrescaling, to be given by the Gumbel distribution with random shift :(1.12) P ( y ) = lim N →∞ P rob ( M N ≥ m N − y ) = E n e − e √ d ( y − z ) o , where the distribution of the random shift variable z depends on details of the behaviour ofcovariance (1.11) for | u − v | ∼ N and | u − v | ∼
1, see the detailed discussion in [15]. Therandom variable Z = e −√ dz is related to the so called derivative martingale associated withthe LCG fields [15] whose distribution is however not known. Recently it has been shownthat the recentering term m N in (1.12) also holds for a randomized model of the Riemannzeta function [2], proved by revealing a special branching structure within the associatedlogarithmic correlations.We see that our conjecture 1.1 for the maximum of characteristic polynomial of large GUEmatrices fully agrees with the predicted structure of the maximum of LCG in dimension d = 1. Note that the expression (1.12) implies that the double-sided Laplace transform ofthe density ρ ( y ) = − ddy P ( y ) for the (shifted) maximum y is related to the density ˜ ρ ( z ) ofthe random variable z as(1.13) E ( e ys ) = Z ρ ( y ) e sy dy = Γ( s + 1) Z ˜ ρ ( z ) e sz dz = Γ( s + 1) E ( e zs )which is in turn equivalent to the Gumbel convolution in eq.(1.6). In fact our formula (1.7)provides the explicit form of the distribution for the derivative martingale of our model, thusgoing considerably beyond the considerations of [15]. rom a quite different perspective, processes similar to (1.8) and (1.9) appeared in thecontext of statistical mechanics of disordered systems when studying extreme values of ran-dom multifractal landscapes supporting spinglass-like thermodynamics [21, 25, 27, 3]. Thelatter link is especially important in the context of the present paper. The idea that it isbeneficial to look at | p N ( θ ) | as a disordered landscape consisting of many peaks and dips,and to think of an associated statistical mechanics problem was put forward in [22, 23]. Itallowed to get quite non-trivial analytical insights into statistics of the maximal value of theCUE polynomial sampled over the full circle θ ∈ [0 , π ], or over its mesoscopic sub-intervals.This was further used to conjecture the associated properties of the modulus of the Riemannzeta-function along the critical line, see some recent advances inpired by that line of researchin [2]. Some relations between between CUE characteristic polynomials and logarithmicallycorrelated processes (in the form of the so-called ”multiplicative chaos” measures introducedby Kahane, see [44] for a review) was recently rigorously verified in [48]. The case of GUEpolynomials however remained outstanding.It is our objective in this paper to provide two separate means of supporting Conjecture1.1. First, we will provide careful and explicit, albeit in part heuristic, analytical arguments.Although our technique is inspired by the approach of [23] it contains new nontrivial fea-tures necessary to overcome challenges arising from the non-uniform eigenvalue density ρ ( x ),reflecting absence of translational invariance for the GUE at the global spectral scale (notee.g. the non-trivial recentering in (1.2)). All this makes actual calculation for the GUE muchmore involved in comparison to the CUE and the limiting random variable u above appearsto be more complicated than its CUE counterpart. Secondly, we will test our Conjecturewith numerical experiments for matrices of size N = 3000 and around 250 ,
000 realizations.This is especially important as part of our analysis is based on very plausible but as yet notfully rigorous considerations. Finally, is natural to expect that the same distribution shouldbe shared by the maximum modulus of characteristic polynomials for Hermitian randommatrices with independent entries taken from the so-called Wigner ensembles, see [19].Before giving the detail of our procedure in the next section we need to quote the followingfundamental asymptotic result obtained by Krasovsky [35] which will be central for ourconsiderations: E k Y j =1 | p N ( x j ) | α j ! = k Y j =1 C ( α j )(1 − x j ) α j / ( N/ α j e (2 x j − − α j N (1.14) × Y ≤ i We are grateful for helpful comments from Christian Webb duringthe preparation of this manuscript. We acknowledge support from EPSRC grant EP/J002763/1“ Insights into Disordered Landscapes via Random Matrix Theory and Statistical Mechanics ”.2. Statistical mechanics approach to the distribution of GUEcharacteristic polynomials. Following the ideas of [23] we recast the problem of computing the value of the globalmaximum of | p N ( x ) | (with an appropriate shift by the mean value) as a statistical mechanicsproblem characterized by the partition function(2.1) Z N ( β ) = N Z − e − βφ N ( x ) ρ ( x ) q dx, β > , q ≥ φ N ( x ) = − | p N ( x ) | − E log | p N ( x ) | )), inverse temperature β > β -independent non-negative parameter q . Specifically, if we define the associated ”freeenergy” as F ( β ) = − β − log Z N ( β ), then(2.2) lim β →∞ F ( β ) = min x ∈ ( − , φ N ( x ) = 2 max x ∈ ( − , [log | p N ( x ) | − E log | p N ( x ) | ] . Note that if compared to a similar partition function for the CUE case the main newfeature in (2.1) is the factor ρ ( x ) q . Although naively the presence of such a factor may seemirrelevant when taking the limit β → ∞ , we will actually see that it plays a very importantrole in supporting our procedure of extracting the free energy for β exceeding some criticalvalue.Now we aim to compute the integer moments of the partition function defined in (2.1):(2.3) E ( Z kN ( β )) = (cid:18) N (cid:19) k Z − . . . Z − E k Y j =1 | p N ( x j ) | β ! k Y j =1 e − β E log | p N ( x j ) | ρ q ( x j ) dx j In the limit N → ∞ the leading asymptotics of the above integral can be extracted byreplacing the factor E (cid:16)Q kj =1 | p N ( x j ) | β (cid:17) with its asymptotics from (1.15). In this way one btains E ( Z kN ( β )) ∼ (cid:18) N (cid:19) β C ( β )(2 /π ) q ! k Z [ − , k k Y j =1 (1 − x j ) β q Y ≤ i One starts with finding a recursionsatisfied by ˜ S k ( a, b, γ ) for integer k which is suitable for the continuation. By writing(2.10) k − Y j =1 Γ( a + b + 2 − ( k + j − γ ) = Γ(2 + a + b − ( k − γ )Γ(2 + a + b − (2 k − γ ) Q kj =1 Γ(2 + a + b − ( k + j − γ )Γ(2 + a + b − (2 k − γ )one sees immediately that(2.11) E ( z kβ ) E ( z k − β ) = Γ( a + 1 − ( k − γ )Γ( b + 1 − ( k − γ )Γ(1 − kγ )Γ(2 + a + b − ( k − γ )Γ(2 + a + b − (2 k − γ )Γ(2 + a + b − (2 k − γ )It is convenient to introduce the moments M β ( s ) of the random variable z β defined for anycomplex s as M β ( s ) = E ( z − sβ ) . We then have E ( z kβ ) = M (1 − k ) , E ( z k − β ) = M (2 − k ) andafter identifying s = 1 − k the recursion (2.11) takes the form(2.12) M β ( s ) M β ( s + 1) = Γ(1 + a + γs )Γ(1 + b + γs )Γ(1 + ( s − γ )Γ(2 + a + b + ( s + 1) γ )Γ(2 + a + b + 2 sγ )Γ(2 + a + b + (2 s + 1) γ )which is now assumed to be valid for any complex s . It is convenient to further use theduplication formula for the Gamma function :(2.13) Γ(2 z ) = 2 z − √ π Γ( z )Γ( z + 1 / s in the denominator. Indeed, we haveΓ(2 + a + b + 2 γs ) = 2 a + b + γ s Γ(1 + γs + ( a + b ) / γs + ( a + b + 3) / / √ π Γ(2 + a + b + γ (2 s + 1))= 2 a + b + γ (2 s +1) Γ(1 + γ ( s + 1 / 2) + ( a + b ) / γ ( s + 1 / 2) + ( a + b + 3) / / √ π so that (2.12) assumes the form M β ( s ) M β ( s + 1) = Γ(1 + a + γs )Γ(1 + b + γs )Γ(1 + γ ( s − a + b + ( s + 1) γ )Γ(1 + γ ( s + 1 / 2) + ( a + b ) / γ ( s + 1 / 2) + ( a + b + 3) / × π a + b )+(4 s +1) γ a + b ) / γs )Γ((3 + a + b ) / γs )Recalling that according to (2.9) in our particular case a = b = q + β we now use theparameterisation a = a + a β , b = b + b β and β -independent constants a , a , b .b . fter this we finally arrive at M β ( s ) M β ( s + 1) =(2.15)Γ(1 + a + β ( s + a ))Γ(1 + b + β ( s + b ))Γ(1 + β ( s − a + b + ( s + 1 + a + b ) β )Γ(1 + β ( s + 1 / a / b / 2) + ( a + b ) / β ( s + (1 + a + b ) / 2) + ( a + b + 3) / × π − a + b ) − (4 s +1+2 a +2 b ) β Γ(1 + ( a + b ) / β ( s + ( a + b ) / a + b ) / β ( s + ( a + b ) / M β ( s ) which satisfies (2.16) for any complex s we follow [25] andintroduce a variant of the Barnes function G β ( x ) which for any ℜ ( x ) > G β ( x ) = x − Q/ 22 ln(2 π ) + Z ∞ dtt e − Q t − e − xt (1 − e − βt )(1 − e − t/β ) + e − t Q/ − x ) + Q/ − xt ! where Q = β + 1 /β . This function satisfies the so-called self-duality relation G β ( x ) = G /β ( x )(2.17)and further posesses a shift property that is central for our studies G β ( x + β ) = β / − βx (2 π ) β − Γ( βx ) G β ( x )(2.18)One can check that G β ( x ) for β = 1 coincides with the standard Barnes function G ( x ) whichis a unique solution of the recursion G ( x + 1) = Γ( x ) G ( x ) satisfying G (1) = 1. Similarly tothe standard Barnes function the general Barnes G β ( x ) has no poles and only zeroes locatedat x = − nβ − m/β , n, m = 0 , , .. . A detailed discussion of properties of functions closelyrelated to G β ( x ) can be found in [40, 41].Let us now define a function M ( G ) β ( s ) of the complex argument s by M ( G ) β ( s ) = π − s B s + B s β β s × G β (cid:16) β (cid:0) s + a + b +12 (cid:1) + a + b β (cid:17) G β (cid:16) β (cid:0) s + a + b (cid:1) + a + b β (cid:17) G β (cid:16) β ( s + a ) + a β (cid:17) G β (cid:16) β ( s + b ) + b β (cid:17) (2.19) × G β (cid:16) β (cid:0) s + a + b +12 (cid:1) + a + b β (cid:17) G β (cid:16) β (cid:0) s + a + b (cid:1) + a + b β (cid:17) G β (cid:16) β ( s + 1 + a + b ) + a + b β (cid:17) G β (cid:16) β ( s − 1) + β (cid:17) where B = 2 β and B = 2( a + b + 1) + β (2 a + 2 b − G β ( β ( s + 1) + c/β ) G β ( βs + c/β ) = (2 π ) β − β / − β s − c Γ( c + β s ) hows that the ratio M ( G ) β ( s ) M ( G ) β ( s +1) reproduces the right-hand side of (2.16) from which we conclude(2.21) M ( G ) β ( s ) M ( G ) β ( s + 1) = M β ( s ) M β ( s + 1)which finally implies that(2.22) M β ( s ) = M Gβ ( s ) M β (1) M ( G ) β (1)where M β (1) ≡ 1. Together with (2.4), (2.5) and the fact that M β (1) = 1, we obtain for β < E ( Z N ( β ) − s ) ∼ (cid:0) N (cid:1) β C ( β )(2 /π ) q Γ(1 − β ) ! − s (1 − s )( β + q +1)+2 β s (1 − s ) M ( G ) β ( s ) M ( G ) β (1)2.2. Duality and the freezing transition. The pair (2.19)-(2.22) solves the problem offinding the complex moments M β ( s ) = E ( z − sβ ) of the random variable z β for any complex s ,and β < 1. Knowledge of such moments can be used to restore the probability distributionof z β , hence of the partition function Z N ( β ), and of its logarithm (the free energy) for large N ≫ 1. Our goal is however to study the limit of the latter as β → ∞ and one thereforeshould have a way of extracting information on the distribution for β > 1. In doing thiswe rely on the freezing transition scenario for logarithmically correlated random landscapes.The background idea of such scenario goes back to [10] and was further advanced and clarifiedin the series of works [21, 25, 26, 27]. In brief, this scenario predicts a phase transition atthe critical value β = 1 and amounts to the following principle: Thermodynamic quantities which for β < are self-dual functions of the inversetemperature β , i.e. functions that remain invariant under the transformation β → β − ,retain for all β > the value they acquired at the point of self-duality β = 1 . Although such a scenario is not yet proven mathematically in full generality and has thestatus of a conjecture supported by physical arguments and available numerics, recently afew nontrivial aspects of freezing were verified within rigorous probabilistic analysis, see e.g.[3, 15, 46] for efforts in this direction.Within that scenario, one of the main outcomes of the analysis performed in [25] is that theself-dual object associated with the distribution of the partition function for logarithmicallycorrelated landscapes is expected to be the appropriately defined Laplace transform: g β ( y ) = E (cid:0) exp (cid:2) − e βy Z N ( β ) / Z eN ( β ) (cid:3)(cid:1) , (2.24)where Z eN ( β ) is a typical scale of the partition function which is extracted from the asymp-totic for the integer moments and in our case can be chosen as(2.25) Z ( e ) N ( β ) = N β [ G ( β + 1)] G (2 β + 1)Γ(1 − β ) (cid:18) π (cid:19) q . Moreover, defining the probability density p β ( y ) by p β ( y ) = − g ′ β ( y ) one can show thatthe double-sided Laplace transform for such a probability density is related to the complex oments ˜ M β ( s ) = E (cid:18) Z N ( β ) Z ( e ) N ( β ) (cid:19) − s of the scaled partition function via the following relation(see eq.(26) of [25]) ln Z ∞−∞ p β ( y ) e ys dy = ln ˜ M β (1 + sβ ) + ln Γ(1 + sβ )(2.26)Actually, as shown in [25] the freezing scenario implies that the variable y whose probabilitydensity is given by p β =1 ( y ) is precisely the fluctuating part of the height of the global min-imum of the random potential which is our main object of interest. Note however that thescale Z ( e ) N ( β ) diverges when approaching the critical point β = 1, and that the associatedfree energy − β log Z ( e ) N ( β ) is self-dual only in the leading order, given by − ( β + β − ) log N .The latter term after freezing at β = 1 yields the leading 2 log N term in our conjectureEq.(1.3) for the maximum, whereas the logarithmically divergent term − β log Γ(1 − β ) af-ter careful re-interpretation results in the second term − log log N , see [27] for the detailedexplanation of that mechanism. The procedure leaves however a certain arbitrariness in theterms of the order of unity in the mean free energy, hence in the overall shift of the positionof the maximum. Let us stress however that apart from such a shift, the shape of the distri-bution function recovered in the framework of the freezing paradigm is completely fixed bythe procedure.Our strategy therefore will be to check if self-duality holds for the right-hand side combi-nation in (2.26) when we substitute our expression for the moments. Before we proceed, itwill be helpful to further expand our expression (2.23). Inserting (2.19) and making use ofthe identity(2.27) 1 G β ( β ( s − 1) + 1 /β ) = Γ(1 + β ( s − G β ( βs + 1 /β ) (2 π ) ( β − / β − / − β ( s − shows that (taking into account all prefactors coming from (2.4), (2.5) and (2.19)) E ( Z N ( β ) − s ) ∼ [ Z ( e ) N ( β )] − s β + B s π − s β β s (2 π ) ( β − / β − / − β ( s − M ( G ) β (1) × Γ(1 + β ( s − G β (cid:16) β (cid:0) s + a + b +12 (cid:1) + a + b β (cid:17) G β (cid:16) β (cid:0) s + a + b (cid:1) + a + b β (cid:17) G β (cid:16) β ( s + a ) + a β (cid:17) G β (cid:16) β ( s + b ) + b β (cid:17) × G β (cid:16) β (cid:0) s + a + b +12 (cid:1) + a + b β (cid:17) G β (cid:16) β (cid:0) s + a + b (cid:1) + a + b β (cid:17) G β (cid:16) β ( s + 1 + a + b ) + a + b β (cid:17) G β (cid:16) βs + β (cid:17) (2.28)A direct inspection makes it clear that the self-duality is only possible if either a = a , b = b or a = b , b = a . For the GUE characteristic polynomials, we have a = b = 1 / a = b = q/ only if q = 1. We therefore have to choose q = 1 to beable to rely upon the freezing scenario allowing to interpret the function p β =1 ( y ) calculatedfrom its Laplace transform via (2.26) as the probability density for the (shifted) global inimum. Using (2.28) with a = a = b = b = we get E ( Z − sN ( β )) ∼ [ Z ( e ) N ( β )] − s (1+ β )( s − × Γ(1 + β ( s − G β (cid:16) β ( s + 1) + β (cid:17) G β (cid:16) β ( s + 1 / 2) + β (cid:17) G β (cid:16) β ( s + 1) + β (cid:17) G β (cid:16) β (cid:0) s + (cid:1) + β (cid:17) G β (cid:16) β ( s + 2) + β (cid:17) G β (cid:16) βs + β (cid:17) c β (2.29)where c β is a constant determined by the condition E ( Z N ( β ) − s ) | s =1 = 1. Inserting (2.19)into the right-hand side of (2.26) (which is now manifestly self-dual) leads to the followingexpression at β = 1: Z ∞−∞ p β =1 ( y ) e ys dy = K s Γ(1 + s ) ˜ M β =1 (1 + s )(2.30)= 1 C K s Γ (1 + s ) G ( s + 7 / G ( s + 3) G ( s + 4) G ( s + 3) G ( s + 6) G ( s + 2) = 1 C K s Γ(1 + s ) G ( s + 7 / Γ( s + 3) G ( s + 1) G ( s + 6) . (2.31)where C = c β =1 and K is a constant which determines the shift in the maximum as discussedbelow (2.26). The value K = 2 π in (1.1) is conjectured from the results of numericalsimulations in the next section. The latter formula (2.31) constitutes our main analyticalresult and finally leads to our Conjecture 1.1.3. Numerical study of the distribution of the maximum modulus of GUEcharacteristic polynomials The purpose of this Section is to provide a numerical test of Conjecture 1.1.3.1. Results. In Figure 3.1 we present a histogram of the recentered and rescaled maximumof the GUE characteristic polynomial, defined by(3.1) y ∗ N := (2 log( N ) − (3 / 2) log(log( N )) − M ∗ N + c ∗ N )(1 + s ∗ N )with M ∗ N defined in (1.2). Here we used the matrix size N = 3000 and 250 , 000 realizationsof the GUE ensemble. The dashed red line is the exact probability density of the randomvariable y defined via its Laplace transform in (1.4). In (3.1) we have recentered and scaledby c ∗ N = 0 . 216 and s ∗ N = 0 . N effects due to the o (1)terms in (1.3). Note that the influence of shift/recentering is already quite small comparedwith the predicted considerably larger (3 / 2) log(log( N )) ∼ . 12 shift. The parameters c ∗ N and s ∗ N were calculated empirically from the mean and variance of y in (1.4) according tothe formula s ∗ N = q Var( y ) / Var( M ∗ N ) − c ∗ N = E ( y ) / ( s ∗ N + 1) − (2 log( N ) − (3 / 2) log(log( N )) − M ∗ N ) , (3.2)as derived by requiring E ( y ∗ N ) = E ( y ) and Var( y ∗ N ) = Var( y ). In Table 1 we display values ofthe parameters c ∗ N and s ∗ N for the studied range of sizes N , as determined empirically fromthe mean and variance of the random variable u . The observed decay with N is certainly 12 −10 −8 −6 −4 −2 0 2 4 600.050.10.150.20.25 Figure 3.1. The centered and scaled maximum as defined by (3.1). Thedashed line is the probability density of the random variable y given in Laplacespace by (1.4). c ∗ N N Figure 3.2. Each triangle represents a value of c ∗ N obtained from (3.2) with250 , 000 realizations. able 1. Finite- N corrections for increasing values of N all with 250 , 000 realizations N c ∗ N s ∗ N 150 0.329 0.331600 0.267 0.2481050 0.244 0.2241500 0.234 0.2121950 0.228 0.2022400 0.221 0.1953000 0.216 0.188 −8 −6 −4 −2 0 2 400.050.10.150.20.250.30.35 Figure 3.3. The inverse Laplace transform of formula (1.7).consistent with asymptotic validity of our Conjecture 1.1, though the convergence to theasymptotic results is too slow to make more definite claims. To resolve further decreaseof the coefficients c ∗ N and s ∗ N would require much larger matrices and is computationallydemanding.Finally, we provide a numerical validation of the decomposition (1.6). In Figure 3.3 weplot the inverse Laplace transform of (1.7) obtained by a direct numerical evaluation of theintegral in the Bromwich inversion formula for the Laplace transform. The positive andnormalized curve clearly corresponds to a bona fide probability density of some real randomvariable y ′ .3.2. Numerical method. The numerical evaluation of the maximum value (1.2) may beconsidered quite a non-trivial problem in its own right, for at least two reasons. Firstly,the characteristic polynomial p N ( x ) having zeros as the eigenvalues of H , displays O ( N )oscillations in the spectral interval [ − , 1] with hugely varying peaks heights. This produces onsiderable clusterings of ‘near-maxima’ which may confuse any naive attempt to find thetrue maximum value. Secondly, the slow changing nature of the correction terms in Conjec-ture 1.1, of order log( N ) and log log( N )) respectively, require one to go to somewhat largematrices to resolve reasonable asymptotic behaviour. The problem is further compoundedby the numerical instability of calculating determinants of such matrices.Our solution to these problems heavily relies on a sparse realization of GUE matrices H originally due to Trotter [47] (see also Dumitriu and Edelman [17]). He discovered thatthe eigenvalues of GUE matrices H have the same joint probability density as those of thefollowing real symmetric tri-diagonal matrix:(3.3) H = 12 √ N N (0 , χ χ N (0 , χ . . . . . . . . . χ N − N (0 , χ N − χ N − N (0 , where N (0 , 2) is a normal random variable with mean 0 and variance 2. The sub-diagonalis composed of random variables χ n having the same density as p χ n where χ n is a χ -square random variable with 2 n degrees of freedom. To compute the maximum value of p N ( x ) = det( xI − H ) = det( xI −H ), we begin by exploiting the known asymptotic behaviour(3.4) 2 E log | p N ( x ) | = N (2 x − − o (1)so that f N ( x ) := 2 log | p N ( x ) | − E log | p N ( x ) | ∼ | det( e − ( x − / − log(2)) ( xI − H )) | (3.5)Further progress is now possible thanks to the fact that determinants of tri-diagonal matricessatisfy a linear recurrence relation. Furthermore, by an appropriate rescaling, the recursioncomputes determinants of all leading principal minors simultaneously, thus computing f j ( x ) for all j = 1 , . . . , N in linear time.Now to find the maximum, we define a mesh M = {− n/ ∆ : n = 0 , . . . , } with∆ ∼ N and evaluate f N ( x ) at each of the points in M . At those points where f N ( x ) ismaximal the Matlab function ‘fminbnd’ is invoked to converge onto the global maximum.Figure 1.2 illustrates the complexity of the problem. Our algorithm is sufficiently precise todistinguish the true maximum (located at x ≈ − . e.g. x ≈ − . References [1] G. W. Anderson, A. Guionnet, and O. Zeitouni. An Introduction to Random Matrices . CambridgeUniversity Press, 2009[2] L.-P. Arguin, D. Belius, and A. J. Harper. Maxima of a randomized Riemann zeta function, and branch-ing random walks. eprint = arXiv:1506.00629 (2015)[3] L.-P. Arguin and O. Zindy. Poisson-Dirichlet Statistics for the Extremes of a Log-Correlated GaussianField. Ann. Appl. Probab. , 1446-1481 (2014)[4] K. Astala, P. Jones, A. Kupiainen, and E. Saksman. Random conformal weldings. Acta Math. ,203–254 (2011)[5] J. Baik, P. Deift and E. Strahov. Products and ratios of characteristic polynomials of random Hermitianmatrices. J. Math. Phys. , 3657-3670 (2003) 6] E. Brezin and S. Hikami. Characteristic Polynomials of Random Matrices. Commun. Math. Phys. ,111-135 (2000)[7] P. Bourgade and J. Kuan. Strong Szeg˝o asymptotics and zeros of the zeta function. Comm. Pure Appl.Math. , 1028-1044 (2014)[8] A. Borodin and E. Strahov. Averages of characteristic polynomials in random matrix theory. Commun.Pure Appl. Math. , 161-253 (2006)[9] D. Bump and A. Gamburd. On the Averages of Characteristic Polynomials from Classical Groups. Commun. Math. Phys. , 227-274 (2006)[10] D. Carpentier and P. Le Doussal. Glass transition of a particle in a random potential, front selection innonlinear renormalization group, and entropic phenomena in Liouville and sinh-Gordon models. Phys.Rev. E , 026110, 33pp (2001)[11] T. Claeys and I. Krasovsky. Toeplitz determinants with merging singularities. eprint = arXiv:1403.3639 (2014)[12] J.B. Conrey, D.W. Farmer, J.P. Keating, M.O. Rubsintein & N.C. Snaith. Integral moments of L -functions. Proc. London. Math. Soc. , 33-104 (2005)[13] J.B. Conrey, P.J. Forrester & N.C. Snaith. Averages of ratios of characteristic polynomials for theclassical compact groups. Int. Math. Res. Not. , 397431 (2005)[14] J.B. Conrey, D.W. Farmer & M.R. Zirnbauer. Autocorrelation of ratios of L -functions. Commun. Num-ber Theory Phys. , 593636 (2008)[15] J. Ding, R. Roy and O. Zeitouni. Convergence of the centered maximum of log-correlated Gaussianfields. eprint = arXiv:1503.04588 (2015)[16] P. Doukhan, G. Oppenheim, and M. Taqqu. Theory and Applications of Long-Range Dependence .Birkhauser, Boston, 2003.[17] I. Dumitriu and A. Edelman. Matrix models for beta ensembles. J. Math. Phys. , 58305847 (2002)[18] B. Duplantier, R. Rhodes, S. Sheffield, and V. Vargas. Log-correlated Gaussian Fields: an overview.eprint = arXiv:1407.5605 (2014)[19] L. Erd¨os, H-T. Yau, and J. Yin. Bulk universality for generalized Wigner matrices. Probab. TheoryRelat. Fields , 341–407 (2012)[20] P. J. Forrester and S. O. Warnaar. The Importance of the Selberg Integral. Bull. Amer. Math. Soc.(N.S.) , 489-534 (2008)[21] Y. V. Fyodorov and J. P. Bouchaud. Freezing and extreme-value statistics in a random energy modelwith logarithmically correlated potential. J. Phys. A: Math. Theor. , no. 37, 372001, 12pp (2008)[22] Y. V. Fyodorov, G.H. Hiary, and J.P. Keating. Freezing Transition, Characteristic Polynomials ofRandom Matrices, and the Riemann Zeta-Function. Phys. Rev. Lett. , 170601, 5pp (2012)[23] Y. V. Fyodorov and J.P. Keating. Freezing Transitions and Extreme Values: Random Matrix Theory, ζ (1 / it ) and Disordered Landscapes. Philos. Trans. R. Soc. Lond. Ser. A Math. Phys. Eng. Sci. ,no. 2007, 20120503, 32pp (2014)[24] Y. V. Fyodorov, B. A. Khoruzhenko, and N. J. Simm. Fractional Brownian Motion with Hurst Index H = 0 and the Gaussian Unitary Ensemble. eprint = arXiv: 1312.0212 (2013)[25] Y. V. Fyodorov, P. Le Doussal, and A. Rosso. Statistical Mechanics of Logarithmic REM: Duality,Freezing and Extreme Value Statistics of 1/f Noises Generated by Gaussian Free Fields. J. Stat. Mech.Theory Exp. , no. 10, P10005, 32 pp (2009)[26] Y. V. Fyodorov , P. Le Doussal, and A. Rosso. Freezing Transition in Decaying Burgers Turbulenceand Random Matrix Dualities. Europhys. Lett. , 60004 , 6 pp (2010)[27] Y. V. Fyodorov, P. Le Doussal, and A. Rosso. Counting Function Fluctuations and Extreme ValueThreshold in Multifractal Patterns: The Case Study of an Ideal 1/f Noise. J. Stat. Phys. , 898–920(2012)[28] Y. V. Fyodorov and E. Strahov. An exact formula for general spectral correlation function of randomHermitian matrices. J. Phys. A: Math. Gen. , no.12, 3203-3213 (2003)[29] S.M. Gonek, C.P. Hughes and J.P. Keating. A Hybrid Euler-Hadamard product for the Riemann zetafunction. Duke Math. J. , 507-549 (2007)[30] C. P. Hughes, J. P. Keating, and N. O’Connell. On the Characteristic Polynomial of a Random UnitaryMatrix. Commun. Math. Phys. , 429–451 (2001) 31] K. Johansson. On fluctuations of eigenvalues of random Hermitian matrices. Duke Math. J. , no. 1,151–204 (1998)[32] J. P. Keating and N. C. Snaith. Random Matrix Theory and ζ (1 / it ). Commun. Math. Phys. ,57–89 (2000)[33] M. Kieburg and T. Guhr. Derivation of determinantal structures for random matrix ensembles in a newway. J. Phys. A: Math. Theor. , 075201, 31pp (2010)[34] N. Kistler. Derrida’s random energy models. From spin glasses to the extremes of correlated randomfields. eprint = arXiv: 1412.0958 (2014)[35] I. V. Krasovsky. Correlations of the characteristic polynomials in the Gaussian unitary ensemble or asingular Hankel determinant. Duke Math. J. , no. 3, 581–619 (2007)[36] M.R. Leadbetter, G. Lindgren & H. Rootzen. Extremes and Related Properties of Random Sequencesand Processes. Springer-Verlag. New York. 1982[37] B. B. Mandelbrot and J. W. van Ness. Fractional Brownian Motions, Fractional Noises and Applications. SIAM Review , no. 4, 422–437 (1968)[38] M. L. Mehta. Random Matrices . Academic Press; 3rd edition, 2004.[39] D. Ostrovsky. Mellin Transform of the Limit Lognormal Distribution. Comm. Math. Phys. , 287-310(2009)[40] D. Ostrovsky. Selberg Integral as a Meromorphic Function. Int. Math. Res. Notes 41 pp (2012)[doi:10.1093/imrn/rns170].[41] D. Ostrovsky. Theory of Barnes Beta Distributions. Electron. Commun. Prob. , no. 59, 116, (2012)[DOI: 10.1214/ECP.v18-2445][42] D. Ostrovsky. On Barnes Beta Distributions, Selberg Integral and Riemann Xi. Forum Mathematicum (2014) [DOI: 10.1515/forum-2013-0149][43] L. Pastur and M. Shcherbina. Eigenvalue Distribution of Large Random Matrices . AMS, 2011.[44] R. Rhodes and V. Vargas. Gaussian Multiplicative Chaos and applications: A review. Probability Surveys , 315-392 (2014) (electronic). DOI: 10.1214/13-PS218.[45] E. Strahov and Y. V. Fyodorov. Universal results for correlations of characteristic polynomials:Riemann-Hilbert approach. Commun. Math. Phys. , 343-382 (2003)[46] E. Subag, O. Zeitouni. Freezing and decorated Poisson point processes. Commun. Math. Phys. ,Issue 1, pp 55-92 (2015)[47] H. Trotter. Eigenvalue distributions of large Hermitian matrices; Wigner’s semicircle law and a theoremof Kac, Murdock, and Szeg¨o. Adv. in Math. , Issue 1, pp 67-82 (1984)[48] C Webb. The Characteristic Polynomial of a Random Unitary Matrix and Gaussian MultiplicativeChaos - the L2-Phase. eprint = arXiv:1410.0939 (2014) Queen Mary University of London, School of Mathematical Sciences, London E1 4NS,United Kingdom e-mails: [email protected] and [email protected]@qmul.ac.uk