Local Spectral Statistics of Gaussian Matrices with Correlated Entries

We prove optimal local law, bulk universality and non-trivial decay for the off-diagonal elements of the resolvent for a class of translation invariant Gaussian random matrix ensembles with correlated entries.


Introduction
Most rigorous works on random matrix ensembles concern either Wigner matrices with independent entries [16,23] (up to the real symmetric or complex hermitian symmetry constraint), or invariant ensembles where the correlation structure of the matrix elements is very specific. Since the existing methods to study Wigner matrices heavily rely on independence, only very few results are available on ensembles with correlated entries, see [10][11][12]19] for the Gaussian case. The global semicircle law in the non Gaussian case with (appropriately) weakly dependent entries has been established via moment method in [22] and via resolvent method in [18]. A similar result for sample covariance matrices was given in [20]. All these works establish limiting spectral density on the macroscopic scale and in models where the dependence is sufficiently weak so that the limiting density of states coincides with that of B László Erdős lerdos@ist.ac.at Oskari H. Ajanki oskari.ajanki@iki.fi Torben Krüger torben.krueger@ist.ac.at 1 IST Austria, Klosterneuburg, Austria the independent case. A more general correlation structure was explored in [5] with a nontrivial limit density, but still only on the global scale, see also [6]. We also mention the very recent proof of the local semicircle law and bulk universality for the adjacency matrix of the d-regular graphs [7,8] which has a completely different specific correlation (due to the requirement that every row contains the same number of ones).
In this paper we consider self-adjoint Gaussian random matrices H with correlated entries. We assume that H is of the form X + X * where the elements of X have a translation invariant correlation structure. Our main result is the optimal local law for H, i.e., we show that the empirical eigenvalue measure of H converges to a deterministic probability density ρ all the way down to the scale N −1 , the typical distance between eigenvalues, as the dimension N of H increases. We also find that the off diagonal elements of the resolvent G(z) := (H − z) −1 with Imz > 0 in the canonical basis are not negligible (unlike in the independent case) and in fact they inherit their decay from the correlation of the matrix elements. As a simple consequence of the local law we get bulk universality. Furthermore, we provide sufficient conditions for the asymptotic eigenvalue density ρ to be supported on a single interval with square root growth at both ends.
The proofs rely on the key observation that the (discrete) Fourier transform H = ( h φθ ) of a translation invariantly correlated self-adjoint random matrix H has independent entries up to an additional symmetry (cf. Lemma 3.2 below). Thus, our recent results [2] on the local law and bulk universality of Wigner type matrices with a general variance matrix can be applied. Some modifications to accommodate this extra symmetry are necessary in the proofs, but they do not influence the final result. The upshot is that in the Fourier space the diagonal elements of G(z) approximately satisfy the equation which constitutes a small perturbation of the quadratic vector equation (QVE), that was extensively analysed in [1,3]. Since the matrix S = (s φθ ) is typically not stochastic, the components m φ (z) of the solution genuinely depend on φ. We establish natural conditions on the correlation structure of H that guarantee that the recently developed theory [1,3] on QVEs is applicable. In particular, the stability of the QVE implies that the solutions of (1.1) and (1.2) are close, i.e., G φφ (z) = m φ (z) + o (1), even for spectral parameters z very close to the real axis, down to the scale Im z N −1 . This yields the local law for the eigenvalue density of H. Moreover, the anisotropic law from [2], applied to H, translates directly into a precise asymptotics for any matrix elements of the resolvent in the canonical basis: The off-diagonal decay of the entries of G(z) thus follows from smoothness properties of m φ (z) in the variable φ. We show that, in turn, this smoothness follows from the decay conditions on the correlation structure of H. Finally, we prove bulk universality of the local spectral statistics of H by using the analogous result from [2] for H and the fact that H and H are isospectral. Gaussian random matrices with translation invariant covariance structure have been analyzed earlier and it has also been realized that the Eq. (1.2) via Fourier transform plays a key role in identifying the limiting density of eigenvalues, see Khorunzhy and Pastur [19,21], Girko [17], as well as Anderson and Zeitouni in [5]. These works, however, were concerned only with the density on macroscopic scales. The off-diagonal decay of the resolvent and the bulk universality require much more detailed information. The current paper in combination with [1,2] presents such a precise analysis.

Set-Up and Main Results
Consider a real symmetric or complex hermitian random matrix, indexed by the large discrete torus of size N , We assume that the matrix is centered, i.e., and that the elements h i j are jointly Gaussian. The covariances of the elements of H are specified by two self-adjoint matrices Here the subtractions in i − k and j − l, etc., are done in the torus T. Let us also denote the graph distance of x ∈ T from the special point 0 ∈ T by |x|. We remark that any random matrix of the form H = X + X * , where X = (x i j ) i, j∈T is centred and translation invariant in the sense that (x i+k, j+l ) i, j∈T has the same law as X for any fixed shift (k, l) ∈ T 2 , has the correlation structure (2.3b).
The following properties of A are needed to state our main results: (D1) Power law decay: There is a positive integer κ, such that In general the solution of the QVE (1.2) specifying the asymptotic density of the states for H may be neither bounded nor stable (cf. [1,Sect. 9]). We will show that certain combinations of the above conditions exclude these issues.
The restrictions on the correlation structure are quantified by the N -independent model parameters ν, κ, ξ 1 , ξ 2 appearing above. We remark that the normalization of (2.4) and (2.5) is chosen for convenience, e.g., we could replace 1 on the right hand side of (2.4) by some finite constant. The set of model parameters depends on our assumptions, e.g., if only (D1) and (R2) are assumed, then κ and ξ 2 are the model parameters. We allow constants appearing in the statements to depend on the model parameters.
For compact statements of our results we define the notion of stochastic domination, introduced in [13,14]. This notion is designed to compare sequences of random variables in the large N limit up to small powers of N on high probability sets.
is a given function, depending only on the model parameters, as well as on an additional tolerance exponent γ ∈ (0, 1). For two sequences, ϕ = (ϕ (N ) ) N and ψ = (ψ (N ) ) N , of non-negative random variables we say that ϕ is stochastically dominated by ψ if for all ε > 0 and D > 0, (2.8) In this case we write ϕ ≺ ψ.
Let us denote the upper complex half plane and the discrete dual torus of T by H := z ∈ C : Im z > 0 , and S := N −1 T, respectively. It was shown in [1] that the quadratic vector equation (QVE) has a unique solution m(z) = (m φ (z)) φ∈S in H S , for every z ∈ H. Our main result is the optimal local law and the decay estimate for the off-diagonal resolvent entries. These are stated in terms of the resolvent G(z) = (G i j (z)) i, j∈T , Theorem 2.2 (Local law for Gaussian matrices with correlated entries) Suppose A is either exponentially decaying (D2) and non-resonant (R1), or decays like a power law (D1) and is strongly non-resonant (R2). Then for any tolerance exponent γ ∈ (0, 1) and uniformly for The vector q(z) = (q x (z)) x∈T inherits the decay type (exponential vs. power law) from A, in the sense that with the constants C > 0 and ν > 0 depending only on the model parameters.
Generally the off-diagonal resolvent entries are not negligible even though (2.13) states only an upper bound. In many cases matching lower bounds can be obtained. For example, for the special model with correlation a xy := e −ν(|x|+|y|) the QVE reduces to a simple scalar equation since a xy factorizes. An elementary calculation shows that in this case as N → ∞, Note that in the general setting of Theorem 2.2 the function π −1 Im q 0 (z) is the harmonic extension of the even probability density to the upper half plane. From (2.11b) it follows that the empirical spectral measure of H approaches the measure with the Lebesgue density ρ as N → ∞. In fact, using a comparison argument (cf. [6, Theorem 1]) this global convergence result extends also to non-Gaussian translation invariant random matrices satisfying (2.3). By applying the general theory for QVEs from [1] we are able to say more about the function q : H → C T , and the associated even probability density ρ : R → [0, ∞).
We remark that there are no explicit conditions on the correlation matrix B in either Theorem 2. Similarly, as in the case of Wigner type matrices the local law implies the bulk universality for Gaussian matrices with correlated entries. However, the q-fullness condition ([2, Definition 1.14]) is replaced by a different non-generacy condition. • H is real symmetric and A is strongly non-resonant (R2); • H is complex hermitian, and there is a constant ξ 3 > 0 such that where b φθ is defined analogously to a φθ in (2.10), and τ + := max{0, τ }, for τ ∈ R.
Then for any parameter ρ 0 > 0 and a smooth compactly supported function F : R n → R, n ∈ N, there exist constants c, C > 0, depending only on ρ 0 , κ, the function F, and either ξ 2 or ξ 3 , such that for any τ ∈ R with ρ(τ ) ≥ ρ 0 the local eigenvalue distribution is universal, Here, E G denotes the expectation with respect to the standard Gaussian ensemble, i.e., with respect to GUE and GOE in the cases of complex hermitian and real symmetric H, respectively, and ρ sc (0) = 1/(2π) is the value of Wigner's semicircle law at the origin.
The following result shows a practical way to construct real symmetric random matrices with translation invariant correlation structure. A similar, but slightly more complicated convolution representation exists for complex hermitian random matrices.

Lemma 2.6 (Linear filtering) Suppose a real symmetric matrix A satisfies the Bochner type condition
for arbitrary matrices W = (w i j ) i, j∈T . Then the random matrix H defined as the convolution,
This lemma is proven at the end of Sect. 5. We introduce the following conventions and notations used throughout this paper.

Convention 2.7 (Constants and comparison relation)
Symbols c, c 1 , c 2 , . . . and C, C 1 , C 2 , . . . denote generic positive and finite constants that depend only on the model parameters. They have a local meaning within a specific proof. For two arbitrary non-negative functions ϕ and ψ defined on some domain U , we write ϕ ψ, or equivalently ψ ϕ, if ϕ(u) ≤ Cψ(u), holds for all u ∈ U . The notation ψ ∼ ϕ is equivalent to both ψ ϕ and ψ ϕ holding at the same time. In this case we say that ψ and ϕ are comparable. In general the relation is called the comparison relation. We also write ψ = ϕ + O(ϑ) if |ψ − ϕ| ϑ.

Structure of the Proof
The proof of Theorem 2.2 splits into three separate parts. In the first part we show how to make H into an almost Wigner type matrix by changing basis. In the second part we describe how the proofs for Wigner type matrices in [2] are modified in order to accommodate some extra dependence in the transformed matrix. In the third part we show that the assumptions on the correlation matrix A imply that the QVE (2.9) has a bounded and sufficiently regular solution m using the general theory developed in [1]. Finally, in the last section we combine the results of the three steps and prove Theorem 2.2.

Mapping H into Wigner Type Matrix by Change of Basis
Since the mapping T → T corresponds to the conjugation by the unitary matrix F = ( f φ,y ) φ∈S,y∈T , with elements the matrices T and T = FTF * have the same spectrum: In the following we analyze random matrices which have independent entries modulo two reflection symmetries.

Definition 3.1 (fourfold correlated ensemble) A random matrix H indexed by a torus is fourfold correlated if h i j and h kl are independent unless
3) The next result shows that the discrete Fourier transform maps Gaussian translation invariant random matrices into Wigner type random matrices with an extra dependence. This connection was first realized by Girko [17] and Khorunzhy and Pastur [19]. It has been later used in [5,6,12].

Lemma 3.2 (Fourier transform) Let H be a (not necessarily Gaussian) random matrix satisfying (2.3). Then the elements of its Fourier transform H satisfy
for every φ, φ , θ, θ ∈ S. If additionally H is Gaussian, then H is fourfold correlated.
We remark that if a xy satisfies the decay estimate (2.4), then a φθ , b φθ N −1 .
Proof The proof of (3.4) is a straightforward computation. We omit further details. From (3.4b) we see that covariances between Re h φθ , Im h φθ and Re h φ θ , Im h φ θ can be non-zero if and only if the condition equivalent to (3.3) holds. Since the covariance matrix captures completely the dependence between the components of a Gaussian random vector the statement about the independence follows trivially.

Local Law for Fourfold Correlation
In this subsection we sketch how to prove a local law for the elements of the Fouriertransformed resolvent by slightly modifying the proof for the Wigner type matrices in [2]. Indeed, the analysis is the same as before, but due to the extra correlation between (φ, θ ) and (−φ, −θ) we have to remove both the rows and columns corresponding to indices φ and −φ from H in order to make it independent of a given row φ. We state a local law for a general self-adjoint random matrix with independent entries apart from a possible correlation of the entries with indices with s i j := E h i j 2 , is uniformly bounded in i and z, and that there exists a constant ε * > 0 to the real axis. Then for any γ > 0 the local law holds uniformly for every z = τ + iη, with η ≥ N γ −1 , and non-random w ∈ C T satisfying max i |w i | ≤ 1: The stochastic domination depends only on γ and ε * , and the constants μ = (μ k ), P, L , p appearing in the estimates contained in the assumptions of Theorem 1.6 from [2]: The extra symmetry condition (3.5) is automatically satisfied by random matrices with the covariance structure (3.4b), but it is generally not needed for the local law to hold (cf. [4] when S is stochastic).
Proof of Theorem 3. 3 We modify slightly the proof of Theorem 1.6 in [2]. The independence of the entries h i j and h −i,− j was used to estimate the off-diagonal resolvent entries and the perturbation d = d(z) of the perturbed QVE satisfied by the diagonal resolvent elements only in the proofs of Lemma 2.1 and Theorem 3.5 in [2]. In order to generalize Lemma 2.1 of [2] we apply the general resolvent identity (2.9) from [2] to replace the entries of G (k) by the corresponding entries of G (k,−k) in the defining formula (2.2) of d k in [2]. This way we obtain a representation for d k as a sum of terms each of which can be individually shown to be small by using either trivial bounds, or by using the large deviation estimates (2.7) similarly as in the proof of Lemma 2.1 in [2]. We will not present these estimates here, since a very similar analysis was carried out in Sect. 5 of [4]. The details for obtaining this representation for d k in the fourfold correlated random matrices are provided in Sect. 5.1 of [4]. We note that (3.9) is equivalent to formula (5.4) in [4] with the symbol ϒ k denoting d k . The off-diagonal resolvent elements are treated similarly by decoupling the dependence between specific rows of H and the entries of G. The treatment of the reflected off-diagonal elements G i,−i is simpler than in Sect. 5.2 of [4] since the extra symmetry (3.5) makes many error terms disappear. Instead of (3.5) another symmetry, h −x,−a = h ax , was assumed in [4]. Since all the factors of the form E h 2 xa in the error terms x in [4] first appeared in the form E h ax h −x,−a , which is zero in our case by (3.5), when following the proof in [4], we may replace the terms E h 2 xa with zeros. The fluctuation averaging ([2, Theorem 3.5]) is extended to fourfold correlated matrices also by slightly modifying the original proof of Theorem 4.7 in [14]. In particular, the arguments do not rely on the stochasticity of S as explained in the proof of Theorem 3.5 in [2]. In order to handle the extra dependence one needs to make a simple modification: The equivalence relation given within the proof needs to be generalized such that for a given index-tuple k = (k 1 , . . . , k p ) ∈ T p , we define r ∼ s to mean that either k r = k s or k r = −k s . This means that for each 'lone index' k one removes the index −k in addition to k from the other resolvent elements within the same monomial. For a more detailed description of the necessary modifications see the proof of Theorem 4.6 in [4].

Anisotropic Local Law for Fourfold Correlation
In order to translate the statements of the local law in Fourier coordinates back to the original coordinates we need an anisotropic local law. Here we consider |z| to be bounded to get simpler estimates. This condition can be easily dropped out if needed. j∈T is a self-adjoint fourfold correlated random matrix with centered entries satisfying the local law at some fixed z, satisfying |z| ≤ 10, where the non-random constant satisfies N −1/2 ≤ ≤ N −ε , for some ε > 0. Then uniformly for all z = τ + iη ∈ H satisfying η ≥ N γ −1 , and all non-random unit vectors u, v ∈ C T : Proof The proof is a straightforward generalization of the method applied in [9] to prove anisotropic local law for random covariance matrices and general Wigner matrices. The proof boils down to showing that for every p ∈ 2N there exists a constant C( p) independent of N such that for every v 2 ≤ 1 the moment bound holds. In the following we will denote generic constants depending only on p by C( p).
Only two minor modifications to the method used in Sect. 7 of [9] are needed. First, since S is not stochastic one needs to take into account that G ii (z) is close to m i (z) instead of an i-independent function such as the Stieltjes transform of the semicircle law, m sc (z). This generalisation was handled in [2] (cf. Theorem 1.12) where the anisotropic local law was proven for Wigner type matrices. As the second modification we need to incorporate the extra dependencies between the matrix elements h ak and h −a,−k into the analysis of [9]. To this end we walk through the key points of the arguments in [9] and point out along the way how the steps are modified to incorporate this extra dependence. The starting point of the argument is to write the right hand side of (3.12) in the form: for an arbitrary even integer p. Let us consider a fixed summand, so that the values of the v-indices b k1 , b k2 are fixed. Here the size of the expectation is naively bounded by p . However, there are 2 p sums over the elements of the 2 -unit vector v. Since v 1 ≤ N 1/2 the naive size of the right hand side of (3.13) is N p/2 p .
The key idea of the proof is to apply recursively the general resolvent identities (cf. [2, (2.9)]) to express the product of resolvent entries in (3.13) as a sum over so-called trivial leaves (cf. [9, Sect. 5.10]) and the sum over C( p) terms (corresponding to the non-trivial leaves in [9]) of the form (3.14) Here We remark that terms of the form (3.14) are coded by expressions A a b ,a b ( ζ ) in the formula (5.45) of [9]. The trivial leaves, exactly as in [9], correspond to products of resolvent entries that remain smaller than p even after summing over the v-indices simply because they contain products of so many off-diagonal resolvent elements of size O(N −ε ) that these together compensate the factor N p/2 originating from the brutal 1 -summation over |v b | (cf. [9,Sect. 5.11]). The classification of the constituents of the product of resolvent entries into the trivial and the nontrivial leaves relies on the concept of maximally expanded resolvent entries [9,Sect. 5.3]. For fourfold correlated matrices we redefine resolvent entries of the type G (B\ab) ab , with a, b ∈ B, as being maximally expanded. Here the set B plays the same role as the black vertices in [9].
From now on we concentrate on the non-trivial leaves of the type (3.14). The key property of these terms is that their expectation factorizes into an expectation over all the entries of H, and the expectation over all the entries of G. From (2.3a) it follows that the expectation over the entries of H can be non-zero only if each entry of H is paired with at least one of the four possible entries of H it is not independent of. As a consequence, either each v-index is paired with at least one other B-index, or there are so many extra entries of H compared to the number of independently summable indices in the products of (3.14) that the small sizes |h ai | ≺ N −1/2 compensate the presence of non-paired v-indices. In order to see that every term of the type (3.14) indeed has these properties one uses the same graph expansion as in [9] to perform the relevant bookkeeping.
The key insight about the combinatorics of the pairing of H-entries is that the number of ways to pair all C( p) of them in (3.14) is bounded by a number only depending on p, say by C( p) C( p) , but not on N . Such a factor can be included in the constant C( p) on the right hand side of (3.12), and is hence harmless. Since h bk may be paired not only with itself but also with h −b,−k , it is now possible that v b gets paired with v −b . However, using |v a | |v −a | ≤ |v a | 2 + |v −a | 2 , (3.15) these terms can be reduced to 2 -norms of v. Let us illustrate the modifications by considering the simplest leading order terms of the type (3.14) when p = 2. Considering the contribution of such terms to the right hand side of (3.13) yields which are stochastically dominated by C( p) 2 [2, (2.10)]) and the local law (3.10). Under the fourfold correlations the pairing produces also terms such as Here the off-diagonal resolvent elements are again stochastically dominated by . The sums over a and c can be bounded using (3.15) and v 2 ≤ 1. Hence, also (3.16) is stochastically dominated by 2 .

Properties of QVE
In this section we show that the assumptions on A in our main theorems guarantee that the induced QVE (2.9) has a sufficiently regular uniformly bounded solution. We show that the quantity q x−y (z) describing the asymptotic value of the off-diagonal resolvent elements If (R0) is assumed we will treat the associated parameters κ, K , ξ 0 as model parameters. By definition (R2) implies a(φ, θ ) ≥ ξ 2 for every φ, θ, and thus (R0) holds with ξ 0 = ξ 2 and D = {[0, 1]}. Assumption (R1) does not imply (R0), but (R1) and (D2) together do ( cf. Lemma 4.2 below).
Instead of directly analyzing the discrete QVE (2.9) we will first establish the correct properties for the solution of the continuous version of (2.9). Afterwards we deduce these properties for the discrete version (2.9) as well. For the transition from the discrete to the continuous version we need certain stability properties of the QVE that were established in [1]. We recall several notations and results from [1]. We will consider QVEs defined on a probability space (X, π) with an operator S in two different setups. When we discuss the discrete QVE (2.9), the setup is In the following, all L p -norms and the scalar products are understood in the appropriate setups (4.5).
The unique solution m to the discrete QVE (2.9) is close to m: Proof We prove (4.6) first. In order to apply the general theory for QVEs we first show that the integral operator A satisfies the conditions A1-A5 from [1]. The qualitative properties A1 and A2 are trivially satisfied. For the property A5 we show that the integral kernel of A K −1 is bounded from below by a constant comparable to one. This follows from (R0) since every element of the (K − 1)-th power of the indecomposable matrix Z is equal to or larger than one (cf. Proposition 4.3 of [1]). For A4 we need to show that the norm A 2→∞ of A as an operator from L 2 [0, 1] to L ∞ [0, 1] is O(1). This follows from (4.1), because (1 + |x|) a xy 1. (4.8) Finally, the normalization A3 of [1] holds if we replace A and m(z; φ) by λ A and λ −1/2 m(λ −1/2 z; φ), respectively, with λ := || A|| −1 ∞→∞ . From (4.8) it follows that λ ∼ 1. Next we show that m(z) is uniformly bounded for z = 0. Indeed, using (4.8) we get From this it follows that Since this implies the condition B1 of [1] (i) of Theorem 4.1 in [1] is applicable in the setup (4.5b). The theorem shows that m(z) ∞ ≤ C(δ) for any |z| ≥ δ with C(δ) depending on δ > 0. The property (R0) is equivalent to property B2 in [1]. Hence by (ii) of Theorem 4.1 in [1] m(z) is uniformly bounded in some neighborhood of z = 0. Combining this with the uniform bound away from z = 0 we get the uniform bound for | m(z; φ)| for all z and φ. In order to bound also the derivative ∂ φ m(z; φ) we differentiate the continuous QVE (4.4) and get (4.9) Using (4.8) and the uniform boundedness of m we finish the proof of (4.6). Im m(τ + iη; φ)dφ, (4.10) analogously to ρ in (2.14), in the continuous setting. With Theorem 1.1 of [1] we see that ρ is a bounded probability density on R. By the continuity (4.8) we also have and hence Theorem 1.9 from [1] yields supp ρ = [− β, β ], for some constant β ∼ 1. Now we prove (4.7) by considering (2.9) as a perturbation of (4.4). Given m we first embed S into [0, 1] canonically, and define the piecewise constant functions (4.12) for φ, θ ∈ [0, 1]. Notice that t (φ, θ ) = a(φ, θ ) and g(z; φ) = m φ (z) when φ, θ ∈ S. Hence it is enough to estimate g − m ∞ . Together with (4.8) this implies In terms of these quantities (2.9) can be written as (4.14) We will now consider g as the solution of the perturbed continuous QVE Using (4.13) we see that the perturbation d is indeed small: Clearly, T 2→∞ ∼ 1 as well. Hence, we know from the general theory (cf. the bound (1.2) of Theorem 1.1 of [1]) that g(z) 2 ≤ 2/ |z|. Using (4.13) we see that for sufficiently large N the operator T is also block fully indecomposable with the same matrix Z and the same partition D as A. Thus we get g(z) 2 1 for all z by (ii) of Theorem 4.1 in [1]. Combining this with (4.16) yields Comparing (4.15) and (4.4) we show that (4.17) implies that the corresponding solutions g and m are close in the sense of (4.7). For this purpose we use the rough stability statement from Theorem 1.10 of [1] to get where c 0 ∼ 1 and λ 1 ∼ 1 are sufficiently small constants left unspecified until the end of the proof. This means that we get stability as long as we stay away from the points ± β. The necessary initial bound inside the indicator function is satisfied for large enough values of |z|, since Here C 1 is a sufficiently large constant. This bound follows from the Stieltjes transform representation of both the solution of the discrete and the continuous QVE (cf. [1, Theorem 1.1]). We use continuity of g and m in z and (4.18) to propagate the initial bound from the regime of large values of |z| to all z ∈ H with dist(z, { β, − β}) ≥ c 0 . In particular, (4.18) remains true even without the indicator function, i.e., It remains to show (4.7) close to the edges by using that the instability at these two points is quadratic. The argument is a simplified version of the one used in a random setting in Sect. 4 of [2]. For the convenience of the reader we show a few details. We restrict to the case |z − β| ≤ c 0 , close to the right edge. The left edge is treated in the same way. For the following analysis we use the stability result, Theorem 4.2 of [2], in the continuous setup (cf. [1, Proposition 8.1]). The theorem yields where the quantity = (z) ≥ 0 is continuous in z and satisfies the cubic inequality Here the constant λ 2 ∼ 1 is independent of c 0 . Note that (4.21) corresponds to (4.10) in [2] and (8.5) in [1], respectively. Combining (4.11), (4.14b) and (4.5d) in [2], the coefficients π k = π k (z) of the cubic equation (4.21) satisfy provided c 0 ∼ 1 is sufficiently small. Since π 1 (z) → 0 as z → β, by decreasing the size c 0 of the neighborhood we are working on, the value of |π 1 | can be made arbitrarily small. This, in turn, implies that the solution of the cubic inequality (4.21) is small, Using this we can make the right hand side of (4.20) smaller than λ 2 /2, say, by decreasing the value of c 0 . Thus, there is a gap in the possible values of the continuous function z → g(z)− m(z) ∞ , in the sense that g − m ∞ / ∈ (λ 2 /2, λ 2 ). Since on the boundary, |z − β| = c 0 , the initial bound, g − m ∞ ≤ λ 2 , holds by (4.19), it propagates to all z with z − β ≤ c 0 . Thus, (4.20) and (4.21) remain true without the indicator functions.
It still remains to bound in (4.20). Since |π 2 | ∼ 1, we may absorb the cubic term in in (4.21). We find that satisfies where := π 1 /(π 2 + ). From this it is easy to see that the bound ≤ N −1/2 can be propagated from the boundary |z − β| = c 0 inside the neighborhood |z − β| ≤ c 0 of the right edge to give N −1/2 everywhere. Using this in (4.20) without the indicator function proves the bound (4.7) at the right edge. A satisfies (D2) and (R1), then also (R0) is satisfied, with the parameters κ, K , ξ 0 depending only on ξ 1 and ν.

Lemma 4.2 If
The part of the proof considering the exponentially decaying correlation matrices relies on the following technical result that is proven in the appendix.

Lemma 4.3 (Jensen-Poisson bound) Suppose f is an analytic function on the complex strip,
then for every ε > 0 there exists δ > 0 depending only on ε, ν, C 1 such that From the exponential decay assumption (D2) it follows that the kernel function a has an analytic extension to the complex strip R ν , where ν > 0 is the exponent from (2.5). Using Lemma 4.3 with f (ζ ) = a(φ, ζ )/ξ 1 for a fixed φ we see that for any ε > 0 there exists δ > 0 depending only on ε such that (4.28) Using (4.28) we now show that A is a block fully indecomposable operator, i.e., (R0) holds. From (4.8) we see that for every φ 1 , φ 2 , θ 1 , θ 2 ∈ [0, 1]. Let K ∈ N be so large that Let us define a partition D = {D k } K k=1 of [0, 1] and a matrix Z = (Z i j ) K i, j=1 , by By the choice of K , we have where e x is the Fourier-basis function. Then we show that q x (z) and q x (z) are so close that (2.13) holds. Let us first assume that A is exponentially decaying, i.e., (D2) holds. Let us periodically extend the kernel function a : [0, 1] 2 → [0, ∞) from (4.1) to all of R 2 . From (2.5) it follows that a can be further analytically extended to the product of complex strips R 2 ν/2 (cf. (4.24)), where ν > 0 is the exponent from (2.5). We will now show that q x (z) decays exponentially in this case. To see this we consider the function (z) : R ν → C, defined by (4.36) In particular, it follows that m(z; φ) = (z; φ) for every φ ∈ [0, 1]. Because a is uniformly continuous and the expression inside the parenthesis on the right hand side of (4.36) is bounded away from zero by a constant comparable to (sup z m(z) ∞ ) −1 when ζ ∈ R, there exists a constant ν < ν such that | (z; ζ )| ≤ C 0 for ζ ∈ R ν . Since a : R 2 ν → C is analytic also (z) : R ν → C is analytic. For any x ∈ Z we thus get by a contour deformation (4.37) where the integrals over the vertical line segments joining ±1 and ±1 − iν cancel each other due to periodicity of the integrand in the horizontal direction. Since x ∈ Z was arbitrary, taking absolute values of (4.37) yields the exponential decay: (4.38) Next we prove that (D1) implies | q x (z)| |x| −κ . To this end let ∂ denote the derivative w.r.t. the variable in [0, 1]. Using e x (φ) = e i2π xφ we get for any k ∈ N: Thus, we need to show that ∂ κ m(z) ∞ 1 uniformly in z. The proof is by induction on the number of derivatives of m. It is based on which follows from (4.9), and the following consequence of (2.4): As the second step of the proof we show that This proves (4.41) for |x| ≤ N 1/2 .
For |x| ≥ N 1/2 we bound q x = q x (z) directly by using the summation of parts where we have dropped two boundary terms of size Thus, estimating each term in the sum over j separately shows that |q x (z)| N −1/2 also in this case.
Next we show that the probability density ρ corresponding to the discrete QVE, via (2.14), is also regular and supported on a single interval.
Proof of Proposition 2.3 Uniform boundedness of m follows from Lemma 4.1. The other statements concerning the density ρ follow by using Theorems 1.1, 1.2, and 1.9 from [1]. As an input for Theorem 1.9 in [1], which shows that the support of ρ is a single interval, we use The first bound follows from N a φθ − a(φ, θ ) ≤ C N −1 . The last bound follows from (4.11). The components q x (z) inherit their analyticity trivially from m(z) since the sum in the definition (2.12) of q x (z) is absolutely summable.

Proofs for Local Law and Bulk Universality
The following is the strongest version of the local law we prove here.
Hence the local law for fourfold correlated matrices, Theorem 3.3, with H and A playing the roles of H and S, applies. In particular, (2.11b) follows. In order to get (2.11a) we use the anisotropic local law (Theorem 3.4). Indeed, fix two arbitrary elements x and y of T and define two unit vectors v and w of C T by setting v φ := N −1/2 e i2π xφ and w θ := N −1/2 e i2π yθ , ∀φ, θ ∈ S. From (3.1) and (2.12) we see that where v · w = i v i w i . Thus the anisotropic local law (3.11) implies (2.11a). The decay estimate (2.13) for q x is already proven in Lemma 4.4.
Next we show that the eigenvalues of H satisfy also the bulk universality provided the elements of h i j contain a small Gaussian GOE/GUE component.
Proof of Corollary 2. 4 We will show that there exists a Gaussian random matrix H (0) and a GOE/GUE matrix U that is independent of H (0) , such that the Gaussian random matrix where ε > 0 equals either ξ 2 or ξ 3 depending on the symmetry class, satisfies (2.3). The matrix H (0) is such that Theorem 2.2 is applicable since the associated correlation matrix A (0) satisfies (D1) and (R2). In particular, the eigenvalues of H (0) satisfy the rigidity estimate (1.33) of Corollary 1.10 in [2]. Here we note that the corollary holds trivially for the fourfold correlated matrix as its proof depends only on the local law and not on the dependence structure of H (0) . The bulk universality is hence proven exactly the same way as Theorem 1.15 in [2], using the method of [15].
Let us first consider the case where H is real symmetric. In this case, the equations (2.3) hold with B = A. Comparing this with the GOE correlation structure, we see that U also satisfies (2.3) with U and in place of H and A = B, respectively. Applying Lemma 3.2 and using φθ = N −1 we obtain the representation, where v φθ are the components of the Fourier transform of some GOE matrix V. Using this representation we can identify the matrix H (0) in (5.1). Namely, we define it in Fouriercoordinates, where the term in the square root is bounded from below by ε/2 by the assumption a φθ ≥ ε/N , and V (0) is a GOE random matrix that is independent of U. Since H (0) and U are independent satisfy (3.4). This immediately yields (2.3) for the matrix (5.1).
Next we consider the case where H is complex self-adjoint. First we remark that for a given pair of hermitian matrices (A, B) there exists a random matrix H satisfying (2.3) if and only if the following hold in the Fourier-space: a φθ ≥ 0, b −φ,−θ = b φθ , and b φθ ≤ a φθ a −φ,−θ , ∀φ, θ ∈ S. (5.4) The necessity of these conditions follow from a φθ = E h φθ 2 , b φθ = E h φθ h −φ,−θ (cf. (3.4b)), and the Cauchy-Schwartz inequality, If θ = −φ, then the identity holds by definition. In order to see that (5.4) is also a sufficient condition for there to exists a random matrix H satisfying (2.3) we consider a fixed index pair (φ, θ ) ∈ S 2 . From the hermitian symmetry and Lemma 3.2 it follows that the two elements h φ,θ and h −φ,−θ determine the four entries of H that may depend on h φ,θ . It is now easily checked that a given 4 ×4-real matrix can be a correlation matrix of the real random vector  (H, A, B). Since φθ = N −1 we see that the Fourier transform of the random matrix (5.1) satisfies (3.4) provided we choose U to be independent of H (0) . This is equivalent to (2.3) and the proof is complete.
From the proof of Corollary 2.4 we read off the convolution representation for symmetric translation invariant random matrices.
Proof of Corollary 2.5 Given the anisotropic local law (3.11) and the uniform boundedness (Lemma 4.1) of the solution m of the QVE (2.9), the delocalization of the eigenvalues is proven exactly the same way as Corollary 1.13 in [2].
Using (6.5) to bound the derivative and writing f := f • ζ 0 we get We will now bound the last integral using the Jensen-Poisson formula, where α j 's are the zeros of f in the unit disk D. The last sum is always non-negative since |α i | ≤ 1 and can be dropped. By splitting the integral into positive and negative parts we obtain an estimate for the integral on the right hand side of (6.7) where we have used (6.6) to get the second inequality. For the last bound we have used f (ω) = f ( ζ 0 (ω)) ≤ C 1 . Plugging this into (6.7) and recalling (6.3) we get φ ∈ [0, 1] : | f (φ)| < δ ≤ 2πC 3 (ζ 0 ) log 2C 1 log(1/δ) .
This finishes the proof as C 3 (ζ 0 ) and C 1 are independent of δ.