25mm25mm25mm26mm

# Out-of-time-ordered correlators for Wigner matrices

###### Abstract.

We consider the time evolution of the _out-of-time-ordered correlator_ (OTOC)
of two general observables $A$ and $B$ in a mean field chaotic quantum system
described by a random Wigner matrix as its Hamiltonian. We rigorously identify
three time regimes separated by the physically relevant _scrambling_ and
_relaxation_ times. The main feature of our analysis is that we express the
error terms in the optimal Schatten (tracial) norms of the observables,
allowing us to track the exact dependence of the errors on their rank. In
particular, for significantly overlapping observables with low rank the OTOC
is shown to exhibit a significant local maximum at the scrambling time, a
feature that may not have been noticed in the physics literature before. Our
main tool is a novel multi-resolvent local law with Schatten norms that
unifies and improves previous local laws involving either the much cruder
operator norm (cf. [10]) or the Hilbert-Schmidt norm (cf. [11]).

Giorgio Cipolloni111Princeton Center for Theoretical Science and Department of
Mathematics, Princeton University, Princeton, NJ 08544, USA.
<EMAIL_ADDRESS>László Erdős222Institute of Science and Technology
Austria, Am Campus 1, 3400 Klosterneuburg, Austria<EMAIL_ADDRESS>Joscha
Henheik222Institute of Science and Technology Austria, Am Campus 1, 3400
Klosterneuburg, Austria<EMAIL_ADDRESS>22footnotetext: Supported by
ERC Advanced Grant “RMTBeyond” No. 101020331.

Key words: Relaxation, Scrambling time, Multi-resolvent local law, Schatten
norm.

2020 Mathematics Subject Classification: 60B20, 82C10.

##  1\. Introduction

A basic feature of a strongly interacting quantum system is that local initial
states become non-local along the unitary time evolution, in particular they
become increasingly harder to distinguish by local observables. The simplest
way to detect this chaotic behavior is to monitor the overlap $\langle
A(t)B\rangle$ of the Heisenberg time evolution
$A(t):=\mathrm{e}^{-\mathrm{i}tH}Ae^{\mathrm{i}tH}$ of an observable $A$ with
another static observable $B$, where $H$ is the Hamiltonian and $A,B$ are
Hermitian operators. Here $\langle M\rangle:=\frac{1}{N}\mathrm{Tr}M$ denotes
the normalized trace of an $N\times N$ matrix, $N$ is the dimension of the
quantum state space. As time goes on, the overlap between two local
observables converges to its stationary value in a process called quantum
thermalisation. Since this stationary value is practically111In a closed
quantum system with finitely many degrees of freedom the initial state is
never fully lost as the stationary value still slightly depends on the
original overlap $\langle AB\rangle$, but it is suppressed by $N$; e.g. it
follows from (1.2) that  $\lim_{t\to\infty}{\mathbb{E}}_{GUE}\langle
A(t)B\rangle=\Big{(}1-\frac{1}{N+1}\Big{)}\langle A\rangle\langle
B\rangle+\frac{1}{N+1}\langle AB\rangle.$ factorized, $\langle A\rangle\langle
B\rangle$, the original observable $A$ becomes hardly detectable from its time
evolution by local observables $B$.

A more refined measure of the dynamically evolving quantum chaos is the out-
of-time-ordered correlator (OTOC) of two observables, defined as222In the
physics literature, the OTOC is usually defined without the factor $1/2$. We
chose it, however, for convenience.

$\mathcal{C}_{A,B}(t):=\frac{1}{2}\big{\langle}\big{|}[A(t),B]\big{|}^{2}\big{\rangle}$
(1.1)

measuring the evolution of the commutator of $A(t)$ and $B$. Starting with
commuting observables, $[A,B]=0$, this quantity initially grows, expressing
how the time evolution $A(t)$ of a local observable spreads (or scrambles) to
non-local degrees of freedom expressed by $B$. The moment when this growth
stops is called333We remark that some papers use slightly different
definition, here we follow the terminology of [27], [17, Section 3.3].  the
scrambling time $t_{*}$. Scrambling is closely related to thermalisation, but
it typically involves non-local observables $B$. Thus in a quantum system with
local interactions, the thermalisation time is smaller than $t_{*}$ and it is
independent of the system size, while the scrambling takes place on a longer
time scale until local information is shared with all degrees of freedom in
the system. Beyond the scrambling time, the OTOC settles to constant value at
a larger time scale called the relaxation time, and then it remains
essentially unchanged.

The fine distinction between thermalisation and scrambling became very popular
in physics about $15$ years ago motivated by the fundamental papers by Hayden
and Preskill [28] and Sekino and Susskind [40] related to the black hole
information paradox. The concept of the OTOC in quantum chaos research was
introduced in Kitaev’s lectures [30] on the connection between the Sachdev-Ye-
Kitaev (SYK) model and black holes. Owing to these fascinating connections,
the physics literature on OTOC in various interacting quantum systems has
become enormous; we refer the reader to the reviews [41, 27] and extensive
references therein. In contrast, OTOC has basically not been considered in the
mathematical literature apart from [34] that studies a very different model
than our current random matrix setup.

Besides the OTOC, quantum chaos has several other signatures: the conventional
one is the spectral statistics of the Hamiltonian. Following E. Wigner’s
groundbreaking observation, in a sufficiently chaotic quantum system the local
eigenvalue statistics are given by the universal Wigner-Dyson distribution
that depends only on the basic symmetries of the system. In the physics
literature the spectral statistics are often described by the spectral form
factor (SFF), or two point spectral correlator, defined as
$r_{2}(t):={\mathbb{E}}|\langle\mathrm{e}^{\mathrm{i}tH}\rangle|^{2}$, where
${\mathbb{E}}$ indicates a statistical averaging over an ensemble of
Hamiltonians. It is well known that the SFF tends to become universal for
large times444See the celebrated slope-dip-ramp-plateau picture, e.g. in
[35]., while it still reflects properties of the actual quantum system
(especially its density of states) for shorter times. A good physics summary
is found in [17], while a recent mathematical analysis of the SFF for general
Wigner matrices was given in [9]; more precise formulas are available for
exactly solvable ensembles [23, 24, 25].

The OTOC is a more refined description of quantum chaos than the SFF, as it
also involves observables. In particular, the SFF misses important features
like the sensitivity of chaos to the locality of the observables or early time
chaos, i.e. the exponential growth of the OTOC for certain strongly
interacting systems like SYK (called fast scramblers [40]) versus the
polynomial growth for slow scramblers like certain weakly chaotic systems
(see, e.g. [26, Section II] and references therein). Note that the SFF can be
recovered from the OTOC by averaging, either over the observables or over the
unitary group in case of unitarily invariant Hamiltonians, see [17]. For
example, if $H$ is a GUE random matrix, then [17, Eq. (57)-(58)])

${\mathbb{E}}_{GUE}\langle A(t)B\rangle=\langle A\rangle\langle
B\rangle+\frac{N^{2}r_{2}(t)-1}{N^{2}-1}\big{[}\langle AB\rangle-\langle
A\rangle\langle B\rangle\big{]},$ (1.2)

with

$r_{2}(t)=\Big{(}\frac{J_{1}(2t)}{t}\Big{)}^{2}+\frac{1}{N}-\frac{1}{N}\Big{(}1-\frac{t}{2N}\Big{)}\mathbf{1}(t\leq
2N),$

where $J_{1}$ is the Bessel function of the first kind of order one. Thus
$r_{2}(t)$ can be expressed from $\langle A(t)B\rangle$. A similar relation
holds between the OTOC and the four point spectral correlator.

The main goal of the current paper is to give a comprehensive mathematical
analysis of the OTOC with general observables $A,B$, when the Hamiltonian $H$
is a Wigner matrix. Wigner matrices represent the Hamiltonian of the most
chaotic quantum systems with matrix elements being independent, identically
distributed (i.i.d.) random variables. In the physics literature, random
matrix theory (RMT) is often used as a test case to see to what extent this
relatively simple model mimics the physics of more complicated systems such as
interacting many-body models (like SYK) or even models with nontrivial spatial
structure (like spin chains). Spectral statistics are remarkably robust,
especially the universality of the large time (so-called plateau) regime of
the SFF has proved to be ubiquitous in many different chaotic quantum systems,
in accordance with the celebrated Bohigas-Giannoni-Schmit conjecture [5]. The
OTOC is a more delicate quantity and, admittedly, its several interesting
features that appear in more realistic strongly coupled systems are not
captured by RMT. For example, the early time exponential growth of the OTOC is
not present in RMT and there is no qualitative difference between the
thermalisation and scrambling times since $H$ is mean field (see [17, Section
3.3] for a detailed analysis). The difference between chaotic and integrable
systems or the effect of their possible coexistence on the OTOC (studied e.g.
in [26]) are also not visible in RMT since Wigner matrices are fully chaotic.
Nevertheless, RMT becomes a good description beyond the scrambling time as
claimed in [17, Section 3.3] and demonstrated in [17, Section 6]. The
calculations are performed under the unitary invariance assumption and without
controlling the error terms.

In the main Theorem 2.2 of the current work, using very different methods, we
rigorously describe the behaviour of $\mathcal{C}_{A,B}(t)$ up to very long
times for general Wigner matrices (no unitary invariance assumed). We mimic
the locality of the observables by considering matrices $A,B$ that are far
from being full rank and track this effect throughout all error terms by using
tracial norms that are sensitive to the rank. We distinguish three time
regimes (see Figure 1 and Section 2.1); for short times (before the scrambling
time $t_{*}$) we find a quadratic growth in $t$; for intermediate times we
find that $\mathcal{C}_{A,B}(t)$ heavily depends on the ranks of $A,B$ and
their overlap $AB$. In particular, we detect a remarkable high peak of
$\mathcal{C}_{A,B}(t)$ when the ranks of $A$ and $B$ are small but their
overlap $\langle AB\rangle$ is still relatively large. To our knowledge, this
observation may be new even in the physics literature. We then identify the
_relaxation time_ , $t_{**}$, when the OTOC saturates, i.e. it becomes
essentially constant (equal its thermal limiting value) with small
oscillations. As expected, in the last regime, after $t_{**}$, our model
behaves universally; a qualitatively similar behaviour has been demonstrated
for several more complicated systems in the physics literature both
theoretically and numerically, see e.g. [17, 26, 41, 36, 27]. While for
technical reasons we cannot consider infinite times, our analysis is valid for
sufficiently long times to see all physically relevant features. For brevity
we carry out the proofs at infinite temperature, but our methods can easily be
extended to any finite temperature and we will give the corresponding formulas
(see Section 2.2 below).

Figure 1. The two curves show the behaviour of $\mathcal{C}_{A,B}(t)$ in two
different scenarios for two commuting traceless observables $A,B$ normalized
to $\langle A^{2}\rangle=\langle B^{2}\rangle=1$. The black curve represents
the case $A=B$ with $\mathrm{rank}(A)=N^{a}$, $a\in[0,1]$ where the OTOC
exhibits a large peak of size $N^{1-a}$ around the scrambling time $t_{*}\sim
1$. Afterwards, it decays to its thermal limiting value (normalized to one)
around the relaxation time $t_{**}\sim N^{\frac{1-a}{3}}$. The red curve
represents the case where $AB=0$. Here, both $t_{*}$ and $t_{**}$ are of order
one, independent of the ranks of $A$ and $B$. For more details see Section
2.1.

From the mathematical point of view, our work is the closest to [8, 10], where
the deterministic leading term of traces of products of observables at
different times, $\langle A_{1}(t_{1})A_{2}(t_{2})A_{3}(t_{3})....\rangle$,
were computed (see also [38, 39], where even the Gaussian fluctuations of such
chains were proven). Clearly the OTOC is a special case of (the difference of
two) such chains. The main novelty is that now we use only Schatten (tracial)
norms555The (normalized) $p$-th Schatten norm of a matrix
$A\in{\mathbb{C}}^{N\times N}$ is defined as $\langle|A|^{p}\rangle^{1/p}$ for
$p\in[1,\infty)$, where $|A|:=(AA^{*})^{1/2}$. of the observables in the
estimates, while [8, 38, 39] used the much cruder operator norm. In
particular, we can extend the time scale for the validity of our description.
More importantly, note that the interesting features of the OTOC are
manifested for small rank observables for which the operator norm is a major
overestimate and conceptually is an overkill.

The main tool is a concentration result, called multi-resolvent local law, for
alternating products of resolvents of random matrices and deterministic
matrices. More precisely, setting $G_{i}:=G(z_{i})$ and considering
deterministic matrices $A_{i}$, the main object of interest is

$G_{1}A_{1}G_{2}A_{2}G_{3}\dots A_{k-1}G_{k}$ (1.3)

for some fixed $k$. We will show that (1.3) concentrates around a
deterministic object and gives an upper bound on the fluctuation. The
interesting regime is the local one, i.e. when $|\mathrm{Im}\,z_{i}|\ll 1$.
Resolvents can then be converted to unitary time evolution
$\mathrm{e}^{\mathrm{i}tH}$ by standard contour integration.

Local laws in general assert that resolvents $G(z)$ tend to become
deterministic (with high probability) in the large $N$ limit even if the
spectral parameter $z$ is very close to the real axis (typically for any
$|\mathrm{Im}\,z|\gg N^{-1}$ in the bulk spectrum). For example, typical
single resolvent local laws for Wigner matrices assert that, for any fixed
$\xi>0$,666Traditionally [19, 31, 4], local laws did not consider arbitrary
test matrix $A$, but only $A=I$ or special rank one projections
$A=\bm{y}\bm{x}^{*}$ leading to the _isotropic local law_ in (1.4). General
$A$ was included later, e.g. in [22].

$\big{|}\langle(G(z)-m(z))A\rangle\big{|}\leq\frac{N^{\xi}\lVert
A\rVert}{N\eta}\,,\qquad\big{|}\langle{\bm{x}},(G(z)-m(z)){\bm{y}}\rangle\big{|}\leq\frac{N^{\xi}\lVert{\bm{x}}\rVert\lVert{\bm{y}}\rVert}{\sqrt{N\eta}}\,,\qquad\text{with}\qquad\eta:=|\mathrm{Im}\,z|\,,$
(1.4)

for a deterministic matrix $A\in{\mathbb{C}}^{N\times N}$ and deterministic
vectors ${\bm{x}},{\bm{y}}\in{\mathbb{C}}^{N}$ with very high probability as
$N$ becomes large. Here, $m(z)$ is the Stieltjes-transform of the Wigner
semicircle law:

$m(z)=m_{\mathrm{sc}}(z):=\int_{\mathbb{R}}\frac{1}{(x-z)}\rho_{\mathrm{sc}}(x)\,\mathrm{d}x,\qquad\quad\rho_{\mathrm{sc}}(x):=\frac{1}{2\pi}\sqrt{[4-x^{2}]_{+}}.$
(1.5)

However, the deterministic limit of a multi-resolvent chain (1.3) is not
simply $m(z_{1})m(z_{2})A_{1}A_{2}\dots$, i.e, one cannot mechanically replace
each $G$ by a scalar $m$, the actual formula is much more complicated, see
(3.2) below. As to the accuracy of this deterministic approximation, besides
$N$ and the imaginary part of the spectral parameter, $\eta=\mathrm{Im}\,z$,
the error term crucially depends on the appropriate norms of $A_{i}$ as well
as on the distinction whether $A_{i}$ is traceless or not. The fact that
traceless observables substantially reduce both the size of the deterministic
limit of (1.3) and of its fluctuation has first been observed and exploited in
[7], see also [10] for a comprehensive analysis of arbitrary long chains. The
results in [10] were optimal both in $N$ and $\eta$, but they all used the
simplest operator norm of $A_{i}$’s in the error term which is far from being
optimal for low rank observables.

Concerning the more accurate norms, only very recent papers [11, 14] started
deviating from the operator norm in the error terms. The main purpose of these
papers was to prove the Eigenstate Thermalisation Hypothesis (ETH) for random
matrices (originally posed by Deutsch in [18]) in its most optimal form,
including low rank observables. Moreover, the key point in [14] was to obtain
ETH also uniformly in the spectrum, including the critical edge regime. This
required to focus on a local law for
$\langle\mathrm{Im}\,GA\mathrm{Im}\,GA\rangle$ and to extract the smallness of
order $\rho^{2}$ at the spectral edge owing to the imaginary part of the
resolvents (here $\rho=\pi^{-1}|\mathrm{Im}\,m_{\mathrm{s}c}|$ is the local
density of states). However, [11, 14] _exclusively_ used the Hilbert-Schmidt
(HS) norm, $\langle|A|^{2}\rangle^{1/2}$, which caused suboptimal
$(N\eta)$-dependence. For the purpose of [11, 14], this suboptimality in
$(N\eta)$ was irrelevant since the proof of ETH relies on local laws in the
$\eta$–regime when $N\eta$ is practically order one.

In contrast to [11, 14], for studying the OTOC at shorter times, we need local
laws in the regime where $N\eta$ is large (since $\eta\sim 1/t$ dictated by
the contour integration and $t\ll N$). In the current paper we prove local
laws that are optimal both in $N$ and $\eta$ and use the optimal Schatten
norms of the observables. Especially, this allows for a more accurate
description of the OTOC in the physically relevant regime of small rank
observables. We, however, do not need to track the $\rho$-dependence or pay
attention to the imaginary parts. Therefore, the current work and [14] are
complementary; they effectively handle very different aspects of the local
law. While the fundamental idea of these two works is similar, both use the
Zigzag strategy described in Section 4, the actual proofs are quite different.
The main focus in [14] was to design and handle contour integral
representations that allowed us to reduce every estimate to resolvent chains
involving only $\mathrm{Im}\,G$’s. In the current paper $\mathrm{Im}\,G$ plays
no role, but we need to track the precise Schatten norms very carefully.

To illustrate the strength of our new result in comparison with the previous
bounds, we present the following three estimates for the simplest case $k=2$,
with $A_{1}=A_{2}=A$, $z_{1}=z_{2}$ in the bulk, and ignoring
$N^{\xi}$-factors for some arbitrarily small $\xi>0$:

$\Big{|}\langle GAGA\rangle-m^{2}\langle
A^{2}\rangle\Big{|}\lesssim\begin{cases}\dfrac{\|A\|^{2}}{N\eta}&\mbox{from
\cite[cite]{[\@@bibref{Number}{multiG}{}{}, Theorem~{}2.5]};}\\\
\sqrt{N\eta}\dfrac{\langle|A|^{2}\rangle}{N\eta}&\mbox{from
\cite[cite]{[\@@bibref{Number}{A2}{}{}, Theorem~{}2.2]}; }\\\
\dfrac{\langle|A|^{2}\rangle}{N\eta}+\dfrac{\langle|A|^{4}\rangle^{1/2}}{N\sqrt{\eta}}&\mbox{from
Theorem~{}\ref{thm:main}.}\end{cases}$ (1.6)

Note that our current result in the last line of (1.6) implies both previous
results since $\langle|A|^{4}\rangle^{1/2}\leq\sqrt{N}\langle|A|^{2}\rangle$
and $\langle|A|^{p}\rangle^{1/p}\leq\|A\|$. While the former bound saturates
for low rank observables, the latter saturates for high rank ones. Therefore,
our new result with Schatten norms optimally interpolates between these two,
at least in the bulk regime.

### Notations

By $\lceil x\rceil:=\min\\{m\in{\mathbb{Z}}\colon m\geq x\\}$ and $\lfloor
x\rfloor:=\max\\{m\in{\mathbb{Z}}\colon m\leq x\\}$ we denote the upper and
lower integer part of a real number $x\in{\mathbb{R}}$. For $k\in{\mathbb{N}}$
we set $[k]:=\\{1,...,k\\}$, and $\langle A\rangle:=d^{-1}\mathrm{Tr}(A)$,
$d\in{\mathbb{N}}$, for the normalised trace of a $d\times d$-matrix $A$,
while $\mathrm{rk}\,A\equiv\mathrm{rank}\,A$ denotes its rank. For positive
quantities $f,g$ we write $f\lesssim g$ resp. $f\gtrsim g$ and mean that
$f\leq Cg$ resp. $f\geq cg$ for some $N$-independent constants $c,C>0$ that
may depend only on the basic control parameters $C_{p}$, see (2.1) in
Assumption 2.1 below. Moreover, we will also write $f\sim g$ in case that
$f\lesssim g$ and $g\lesssim f$.

We denote vectors by bold-faced lower case Roman letters
$\bm{x},\bm{y}\in{\mathbb{C}}^{N}$, for some $N\in{\mathbb{N}}$, and define

$\langle\bm{x},\bm{y}\rangle:=\sum_{i}\bar{x}_{i}y_{i}\,,\qquad
A_{\bm{x}\bm{y}}:=\langle\bm{x},A\bm{y}\rangle\,.$

Matrix entries are indexed by lower case Roman letters $a,b,c,...,i,j,k,...$
from the beginning or the middle of the alphabet and unrestricted sums over
those are always understood to be over $\\{1,...,N\\}$.

We will use the concept _’with very high probability’_ , meaning that for any
fixed $D>0$, the probability of an $N$-dependent event is bigger than
$1-N^{-D}$ for all $N\geq N_{0}(D)$. Also, we will use the convention that
$\xi>0$ denotes an arbitrarily small positive exponent, independent of $N$.
Moreover, we introduce the common notion of _stochastic domination_ (see,
e.g., [20]): For two families

$X=\left(X^{(N)}(u)\mid N\in{\mathbb{N}},u\in
U^{(N)}\right)\quad\text{and}\quad Y=\left(Y^{(N)}(u)\mid
N\in{\mathbb{N}},u\in U^{(N)}\right)$

of non-negative random variables indexed by $N$, and possibly an additional
parameter $u$ from a parameter space $U^{(N)}$, we say that $X$ is
stochastically dominated by $Y$, if for all $\epsilon,D>0$ we have

$\sup_{u\in
U^{(N)}}\mathbf{P}\left[X^{(N)}(u)>N^{\epsilon}Y^{(N)}(u)\right]\leq N^{-D}$

for large enough $N\geq N_{0}(\epsilon,D)$. In this case we write $X\prec Y$.
If for some complex family of random variables we have $|X|\prec Y$, we also
write $X=O_{\prec}(Y)$.

##  2\. Main results

We consider $N\times N$ Wigner matrices $W$, i.e. $W$ is a random real
symmetric or complex Hermitian matrix $W=W^{*}$ with independent entries (up
to the Hermitian symmetry) and with single entry distributions
$w_{aa}\stackrel{{\scriptstyle\mathrm{d}}}{{=}}N^{-1/2}\chi_{\mathrm{d}}$, and
$w_{ab}\stackrel{{\scriptstyle\mathrm{d}}}{{=}}N^{-1/2}\chi_{\mathrm{od}}$,
for $a>b$. The random variables $\chi_{\mathrm{d}},\chi_{\mathrm{od}}$ satisfy
the following assumptions.777By inspecting our proof, it is easy to see that
actually we do not need to assume that the off-diagonal entries of $W$ are
identically distributed. We only need that they all have the same second
moments, but higher moments can be different.

######  Assumption 2.1.

The off-diagonal distribution $\chi_{\mathrm{od}}$ is a real or complex
centered random variable, ${\mathbb{E}}\chi_{\mathrm{od}}=0$, satisfying
${\mathbb{E}}|\chi_{\mathrm{od}}|^{2}=1$. The diagonal distribution is a real
centered random variable, ${\mathbb{E}}\chi_{\mathrm{d}}=0$. Furthermore, we
assume the existence of high moments, i.e. for any $p\in{\mathbb{N}}$ there
exists $C_{p}>0$ such that

${\mathbb{E}}\big{[}|\chi_{\mathrm{d}}|^{p}+|\chi_{\mathrm{od}}|^{p}\big{]}\leq
C_{p}\,.$ (2.1)

Our main result, Theorem 2.2 below, concerns the Heisenberg time evolution
$A(t):=\mathrm{e}^{-\mathrm{i}Wt}A\mathrm{e}^{\mathrm{i}Wt}$ of a fixed
deterministic self-adjoint observable $A=A^{*}\in{\mathbb{C}}^{N\times N}$
governed by the Wigner matrix $W$. More precisely, we consider the _out-of-
time-ordered correlator (OTOC)_

$\mathcal{C}_{A,B}(t):=\frac{1}{2}\big{\langle}\big{|}[A(t),B]\big{|}^{2}\big{\rangle}=\big{\langle}A(t)^{2}B^{2}\big{\rangle}-\big{\langle}A(t)BA(t)B\big{\rangle}$
(2.2)

with another self-adjoint observable $B=B^{*}\in{\mathbb{C}}^{N\times N}$,
consisting of a two-point and a four-point888We remark that some papers (see,
e.g., [37]) refer to the four-point part $\mathcal{F}_{A,B}(t)$ alone as the
OTOC.  part,

$\mathcal{D}_{A,B}(t):=\big{\langle}A(t)^{2}B^{2}\big{\rangle}\qquad\text{and}\qquad\mathcal{F}_{A,B}(t):=\big{\langle}A(t)BA(t)B\big{\rangle}\,,\quad\text{respectively.}$
(2.3)

In the formulation of Theorem 2.2, a key role is played by the Fourier
transform of the semicircular density (1.5),

$\varphi(t):=\widehat{\rho_{\mathrm{sc}}}(t)=\int_{-2}^{2}\mathrm{d}x\rho_{\mathrm{sc}}(x)\mathrm{e}^{\mathrm{i}xt}=\frac{J_{1}(2t)}{t}\,,\quad
t\in\mathbb{R},$ (2.4)

where we recall that $J_{1}$ is the Bessel function of the first kind of order
one. We note that by standard asymptotics of the Bessel functions on the real
line, it holds that

$J_{1}(s)=\begin{cases}-\mathrm{sgn}(s)\cos\left(|s|+\dfrac{\pi}{4}\right)\sqrt{\dfrac{2}{\pi|s|}}+\mathcal{O}\left(\dfrac{1}{|s|^{3/2}}\right)\quad&s\to\pm\infty\,,\\\\[11.38109pt]
\dfrac{s}{2}-\dfrac{1}{2}\left(\dfrac{s}{2}\right)^{3}+\mathcal{O}(|s|^{5})\quad&s\to
0\,.\end{cases}$ (2.5)

For simplicity, we formulate our main result only for traceless matrices,
$\langle A\rangle=\langle B\rangle=0$ and at infinite temperature. For general
observables, see Remark 2.3 and for finite temperature, see Section 2.2.

######  Theorem 2.2 (OTOC for Wigner matrices).

Let $W$ be a Wigner matrix satisfying Assumption 2.1 and let
$A,B\in{\mathbb{C}}^{N\times N}$ be self-adjoint deterministic matrices which
are traceless, $\langle A\rangle=\langle B\rangle=0$. Fix any $\epsilon>0$.
Then, the OTOC (2.2) satisfies

$\begin{split}\mathcal{C}_{A,B}(t)=&\langle A^{2}\rangle\langle
B^{2}\rangle\big{[}1-\varphi(t)^{2}\big{]}+2\langle
AB\rangle^{2}\varphi(t)^{2}\,\big{[}\varphi(2t)-\varphi(t)^{2}\big{]}\\\\[5.69054pt]
&+\langle A^{2}B^{2}\rangle\,\varphi(t)^{2}-\langle
ABAB\rangle\varphi(t)^{4}+\mathcal{O}_{\prec}\big{(}\mathcal{E}_{A,B}(t,N)\big{)}\end{split}$
(2.6)

with an error term

$\mathcal{E}=\mathcal{E}_{A,B}(t,N):=\mathrm{e}^{t/N^{1/2-\epsilon}}\left(\frac{|t|^{4}\langle
A^{2}\rangle^{2}}{N}+\frac{|t|\langle
A^{8}\rangle^{1/2}}{N}\right)^{1/2}\left(\frac{|t|^{4}\langle
B^{2}\rangle^{2}}{N}+\frac{|t|\langle B^{8}\rangle^{1/2}}{N}\right)^{1/2}.$
(2.7)

The proof of Theorem 2.2 is given in Section 3.3 below. It is based on a novel
_multi-resolvent local law_ with error terms involving optimal Schatten norms,
see Theorem 3.3.

######  Remark 2.3.

We have several comments on Theorem 2.2:

  * (i)

_[Non-traceless observables]_ For general observables $A,B$, we can decompose
them into a tracial and traceless part, $A=:\langle A\rangle+\mathring{A}$,
and similarly for $B$. The tracial parts, $\langle A\rangle$ and $\langle
B\rangle$, then commute with the unitary time evolution
$\mathrm{e}^{\mathrm{i}tW}$ and one straightforwardly obtains a result similar
to Theorem 2.2 (see also Remark 3.5).

  * (ii)

_[Variance of fluctuations]_ The size of the fluctuations around the
deterministic leading term in (2.6), i.e. the variance of
$\mathcal{C}_{A,B}(t)$, is explicitly computable, following the arguments
leading to [38, Lemma 2.5]. The result is expressible purely in terms of
Schatten norms of $A$ and $B$ (cf. [38, Lemma 2.5 and Definition 3.4]),
however in [38, 39] the error terms are still in terms of crude operator
norms.

  * (iii)

_[Gaussianity]_ It is also possible to prove Gaussianity of the fluctuations
of $\mathcal{C}_{A,B}(t)$ (cf. [38, Theorem 2.7 and Corollary 2.12]) by
showing an approximate Wick theorem for resolvent chains, similarly to [38,
Theorem 3.6]. However, we refrain from doing this for brevity of the current
paper.

###  2.1. Physical interpretation of Theorem 2.2 by two examples

We will now discuss the behavior of $\mathcal{C}_{A,B}(t)$ in two exemplary
and extreme situations of observables $A,B$. In the first example we will set
the two observables identical, in the second we will assume that their product
vanishes. More concretely, we define

$\begin{split}A&=N^{\frac{1-a}{2}}\mathrm{diag}(1,-1,...,1,-1,0,...,0)\,,\qquad
a\in[0,1)\,,\\\
B&=N^{\frac{1-b}{2}}\mathrm{diag}(0,...,0,1,-1,...,1,-1)\,,\qquad
b\in[0,1)\,,\end{split}$ (2.8)

in such a way that $AB=BA=0$, and $\langle A\rangle=\langle B\rangle=0$ as
well as $\langle A^{2}\rangle=\langle B^{2}\rangle=1$, i.e. $A$ (resp. $B$)
contains $N^{a}$-many (resp. $N^{b}$-many) non-zero entries on the diagonal.

Example 1. For the first example, we have

$\mathcal{C}_{A,A}(t)=\langle
A^{2}\rangle^{2}\big{[}1-\varphi(t)^{2}\big{\\{}1-2\varphi(2t)+2\varphi(t)^{2}\\}\big{]}+\langle
A^{4}\rangle\,\big{[}\varphi(t)^{2}-\varphi(t)^{4}\big{]}+\mathcal{O}_{\prec}\big{(}\mathcal{E}_{A,A}(t,N)\big{)}\,.$
(2.9)

Example 2. For the second example, we have

$\mathcal{C}_{A,B}(t)=\langle A^{2}\rangle\langle
B^{2}\rangle\big{[}1-\varphi(t)^{2}\big{]}+\mathcal{O}_{\prec}\big{(}\mathcal{E}_{A,B}(t,N)\big{)}\,.$
(2.10)

The key features of these two examples (2.9)–(2.10) are briefly summarized in
Table 1. Ignoring the respective error terms, $\mathcal{C}_{A,A}(t)$ and
$\mathcal{C}_{A,B}(t)$ are schematically depicted in Figure 1. Now we comment
on each time regime.

  * (i)

[Short-time regime] By the asymptotics in (2.5), we have the short-time
asymptotic

$\displaystyle\mathcal{C}_{A,A}(t)$ $\displaystyle=t^{2}\,\big{(}\langle
A^{4}\rangle-\langle A^{2}\rangle^{2}\big{)}+\mathcal{O}\big{(}|t|^{4}\langle
A^{4}\rangle\big{)}+\mathcal{O}_{\prec}\big{(}\mathcal{E}\big{)}\qquad\qquad$
$\displaystyle\mbox{(Example 1)}\,,$ $\displaystyle\mathcal{C}_{A,B}(t)$
$\displaystyle=t^{2}\,\langle A^{2}\rangle\langle
B^{2}\rangle+\mathcal{O}\big{(}|t|^{4}\langle A^{2}\rangle\langle
B^{2}\rangle\big{)}+\mathcal{O}_{\prec}\big{(}\mathcal{E}\big{)}\qquad\qquad$
$\displaystyle\mbox{(Example 2)}\,.$

This shows, that the OTOC for Wigner matrices does _not_ exhibit the
exponential increase $\sim\mathrm{e}^{2\lambda t}$, observed for quantum
systems with a _classically_ chaotic analogue, where $\lambda$ is the Lyapunov
exponent of the classical system. Instead, the OTOC (2.2) behaves
polynomially, as expected for quantum chaotic systems without a classical
analogue (e.g. for certain spin$-\tfrac{1}{2}$ chains [26]).

  * (ii)

[Scrambling time] The monotonous growth of both $\mathcal{C}_{A,A}(t)$ and
$\mathcal{C}_{A,B}(t)$ stops at a time of order one (using elementary
properties of $\varphi$ from (2.4)), hence the scrambling time is $t_{*}\sim
1$. However, the maximally attained value strongly differs for the two
examples: While $\mathcal{C}_{A,A}(t_{*})\sim N^{1-a}$ heavily depends on $a$
(i.e. the rank of $A$), the peak of $\mathcal{C}_{A,B}(t_{*})\sim 1$ is
independent of the ranks of $A$ and $B$.

  * (iii)

[Intermediate regime up to the relaxation time] The following interval of
intermediate times is characterized by a decay of the OTOC (2.6) towards its
thermal limiting value $\langle A^{2}\rangle\langle B^{2}\rangle=1$ up to the
relaxation time $t_{**}$. This regime is also quite different for the two
examples (2.9)–(2.10): While for Example 1 the interval of intermediate times
is given by $t\in[t_{*},t_{**}]$ where $t_{**}\sim N^{(1-a)/3}$, in Example 2
the relaxation time $t_{**}\sim 1$ is comparable with the scrambling time.
However, for technical reasons, the _entire_ interval of intermediate times is
only accessible if $a>5/11$ for Example 1, and $a+b>2/3$ for Example 2,
otherwise the leading terms in (2.9)–(2.10) become smaller than their
respective error term. In comparison, computing the OTOC with the operator
norm in the error terms [10, Corollary 2.7], would lead to the (more
restrictive) conditions $a>5/8$ for Example 1, and $a+b>1$ for Example 2.

  * (iv)

[Long-time regime] In the consecutive long-time regime, i.e. $t\gg
N^{(1-a)/3}$ for Example 1 and $t\gg 1$ for Example 2, we find the OTOC (2.6)
to concentrate around its thermal limiting value $\langle A^{2}\rangle\langle
B^{2}\rangle=1$ with small oscillations. This confirms the expectation, that
the OTOC in strongly chaotic systems exhibits only small fluctuations for long
times. These are accessible up to

$t\leq
N^{\min\left\\{\frac{1}{4},\,\frac{3a-1}{2}\right\\}-\epsilon}\qquad\text{and}\qquad
t\leq
N^{\min\left\\{\frac{1}{4},\,\frac{3a+1}{10},\,\frac{3b+1}{10},\,\frac{3(a+b)-2}{4}\right\\}-\epsilon}$

for Example 1 and 2, respectively. Again, in comparison with Theorem 2.2, the
operator norm error terms from [10, Cor. 2.7] would lead to the constraints
$t\leq N^{\frac{2a-1}{2}-\epsilon}$ for Example 1, and $t\leq
N^{\frac{a+b-1}{2}-\epsilon}$ for Example 2.

| Ex. 1 (2.9): $A=B$, $\mathrm{rk}\,A=N^{a}$ | Ex. 2 (2.10): $AB=0$, $\mathrm{rk}\,A=N^{a}$, $\mathrm{rk}\,B=N^{b}$
---|---|---
short times | $\mathcal{C}_{A,A}(t)\sim t^{2}\,\big{(}\langle A^{4}\rangle-\langle A^{2}\rangle^{2}\big{)}$ | $\mathcal{C}_{A,B}(t)\sim t^{2}\,\langle A^{2}\rangle\langle B^{2}\rangle$
scrambling time |  $t_{*}\sim 1$ and $\mathcal{C}_{A,A}(t_{*})\sim N^{1-a}$ |  $t_{*}\sim 1$ and $\mathcal{C}_{A,B}(t_{*})\sim 1$
intermediate times | full access if $a>5/11$ | full access if $a+b>2/3$
relaxation time |  $t_{**}\sim N^{\frac{1-a}{3}}$ and $\mathcal{C}_{A,A}(t_{**})\sim\langle A^{2}\rangle^{2}$ |  $t_{**}\sim 1$ and $\mathcal{C}_{A,B}(t_{**})\sim\langle A^{2}\rangle\langle B^{2}\rangle$
long times | up to $t\leq N^{\min\left\\{\frac{1}{4},\,\frac{3a-1}{2}\right\\}-\epsilon}$ | up to $t\leq N^{\min\left\\{\frac{1}{4},\,\frac{3a+1}{10},\,\frac{3b+1}{10},\,\frac{3(a+b)-2}{4}\right\\}-\epsilon}$
Table 1. Overview of the two examples (2.9) and (2.10).

###  2.2. Finite temperature case

Theorem 2.2 can easily be extended to the case of finite temperature,
$\beta=1/T>0$. The OTOC (2.2) now is given by

$\mathcal{C}^{(\beta)}_{A,B}(t):=\frac{1}{2}\frac{\mathrm{Tr}\big{[}|[A(t),B]|^{2}\mathrm{e}^{-\beta
W}\big{]}}{Z}$ (2.11)

with partition function $Z=\mathrm{Tr}\big{[}\mathrm{e}^{-\beta W}\big{]}$. In
this case, the analog of (2.6) in the regime999We restrict to this regime for
simplicity, as all the error terms can be absorbed into $\prec$. $\beta\ll\log
N$ reads

$\begin{split}\mathcal{C}^{(\beta)}_{A,B}(t)=&\langle A^{2}\rangle\langle
B^{2}\rangle\big{[}1-\varphi(t)^{2}\big{]}+\langle
AB\rangle^{2}\frac{\varphi(t)}{\varphi(\mathrm{i}\beta)}\,\mathrm{Re}\,\bigg{[}\varphi(2t)\varphi(t+\mathrm{i}\beta)+\varphi(t)\varphi(2t+\mathrm{i}\beta)-2\varphi(t)^{2}\varphi(t+\mathrm{i}\beta)\bigg{]}\\\\[5.69054pt]
&+\langle A^{2}B^{2}\rangle\,\varphi(t)^{2}-\langle
ABAB\rangle\frac{\varphi(t)^{3}\mathrm{Re}\,[\varphi(t+\mathrm{i}\beta)]}{\varphi(\mathrm{i}\beta)}+\mathcal{O}_{\prec}\left(\mathcal{E}_{A,B}(t,N)\right)\,,\end{split}$
(2.12)

where $\mathcal{E}$ is from (2.7). Here $\varphi(z)$ is the complex extension
of $\varphi(t)$ for $z\in{\mathbb{C}}$; note that $\varphi(z)$ is generically
complex but $\varphi(\mathrm{i}\beta)$ is real.

Figure 2. Depicted are four curves illustrating the influence of $\beta=1/T$
on the OTOC $\mathcal{C}_{A,A}^{(\beta)}(t)$ up to intermediate times for
Example 1 from Section 2.1 (i.e. normalized to $\langle A^{2}\rangle=1$ with
$\mathrm{rank}\,A=N^{\frac{1-a}{2}}$ and $a\in[0,1]$). As $\beta$ increases,
the characteristic rank-dependent peak of size $\sim N^{1-a}$ around the
scrambling time $t_{*}\sim 1$ becomes more pronounced and very slightly
shifted to the left.

Moreover, using the asymptotics of the Bessel function in the complex plane,
we have

$\varphi(z)=-\sqrt{\frac{1}{\pi
z^{3}}}\left(\cos\left(2z+\frac{\pi}{4}\right)+\mathrm{e}^{2|\mathrm{Im}\,z|}\mathcal{O}\left(|z|^{-1}\right)\right)\qquad\text{for}\qquad|z|\to\infty\quad\text{with}\quad|\arg
z|<\pi\,.$

In particular, the thermal limiting value of $\mathcal{C}^{(\beta)}_{A,B}(t)$
is _independent_ of $\beta$ at least in our regime $\beta\ll\log N$. Note that
for much larger $\beta\gtrsim\sqrt{N}$ physics calculations predict a
temperature dependence of the thermal limiting value of the OTOC , see [36,
Eqs. (3.8)–(3.9)].

However, before the long-time regime and neglecting the error term in (2.12),
we find a strong dependence of the OTOC on temperature for Example 1 from
Section 2.1 (cf. (2.8)–(2.9)) as illustrated in Figure 2. In contrast to that,
for Example 2 from Section 2.1 (cf. (2.8) and (2.10)), the whole OTOC curve
(as depicted in Figure 1) is independent of temperature at any time. This is
because all the $\beta$ dependent terms in (2.12) drop out since $AB=BA=0$.

##  3\. Proof of Theorem 2.2: Multi-resolvent local law with Schatten norms

Theorem 2.2 relies on a new multi-resolvent local law for alternating chains
(1.3) with deterministic matrices via simple contour integration (see Section
3.3). Its main novelty, using Schatten norms for the observables and still
keeping optimality in $N$ and $\eta$, has already been explained in the
Introduction. After collecting some preliminary information, in Section 3.2 we
present our new local law, Theorem 3.3, then in Section 3.3 we quickly
complete the proof of Theorem 2.2. Starting from Section 4 we will focus on
the proof of Theorem 3.3.

###  3.1. Preliminaries on the deterministic approximation

Before stating our main technical result, we introduce some additional
notation. Given a non-crossing partition $\pi$ of the set
$[k]:=\set{1,\ldots,k}$ arranged in cyclic order, the partial trace
$\mathrm{pTr}_{\pi}$ is defined as

$\mathrm{pTr}_{\pi}(B_{1},\ldots,B_{k-1}):=\prod_{S\in\pi\setminus\mathfrak{B}(k)}\left\langle\prod_{j\in
S}B_{j}\right\rangle\prod_{j\in\mathfrak{B}(k)\setminus\set{k}}B_{j},$ (3.1)

with $\mathfrak{B}(k)\in\pi$ denoting the unique block containing $k$. Then,
for generic $B_{i}$’s, the deterministic approximation of (1.3) is given by
(see [8, Theorem 3.4])

$M_{[1,k]}=M(z_{1},B_{1},\ldots,B_{k-1},z_{k}):=\sum_{\pi\in\mathrm{NC}([k])}\mathrm{pTr}_{K(\pi)}(B_{1},\ldots,B_{k-1})\prod_{S\in\pi}m_{\circ}[S],$
(3.2)

where $\mathrm{NC}([k])$ denotes the non-crossing partitions of the set $[k]$,
and $K(\pi)$ denotes the Kreweras complement of $\pi$ (see [8, Definition 2.4]
and [32]). Furthermore, for any subset $S\subset[k]$ we define
$m[S]:=m_{\mathrm{sc}}[\bm{z}_{S}]$ as the iterated divided difference of
$m_{\mathrm{sc}}$ evaluated in $\bm{z}_{S}:=\\{z_{i}:i\in S\\}$, and by
$m_{\circ}[\cdot]$ denote the free-cumulant transform of $m[\cdot]$ which is
uniquely defined implicitly by the relation

$m[S]=\sum_{\pi\in\mathrm{NC}(S)}\prod_{S^{\prime}\in\pi}m_{\circ}[S^{\prime}],\qquad\forall
S\subset[k],$ (3.3)

e.g. $m_{\circ}[i,j]=m[\set{i,j}]-m[\set{i}]m[\set{j}]$. Throughout the paper,
we will often use the fact that $m_{\mathrm{sc}}[\bm{z}_{S}]$ can be written
as follows

$m_{\mathrm{sc}}[\set{z_{i}:i\in
S}]=\int_{-2}^{2}\rho_{\mathrm{sc}}(x)\prod_{i\in
S}\frac{1}{(x-z_{i})}\mathrm{d}x.$ (3.4)

In order to formulate bounds on the deterministic approximation $M$ as well as
the local law bounds in a concise way, we introduce the following weighted
Schatten norms.

######  Definition 3.1 ($\ell$-weighted Schatten norms).

For $B\in{\mathbb{C}}^{N\times N}$, $\ell>0$ and $p\in[2,\infty]$, we define
the _$\ell$ -weighted $(2,p)$-Schatten norm_ as

${\left|\kern-1.29167pt\left|\kern-1.29167pt\left|B\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{p,\ell}:=\frac{\langle|B|^{2}\rangle^{1/2}}{\ell^{1/2}}+\frac{\langle|B|^{p}\rangle^{1/p}}{\ell^{1/p}}\quad\text{for}\quad
p<\infty\qquad\text{and}\qquad{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|B\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{\infty,\ell}:=\frac{\langle|B|^{2}\rangle^{1/2}}{\ell^{1/2}}+\|B\|\,.$
(3.5)

By elementary inequalities, we have

${\left|\kern-1.29167pt\left|\kern-1.29167pt\left|B\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{p,\ell}\lesssim{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|B\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{q,\ell}\qquad\text{and}\qquad{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|B\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{q,\ell}\lesssim\left[1+(N\ell)^{\frac{q-p}{pq}}\right]\,{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|B\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{p,\ell}\,,\quad\mbox{$2\leq
p\leq q\leq\infty$.}$ (3.6)

The next lemma gives a bound on the deterministic approximation $M_{[1,k]}$
with traceless observables; its proof is given in Appendix A. We will
henceforth follow the convention that the letter $B$ is used for generic
matrices, while $A$ is reserved for traceless ones.

######  Lemma 3.2 ($M$-bound).

Fix $k\geq 1$. Consider spectral parameters
$z_{1},...,z_{k}\in{\mathbb{C}}\setminus{\mathbb{R}}$ and traceless matrices
$A_{1},...,A_{k}\in{\mathbb{C}}^{N\times N}$. Denoting
$\ell:=\min_{j\in[k]}\eta_{j}\rho_{j}$ with $\eta_{j}:=|\mathrm{Im}\,z_{j}|$
and $\rho_{j}:=\pi^{-1}|\mathrm{Im}\,m(z_{j})|$, it holds that101010We point
out that the bound can be improved to
$\ell\prod_{i\in[k]}\left(\ell^{(1-k)/(2k)}\langle|A_{i}|^{3}\rangle^{1/k}\langle|A_{i}|^{2}\rangle^{(k-3)/(2k)}+\ell^{-1/k}\langle|A_{i}|^{k}\rangle^{1/k}\right)$
in case of odd $k\geq 3$, but we do not follow this improvement for brevity
and ease of notation.

$\left|\langle{M}_{[1,k]}(z_{1},A_{1},...,z_{k})A_{k}\rangle\right|\lesssim\mathbf{1}(k\geq
2)\
\ell\prod_{i\in[k]}{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{k,\ell}\,,$
(3.7)

where
${\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{k,\ell}$
was defined in (3.5).

It is easy to see that the bound (3.7) is optimal for $k$ even (cf. Footnote
10) and in case that $A_{i}=A$ for all $i\in[k]$ and $z_{i}\in\\{z,\bar{z}\\}$
for some $z\in{\mathbb{C}}\setminus{\mathbb{R}}$ in the _bulk_ , i.e.
$\mathrm{Re}\,z\in[-2+\delta,2-\delta]$ for some $\delta>0$.

###  3.2. Multi-resolvent local law

We now formulate our main technical result, an _averaged multi-resolvent local
law with Schatten norms_ , in Theorem 3.3; its proof is given in Section 4.
The corresponding _isotropic multi-resolvent local law_ will be formulated in
Theorem 6.2 later.

######  Theorem 3.3 (Averaged multi-resolvent local law with Schatten norms).

Let $W$ be a Wigner matrix satisfying Assumption 2.1, and fix
$k\in{\mathbb{N}}$. Consider spectral parameters
$z_{1},...,z_{k}\in{\mathbb{C}}\setminus{\mathbb{R}}$, the associated
resolvents $G_{j}=G(z_{j}):=(W-z_{j})^{-1}$, and traceless deterministic
matrices $A_{1},...,A_{k}\in{\mathbb{C}}^{N\times N}$. Denote

$\eta_{j}:=|\mathrm{Im}\,z_{j}|\,,\qquad\rho_{j}:=\pi^{-1}|\mathrm{Im}\,m(z_{j})|\,,\qquad\ell:=\min_{j\in[k]}\eta_{j}\rho_{j},$

and let ${M}_{[1,k]}$ be given in (3.2). Then, for any fixed $\epsilon>0$ and
recalling (3.5), we have

$\left|\langle{G}_{1}A_{1}...{G}_{k}A_{k}\rangle-\langle{M}_{[1,k]}A_{k}\rangle\right|\prec\frac{1}{N}\prod_{i\in[k]}{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{2k,\ell}\,,$
(3.8)

uniformly in spectral parameters satisfying $N\ell\geq N^{\epsilon}$ and
$\max_{j}|z_{j}|\leq N^{1/\epsilon}$.

Notice that, in the bulk regime, $\ell$ essentially agrees with
$\eta:=\min_{i\in[k]}|\mathrm{Im}\,z_{i}|$, since $\rho_{i}\sim 1$. However,
for (3.8) to be valid in any regime, the standard condition $N\eta\gg 1$ for
the local law in the bulk needs to be replaced by $N\ell\gg 1$. This ensures
that we are at the _mesoscopic scale_ , i.e. there are many eigenvalues in a
local window of size $\eta_{i}$ around each $\mathrm{Re}\,z_{i}$.

As already mentioned in the Introduction around (1.6), Theorem 3.3 unifies and
improves the previous local laws (in the bulk spectrum) with an error term
involving either the operator norm (see [10, Eq. (2.11a) in Theorem 2.5]) or
Hilbert-Schmidt norm (see [11, Theorem 2.2]). This follows by estimating (in
the relevant $\ell\lesssim 1$ regime)

${\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{2k,\ell}\lesssim\begin{cases}\dfrac{\|A_{i}\|}{\ell^{1/2}}\qquad&\text{for
\ \cite[cite]{[\@@bibref{Number}{multiG}{}{}, Eq.~{}(2.11a) in Theorem
2.5]}}\\\\[5.69054pt]
(N\ell)^{\frac{k-1}{2k}}\dfrac{\langle|A_{i}|^{2}\rangle^{1/2}}{\ell^{1/2}}\qquad&\text{for
\ \cite[cite]{[\@@bibref{Number}{A2}{}{}, Theorem~{}2.2]}}\end{cases}$ (3.9)

in (3.8) for $N\ell>1$ and every $i\in[k]$ by means of elementary
inequalities.

######  Remark 3.4 (Optimality).

The bound (3.8) is optimal in case that $A_{i}=A$ for all $i\in[k]$ and
$z_{i}\in\\{z,\bar{z}\\}$ for some $z\in{\mathbb{C}}\setminus{\mathbb{R}}$
with $|\mathrm{Im}\,z|\geq N^{-1+\epsilon}$ in the bulk. For $W$ being GUE,
this can easily be checked by spectral decomposition of the resolvents and
using Weingarten calculus [16] for the Haar-distributed eigenvectors. The so-
called “ladder diagram” gives the (first) Hilbert-Schmidt term in the estimate

$\left|\langle(GA)^{k}\rangle-\langle
M_{[1,k]}A\rangle\right|\prec\frac{{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{2k,\ell}^{k}}{N}\sim\frac{\langle|A|^{2}\rangle^{k/2}}{N\ell^{k/2}}+\frac{\langle|A|^{2k}\rangle^{1/2}}{N\ell^{1/2}}\,.$
(3.10)

Interestingly, in some regimes the second term in (3.10) (non ladder diagram)
gives the main contribution, defying the general belief that always the ladder
diagrams are the leading terms.

######  Remark 3.5 (Extensions).

In Theorem 3.3, each $G$ may also be replaced by a product of $G$’s or even
$|G|$’s (absolute value). We refrain from stating these results explicitly, as
they can be easily obtained from appropriate integral representations (see
(5.24) and (5.28) below). We formulate only the following example for $k=2$
and identical observables $B_{1}=B_{2}$ for illustration. Let
$B\in{\mathbb{C}}^{N\times N}$ be an arbitrary (i.e. not necessarily
traceless) deterministic matrix. Then, decomposing $B=\mathring{B}+\langle
B\rangle$, we have

$\begin{split}\langle G_{1}BG_{2}B\rangle&=\langle B\rangle^{2}\langle
G_{1}G_{2}\rangle+\langle B\rangle\langle
G_{1}G_{2}\mathring{B}\rangle+\langle B\rangle\langle
G_{2}G_{1}\mathring{B}\rangle+\langle
G_{1}\mathring{B}_{1}G_{2}\mathring{B}_{2}\rangle\\\ &=\frac{m_{1}m_{2}\langle
B\rangle^{2}}{1-m_{1}m_{2}}+m_{1}m_{2}\langle\mathring{B}^{2}\rangle+\mathcal{O}_{\prec}\left(\frac{|\langle
B\rangle|^{2}}{N\ell^{2}}+\frac{|\langle
B\rangle|\langle|\mathring{B}|^{2}\rangle^{1/2}}{N\ell^{3/2}}+\frac{\langle|\mathring{B}|^{2}\rangle}{N\ell}+\frac{\langle|\mathring{B}|^{4}\rangle^{1/2}}{N\ell^{1/2}}\right)\,.\end{split}$
(3.11)

The statement for different observables $B_{1},B_{2}$ is analogous.

###  3.3. Out-of-time-ordered correlators: Proof of Theorem 2.2

In order to prove Theorem 2.2, we distinguish the following three time
regimes,

$\mathrm{(i)}\ |t|<1\,,\qquad\mathrm{(ii)}\ 1\leq|t|\leq
N^{(1-\epsilon)/2}\,,\quad\text{and}\quad\mathrm{(iii)}\
|t|>N^{(1-\epsilon)/2}$ (3.12)

for some small fixed $\epsilon>0$ from Theorem 2.2. In the following, we focus
on the most complicated case (ii) and discuss the other two cases briefly at
the end of this section. Since the arguments in this section are fairly
standard, we will leave some irrelevant technical details to the reader.

For $1\leq|t|\leq N^{(1-\epsilon)/2}$, we employ the integral representation

$\mathrm{e}^{\pm\mathrm{i}tW}=\frac{1}{2\pi\mathrm{i}}\oint_{\Gamma}\mathrm{e}^{\pm\mathrm{i}tz}G(z)\mathrm{d}z$
(3.13)

with the contour

$\Gamma\equiv\Gamma_{t,R}:=\partial\big{(}[-R,R]\times\mathrm{i}[-|t|^{-1},|t|^{-1}]\big{)}$
(3.14)

parametrized counterclockwise and the parameter $R$ chosen as $R=N^{\kappa}$
for some arbitrarily small (but fixed) $\kappa>0$. In this way, we can write
the four-point part $\mathcal{F}_{A,B}(t)$ in (2.3) of the OTOC (2.2) as

$\mathcal{F}_{A,B}(t)=\left(\prod_{i\in[4]}\oint_{\Gamma}\frac{\mathrm{d}z_{i}}{2\pi\mathrm{i}}\right)\mathrm{e}^{\mathrm{i}t(z_{1}+z_{3}-z_{2}-z_{4})}\langle
G(z_{1})AG(z_{2})BG(z_{3})AG(z_{4})B\rangle\,.$ (3.15)

For the part, where _all_ $z_{i}$’s, $i\in[4]$, run on the horizontal parts of
the contour $\Gamma$, we can replace $\langle G(z_{1})A...G(z_{4})B\rangle$ by
$\langle M(z_{1},A,...,z_{4})B\rangle$ at the price of an error
$\mathcal{O}_{\prec}(\mathcal{E})$ with $\mathcal{E}=\mathcal{E}_{A,B}(t)$
from (2.7) by means of Theorem 3.3 for $k=4$. Here, we used $\ell\gtrsim
N^{-2\kappa}|t|^{-2}$, which follows by
$\rho(z_{i})\gtrsim|\mathrm{Im}\,z_{i}|(1+\mathrm{dist}(z_{i},[-2,2])^{2})^{-1}$,
and that $\kappa>0$ is arbitrarily small and hence $N^{\kappa}$ can be
absorbed into $\prec$.

If _all_ the $z_{i}$’s run on the vertical parts of $\Gamma$, we employ the
global law version of our main result, Proposition 4.1 below, to replace
$\langle G(z_{1})A...G(z_{4})B\rangle$ by $\langle
M(z_{1},A,...,z_{4})B\rangle$, now at the price of an error
$\mathcal{O}_{\prec}\big{(}\langle A^{8}\rangle^{1/4}\langle
B^{8}\rangle^{1/4}/N\big{)}$, which can easily be included in
$\mathcal{O}_{\prec}(\mathcal{E})$. We point out that Theorem 3.3 cannot be
used in this regime, since $\ell=0$ when one of the $z_{i}$’s crosses the real
axis.

In the remaining cases, where some of the $z_{i}$’s run on the horizontal
parts and others run on the vertical parts, we can no longer cleanly apply
either the local or global law (Theorem 3.3 and Proposition 4.1,
respectively). Instead, we treat this situation by expanding around the case
where all $z_{i}$ are on the horizontal parts. More precisely, in case that,
say, $z_{1}$ runs on the right vertical part of the contour $\Gamma$ and
$z_{2},....,z_{4}$ are fixed on the horizontal parts, we employ analyticity of
$G(z_{1})$ away from $[-3,3]$ (since $\|W\|\leq 2+\epsilon$ with very high
probability). This enables us to write $g(z_{1}):=\langle
G(z_{1})A...G(z_{4})B\rangle$ at $z_{1}=R+\mathrm{i}q$ with
$q\in[-|t|^{-1},|t|^{-1}]$, by Taylor expansion around
$\zeta:=R+\mathrm{i}|t|^{-1}$, as

$g(R+\mathrm{i}q)=\sum_{j=0}^{K}\frac{(q-|t|^{-1})^{j}}{2\pi\mathrm{i}}\oint_{C_{1/(2|t|)}}\frac{g(w)}{(w-\zeta)^{j+1}}\mathrm{d}w+\mathcal{O}\big{(}\widetilde{\mathcal{E}}\big{)}\,,\quad\widetilde{\mathcal{E}}:=\frac{\langle
A^{4}\rangle^{1/2}\langle B^{4}\rangle^{1/2}}{R^{K+1}|t|^{K-2}}$ (3.16)

for any $K\in{\mathbb{N}}$ to be chosen below, where $C_{1/(2|t|)}$ is the
circle of radius $1/(2|t|)$ centered around $\zeta$. In (3.16), we used that
the $(K+1)^{\rm th}$ derivative of $g(z_{1})$ on the vertical part of the
contour is (deterministically) bounded as

$\big{|}g^{(K+1)}(z_{1})\big{|}\lesssim\big{|}\big{\langle}\big{(}G(z_{1})\big{)}^{K+1}AG(z_{2})BG(z_{3})AG(z_{4})B\big{\rangle}\big{|}\lesssim\frac{\langle
A^{4}\rangle^{1/2}\langle B^{4}\rangle^{1/2}|t|^{3}}{R^{K+1}},$

since $\|G(z_{i})\|\leq|t|$ for $i\in\\{2,3,4\\}$. Using the representation
(3.16), we can now replace $g(w)$ by $h(w):=\langle M(w,A,...,z_{4})B\rangle$,
again at the expense of an error $\mathcal{O}_{\prec}(\mathcal{E})$ by means
of Theorem 3.3. As one can easily see that $h$ is also analytic around
$R+\mathrm{i}q$ and satisfies the same relation as $g$ in (3.16), we find that
$|g(R+\mathrm{i}q)-h(R+\mathrm{i}q)|\prec\mathcal{E}+\widetilde{\mathcal{E}}$.
Hence, in order to absorb $\widetilde{\mathcal{E}}$ into $\mathcal{E}$, it
remains to choose $K$ in (3.16) as $K:=\lceil\kappa^{-1}\rceil$.

Therefore, since integrating the $\mathcal{O}_{\prec}(\mathcal{E})$-bound in
(3.15) only adds some $N^{\kappa}$-factors from the length of $\Gamma$ (note
that $|\mathrm{e}^{\mathrm{i}t(z_{1}+z_{3}-z_{2}-z_{4})}|\lesssim 1$), which
can easily be absorbed into $\prec$, we find that

$\mathcal{F}_{A,B}(t)=\left(\prod_{i\in[4]}\oint_{\Gamma}\frac{\mathrm{d}z_{i}}{2\pi\mathrm{i}}\right)\mathrm{e}^{\mathrm{i}t(z_{1}+z_{3}-z_{2}-z_{4})}\langle
M(z_{1},A,z_{2},B,z_{3},A,z_{4})B\rangle+\mathcal{O}_{\prec}\big{(}\mathcal{E}\big{)}\,.$
(3.17)

Similarly, decomposing $A^{2}=(A^{2}-\langle A^{2}\rangle)+\langle
A^{2}\rangle$ and analogously for $B$, the two-point part
$\mathcal{D}_{A,B}(t)$ from (2.3) is given by, again for times $1\leq|t|\leq
N^{(1-\epsilon)/2}$,

$\mathcal{D}_{A,B}(t)=\langle A^{2}\rangle\langle
B^{2}\rangle+\left(\prod_{i\in[2]}\oint_{\Gamma}\frac{\mathrm{d}z_{i}}{2\pi\mathrm{i}}\right)\mathrm{e}^{\mathrm{i}t(z_{1}-z_{2})}\langle
M(z_{1},(A^{2}-\langle A^{2}\rangle),z_{2})(B^{2}-\langle
B^{2}\rangle)\rangle+\mathcal{O}_{\prec}(\mathcal{E})\,.$ (3.18)

In the other two regimes, $|t|<1$ and $|t|>N^{(1-\epsilon)/2}$, we follow the
above arguments for $1\leq|t|\leq N^{1/2-\epsilon}$, but replace the contour
$\Gamma$ from (3.14) by111111We choose the first contour only for $t\neq 0$.
If $t=0$, the lhs. of (2.6) carries no randomness, and we find
$\mathcal{C}_{A,B}(0)=\langle A^{2}B^{2}\rangle-\langle ABAB\rangle$.

$\\{z\in{\mathbb{C}}:\mathrm{dist}(z,[-2,2])=|t|^{-1}\\}\quad\text{and}\quad\partial\big{(}[-R,R]\times\mathrm{i}[-N^{(1-\epsilon)/2},N^{(1-\epsilon)/2}]\big{)}\,,$

respectively. This results in error terms in identities analogous to
(3.17)–(3.18), that can easily be seen to be bounded by
$\mathcal{O}_{\prec}(\mathcal{E})$, just as before, by application of
Proposition 4.1 and Theorem 3.3, respectively.

It remains to explicitly compute the deterministic terms in (3.17)–(3.18),
which, by using the formulas (3.1)–(3.4), straightforwardly results in the
expression given in (2.6). This concludes the proof of Theorem 2.2. ∎

##  4\. Zigzag strategy: Proof of Theorem 3.3

In this section we prove our main technical result from Section 3, the multi-
resolvent local law in Theorem 3.3. Its proof is conducted via the
_characteristic flow method_ [2, 6, 29, 1, 33, 3] followed by a _Green
function comparison_ (GFT) argument. A combination of these tools, which we
coin the _Zigzag strategy_ , has first been used in [12, 14, 15]. It consists
of the following three steps:

  * 1.

Global law. Proof of a multi-resolvent _global law_ , i.e. for spectral
parameters “far away” from the spectrum,
$\min_{j}\mathrm{dist}(z_{j},[-2,2])\geq\delta$ for some small $\delta>0$ (see
Proposition 4.1 below).

  * 2.

Characteristic flow. Propagate the global law to a _local law_ by considering
the evolution of the Wigner matrix $W$ along the Ornstein-Uhlenbeck flow

$\mathrm{d}W_{t}=-\frac{1}{2}W_{t}\mathrm{d}t+\frac{\mathrm{d}B_{t}}{\sqrt{N}},\qquad
W_{0}=W\,,$ (4.1)

with $B_{t}$ a standard real symmetric or complex Hermitian Brownian motion,
thereby introducing an order one Gaussian component (see Proposition 4.3). The
spectral parameters of the resolvents evolve from the global regime to the
local regime according to the _characteristic (semicircular) flow_

$\partial_{t}z_{i,t}=-m(z_{i,t})-\frac{z_{i,t}}{2}\,.$ (4.2)

The simultaneous effect of these two evolutions is a key cancellation of two
large terms.

  * 3.

Green function comparison. Remove the Gaussian component by a Green function
comparison (GFT) argument (see Proposition 4.4).

As the first step, we have the following global law, the proof of which is
completely analogous to the proofs presented in [10, Appendix B], [11,
Appendix A], and [14, Appendix A], and so omitted. Proposition 4.1 is stated
for a general deterministic matrix $B$ since the traceless condition plays no
role in this case.

######  Proposition 4.1 (Step 1: Global law).

Let $W$ be a Wigner matrix satisfying Assumption 2.1, and fix
$k\in{\mathbb{N}}$ and $\delta>0$. Consider spectral parameters
$z_{1},...,z_{k}\in{\mathbb{C}}\setminus{\mathbb{R}}$, the associated
resolvents $G_{j}=G(z_{j}):=(W-z_{j})^{-1}$, and deterministic matrices
$B_{1},...,B_{k}\in{\mathbb{C}}^{N\times N}$. Then, uniformly in spectral
parameters satisfying $d:=\min_{j\in[k]}\mathrm{dist}(z_{j},[-2,2])\geq\delta$
and deterministic matrices $B_{1},...,B_{k}$, it holds that

$\left|\langle
G_{1}B_{1}...G_{k}B_{k}\rangle-\langle{M}_{[1,k]}B_{k}\rangle\right|\prec\frac{1}{N}\prod_{i\in[k]}{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|B_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{2k,d}\,.$
(4.3)

Next, using Proposition 4.1 as an input, we derive Theorem 3.3 for Wigner
matrices which have an order one Gaussian component, as formulated in
Proposition 4.3. For this purpose we consider the evolution of the Wigner
matrix $W$ along the Ornstein-Uhlenbeck flow (4.1) and define its resolvent
$G_{t}(z):=(W_{t}-z)^{-1}$ with $z_{i}\in{\mathbb{C}}\setminus{\mathbb{R}}$.
Even if not stated explicitly we will always consider this flow only for short
times, i.e. for $0\leq t\leq T$, where the maximal time $T$ is smaller than
some ${\gamma}>0$. Note that the first two moments of $W_{t}$ are preserved
along the flow (4.1), and hence the self-consistent density of states of
$W_{t}$ is unchanged; it remains the standard semicircle law. We now want to
compute the deterministic approximation to an alternating product of
resolvents and deterministic matrices $A_{1},A_{2},...$ with trace zero

${G}_{t}(z_{1})A_{1}{G}_{t}(z_{2})A_{2}{G}_{t}(z_{3})A_{3}\dots,$ (4.4)

and have a very precise estimate of the error term.

In fact, we will let the spectral parameters evolve with time with a carefully
chosen equation that will conveniently cancel some error terms in the time
evolution of (4.4). The corresponding equation will be the characteristic
equation for the semicircular flow, i.e. given by the first order ODE (4.2)
(see [14, Figure 1] for an illustration of this flow). We remark that, along
the characteristics we have

$\partial_{t}m(z_{i,t})=-\partial_{z}m(z_{i,t})\left(m(z_{i,t})+\frac{z_{i,t}}{2}\right)=-\partial_{z}m(z_{i,t})\left(-\frac{1}{2m(z_{i,t})}+\frac{m(z_{i,t})}{2}\right)=\frac{m(z_{i,t})}{2},$
(4.5)

where in the last two equalities we used the defining equation
$m(z)^{2}+zm(z)+1=0$ of the Stieltjes transform of the semicircular law. In
particular, this implies that $\rho_{i,s}\sim\rho_{i,t}$ for any $0\leq s\leq
t$, where we denoted $\rho_{i,t}:=\pi^{-1}|\mathrm{Im}\,m(z_{i,t})|$. In
contrast to that, the behavior of the $\eta_{i,t}:=|\mathrm{Im}\,z_{i,t}|$
depends on the regime: in the bulk $\eta_{i,t}$ decreases linearly in time
with a speed of order one, at the edge the decay is still linear but with a
speed depending on the size of the local density of states. By standard ODE
theory we obtain the following lemma:

######  Lemma 4.2 (see Lemma 3.2 in [14]).

Fix an $N$–independent $\gamma>0$, fix $0<T<\gamma\normalcolor$, and pick
$z\in{\mathbb{C}}\setminus{\mathbb{R}}$. Then there exists an initial
condition $z_{0}$ such that the solution $z_{t}$ of (4.2) with this initial
condition $z_{0}$ satisfies $z_{T}=z$. Furthermore, there exists a constant
$C>0$ such that $\mathrm{dist}(z_{0},[-2,2])\geq C\normalcolor T$.

The spectral parameters evolving by (4.1) will have the property that

${G}_{t}(z_{1,t})A_{1}\dots
A_{k-1}{G}_{t}(z_{k,t})-M_{[1,k],t}\approx{G}_{0}(z_{1,0})A_{1}\dots
A_{k-1}{G}_{0}(z_{k,0})-M_{[1,k],0},$ (4.6)

with $M_{[1,k],t}:=M(z_{1,t},A_{1},\dots,A_{k-1},z_{k,t})$, for any $0\leq
t\leq T$. Note that the deterministic approximation $M_{[1,k],t}$ depends on
time only through the time dependence of the spectral parameters, the
deterministic approximation of (4.4) with fixed spectral parameters does not
depend on time, i.e. it is unchanged along the whole flow (4.1).

######  Proposition 4.3 (Step 2: Characteristic flow).

Fix any $\epsilon,\gamma>0$, a time $0\leq T\leq\gamma$, and
$K\in{\mathbb{N}}$. Consider
$z_{1,0},\dots,z_{K,0}\in{\mathbb{C}}\setminus{\mathbb{R}}$ as initial
conditions to the solution $z_{j,t}$ of (4.2) for $0\leq t\leq T$ and define
$G_{j,t}:=G_{t}(z_{j,t})$ as well as $\eta_{j,t}:=|\mathrm{Im}\,z_{j,t}|$ and
$\rho_{j,t}:=\pi^{-1}|\mathrm{Im}\,m(z_{j,t})|$.

Let $k\leq K$, define $\ell_{t}:=\min_{j\in[k]}\eta_{j,t}\rho_{j,t}$ and
recall (3.5). Then, assuming that

$\left|\langle{G}_{1,0}A_{1}...{G}_{k,0}A_{k}\rangle-\langle{M}_{[1,k],0}A_{k}\rangle\right|\prec\frac{1}{N}\prod_{i\in[k]}{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{2k,\ell_{0}}\,,$
(4.7)

holds uniformly for any $k\leq K$, any choice of deterministic traceless
$A_{1},...,A_{k}$ and any choice of $z_{i,0}$’s such that $N\ell_{0}\geq
N^{\epsilon}$ and $\max_{j\in[k]}|z_{j,0}|\leq N^{1/\epsilon}$, then we have

$\left|\langle{G}_{1,T}A_{1}...{G}_{k,T}A_{k}\rangle-\langle{M}_{[1,k],T}A_{k}\rangle\right|\prec\frac{1}{N}\prod_{i\in[k]}{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{2k,\ell_{T}}$
(4.8)

for any $k\leq K$, again uniformly in traceless matrices $A_{i}$ and in
spectral parameters satisfying $N\ell_{T}\geq N^{\epsilon}$ and
$\max_{j\in[k]}|z_{j,T}|\leq N^{1/\epsilon}$.

The proof of Proposition 4.3 is given in Section 5. As the third and final
step, we show that the additional Gaussian component introduced in Proposition
4.3 can be removed using a Green function comparison (GFT) argument at the
price of a negligible error. The proof of this proposition is presented in
Section 6.

######  Proposition 4.4 (Step 3: Green function comparison).

Let $H^{(\bm{v})}$ and $H^{(\bm{w})}$ be two $N\times N$ Wigner matrices with
matrix elements given by the random variables $v_{ab}$ and $w_{ab}$,
respectively, both satisfying Assumption 2.1 and having matching moments up to
third order,121212This condition can easily be relaxed to being matching up to
an error of size $N^{-2}$ as done, e.g., in [21, Theorem 16.1]. i.e.

${\mathbb{E}}\bar{v}_{ab}^{u}v_{ab}^{s-u}={\mathbb{E}}\bar{w}_{ab}^{u}w_{ab}^{s-u}\,,\quad
s\in\\{0,1,2,3\\}\,,\quad u\in\\{0,...,s\\}\,.$ (4.9)

Fix $K\in{\mathbb{N}}$ and consider spectral parameters
$z_{1},...,z_{K}\in{\mathbb{C}}\setminus{\mathbb{R}}$ satisfying
$\min_{j}N\eta_{j}\rho_{j}\geq N^{\epsilon}$ and $\max_{j}|z_{j}|\leq
N^{1/\epsilon}$ for some $\epsilon>0$ and the associated resolvents
$G_{j}^{(\\#)}=G^{(\\#)}(z_{j}):=(H^{(\\#)}-z_{j})^{-1}$ with
$\\#=\bm{v},\bm{w}$. Pick traceless matrices
$A_{1},...,A_{K}\in{\mathbb{C}}^{N\times N}$.

Assume that, for $H^{(\bm{v})}$, we have the following bounds (writing
$G_{j}\equiv G_{j}^{(\\#)}$ for brevity): For any $k\leq K$, consider any
subset of cardinality $k$ of the $K$ spectral parameters, and similarly,
consider any subset of cardinality $k$ of the deterministic matrices.
Relabeling both of them by $[k]$, setting
$\ell:=\min_{j\in[k]}\eta_{j}\rho_{j}$ and recalling (3.5), we have that

$\left|\langle{G}_{1}A_{1}...{G}_{k}A_{k}\rangle-\langle{M}_{[1,k]}A_{k}\rangle\right|\prec\frac{1}{N}\prod_{i\in[k]}{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{2k,\ell}$
(4.10)

uniformly in all choices of subsets of $z$’s and $A$’s.

Then, (4.10) also holds for the ensemble $H^{(\bm{w})}$, uniformly all choices
of subsets of $z$’s and $A$’s.

We are now ready to finally conclude the proof of Theorem 3.3. Fix $T>0$, and
fix $z_{1},\dots,z_{k}\in{\mathbb{C}}\setminus{\mathbb{R}}$ such that
$\min_{i}N\eta_{i}\rho_{i}\geq N^{\epsilon}$, and let $z_{i,0}$ be the initial
conditions of the characteristics (4.2) chosen so that $z_{i,T}=z_{i}$ (this
is possible thanks to Lemma 4.2). Then, the assumption (4.7) of Proposition
4.3 is satisfied for those $z_{i,0}$ by Proposition 4.1 with $\delta=CT$,
since $d\gtrsim\ell_{0}$ and where $C>0$ is the constant from Lemma 4.2. We
can thus use Proposition 4.3 to show that (4.8) holds. Finally, the Gaussian
component added in Proposition 4.3 is removed using Proposition 4.4 with the
aid of a complex version (see [14, Lemma A.2]) of the standard moment-matching
lemma [21, Lemma 16.2]. ∎

##  5\. Characteristic flow: Proof of Proposition 4.3

In this section, we give the proof of Proposition 4.3. The argument is divided
into three parts.

  * (i)

We begin by introducing several deterministic and stochastic control
quantities, which play a fundamental role throughout the rest of this paper.
The stochastic control quantities (some normalized differences between a
resolvent chain and the corresponding $M$-term, see (5.12)–(5.13) below)
satisfy a system of _master inequalities_ (see Propositions 5.2–5.3).

  * (ii)

Taking these master inequalities together with Proposition 4.1 as inputs, we
will prove Proposition 4.3 in Section 5.1, thus concluding Step 2 of the
argument in Section 4.

  * (iii)

Afterwards, based on Proposition 4.1 and several relations among the
_deterministic_ control quantities (see Lemma 5.1), we will prove the master
inequalities in Section 5.2.

To keep the presentation simpler, within this section we assume that
${\mathbb{E}}\chi^{2}_{\mathrm{od}}=0$ and that
${\mathbb{E}}\chi^{2}_{\mathrm{d}}=1$, i.e. we consider only complex Wigner
matrices. The modifications for the general case are analogous to [14, Section
4], and so omitted.

We recall our choice of the characteristics

$\partial_{t}z_{i,t}=-m(z_{i,t})-\frac{z_{i,t}}{2},$ (5.1)

note that $t\mapsto\eta_{i,t}$ is decreasing and that
$\rho_{i,t}\sim\rho_{i,s}$ for any $0\leq s\leq t$ (see (4.5) and the
paragraph below it for more details). Additionally, we record the following,
trivially checkable, simple integration rule which will be used many times in
this section:

$\int_{0}^{t}\frac{1}{\eta_{s}^{\alpha}}\mathrm{d}s\lesssim\frac{\log
N}{\eta_{t}^{\alpha-2}\ell_{t}}\quad\text{with}\quad\eta_{s}:=\min_{i}\eta_{i,s}\,,\quad\ell_{s}:=\min_{i}\eta_{i,s}\rho_{i,s}\,.$
(5.2)

Along the characteristics, using the short–hand notation
$G_{i,t}=(W_{t}-z_{i,t})^{-1}$, with $W_{t}$ being the solution of (4.1), by
Itô’s formula, we have131313We point out that (5.3) holds for any matrix
$A_{i}\in{\mathbb{C}}^{N\times N}$, i.e. we did not use that the $A_{i}$’s are
traceless.

$\begin{split}\mathrm{d}\langle(G_{[1,k],t}-M_{{[1,k]},t})A\rangle&=\frac{1}{\sqrt{N}}\sum_{a,b=1}^{N}\partial_{ab}\langle{G}_{[1,k],t}A\rangle\mathrm{d}B_{ab,t}+\frac{k}{2}\langle{G}_{[1,k],t}A\rangle\mathrm{d}t+\sum_{i,j=1\atop
i<j}^{k}\langle{G}_{[i,j],t}\rangle\langle{G}_{[j,i],t}\rangle\mathrm{d}t\\\
&\quad+\sum_{i=1}^{k}\langle
G_{i,t}-m_{i,t}\rangle\langle{G}^{(i)}_{[1,k],t}A\rangle\mathrm{d}t-\partial_{t}\langle
M_{[1,k],t}A\rangle\mathrm{d}t,\end{split}$ (5.3)

where $\partial_{ab}$ denotes the directional derivative $\partial_{w_{ab}}$.
We also set

${G}_{[i,j],t}:=\begin{cases}G_{i,t}A_{i}\dots A_{j-1}G_{j,t}&\mathrm{if}\quad
i<j\\\ G_{i,t}&\mathrm{if}\quad i=j\\\ G_{i,t}A_{i,t}\dots
G_{k,t}A_{k}G_{1,t}A_{1}\dots A_{j-1}G_{j,t}&\mathrm{if}\quad i>j,\end{cases}$
(5.4)

and defined ${G}^{(l)}_{[i,j],t}$ as $G_{[i,j],t}$ but with the $l$–th
resolvent substituted by $G_{l,t}^{2}$. The last definition for $i>j$ reflects
the cyclicity of the trace, since ${G}_{[i,j],t}$ will be needed in a tracial
situation. For any $i,j\in[k]$, we denote the deterministic approximation of
$G_{[i,j],t}$ by $M_{[i,j],t}$, with $M_{[i,j],t}$ being defined as in (3.2)
with $[1,k]$ replaced by $[i,j]$ if $i<j$ and by $[1,j]\cup[i,k]$ if $i>j$.
Deterministic control quantities: mean and standard deviation size. We now
introduce two deterministic control quantities that measure the size of the
mean of long chains $\langle G_{1}A_{1}...G_{k}A_{k}\rangle$ and their
standard deviation, respectively, in terms of the $\ell$-weighted Schatten
norms of $A_{1},...,A_{k}$ (recall Definition 3.1) and the spectral parameters
of the resolvents $G_{i}=G(z_{i})$. These are given by

$\begin{split}{\mathfrak{m}}_{k}(\ell;\bm{A}_{J_{k}}):=\mathbf{1}(k\geq
2)\prod_{\alpha\in
J_{k}}{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{\alpha}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{k,\ell}\qquad\text{and}\qquad{\mathfrak{s}}_{k}(\ell;\bm{A}_{J_{k}}):=\prod_{\alpha\in
J_{k}}{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{\alpha}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{2k,\ell}\,,\end{split}$
(5.5)

and will be called the _mean size_ and the _standard deviation size_ ,
respectively. They are functions of a positive number $\ell$ (usually given by
$\ell=\min_{i}\eta_{i}\rho_{i}$) and a multiset
$\bm{A}_{J_{k}}=(A_{\alpha})_{\alpha\in J_{k}}$ of deterministic matrices of
cardinality $|J_{k}|=k$.

In the following, we shall frequently drop the symbol $\bm{A}$ from the
definitions in (5.5), i.e. write ${\mathfrak{m}}_{k}(\ell;J_{k})$ and
${\mathfrak{s}}_{k}(\ell;J_{k})$ instead of
${\mathfrak{m}}_{k}(\ell;\bm{A}_{J_{k}})$ and
${\mathfrak{s}}_{k}(\ell;\bm{A}_{J_{k}})$, respectively. Moreover, for time
dependent spectral parameters, we will also use the notation
${\mathfrak{m}}_{k}(t)={\mathfrak{m}}_{k}(t;J_{k})={\mathfrak{m}}_{k}(\ell_{t};\bm{A}_{J_{k}})$,
with $\ell_{t}=\min_{i}\eta_{i,t}\rho_{i,t}$, and a similar notation for the
standard deviation size. In particular, we may often omit the dependence on
the deterministic matrix $A$. (More generally, we will omit every argument of
(5.5), whenever it is clear what they are and it does not lead to confusion.)
For example, for $i<j$, we have, with $\eta_{t}:=\min_{i}\eta_{i,t}$,

$|\langle
M_{[i,j],t}A_{j}\rangle|\lesssim\ell_{t}\,\mathfrak{m}_{j-i+1}(t),\qquad\quad|\langle
M_{[i,j],t}\rangle|\lesssim(\ell_{t}/\eta_{t})\mathfrak{m}_{j-i}(t).$ (5.6)

The first bound follows from (3.7); the second bound in (5.6) can easily be
obtained by arguments entirely analogous to the ones leading to Lemma A.1 in
Appendix A. The bounds (5.6) justify the terminology that $\mathfrak{m}_{k}$
is the mean size of a resolvent chain $\langle
G_{1}A_{1}...G_{k}A_{k}\rangle$, which is well approximated by $\langle
M_{[1,k]}A_{k}\rangle$. The additional $\ell_{t}$ factor was introduced to
simplify formulas later.

In the remainder of this section we will often use the following lemma about
some relations among the sizes $\mathfrak{m}_{k}(t)$, ${\mathfrak{s}}(t)$ from
(5.5). The proofs are immediate consequences of (3.6) and
$\ell_{s}\gtrsim\ell_{t}$ for $s\leq t$ (recall the discussion below (4.5))
and hence omitted.

######  Lemma 5.1 ($\mathfrak{m}/{\mathfrak{s}}$-relations).

Let $k\geq 1$ and consider time-dependent spectral parameters
$z_{1,t},...,z_{k,t}$ and a multiset of traceless matrices $\bm{A}_{J_{k}}$ as
above. Then, the mean size $\mathfrak{m}_{k}$ and the standard variation size
${\mathfrak{s}}_{k}$ satisfy the following inequalities.

  * (i)

Super-multiplicativity, i.e. for any $1\leq j\leq k-1$ we have

$\mathfrak{m}_{j}(t;J_{j})\,\mathfrak{m}_{k-j}(t;J_{k-j})\lesssim\,\mathfrak{m}_{k}(t;J_{k})\qquad\text{and}\qquad{\mathfrak{s}}_{j}(t;J_{j})\,{\mathfrak{s}}_{k-j}(t;J_{k-j})\lesssim\,{\mathfrak{s}}_{k}(t;J_{k})$
(5.7)

for all disjoint decompositions $J_{k}=J_{j}\dot{\cup}J_{k-j}$.

  * (ii)

The mean size can be upper bounded by the standard deviation size as

${\mathfrak{m}}_{k}(t;J_{k})\lesssim\,{\mathfrak{s}}_{k}(t;J_{k})\qquad\text{and}\qquad\sqrt{{\mathfrak{m}}_{2k}(t;J_{k}\cup
J_{k})}\lesssim{\mathfrak{s}}_{k}(t;J_{k})\,,$ (5.8)

where $J_{k}\cup J_{k}$ denotes the union of the multiset $J_{k}$ with itself.

  * (iii)

The standard deviation size satisfies the doubling inequality (recall
$\ell_{t}=\min_{j\in[k]}\eta_{j,t}\rho_{j,t}$)

${\mathfrak{s}}_{2k}(t;J_{k}\cup
J_{k})\lesssim[1+\sqrt{N\ell_{t}}]\big{(}{\mathfrak{s}}_{k}(t;J_{k})\big{)}^{2}\,.$
(5.9)

  * (iv)

Monotonicity in time: for any $0\leq s\leq t$ we have

$\mathfrak{m}_{k}(s;J_{k})\,\ell_{s}^{\alpha}\lesssim\mathfrak{m}_{k}(t;J_{k})\,\ell_{t}^{\alpha}\qquad{\mathfrak{s}}_{k}(s;J_{k})\,\ell_{s}^{\beta}\lesssim{\mathfrak{s}}_{k}(t;J_{k})\,\ell_{t}^{\beta}$
(5.10)

for all $\alpha\in[0,1]$ and $\beta\in[0,1/2]$.

Stochastic control quantities: Master and reduction inequalities. Using the
notation from (5.5), the goal is thus to prove that

$\langle G_{[1,k],T}A_{k}\rangle-\langle M_{[1,k],T}A_{k}\rangle=\langle
G_{[1,k],0}A_{k}\rangle-\langle
M_{[1,k],0}A_{k}\rangle+\mathcal{O}_{\prec}\Big{(}\frac{{\mathfrak{s}}_{k}(T)}{N}\Big{)}\,,$
(5.11)

uniformly in the spectrum and uniformly in traceless deterministic matrices
$A_{i}$, for some fixed $T\lesssim 1$. We may henceforth assume that all the
$A_{i}$’s are Hermitian; the general case follows by multilinearity.

For the purpose of proving (5.11), recall the notation
$\ell_{t}:=\min_{i}\eta_{i,t}\rho_{i,t}$ from (5.2) and define

$\Phi_{k}(t)=\Phi_{k}(t;\bm{z}_{[1,k]};\bm{A}_{[1,k]}):=\frac{N}{{\mathfrak{s}}_{k}(\ell_{t};\bm{A}_{[1,k]})}\big{|}\langle(G_{[1,k],t}-M_{[1,k],t})A_{k}\rangle\big{|}\quad\text{for
all}\quad k\geq 1\,.$ (5.12)

Here, $\bm{A}_{[1,k]}=(A_{1},\dots,A_{k})$ and
$\bm{z}_{[1,k]}=(z_{1},...,z_{k})$ with $z_{i}=z_{i,0}$ (initial condition)
denote multisets of deterministic matrices and spectral parameters,
respectively. We now briefly comment on the definition (5.12). We chose the
pre–factor in the definition of $\Phi_{k}(t)$ so that eventually $\Phi_{k}(t)$
will be of order one with high probability (cf. (5.11)). However, we will not
be able to prove this directly, we first prove that
$\Phi_{k}(t)\prec\sqrt{N\ell_{t}}$ and then, using this bound as an input, we
prove the desired bound $\Phi_{k}(t)\prec 1$. To implement this technically we
introduce another quantity

$\Psi_{k}(t)=\Psi_{k}(t;\bm{z}_{[1,k]};\bm{A}_{[1,k]}):=\frac{N}{{\mathfrak{s}}_{k}(\ell_{t};\bm{A}_{[1,k]})\sqrt{N\ell_{t}}}\big{|}\langle(G_{[1,k],t}-M_{[1,k],t})A_{k}\rangle\big{|}\quad\text{for
all}\quad k\geq 1\,,$ (5.13)

which we will show to be bounded by one (note that $\Psi_{k}(t)\prec 1$
implies $\Phi_{k}(t)\prec\sqrt{N\ell_{t}}$), and then show that
$\Psi_{k}(t)\prec 1$ in fact implies $\Phi_{k}(t)\prec 1$.

We start considering $\Psi_{k}(t)$; note that by (4.7) it follows

$\Psi_{k}(0)\lesssim\Phi_{k}(0)\prec 1,$ (5.14)

for any $k\geq 1$. For this purpose we will derive a series of master
inequalities for these quantities with the following structure. We assume that

$\Psi_{k}(t)\prec\psi_{k}$ (5.15)

holds for some deterministic control parameter $\psi_{k}\geq 1$, uniformly in
deterministic matrices $A_{i}$, times $0\leq t\leq T$ and in spectral
parameters with $\ell_{t}\geq N^{-1+\epsilon}$ for some small fixed
$\epsilon>0$ (we stress that the $\psi_{k}$’s depend neither on time, nor on
the spectral parameters $z_{i,t}$, nor on the deterministic matrix $A$).
Starting from (5.15) we derive an improved upper bound for $\Psi_{k}(t)$ and
show that, by _iterating_ these inequalities, we indeed obtain the desired
$\Psi_{k}(t)\prec 1$. The main inputs to prove this fact are the _master
inequalities_ in (5.16), which informally states that if
$\Psi_{l}(t)\prec\psi_{l}$ for any $l=1,\dots,k+\mathbf{1}(k\,\mathrm{odd})$,
then this actually implies the improved bound (5.16).

######  Proposition 5.2 (Master inequalities).

Fix $k\in{\mathbb{N}}$. Assume that $\Psi_{l}(t)\prec\psi_{l}$ for some
deterministic control parameters $\psi_{l}$, for any $1\leq l\leq
k+\bm{1}(k\,\mathrm{odd})$, uniformly in $t\in[0,T]$. Then

$\Psi_{k}(t)\prec
1+\sum_{j=1}^{k-2}\psi_{j}+\frac{1}{(N\ell_{T})^{1/4}}\sum_{j=1}^{k}\widetilde{\psi}_{j}\widetilde{\psi}_{k-j},$
(5.16)

uniformly in $t\in[0,T]$. Here we set
$\widetilde{\psi}_{l}:=\psi_{l}+\bm{1}(l\,\mathrm{odd})\sqrt{\psi_{l-1}\psi_{l+1}}$
and $\psi_{0}:=1$.

Using (5.16), in the next section we show that $\Psi_{k}(t)\prec 1$ by an
iterative procedure. Then, to conclude $\Phi_{k}(t)\prec 1$, we rely on the
following proposition, which will eventually prove Proposition 4.3.

######  Proposition 5.3.

Fix $k\in{\mathbb{N}}$, assume that $\Phi_{l}(t)\prec 1$, for $l\leq k-2$
uniformly in $t\in[0,T]$, and that $\Phi_{l}(t)\prec\sqrt{N\ell_{t}}$ for
$1\leq l\leq 2k$ and $t\in[0,T]$. Then

$\Phi_{k}(t)\prec 1,$ (5.17)

uniformly in $t\in[0,T]$.

###  5.1. Closing the hierarchy: Proof of Proposition 4.3

To show that Proposition 5.2 in fact implies $\Psi_{k}(t)\prec 1$ we rely on
the following procedure, which we refer to as _iteration_ (see, e.g., [13,
Lemma 4.1]).

######  Lemma 5.4 (Iteration).

Fix $k\in{\mathbb{N}}$, $T>0$, and $N$-independent constants
$\epsilon,\delta>0$, $\alpha\in(0,1)$ and $D>0$. Let $X$ be a random variable
depending on $k$ time dependent spectral parameters $z_{1,t},...,z_{k,t}$,
$t\in[0,T]$, and recall that $\ell_{t}=\min_{j\in[k]}\eta_{j,t}\rho_{j,t}$.
Assume that the a-priori bound $X\prec N^{D}$ holds uniformly in $t\in[0,T]$,
$\ell_{t}\geq N^{-1+\epsilon}$. Suppose that there is a deterministic quantity
$x$ (may depend on $\ell_{T}$ and $N$) such that for any fixed
$l\in{\mathbb{N}}$ the fact that $X\prec x$ uniformly for $t\in[0,T]$ and
$\ell_{t}\geq N^{-1+l\epsilon}$ implies141414Here the scalar $A>0$ should not
be confused with the matrices $A_{i}$ which appear throughout the proof.

$X\prec A+\frac{x}{B}+x^{1-\alpha}C^{\alpha},$ (5.18)

uniformly for $t\in[0,T]$ and $\ell_{t}\geq N^{-1+(l+l^{\prime})\epsilon}$,
for some constants151515The constants $A,B,C$ may depend on $N$, while $l$ and
$l^{\prime}$ are independent of $N$. $l^{\prime}\in{\mathbb{N}}$,
$B\geq\delta>0$, $A,C>0$. Then, iterating (5.18) finitely many times, we
obtain

$X\prec A+C,$

uniformly for $t\in[0,T]$ and $\ell_{t}\geq N^{-1+(1+Kl^{\prime})\epsilon}$,
for some $K=K(\alpha,D,\delta)$.

We are now ready to present the proof of Proposition 4.3.

###### Proof of Proposition 4.3.

The proof of this proposition is divided into two steps: we first prove that
$\Psi_{k}(t)\prec 1$ for any $k$ and then use this information as an input to
conclude that $\Phi_{k}(t)\prec 1$ for any $k$. These two bounds are quite
different in spirit. The first one is an a-priori bound so the iterative
procedure behind its proof must be self-improving. This is reflected in the
triangular structure of (5.16): the right hand side contains quantities with
index at most $k$ (or $k+1$ when $k$ is odd) and terms with the highest index
come with a small prefactor (recall that $N\ell_{T}$ is large). This makes the
system (5.16) closable. The iteration is done essentially for each fixed $k$
and then we use an induction on $k$. Owing to the parity issue in the
definition of $\widetilde{\psi}$, we use a step-two induction, but this is
just a small technicality. The second bound $\Phi_{k}(t)\prec 1$ is quite
different, since its proof relies on the a-priori bound obtained in the first
step. The key point is that in order to prove $\Phi_{k}(t)\prec 1$ for some
fixed $k$, we need to know the a-priori bound
$\Phi_{l}(t)\prec\sqrt{N\ell_{t}}$ for any $l\leq 2k$, i.e. without the
a-priori bound the system of inequalities behind the proof of
$\Phi_{k}(t)\prec 1$ would not be closable. This explains why we need to
proceed in two stages.

We now start we the proof of $\Psi_{k}(t)\prec 1$. We first prove this for
$k=1,2$ and then, using an inductive argument, we show that the same bound
holds for any $k\geq 3$. By (5.16) for $k=1,2$ we have

$\Psi_{1}(t)\prec
1+\frac{\psi_{1}+\sqrt{\psi_{2}}}{(N\ell_{T})^{1/4}},\qquad\quad\Psi_{2}(t)\prec
1+\frac{\psi_{1}^{2}+\psi_{2}}{(N\ell_{T})^{1/4}}$ (5.19)

Using iteration we then obtain (all estimates are uniform in $t\in[0,T]$)

$\Psi_{1}(t)\prec
1+\frac{\sqrt{\psi_{2}}}{(N\ell_{T})^{1/4}},\qquad\quad\Psi_{2}(t)\prec
1+\frac{\psi_{1}^{2}}{(N\ell_{T})^{1/4}}.$ (5.20)

Plugging the first inequality into the second one and using iteration, this
immediately gives that $\Psi_{1}(t)+\Psi_{2}(t)\prec 1$.

Next, we proceed with the induction step. Fix an even $k\in{\mathbb{N}}$, and
assume that $\Psi_{l}(t)\prec 1$, for $l\leq k-2$, then by (5.16) we have

$\Psi_{k-1}(t)\prec
1+\frac{\psi_{k-1}+\sqrt{\psi_{k}}}{(N\ell_{T})^{1/4}},\qquad\quad\Psi_{2}(t)\prec
1+\psi_{k-1}+\frac{\psi_{k}+\sqrt{\psi_{k}}}{(N\ell_{T})^{1/4}}$ (5.21)

By iteration, we have

$\Psi_{k-1}(t)\prec
1+\frac{\sqrt{\psi_{k}}}{(N\ell_{T})^{1/4}},\qquad\quad\Psi_{2}(t)\prec
1+\psi_{k-1},$ (5.22)

which concludes $\Psi_{l}(t)\prec 1$, for $l\leq k$, by plugging the first
inequality into the second one and using iteration once again.

Finally, to conclude $\Phi_{k}(t)\prec 1$, we proceed by a step-two induction
on $k$. For $k=1,2$, the assumption $\Phi_{l}(t)\prec\sqrt{N\ell_{t}}$, for
$l\leq 4$, of Proposition 5.3 is satisfied, and so we have
$\Phi_{1}(t)+\Phi_{2}(t)\prec 1$. Then we proceed with the induction step,
i.e. for a fixed even $k\in{\mathbb{N}}$ we assume that
$\Phi_{l}(t)\prec\sqrt{N\ell_{t}}$, for $l\leq 2k$, and $\Phi_{l}(t)\prec 1$,
for $l\leq k-2$. Then, by Proposition 5.3, we have $\Phi_{l}(t)\prec 1$ for
$l\leq k$; this concludes the induction step, hence the proof. ∎

###  5.2. Master and reduction inequalities: Proofs of Propositions 5.2–5.3

We first present the proof of Proposition 5.2 in detail and then at the end of
this section we explain the minor changes to obtain Proposition 5.3.

###### Proof of Proposition 5.2.

To keep the notation simple, from now on we often omit the $t$–dependence from
the resolvents $G_{i}=G_{i,t}$ and from their deterministic approximations
$M_{[i,j]}=M_{[i,j],t}$, but we will still keep the $t$–dependence in the
spectral parameters, in $\Psi_{k}(t)$ and in the quantities (5.5);
additionally, we stress that $\psi_{k}$ does not depend on time. All estimates
are uniform in the time parameters $t,s\in[0,T]$.

Note that in (5.13) we defined $\Psi_{k}(t)$ only for alternating chains of
single resolvents and deterministic $A$’s, i.e. no $G_{i}^{k}$ appears.
However, in (5.3) we naturally get chains involving $G_{i}^{2}$. For these
terms we use the estimate (recall the estimate for $\langle M_{[i,j]}\rangle$
from (5.6)):

$\big{|}\langle
G_{[i,j]}-M_{[i,j]}\rangle\big{|}\prec\frac{\psi_{j-i}\mathfrak{s}_{j-i}(t)}{N\eta_{t}}.$
(5.23)

For $\mathrm{Im}\,z_{i,t}\mathrm{Im}\,z_{j,t}>0$ this trivially follows by
integral representation

$G_{j}G_{i}=\frac{1}{2\pi\mathrm{i}}\int_{\Gamma}\frac{G(z)}{(z-z_{i,t})(z-z_{j,t})}\,\mathrm{d}z,$
(5.24)

for ”linearizing” the product of the first and last $G$’s in $\langle
G_{[i,j]}\rangle=\langle G_{i}A_{i}\dots A_{j-i}G_{j}\rangle$ after using
cyclicity of the trace. Here $\Gamma$ is a contour which lies in the region
$\\{z\in{\mathbb{C}}:|\mathrm{Im}\,z|\rho(z)\geq\ell_{t}/2\\}$. For
$\mathrm{Im}\,z_{i,t}\mathrm{Im}\,z_{j,t}<0$ we ”linearize” by the resolvent
identity. Using (5.24) will change the value of the imaginary part of the
spectral parameters, so the domain on which the inequalities below hold
(characterized by $N\ell\geq N^{\epsilon}$) may change from time to time,
i.e., say $\epsilon\to\epsilon/2$. However, this can happen only finitely many
times as it does not affect the $\prec$–bound (see [13, Figure 2] for a
detailed discussion of this minor technicality).

To describe the evolution of $M_{[1,k],t}$ in (5.3) we rely on the following
lemma.

######  Lemma 5.5 (Lemma 4.8 of [14]).

For any deterministic matrices $A_{i}\in{\mathbb{C}}^{N\times N}$ (i.e. not
necessarily traceless), we have

$\partial_{t}\langle M_{[1,k],t}A_{k}\rangle=\frac{k}{2}\langle
M_{[1,k],t}A_{k}\rangle+\sum_{i,j=1,\atop i<j}^{k}\langle
M_{[i,j],t}\rangle\langle M_{[j,i],t}\rangle.$ (5.25)

Adding and subtracting the deterministic approximation of each term in (5.3)
and using Lemma 5.5 we thus obtain

$\begin{split}\mathrm{d}\langle(G_{[1,k],t}-M_{{[1,k]},t})A_{k}\rangle&=\frac{1}{\sqrt{N}}\sum_{a,b=1}^{N}\partial_{ab}\langle{G}_{[1,k],t}A_{k}\rangle\mathrm{d}B_{ab,t}+\frac{k}{2}\langle{G}_{[1,k],t}A_{k}\rangle\mathrm{d}t+\sum_{i,j=1\atop
i<j}^{k}\langle{G}_{[i,j],t}-M_{[i,j],t}\rangle\langle
M_{[j,i],t}\rangle\mathrm{d}t\\\ &\quad+\sum_{i,j=1\atop i<j}^{k}\langle
M_{[i,j],t}\rangle\langle
G_{[j,i],t}-M_{[j,i],t}\rangle\mathrm{d}t+\sum_{i,j=1\atop
i<j}^{k}\langle{G}_{[i,j],t}-M_{[i,j],t}\rangle\langle
G_{[j,i],t}-M_{[j,i],t}\rangle\mathrm{d}t\\\ &\quad+\sum_{i=1}^{k}\langle
G_{i,t}-m_{i,t}\rangle\langle{G}^{(i)}_{[1,k],t}A_{k}\rangle\mathrm{d}t.\end{split}$
(5.26)

We point out that, using a simple change of variables, the term
$\langle{G}_{[1,k]}A_{k}\rangle$ amounts to a simple rescaling $e^{kt/2}$, so
we will ignore it. The quadratic variation of the martingale term in (5.26) is
given by

$\frac{1}{N}\sum_{a,b=1}^{N}\big{|}\partial_{ab}\langle
G_{[1,k]}A_{k}\rangle\big{|}^{2}\lesssim\sum_{i=1}^{k}\frac{\langle\mathrm{Im}\,G_{i}(A_{i}\mathrm{Im}\,G_{i+1}\dots
A_{i-1})\mathrm{Im}\,G_{i}(A_{i}\mathrm{Im}\,G_{i+1}\dots
A_{i-1})^{*}\rangle}{N^{2}\eta_{i,t}^{2}},$ (5.27)

where we also used the Ward-identity
$G_{i}G_{i}^{*}=\eta^{-1}_{i}\mathrm{Im}\,G_{i}$. Notice that the rhs. of
(5.27) naturally contains chains of $2k$ resolvents. However, to have a closed
system of master inequalities for products of resolvents of length $k$, we
split the chain of length $2k$ into the product of two chains of length $k$.
For this purpose we use the following _reduction inequality_ which will be
proven at the end of this section. For any fixed matrices
$R,Q\in{\mathbb{C}}^{N\times N}$, and spectral parameters
$z,w\in{\mathbb{C}}\setminus{\mathbb{R}}$, we have161616We point out that the
bound $\Psi_{k}(t)\leq\psi_{k}$ holds also for chains when some resolvents
$G$’s are replaced by their absolute value $|G|$. This can be easily seen
using the integral representation
$|G(E+\mathrm{i}\eta)|=\frac{2}{\pi}\int_{0}^{\infty}\frac{\mathrm{Im}\,G(E+\mathrm{i}\sqrt{\eta^{2}+v^{2}})}{\sqrt{\eta^{2}+v^{2}}}\,\mathrm{d}v\,,$
(5.28) together with the bound (5.13) for chains containing only resolvents
(see e.g. [10, Lemma 5.1]).

$\langle\mathrm{Im}\,G(z)QG(w)R\mathrm{Im}\,G(z)R^{*}G(w)^{*}Q^{*}\rangle\lesssim
N\langle\mathrm{Im}\,G(z)Q|G(w)|Q^{*}\rangle\langle\mathrm{Im}\,G(z)R^{*}|G(w)|R\rangle.$
(5.29)

We now focus on the case $k$ being even for notational simplicity. Combining
(5.27) with (5.29) used for $z=z_{1},w=z_{k/2}$ and

$Q=A_{1}G_{2}\dots G_{k/2-1}A_{k/2-1}\,,\qquad\qquad\quad
R=A_{k/2}G_{k/2+1}\dots G_{k}A_{k},$

together with $\Psi_{k}(t)\prec\psi_{k}$, we obtain the following bound for
the quadratic variation

$\frac{1}{N}\sum_{a,b=1}^{N}\big{|}\partial_{ab}\langle
G_{[1,k],s}A_{k}\rangle\big{|}^{2}\prec\frac{\mathfrak{m}_{k}(s)^{2}}{N\eta_{s}^{2}}+\frac{N\ell_{t}\mathfrak{s}(s)^{2}\psi_{k}^{2}}{N^{3}\eta_{s}^{2}}.$
(5.30)

Then, using the Burkholder–Davis–Gundy (BDG) inequality, we conclude that the
martingale term in (5.26) is bounded by

$\begin{split}\frac{N}{{\mathfrak{s}_{k}(t)}\sqrt{N\ell_{t}}}\left[\int_{0}^{t}\left(\frac{\ell_{s}\mathfrak{m}_{k}(s)^{2}}{N\eta_{s}^{2}}+\frac{N\ell_{s}\mathfrak{s}_{k}(s)^{2}\psi_{k}^{2}}{N^{3}\eta_{s}^{2}}\right)\,\mathrm{d}s\right]^{1/2}&\lesssim\frac{\mathfrak{m}_{k}(t)}{\mathfrak{s}_{k}(t)}+\frac{\psi_{k}}{\sqrt{N}}\left(\int_{0}^{t}\frac{1}{\eta_{s}^{2}}\,\mathrm{d}s\right)^{1/2}\\\
&\lesssim 1+\frac{\psi_{k}}{\sqrt{N\ell_{t}}},\end{split}$ (5.31)

with very high probability. Here in the first inequality we used (5.10), and
in the second inequality we used (5.8) for the first term and (5.2) for the
second term.

Similarly, for odd $k$, using (5.29) for $z=z_{1},w=z_{(k-1)/2}$ and

$Q=A_{1}G_{2}\dots G_{(k-3)/2}A_{(k-3)/2}\,,\qquad\qquad\quad
R=A_{(k-1)/2}G_{(k+1)/2}\dots G_{k}A_{k}$

we get a bound $1+\sqrt{\psi_{k-1}\psi_{k+1}}/(N\ell_{t})^{1/4}$ for the
martingale term.

Next, using (5.23), we estimate the contribution of the last term in (5.26)

$\begin{split}\frac{N}{{\mathfrak{s}_{k}(t)}\sqrt{N\ell_{t}}}\int_{0}^{t}\langle
G_{i,s}-m_{i,s}\rangle\langle
G_{[1,k],s}^{(i)}A_{k}\rangle\,\mathrm{d}s&\prec\frac{N}{{\mathfrak{s}_{k}(t)}\sqrt{N\ell_{t}}}\int_{0}^{t}\frac{1}{N\eta_{s}^{2}}\left(\ell_{s}\mathfrak{m}_{k}(s)+\frac{\psi_{k}\mathfrak{s}_{k}(s)\sqrt{N\ell_{s}}}{N}\right)\,\mathrm{d}s\\\
&\lesssim\frac{1}{\sqrt{N\ell_{t}}}+\frac{\psi_{k}}{N\ell_{t}},\end{split}$
(5.32)

where in the first inequality we used $|\langle
G_{i,s}-m_{i,s}\rangle|\prec(N\eta_{s})^{-1}$ from (1.4), and in the second
inequality we used (5.8), (5.10) and (5.2).

Using again (5.23), we estimate the last term in the first line and the first
term in the second line of (5.26) as follows:

$\begin{split}\frac{N}{{\mathfrak{s}_{k}(t)}\sqrt{N\ell_{t}}}&\int_{0}^{t}\sum_{i,j=1\atop
i<j}^{k-1}\langle M_{[i,j],s}\rangle\langle
G_{[j,i],s}-M_{[j,i],s}\rangle\,\mathrm{d}s\\\
&\prec\frac{N}{{\mathfrak{s}_{k}(t)}\sqrt{N\ell_{t}}}\int_{0}^{t}\sum_{i,j=1\atop
i+1<j}^{k-1}\frac{\ell_{s}\mathfrak{m}_{j-i}(s)}{\eta_{s}}\frac{\psi_{k-(j-i)}\mathfrak{s}_{k-j+i}(s)\sqrt{N\ell_{s}}}{N\eta_{s}}\,\mathrm{d}s\lesssim\sum_{j=1}^{k-2}\psi_{j},\end{split}$
(5.33)

where in the first inequality we used that $\langle M_{[i,i+1],s}\rangle=0$
and (5.6) to estimate $\langle M_{[i,j],s}\rangle$, and in the second
inequality we used (5.7)–(5.8), (5.10) and (5.2).

Finally, for the last term in the second line of (5.26) we estimate:

$\begin{split}&\frac{N}{{\mathfrak{s}_{k}(t)}\sqrt{N\ell_{t}}}\int_{0}^{t}\sum_{i,j=1\atop
i<j}^{k-1}\langle G_{[i,j],s}-M_{[i,j],s}\rangle\langle
G_{[j,i],s}-M_{[j,i],s}\rangle\,\mathrm{d}s\\\
&\qquad\qquad\qquad\quad\prec\frac{N}{{\mathfrak{s}_{k}(t)}\sqrt{N\ell_{t}}}\sum_{i,j=1\atop
i<j}^{k-1}\int_{0}^{t}\frac{\psi_{j-i}\mathfrak{s}_{j-i}(s)\sqrt{N\ell_{s}}}{N\eta_{s}}\frac{\psi_{k-(j-i)}\mathfrak{s}_{k-j+i}(s)\sqrt{N\ell_{s}}}{N\eta_{s}}\,\mathrm{d}s\\\
&\qquad\qquad\qquad\quad\lesssim\frac{1}{\sqrt{N\ell_{t}}}\sum_{j=1}^{k-1}\psi_{j}\psi_{k-j},\end{split}$
(5.34)

where in the second inequality we used (5.7), (5.10) and (5.2). Collecting all
these bounds, using that $\ell_{t}\gtrsim\ell_{T}$ and $N\ell_{T}\geq 1$, we
obtain

$\Psi_{k}(t)\prec
1+\frac{\psi_{k}+\bm{1}(k\,\mathrm{odd})\sqrt{\psi_{k-1}\psi_{k+1}}}{(N\ell_{T})^{1/4}}+\sum_{j=1}^{k-2}\psi_{j}+\frac{1}{(N\ell_{T})^{1/2}}\sum_{j=1}^{k-1}\psi_{j}\psi_{k-j}.$

Finally, using that $\psi_{l}\geq 1$ and $N\ell_{T}\geq 1$, we conclude
(5.16). ∎

###### Proof of Proposition 5.3.

The proof of (5.17) is very similar to the one of (5.16), for this reason we
only explain the minor differences. All the terms in (5.26) are estimated
exactly in the same way as in the proof of Proposition 5.2 with the exception
of the martingale term. In fact, this is the only step where the estimate for
$\Psi_{k}$ (first step) differs from the estimate on $\Phi_{k}$ (second step).
Estimating a longer chain by two smaller ones loses a certain $N\ell$ factor;
this loss is unavoidable in the first step, but it can be avoided in the
second one, once the a-priori bound is available.

More precisely, the estimates in (5.32), (LABEL:eq:prob3aaa), after
multiplying them by a factor $\sqrt{N\ell_{t}}$, are the same with the only
minor difference that we can now use $\psi_{l}=1$. The estimate (5.33) becomes
(recall that $\langle M_{[i,i+1],s}\rangle=0$)

$\frac{N}{{\mathfrak{s}_{k}(t)}}\int_{0}^{t}\sum_{i,j=1\atop
i+1<j}^{k-1}\langle M_{[i,j],s}\rangle\langle
G_{[j,i],s}-M_{[j,i],s}\rangle\,\mathrm{d}s\prec\frac{N}{{\mathfrak{s}_{k}(t)}}\int_{0}^{t}\sum_{i,j=1\atop
i+1<j}^{k-1}\frac{\ell_{s}\mathfrak{m}_{j-i}(s)}{\eta_{s}}\frac{\Phi_{k-(j-i)}\mathfrak{s}_{k-j+i}(s)}{N\eta_{s}}\,\mathrm{d}s\lesssim
1,$ (5.35)

where the only difference with (5.33) is that now in the second inequality we
used that $\Phi_{l}(t)\prec 1$, for $l\leq k-2$, by the induction assumption.
Instead of the bound for the quadratic variation (5.27), using that
$\Phi_{2k}(t)\prec\sqrt{N\ell_{t}}$, the estimate in (5.31) is replaced by

$\begin{split}\frac{N}{\mathfrak{s}_{k}(t)}\left[\int_{0}^{t}\frac{\langle\mathrm{Im}\,G_{1,s}(A_{1}\dots
A_{k})\mathrm{Im}\,G_{1,s}(A_{k}\dots
A_{1})^{*}\rangle}{N^{2}\eta_{s}^{2}}\,\mathrm{d}s\right]^{1/2}&\prec\frac{N}{\mathfrak{s}_{k}(t)}\left[\int_{0}^{t}\left(\frac{\mathfrak{m}_{2k}(s)\ell_{s}}{N^{2}\eta_{s}^{2}}+\frac{\Phi_{2k}(s)\mathfrak{s}_{2k}(s)}{N^{3}\eta_{s}^{2}}\right)\,\mathrm{d}s\right]^{1/2}\\\
&\lesssim\left[\int_{0}^{t}\left(\frac{\ell_{t}}{\eta_{s}^{2}}+\frac{\sqrt{N\ell_{t}}}{N\eta_{s}^{2}}\right)\,\mathrm{d}s\right]^{1/2}\lesssim
1\,,\end{split}$ (5.36)

where in the second inequality we used (5.8), (5.10) and in the last
inequality we used (5.2). This concludes the proof of Proposition 5.3. ∎

We conclude this section with the proof of (5.29).

###### Proof of (5.29).

By spectral decomposition of $W$ ($\lambda_{i},{\bm{u}}_{i}$ being its
eigenvalues and eigenvectors) we have

$\begin{split}\langle\mathrm{Im}\,G(z)QG(w)R\mathrm{Im}\,G(z)R^{*}G(w)^{*}Q^{*}\rangle&=\frac{|\mathrm{Im}\,z|^{2}}{N}\sum_{i,j,k,l=1}^{N}\frac{\langle{\bm{u}}_{i},Q{\bm{u}}_{j}\rangle\langle{\bm{u}}_{j},R{\bm{u}}_{k}\rangle\langle{\bm{u}}_{k},R^{*}{\bm{u}}_{l}\rangle\langle{\bm{u}}_{l},Q^{*}{\bm{u}}_{i}\rangle}{|\lambda_{i}-z|^{2}(\lambda_{j}-w)|\lambda_{k}-z|^{2}(\lambda_{l}-\overline{w})}\\\
&\lesssim
N\langle\mathrm{Im}\,G(z)Q|G(w)|Q^{*}\rangle\langle\mathrm{Im}\,G(z)R^{*}|G(w)|R\rangle,\end{split}$
(5.37)

where in the last step we used Schwarz inequality

$\big{|}\langle{\bm{u}}_{i},Q{\bm{u}}_{j}\rangle\langle{\bm{u}}_{j},R{\bm{u}}_{k}\rangle\langle{\bm{u}}_{k},R^{*}{\bm{u}}_{l}\rangle\langle{\bm{u}}_{l},Q^{*}{\bm{u}}_{i}\rangle\big{|}\lesssim\big{|}\langle{\bm{u}}_{i},Q{\bm{u}}_{j}\rangle\langle{\bm{u}}_{k},R^{*}{\bm{u}}_{l}\rangle\big{|}^{2}+\big{|}\langle{\bm{u}}_{j},R{\bm{u}}_{k}\rangle\langle{\bm{u}}_{l},Q^{*}{\bm{u}}_{i}\rangle\big{|}^{2}\,.\qed$

##  6\. Green function comparison: Proof of Proposition 4.4

In this section, we remove the Gaussian component introduced in Proposition
4.3 using a Green function comparison (GFT) argument, i.e. we prove
Proposition 4.4. The basic idea is the same as in [14, Section 5]: We perform
a _self-consistent_ GFT (i.e., given a local law for one ensemble, we aim to
prove it for a different one) using an entry-by-entry Lindeberg replacement
strategy in $O(N^{2})$ many steps. Note that, unlike here, in typical
applications of GFT to answer universality questions, the local law is given
as an a-priori input. Prior to [14] the GFT has been used in a similar spirit
by Knowles and Yin [31] in order to prove a single resolvent local law for
ensembles, where the deterministic approximation $M$ to $G$ is no longer a
multiple of the identity (e.g. deformed Wigner matrices). In contrast to our
approach, they used a continuous interpolation between ensembles, but we stick
with the entrywise Lindeberg replacement, which is easier to adjust to
multiple resolvents, similarly as in [14].

A characteristic property of the Lindeberg strategy, is that along the
replacement procedure _isotropic resolvent chains_ naturally arise. In
particular, we have to consider the isotropic analog of Theorem 3.3, the
average local law, as well, and need to show that also

$\big{(}G_{1}A_{1}...G_{k}A_{k}G_{k+1}\big{)}_{\bm{x}\bm{y}}$ (6.1)

concentrates around a deterministic value
$\big{(}M(z_{1},A_{1},...,z_{k},A_{k},z_{k+1})\big{)}_{\bm{x}\bm{y}}$ with $M$
given by (3.2). This will also be done via the Zigzag strategy (recall the
outline in the beginning of Section 4) in Section 6.1.

First, analogously to Lemma 3.2, we have the following bound on the
deterministic approximation of (6.1), the proof of which is deferred to
Appendix A.

######  Lemma 6.1 (Isotropic $M$-bounds).

Assume the setting of Lemma 3.2 but with $k+1$ (instead of $k$) spectral
parameters $z_{1},...,z_{k+1}\in{\mathbb{C}}\setminus{\mathbb{R}}$. Then, for
deterministic vectors $\bm{x},\bm{y}\in{\mathbb{C}}^{N}$ with
$\|\bm{x}\|,\|\bm{y}\|\lesssim 1$, it holds that171717Analogously to Footnote
10, we point out that the case of $k$ odd admits an improved bound by
$\prod_{i\in[k]}\|A_{i}\|^{\frac{1}{k}}\,{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{\infty,\ell}^{1-\frac{1}{k}}$,
but we do not follow this improvement for simplicity.

$\left|\langle\bm{x},{M}_{[1,k+1]}(z_{1},A_{1},...,A_{k},z_{k+1})\bm{y}\rangle\right|\lesssim\prod_{i\in[k]}{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{\infty,\ell}\,,$
(6.2)

where
${\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{\infty,\ell}$
has been introduced in Definition 3.1.

The analog of Theorem 3.3 is the following _isotropic multi-resolvent local
law_. The proof is given in Section 6.1.

######  Theorem 6.2 (Isotropic multi-resolvent local laws).

Assume the setting of Theorem 3.3 but with $k+1$ (instead of $k$) spectral
parameters $z_{1},...,z_{k+1}\in{\mathbb{C}}\setminus{\mathbb{R}}$ and let
$\bm{x},\bm{y}\in{\mathbb{C}}^{N}$ be deterministic vectors with
$\|\bm{x}\|,\|\bm{y}\|\lesssim 1$. Then it holds that

$\begin{split}\left|\langle\bm{x},{G}_{1}A_{1}...A_{k}{G}_{k+1}\bm{y}\rangle-\langle\bm{x},{M}_{[1,k+1]}\bm{y}\rangle\right|\prec\frac{1}{\sqrt{N\ell}}\prod_{i\in[k]}{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{\infty,\ell}\,,\end{split}$
(6.3)

where
${\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{\infty,\ell}$
has been introduced in Definition 3.1.

Analogously to Theorem 3.3 and (3.9), this unifies and improves the previous
local laws with operator norm (see [10, Eq. (2.11b) in Theorem 2.5]) and
Hilbert-Schmidt norm (see [11, Corollary 2.4]) in the bulk of the spectrum.
This follows by estimating (in the relevant $\ell\lesssim 1$ regime)

${\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{\infty,\ell}\lesssim\begin{cases}\dfrac{\|A_{i}\|}{\ell^{1/2}}\qquad&\text{for
\ \cite[cite]{[\@@bibref{Number}{multiG}{}{}, Eq.~{}(2.11b) in Theorem
2.5]}}\\\\[5.69054pt]
(N\ell)^{1/2}\dfrac{\langle|A_{i}|^{2}\rangle^{1/2}}{\ell^{1/2}}\qquad&\text{for
\ \cite[cite]{[\@@bibref{Number}{A2}{}{}, Corollary~{}2.4]}}\end{cases}$ (6.4)

in (6.3) for $N\ell>1$ and every $i\in[k]$ by means of elementary
inequalities.

Note that in the isotropic law (6.3) we use the
${\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{\infty,\ell}$
norm instead of the
${\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{2k,\ell}$
norm in the corresponding averaged law (3.8). By taking
$A_{k+1}:=N\bm{y}\bm{x}^{*}$ (assume for simplicity that
$\langle\bm{x},\bm{y}\rangle=0$) in (3.8) for $k\to k+1$, it would in fact be
possible to obtain an isotropic law immediately from Theorem 3.3. However, the
bound provided in (6.3) is stronger than that, as can be seen by means of
(3.6) and
${\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{k+1}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{2(k+1),\ell}\sim
N/(N\ell)^{\frac{1}{2(k+1)}}$, which yield

$\frac{1}{\sqrt{N\ell}}\prod_{i\in[k]}{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{\infty,\ell}\lesssim\frac{1}{N}\prod_{i\in[k+1]}{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{2(k+1),\ell}\,.$
(6.5)

In case that all $A_{i}$ for $i\in[k]$ have large rank, the lhs. of (6.5) is
in fact much smaller (by some inverse $(N\ell)$-power) than the rhs. of (6.5).

As for Theorem 3.3, we now give a concrete example (6.6) how to use Theorem
6.2 for general (i.e. not necessarily traceless) matrices.

######  Example 6.3.

For $k=1$, by (6.3), we have

$\left|\langle\bm{x},G_{1}BG_{2}\bm{y}\rangle-
m_{1}m_{2}B_{\bm{x}\bm{y}}-\frac{m_{1}m_{2}\langle
B\rangle\langle\bm{x},\bm{y}\rangle}{1-m_{1}m_{2}}\right|\prec\,\frac{|\langle
B\rangle|}{\sqrt{N\ell^{3}}}+\frac{\langle|\mathring{B}|^{2}\rangle^{1/2}}{\sqrt{N\ell^{2}}}+\frac{\|\mathring{B}\|}{\sqrt{N\ell}}\,,$
(6.6)

for a general matrix $B=\mathring{B}+\langle B\rangle$, completely analogously
to (3.11).

###  6.1. Isotropic law: Proof of Theorem 6.2

Analogously to the proof of the averaged law, Theorem 3.3, the proof of the
isotropic law, Theorem 6.2, is also conducted via the Zigzag strategy with
natural modifications, so we will be very brief. The initial step, the global
law, has already been proven in [10, Theorem 2.5] (see Eq. (2.11b) in the
$d\geq 1$ regime, which is even stronger than (6.7) below).

######  Proposition 6.4 (Step 1: Isotropic global law).

Let $W$ be a Wigner matrix satisfying Assumption 2.1, and fix
$k\in{\mathbb{N}}$ and $\delta>0$. Consider spectral parameters
$z_{1},...,z_{k+1}\in{\mathbb{C}}\setminus{\mathbb{R}}$, the associated
resolvents $G_{j}=G(z_{j}):=(W-z_{j})^{-1}$, and deterministic matrices
$B_{1},...,B_{k}\in{\mathbb{C}}^{N\times N}$. Then, uniformly in spectral
parameters satisfying
$d:=\min_{j\in[k+1]}\mathrm{dist}(z_{j},[-2,2])\geq\delta$, deterministic
matrices $B_{1},...,B_{k}$ and deterministic vectors $\bm{x},\bm{y}$ with
$\|\bm{x}\|,\|\bm{y}\|\lesssim 1$, it holds that

$\left|\big{(}G_{1}B_{1}...B_{k}G_{k+1}\big{)}_{\bm{x}\bm{y}}-\big{(}{M}_{[1,k+1]}\big{)}_{\bm{x}\bm{y}}\right|\prec\frac{\prod_{i\in[k]}{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|B_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{\infty,d}}{\sqrt{Nd}}\,.$
(6.7)

In the second step, the global law is propagated to a local law through the
characteristic flow (4.1)–(4.2), thereby introducing an order one Gaussian
component. The proof of Proposition 6.5 is postponed to Appendix A.2.1.

######  Proposition 6.5 (Step 2: Isotropic characteristic flow).

Fix any $\epsilon,\gamma>0$, a time $0\leq T\leq\gamma$, and
$K\in{\mathbb{N}}$. Consider
$z_{1,0},\dots,z_{K+1,0}\in{\mathbb{C}}\setminus{\mathbb{R}}$ as initial
conditions to the solution $z_{j,t}$ of (4.2) for $0\leq t\leq T$ and define
$G_{j,t}:=G_{t}(z_{j,t})$ as well as $\eta_{j,t}:=|\mathrm{Im}\,z_{j,t}|$ and
$\rho_{j,t}:=\pi^{-1}|\mathrm{Im}\,m(z_{j,t})|$.

Let $k\leq K$, define $\ell_{t}:=\min_{j\in[k+1]}\eta_{j,t}\rho_{j,t}$ and
recall (3.5). Assuming that

$\left|\big{(}{G}_{1,0}A_{1}...A_{k}{G}_{k+1,0}\big{)}_{\bm{x}\bm{y}}-\big{(}{M}_{[1,k+1],0}\big{)}_{\bm{x}\bm{y}}\right|\prec\frac{1}{\sqrt{N\ell_{0}}}\prod_{i\in[k]}{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{\infty,\ell_{0}}$
(6.8)

holds uniformly for any $k\leq K$, any choice of deterministic traceless
$A_{1},...,A_{k}$, any choice of $z_{i,0}$’s such that $N\ell_{0}\geq
N^{\epsilon}$ and $\max_{j\in[k]}|z_{j,0}|\leq N^{1/\epsilon}$, and all
deterministic vectors $\|\bm{x}\|,\|\bm{y}\|\lesssim 1$, then we have

$\left|\big{(}{G}_{1,T}A_{1}...A_{k}{G}_{k+1,T}\big{)}_{\bm{x}\bm{y}}-\big{(}{M}_{[1,k+1],T}\big{)}_{\bm{x}\bm{y}}\right|\prec\frac{1}{\sqrt{N\ell_{T}}}\prod_{i\in[k]}{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{\infty,\ell_{T}}$
(6.9)

for any $k\leq K$, again uniformly in traceless matrices $A_{i}$, in
deterministic vectors $\bm{x},\bm{y}$ with $\|\bm{x}\|,\|\bm{y}\|\lesssim 1$
and in spectral parameters satisfying $N\ell_{T}\geq N^{\epsilon}$ and
$\max_{j\in[k]}|z_{j,T}|\leq N^{1/\epsilon}$.

In the third and final step, we remove the Gaussian component introduced in
Proposition 6.5 by a GFT argument. The proof of Proposition 6.6 is given in
Section 6.2.1 below.

######  Proposition 6.6 (Step 3: Isotropic Green function comparison).

Let $H^{(\bm{v})}$ and $H^{(\bm{w})}$ be two $N\times N$ Wigner matrices with
matrix elements given by the random variables $v_{ab}$ and $w_{ab}$,
respectively, both satisfying Assumption 2.1 and having matching moments up to
third order, i.e.

${\mathbb{E}}\bar{v}_{ab}^{u}v_{ab}^{s-u}={\mathbb{E}}\bar{w}_{ab}^{u}w_{ab}^{s-u}\,,\quad
s\in\\{0,1,2,3\\}\,,\quad u\in\\{0,...,s\\}\,.$ (6.10)

Fix $K\in{\mathbb{N}}$ and consider spectral parameters
$z_{1},...,z_{K+1}\in{\mathbb{C}}\setminus{\mathbb{R}}$ satisfying
$\min_{j}N\eta_{j}\rho_{j}\geq N^{\epsilon}$ and $\max_{j}|z_{j}|\leq
N^{1/\epsilon}$ for some $\epsilon>0$ and the associated resolvents
$G_{j}^{(\\#)}=G^{(\\#)}(z_{j}):=(H^{(\\#)}-z_{j})^{-1}$ with
$\\#=\bm{v},\bm{w}$. Pick traceless matrices
$A_{1},...,A_{K}\in{\mathbb{C}}^{N\times N}$.

For any $k\leq K$, consider any subset of cardinality $k+1$ of the $K+1$
spectral parameters, and similarly, consider any subset of cardinality $k$ of
the deterministic matrices. Relabeling them by $[k+1]$ and $[k]$,
respectively, setting $\ell:=\min_{j\in[k+1]}\eta_{j}\rho_{j}$ and recalling
(3.5), we assume that for $G_{j}=G_{j}^{(\bm{v})}$ we have

$\left|\big{(}{G}_{1}A_{1}...A_{k}{G}_{k+1}\big{)}_{\bm{x}\bm{y}}-\big{(}{M}_{[1,k+1]}\big{)}_{\bm{x}\bm{y}}\right|\prec\frac{1}{\sqrt{N\ell}}\prod_{i\in[k]}{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{\infty,\ell}\,,$
(6.11)

uniformly in all choices of subsets of $z$’s and $A$’s and in bounded
deterministic vectors $\|\bm{x}\|,\|\bm{y}\|\lesssim 1$.

Then, (6.11) also holds for the ensemble $H^{(\bm{w})}$, i.e. for
$G_{j}=G_{j}^{(\bm{w})}$, uniformly in all choices of subsets of $z$’s and
$A$’s and in bounded deterministic vectors $\|\bm{x}\|,\|\bm{y}\|\lesssim 1$.

Based on Propositions 6.4–6.6, the proof of Theorem 6.2 is concluded in the
same way as the proof of Theorem 3.3 in the end of Section 4. ∎

###  6.2. GFT argument: Proof of Propositions 4.4 and 6.6

The principal idea of the GFT argument is the same as in [14, Section 5] and
even the detailed argument almost directly translates to our case. The main
difference is that we now use the conceptual _mean and standard deviation
sizes_ ${\mathfrak{m}}/{\mathfrak{s}}$ (recall (5.5)) and ${\mathfrak{m}}^{\rm
iso}/{\mathfrak{s}}^{\rm iso}$ (see (6.13) and (A.5) below) as basic
deterministic control quantities; while in [14, Section 5] they were not
introduced explicitly as they essentially boiled down to simple $N$-powers.
More precisely, in view of the bounds [14, Eqs. (2.17)–(2.18) and
(2.20)–(2.21)], we simply replace (ignoring all the $\rho$-factors and using
the normalization $\langle|A_{i}|^{2}\rangle=1$ for the just mentioned terms
in [14])

$\frac{N^{k/2}}{N\ell}\longrightarrow{\mathfrak{m}}_{k}\,,\qquad\frac{N^{k/2}}{(N\ell)^{1/2}}\longrightarrow{\mathfrak{s}}_{k}\,,\qquad
N^{k/2}\longrightarrow{\mathfrak{m}}_{k}^{\rm
iso}\,,\qquad\text{and}\qquad\frac{N^{k/2}}{\ell^{1/2}}\longrightarrow{\mathfrak{s}}_{k}^{\rm
iso}\,.$ (6.12)

Given this similarity (6.12) we will henceforth be very brief and only point
out a few minor adaptions of the proof in [14, Section 5] to our new
conceptual notations in Sections 6.2.1–6.2.2, discussing the isotropic and
averaged case, respectively. For further simplicity of notation, we shall drop
all irrelevant sub-scripts of resolvents and deterministic matrices, i.e.
write $G=G_{j}$ and $A=A_{j}$. Additionally, we will also drop the
$\widecheck{G}$\- and $G^{(\gamma)}$-notations (see [14, Eq. (5.10)]),
indicating the precise step in the replacement procedure, as it will be
irrelevant for the modifications discussed below.

####  6.2.1. Part (a): Proof of the isotropic law, Proposition 6.6 (cf.
Section 5.2 in [14])

We start with some preliminaries. In order to express the bounds in (6.2) and
(6.3) concisely, we employ the ${\mathfrak{m}}^{\rm iso}/{\mathfrak{s}}^{\rm
iso}$-notation (the _isotropic mean and standard deviation sizes_ ,
analogously to (5.5) in the average case) to write (6.2) and (6.3) as

$\left|({M}_{[1,k+1]})_{\bm{x}\bm{y}}\right|\lesssim({\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{\infty,\ell})^{k}=:{\mathfrak{m}}_{k}^{\rm
iso}\quad\text{and}\quad\left|\big{(}(GA)^{k}G\big{)}_{\bm{x}\bm{y}}-\big{(}{M}_{[1,k+1]}\big{)}_{\bm{x}\bm{y}}\right|\prec\frac{({\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{\infty,\ell})^{k}}{\sqrt{N\ell}}=:\frac{{\mathfrak{s}}_{k}^{\rm
iso}}{N^{1/2}}\,,$ (6.13)

respectively (the general definition of ${\mathfrak{m}}^{\rm
iso}/{\mathfrak{s}}^{\rm iso}$ is given in (A.5)). Moreover, we also write
$M_{j-i+1}\equiv M_{[i,j]}$ for all $1\leq i<j\leq k+1$ with a slight abuse of
notation (see [14, Eq. (5.33)]). Lastly, we will heavily use the following
relations, proven in parts (i)–(ii) of Lemma A.2,

$\mathfrak{m}^{\rm iso}_{j}\,\mathfrak{m}_{k-j}^{\rm
iso}\lesssim\mathfrak{m}_{k}^{\rm iso}\,,\qquad{\mathfrak{s}}^{\rm
iso}_{j}\,{\mathfrak{s}}^{\rm iso}_{k-j}\lesssim\ell^{-1/2}{\mathfrak{s}}^{\rm
iso}_{k}\,,\quad\text{and}\quad\mathfrak{m}_{k}^{\rm
iso}\lesssim\ell^{1/2}\,{\mathfrak{s}}_{k}^{\rm iso}\,$ (6.14)

for all $k\in{\mathbb{N}}$ and $0\leq j\leq k$ using the conventions
${\mathfrak{m}}_{0}^{\rm iso}:=1$ and ${\mathfrak{s}}_{0}^{\rm
iso}:=\ell^{-1/2}$ (recall $\ell=\min_{i}\eta_{i}\rho_{i}$). With
(6.13)–(6.14) at hand, we can now discuss the two bits of the argument in [14,
Section 5.2], which are not completely straightforward to adapt to our current
setting (cf. Case (i) and Case (ii) in [14, Section 5.2.2]). We will refer to
explicit equation numbers within a longer proof in [14] that we do not repeat
here, hence the reader needs to be familiar with [14, Section 5.2] to follow
the details. However, to facilitate a high level understanding without going
into details, we recall the novel idea in [14]. Traditional GFT proofs via the
Lindeberg strategy estimated the change of, say,
${\mathbb{E}}\big{|}\big{(}(GA)^{k}G\big{)}_{\bm{x}\bm{y}}\big{|}^{p}$ after
each replacement ${\bm{v}}_{ab}\to{\bm{w}}_{ab}$ and showed that it is bounded
by $o(N^{-2})$-times its natural target size, so $N^{2}$ replacements were
affordable. In other words the telescopic summation over $(a,b)\in[N]^{2}$ was
done trivially. This does not hold for the self-consistent GFT argument in
[14]: the replacement for some $(a,b)$ pairs may be too large, but their sum
is still affordable. We will refer to it as the gain from the summation idea.
This is related to our efforts to control the observables in tracial norms;
the simplest toy example is the estimate

$\frac{1}{N^{2}}\sum_{ab}|A_{ab}|\leq N^{-1/2}\langle|A|^{2}\rangle^{1/2}$
(6.15)

which is optimal. However, if we use the best available bound
$|A_{ab}|\leq\|A\|$ for each summand individually, then the lhs. of (6.15) is
overestimated by $\|A\|$ which is much bigger than the rhs. With keeping this
idea in mind, we now return to the detailed modifications in [14, Section
5.2.2].

The terms considered in Case (i) are the ones for which there are enough
$(G-M)$-type terms to balance the $N$-prefactor (cf. (6.16) and (6.17)) by
solely using (6.14). The terms in Case (ii) require to _gain from the
summation_ over all steps in the replacement procedure.

  * Case (i):

First, the analog of [14, Eq. (5.34)] becomes (neglecting the irrelevant
$|\Psi_{k}|^{p-1}$\- and $N^{\xi}$-factors)

$\begin{split}&\frac{N^{1/2}}{{\mathfrak{s}}_{k}^{\rm
iso}}\sum_{\begin{subarray}{c}0\leq k_{l}\leq k-1:\\\
\sum_{l}k_{l}=k\end{subarray}}{\mathbb{E}}\Big{[}\big{|}\big{(}({G}A)^{k_{1}}{G}-M_{k_{1}+1}\big{)}_{\bm{x}\bm{e}_{i}}\big{(}M_{k_{2}+1}\big{)}_{\bm{e}_{j}\bm{e}_{j}}\big{(}M_{k_{3}+1}\big{)}_{\bm{e}_{i}\bm{e}_{i}}\big{(}M_{k_{4}+1}\big{)}_{\bm{e}_{j}\bm{e}_{j}}\big{(}M_{k_{5}+1}\big{)}_{\bm{e}_{i}\bm{y}}\big{|}\\\
&\quad+\big{|}\big{(}({G}A)^{k_{1}}{G}-M_{k_{1}+1}\big{)}_{\bm{x}\bm{e}_{i}}\big{(}({G}A)^{k_{2}}{G}-M_{k_{2}+1}\big{)}_{\bm{e}_{j}\bm{e}_{j}}\big{(}M_{k_{3}+1}\big{)}_{\bm{e}_{i}\bm{e}_{i}}\big{(}M_{k_{4}+1}\big{)}_{\bm{e}_{j}\bm{e}_{j}}\big{(}M_{k_{5}+1}\big{)}_{\bm{e}_{i}\bm{y}}\big{|}+...\Big{]}\\\
\lesssim\,&\frac{N^{1/2}}{{\mathfrak{s}}_{k}^{\rm
iso}}\sum_{\begin{subarray}{c}0\leq k_{l}\leq k-1:\\\
\sum_{l}k_{l}=k\end{subarray}}\left[\frac{{\mathfrak{s}}_{k_{1}}^{\rm
iso}}{N^{1/2}}\left(\prod_{l=2}^{5}{\mathfrak{m}}_{k_{l}}^{\rm
iso}\right)+\frac{{\mathfrak{s}}_{k_{1}}^{\rm iso}{\mathfrak{s}}_{k_{2}}^{\rm
iso}}{N}\left(\prod_{l=3}^{5}{\mathfrak{m}}_{k_{l}}^{\rm
iso}\right)+...\right]\lesssim\left[1+\frac{1}{(N\ell)^{1/2}}+...\right]\lesssim
1\,.\end{split}$ (6.16)

In the penultimate step we used the relations in (6.14) multiple times to
estimate ${\mathfrak{s}}_{k_{1}}^{\rm
iso}\prod_{l=2}^{5}{\mathfrak{m}}_{k_{l}}^{\rm
iso}\lesssim{\mathfrak{s}}_{k}^{\rm iso}$ and ${\mathfrak{s}}_{k_{1}}^{\rm
iso}{\mathfrak{s}}_{k_{l}}^{\rm iso}\prod_{l=3}^{5}{\mathfrak{m}}_{k_{l}}^{\rm
iso}\lesssim\ell^{-1/2}{\mathfrak{s}}_{k}^{\rm iso}$, and similarly for the
other analogous terms indicated by dots. The last step in (6.16) is due to
$N\ell>1$.

  * Case (ii):

Second, the key trick in [14, Section 5.2.2], the _gain from summations_
described in [14, Example 5.7], turns into the following. The _trivial
estimate_ , cf. [14, Eq. (5.36)], reads (again neglecting the irrelevant
$|\Psi_{k}|^{p-1}$\- and $N^{\xi}$-factors)

$\begin{split}&\frac{N^{1/2}}{{\mathfrak{s}}_{k}^{\rm
iso}}\sum_{\begin{subarray}{c}0\leq k_{l}\leq k:\\\
\sum_{l}k_{l}=k\end{subarray}}\left[\big{|}\big{(}M_{k_{1}+1}\big{)}_{\bm{x}\bm{e}_{i}}\big{(}M_{k_{2}+1}\big{)}_{\bm{e}_{j}\bm{e}_{j}}\big{(}M_{k_{3}+1}\big{)}_{\bm{e}_{i}\bm{e}_{i}}\big{(}M_{k_{4}+1}\big{)}_{\bm{e}_{j}\bm{e}_{j}}\big{(}M_{k_{5}+1}\big{)}_{\bm{e}_{i}\bm{y}}\big{|}+...\right]\\\
\lesssim\,&\frac{N^{1/2}}{{\mathfrak{s}}_{k}^{\rm
iso}}\sum_{\begin{subarray}{c}0\leq k_{l}\leq k:\\\
\sum_{l}k_{l}=k\end{subarray}}\left[\left(\prod_{l=1}^{5}{\mathfrak{m}}_{k_{l}}^{\rm
iso}\right)+...\right]\lesssim\left(N\ell\right)^{1/2}\,,\end{split}$ (6.17)

where we again used (6.14) multiple times. Analogously to [14, Eq. (5.37)],
this can be improved on by averaging over all replacement steps: Fixing one
constellation of $k_{l}$’s in (6.17), we find

$\begin{split}&\frac{N^{1/2}}{{\mathfrak{s}}_{k}^{\rm
iso}}\frac{1}{N^{2}}\sum_{i,j}\left[\big{|}\big{(}M_{k_{1}+1}\big{)}_{\bm{x}\bm{e}_{i}}\big{(}M_{k_{2}+1}\big{)}_{\bm{e}_{j}\bm{e}_{j}}\big{(}M_{k_{3}+1}\big{)}_{\bm{e}_{i}\bm{e}_{i}}\big{(}M_{k_{4}+1}\big{)}_{\bm{e}_{j}\bm{e}_{j}}\big{(}M_{k_{5}+1}\big{)}_{\bm{e}_{i}\bm{y}}\big{|}+...\right]\\\
&\qquad\qquad\lesssim\,(N\ell)^{1/2}\frac{1}{N^{2}}\sum_{i,j}\left[\frac{\big{|}(M_{k_{1}+1})_{\bm{x}\bm{e}_{i}}\big{|}}{{\mathfrak{m}}_{k_{1}}^{\rm
iso}}+...\right]\lesssim(N\ell)^{1/2}\frac{1}{N^{1/2}}\lesssim 1\end{split}$
(6.18)

in the relevant $\ell\leq 1$ regime (the opposite regime being covered by
Proposition 6.4). To go to the second line, we employed (6.14). In the
penultimate step, we used

$\sum_{i}\big{|}(M_{k_{1}+1})_{\bm{x}\bm{e}_{i}}\big{|}\leq\sqrt{N}\sqrt{\big{(}|M_{k_{1}+1}|^{2}\big{)}_{\bm{x}\bm{x}}}\lesssim\sqrt{N}\,{\mathfrak{m}}_{k_{1}}^{\rm
iso}\,,$ (6.19)

which follows by a Schwarz inequality, completely analogously to [14, Eq.
(5.38)–(5.39)].

With these two slight adjustments, which rest on the estimates in (6.14), the
proof in [14, Section 5.2] can entirely be translated to our current setting
in a straightforward way, so we omit further details. This concludes the proof
of Proposition 6.6. ∎

####  6.2.2. Part (b): Proof of the averaged law, Proposition 4.4 (cf. Section
5.3 in [14])

For the averaged law, in addition to the bounds in (6.13)–(6.14), we will also
use the relations

${\mathfrak{m}}_{k}\lesssim{\mathfrak{m}}_{k}^{\rm
iso}\,,\qquad{\mathfrak{s}}_{k}\lesssim\ell^{1/2}\,{\mathfrak{s}}_{k}^{\rm
iso}\,,\quad\text{and}\quad{\mathfrak{s}}_{k}^{\rm
iso}\lesssim\,N\,\big{[}1+(N\ell)^{-1/2}\big{]}\,{\mathfrak{s}}_{k}$ (6.20)

from (A.11)–(A.12) in Lemma A.3 for all $k\in{\mathbb{N}}$. Just as in the
case of Section 6.2.1, there are two bits of the argument in [14, Section
5.3], that are not entirely straightforward to adjust to the setting of the
current paper (cf. Case (i) and Case (ii) in [14, Section 5.3]).

Analogously to the isotropic case in Section 6.2.1, while the terms considered
in Case (i) solely use (6.14) and (6.20), the terms in Case (ii) require to
_gain from the summation_ over all steps in the replacement procedure.

  * Case (i):

First, consider [14, Eq. (5.54)] with $d=2$ (for concreteness). Then, the
estimate turns into (neglecting the irrelevant $N^{\xi}$-factor)

$\left(\frac{N}{{\mathfrak{s}}_{k}}\right)^{2}\frac{1}{N^{2}}\sum_{\begin{subarray}{c}k_{1}+k_{2}=k\\\
k_{1}^{\prime}+k_{2}^{\prime}=k\end{subarray}}\left[{\mathfrak{m}}_{k_{1}}^{\rm
iso}\frac{{\mathfrak{s}}_{k_{2}}^{\rm
iso}}{N^{1/2}}{\mathfrak{m}}_{k_{1}^{\prime}}^{\rm
iso}\frac{{\mathfrak{s}}_{k_{2}^{\prime}}^{\rm
iso}}{N^{1/2}}+{\mathfrak{m}}_{k_{1}}^{\rm iso}{\mathfrak{m}}_{k_{2}}^{\rm
iso}\frac{{\mathfrak{s}}_{k_{1}^{\prime}}^{\rm
iso}}{N^{1/2}}\frac{{\mathfrak{s}}_{k_{2}^{\prime}}^{\rm
iso}}{N^{1/2}}+...\right]\lesssim\left[1+\frac{1}{N\ell}\right]\lesssim 1\,,$

where we used (6.14) and (6.20) to bound ${\mathfrak{m}}_{k_{1}}^{\rm
iso}{\mathfrak{s}}_{k_{2}}^{\rm
iso}\lesssim\ell^{-1/2}[1+(N\ell)^{1/2}]{\mathfrak{s}}_{k}$ as well as
${\mathfrak{m}}_{k_{1}}^{\rm iso}{\mathfrak{m}}_{k_{2}}^{\rm
iso}\lesssim[1+(N\ell)^{1/2}]{\mathfrak{s}}_{k}$ and
${\mathfrak{s}}_{k_{1}^{\prime}}^{\rm iso}{\mathfrak{s}}_{k_{2}^{\prime}}^{\rm
iso}\lesssim\ell^{-1}[1+(N\ell)^{1/2}]{\mathfrak{s}}_{k}$ in the penultimate
step, and the last step used $N\ell>1$.

  * Case (ii):

Second, we again discuss how to _gain from the summation_ , as originally
explained in [14, Example 5.11]. The _trivial estimate_ from [14, Eq. (5.55)],
again neglecting the irrelevant $|\Psi_{k}|^{p-1}$ factor and fixing one
constellation of $k_{l}$’s summing up to $k$, becomes (neglecting the
irrelevant $N^{\xi}$-factor)

$\begin{split}&({\mathfrak{s}}_{k})^{-1}\left[\big{|}\big{(}M_{k_{1}+1}\big{)}_{\bm{e}_{i}\bm{e}_{i}}\big{(}M_{k_{2}+1}\big{)}_{\bm{e}_{j}\bm{e}_{j}}\big{(}M_{k_{3}+1}\big{)}_{\bm{e}_{i}\bm{e}_{i}}\big{(}M_{k_{4}+1}\big{)}_{\bm{e}_{j}\bm{e}_{j}}\big{|}+...\right]\\\
&\qquad\qquad\qquad\lesssim\,({\mathfrak{s}}_{k})^{-1}\left[\prod_{l}{\mathfrak{m}}_{k_{l}}^{\rm
iso}+...\right]\lesssim(N\ell)^{1/2}\,,\end{split}$ (6.21)

where we used (6.14) and (6.20). Analogously to [14, Eq. (5.59)], we can again
improve upon this by averaging over all replacement positions. The key for
this _gain from the summation_ is the bound

$\langle|M_{k+1}|^{2}\rangle\lesssim({\mathfrak{s}}_{k})^{2}\qquad\forall
k\in{\mathbb{N}}\,,$ (6.22)

which can be obtained completely analogously to [14, Lemma 5.12]. Armed with
(6.22), the analog of the improved bound [14, Eq. (5.59)] now reads (again
neglecting the irrelevant $N^{\xi}$-factor)

$\begin{split}&({\mathfrak{s}}_{k})^{-1}\frac{1}{N^{2}}\sum_{i,j}\left[\big{|}\big{(}M_{k_{1}+1}\big{)}_{\bm{e}_{i}\bm{e}_{i}}\big{(}M_{k_{2}+1}\big{)}_{\bm{e}_{j}\bm{e}_{j}}\big{(}M_{k_{3}+1}\big{)}_{\bm{e}_{i}\bm{e}_{i}}\big{(}M_{k_{4}+1}\big{)}_{\bm{e}_{j}\bm{e}_{j}}\big{|}+...\right]\\\
\lesssim&({\mathfrak{s}}_{k})^{-1}\left[\prod_{l\in[4]}\left(\frac{1}{N}\sum_{i}\big{|}\big{(}M_{k_{l}+1}\big{)}_{\bm{e}_{i}\bm{e}_{i}}\big{|}^{2}\right)^{1/2}+...\right]\lesssim({\mathfrak{s}}_{k})^{-1}\left[\prod_{l\in[4]}{\mathfrak{s}}_{k_{l}}+...\right]\lesssim
1\,,\end{split}$ (6.23)

where in the first step we employed a trivial Schwarz inequality.

With the above two slight adjustments at hand, which basically rest on the
estimates in (6.14), (6.20) and (6.22), we can straightforwardly translate the
entire proof of the averaged part in [14, Section 5.3] to our current setting.
This finishes our discussion of the adjustments of the arguments from [14,
Section 5.2] and thus the proof of Proposition 4.4. ∎

##  Appendix A Additional technical results

In this section we prove several additional technical results which were used
in the main sections.

###  A.1. Bound on the deterministic approximation: Proofs of Lemmas 3.2 and
6.1

We will deduce Lemmas 3.2 and 6.1 from the following stronger bound. Part (a)
of Lemma A.1 is already proven in [14, Lemma A.1] applied to the special case
$\mathfrak{I}_{k}=\emptyset$. The proof of part (b) is analogous to that and
hence omitted.

######  Lemma A.1 (cf. Lemma A.1 in [14]).

Fix $k\geq 1$ and let $A_{1},...,A_{k}\in{\mathbb{C}}^{N\times N}$ be
traceless deterministic matrices.

  * (a)

Consider spectral parameters
$z_{1},...,z_{k}\in{\mathbb{C}}\setminus{\mathbb{R}}$ and define
$\eta:=\min_{j\in[k]}|\mathrm{Im}\,z_{j}|$. Then, for every $1\leq
s\leq\lfloor k/2\rfloor$ and $\pi\in\mathrm{NC}([k])$ with $|\pi|=k+1-s$, it
holds that

$\left|\langle\mathrm{pTr}_{K(\pi)}(A_{1},\ldots,A_{k-1})A_{k}\rangle\prod_{S\in\pi}m_{\circ}[S]\right|\lesssim\frac{1}{\eta^{s-1}}\left(\prod_{\begin{subarray}{c}S\in
K(\pi)\\\ |S|\geq 2\end{subarray}}\prod_{j\in
S}\langle|A_{j}|^{|S|}\rangle^{\frac{1}{|S|}}\right)\,.$ (A.1)

with $m_{\circ}[S]$ being defined in (3.3). For $s>\lfloor k/2\rfloor$ the
lhs. of (A.1) equals zero.

  * (b)

Consider spectral parameters
$z_{1},...,z_{k+1}\in{\mathbb{C}}\setminus{\mathbb{R}}$ and define
$\eta:=\min_{j\in[k+1]}|\mathrm{Im}\,z_{j}|$. Then, for every $1\leq
s\leq\lceil(k+1)/2\rceil$ and $\pi\in\mathrm{NC}([k+1])$ with $|\pi|=k+1-s$,
it holds that

$\left\|\mathrm{pTr}_{K(\pi)}(A_{1},\ldots,A_{k})\prod_{S\in\pi}m_{\circ}[S]\right\|\lesssim\frac{1}{\eta^{s-1}}\left(\prod_{\begin{subarray}{c}S\in
K(\pi)\setminus\mathfrak{B}:\\\ |S|\geq 2\end{subarray}}\,\prod_{j\in
S}\langle|A_{j}|^{|S|}\rangle^{\frac{1}{|S|}}\right)\left(\prod_{j\in\mathfrak{B}\setminus\\{k+1\\}}\|A_{j}\|\right)\,.$
(A.2)

where $\mathfrak{B}\equiv\mathfrak{B}(k+1)\in K(\pi)$ being the unique block
containing $k+1$ (recall (3.1)). For $s>\lceil(k+1)/2\rceil$ the lhs. of (A.2)
equals zero.

Given Lemma A.1, analogously to [14, Appendix A.1], the proofs of Lemmas 3.2
and 6.1 immediately follow after realizing that $\eta\gtrsim\ell$ and applying
Hölder’s and Young’s inequality. To show this mechanism, consider (A.1) for an
example where $k=6$ and $s=2$, and focus on the non-crossing partition
$\pi=\\{15|2|3|4|6\\}$ for which $K(\pi)=\\{1234|56\\}$. In this case, the
rhs. of (A.1) can be estimated as

$\begin{split}\eta^{-1}\left(\prod_{i=1}^{4}\langle|A_{i}|^{4}\rangle^{1/4}\right)\left(\prod_{i=5}^{6}\langle|A_{i}|^{2}\rangle^{1/2}\right)&\lesssim\ell\left(\prod_{i=1}^{4}\big{[}\langle|A_{i}|^{6}\rangle^{1/8}\ell^{-1/8}\big{]}\,\big{[}\langle|A_{i}|^{2}\rangle^{1/8}\ell^{-1/8}\big{]}\right)\left(\prod_{i=5}^{6}\langle|A_{i}|^{2}\rangle^{1/2}\ell^{-1/2}\right)\\\
&\lesssim\ell\prod_{i\in[6]}\left(\frac{\langle|A_{i}|^{6}\rangle^{1/6}}{\ell^{1/6}}+\frac{\langle|A_{i}|^{2}\rangle^{1/2}}{\ell^{1/2}}\right)\,.\end{split}$

In the first step, we employed $\eta\gtrsim\ell$ and Hölder’s inequality in
the form
$\langle|A_{i}|^{4}\rangle^{1/4}\leq\langle|A_{i}|^{6}\rangle^{1/8}\langle|A_{i}|^{2}\rangle^{1/8}$
for $i\in[4]$. In the second step, we then used Young’s inequality for
$i\in[4]$ and added $\langle|A_{i}|^{6}\rangle^{1/6}\ell^{-1/6}$ for $i=5,6$
to complete the
${\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{i}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{6,\ell}$
norm for all $i\in[6]$ (recall (3.5) and (3.7)). ∎

###  A.2. GFT and isotropic local law: Additional proofs for Section 6

####  A.2.1. Isotropic law: Proof of Theorem 6.2

In this section we want to study chains of the form

$(G_{1}A_{1}\dots
A_{k}G_{k+1})_{{\bm{x}}{\bm{y}}}:=\langle{\bm{x}},G_{1}A_{1}\dots
A_{k}G_{k+1}{\bm{y}}\rangle,$ (A.3)

for unit deterministic vectors ${\bm{x}},{\bm{y}}$.

Following the notation (5.4), by $G_{[1,k+1],t}$ we denote the evolution of
the quantity in (A.3) along the Ornstein-Uhlenbeck flow (4.1) with the
characteristic equation (4.2). Then, by (5.3) and Lemma 5.5 (used for $k\to
k+1$), choosing $A_{k+1}=N{\bm{y}}{\bm{x}}^{*}$, we obtain the flow

$\begin{split}\mathrm{d}(G_{[1,k+1],t}-M_{{[1,k+1]},t})_{{\bm{x}}{\bm{y}}}&=\frac{1}{\sqrt{N}}\sum_{a,b=1}^{N}\partial_{ab}({G}_{[1,k+1],t})_{{\bm{x}}{\bm{y}}}\mathrm{d}B_{ab,t}+\frac{k}{2}(G_{[1,k+1],t}-M_{{[1,k+1]},t})_{{\bm{x}}{\bm{y}}}\mathrm{d}t\\\
&\quad+\sum_{i,j=1\atop
i<j}^{k+1}\langle{G}_{[i,j],t}-M_{[i,j],t}\rangle({M}_{[1,i]\cup[j,k+1],t})_{{\bm{x}}{\bm{y}}}\mathrm{d}t\\\
&\quad+\sum_{i,j=1\atop
i<j}^{k+1}\langle{M}_{[i,j],t}\rangle({G}_{[1,i]\cup[j,k+1],t}-M_{[1,i]\cup[j,k+1],t})_{{\bm{x}}{\bm{y}}}\mathrm{d}t\\\
&\quad+\sum_{i,j=1\atop i<j}^{k+1}\langle
G_{[i,j],t}-{M}_{[i,j],t}\rangle({G}_{[1,i]\cup[j,k+1],t}-M_{[1,i]\cup[j,k+1],t})_{{\bm{x}}{\bm{y}}}\mathrm{d}t\\\
&\quad+\sum_{i=1}^{k+1}\langle
G_{i,t}-m_{i,t}\rangle({G}^{(i)}_{[1,k+1],t})_{{\bm{x}}{\bm{y}}}\mathrm{d}t,\end{split}$
(A.4)

where we used the notation

$({G}_{[1,i]\cup[j,k+1],t})_{{\bm{x}}{\bm{y}}}:=(G_{1,t}A_{1}\dots
A_{i-1}G_{i,t}G_{j,t}A_{j}\dots A_{k-1}G_{k+1,t})_{{\bm{x}}{\bm{y}}}$

and similarly for $M_{[1,i]\cup[j,k+1],t}$.
Deterministic control quantities. Next, we introduce the new _isotropic_
control quantities, the analogues of the _averaged_ quantities defined in
(5.5) (recall (3.5) for the definition of
${\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{\alpha}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{\infty,\ell}$):

$\begin{split}\mathfrak{m}_{k}^{\rm
iso}(\ell;\bm{A}_{J_{k}}):=\prod_{\alpha\in
J_{k}}{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{\alpha}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{\infty,\ell}\qquad\text{and}\qquad{\mathfrak{s}}_{k}^{\rm
iso}(\ell;\bm{A}_{J_{k}}):=\ell^{-1/2}\prod_{\alpha\in
J_{k}}{\left|\kern-1.29167pt\left|\kern-1.29167pt\left|A_{\alpha}\right|\kern-1.29167pt\right|\kern-1.29167pt\right|}_{\infty,\ell}\,,\end{split}$
(A.5)

which will be called the _isotropic mean size_ and _isotropic standard
deviation size_ , respectively. They are functions of a positive number $\ell$
(usually given by $\ell=\min_{i}\eta_{i}\rho_{i}$ for some spectral parameters
$z_{i}$) and a multiset $\bm{A}_{J_{k}}=(A_{\alpha})_{\alpha\in J_{k}}$ of
deterministic matrices of cardinality $|J_{k}|=k$.

Similarly to the paragraph below (5.5), we may often omit the arguments of
$\mathfrak{m}_{k}^{\rm iso}$ and ${\mathfrak{s}}_{k}^{\rm iso}$, write $J_{k}$
instead of $\bm{A}_{J_{k}}$, and for time dependent spectral parameters we use
the short–hand notations $\mathfrak{m}_{k}^{\rm
iso}(t)={\mathfrak{m}}_{k}^{\rm iso}(t;J_{k})={\mathfrak{m}}_{k}^{\rm
iso}(\ell_{t};\bm{A}_{J_{k}})$ and ${\mathfrak{s}}_{k}^{\rm
iso}(t)={\mathfrak{s}}_{k}^{\rm iso}(t;J_{k})={\mathfrak{s}}_{k}^{\rm
iso}(\ell_{t};\bm{A}_{J_{k}})$. Note that with this definitions, by (6.2), we
have (analogously to (5.6)), with $\eta_{t}:=\min_{i}\eta_{i,t}$,

$\big{|}(M_{[i,j],t})_{{\bm{x}}{\bm{y}}}\big{|}\lesssim\mathfrak{m}_{j-i}^{\rm
iso}(t)\qquad\text{and}\qquad\big{|}(M_{[1,i]\cup[j,k+1],t})_{\bm{x}\bm{y}}\big{|}\lesssim\eta_{t}^{-1}{\mathfrak{m}}_{k-(j-i)}^{\rm
iso}(t)\,.$ (A.6)

We now record the following relations about $\mathfrak{m}_{k}^{\rm
iso}/{\mathfrak{s}}_{k}^{\rm iso}$, analogously to Lemma 5.1. The proof is
again an immediate consequence of (3.6) and $\ell_{s}\gtrsim\ell_{t}$ for
$s\leq t$ (recall the discussion below (4.5)) and hence omitted.

######  Lemma A.2 ($\mathfrak{m}^{\rm iso}/{\mathfrak{s}}^{\rm
iso}$-relations).

Let $k\geq 1$ and consider time-dependent spectral parameters
$z_{1,t},...,z_{k+1,t}$ and a multiset of traceless matrices $\bm{A}_{J_{k}}$
as above. Set $\ell_{t}:=\min_{j\in[k+1]}\eta_{j,t}\rho_{j,t}$. Then, the mean
size $\mathfrak{m}^{\rm iso}_{k}$ and the standard deviation size
$\mathfrak{s}_{k}^{\rm iso}$ satisfy the following inequalities (using the
conventions ${\mathfrak{m}}_{0}^{\rm iso}(t):=1$ and ${\mathfrak{s}}_{0}^{\rm
iso}(t):=\ell_{t}^{-1/2}$).

  * (i)

Super-multiplicativity; i.e. for any $0\leq j\leq k$ it holds that

$\begin{split}\mathfrak{m}^{\rm iso}_{j}({t;J_{j}})\,\mathfrak{m}_{k-j}^{\rm
iso}({t;J_{k-j}})\lesssim\mathfrak{m}_{k}^{\rm
iso}({t;J_{k}})\quad\text{and}\quad{\mathfrak{s}}^{\rm
iso}_{j}(t;{J_{j}})\,{\mathfrak{s}}^{\rm
iso}_{k-j}(t;{J_{k-j}})\lesssim\ell_{t}^{-1/2}{\mathfrak{s}}^{\rm
iso}_{k}(t;{J_{k}})\end{split}$ (A.7)

for all disjoint decompositions ${J_{k}}={J_{j}}\dot{\cup}{J_{k-j}}$.

  * (ii)

The mean size can be upper bounded by the standard deviation size as

$\mathfrak{m}_{k}^{\rm
iso}(t;{J_{k}})\lesssim\ell_{t}^{1/2}\,{\mathfrak{s}}_{k}^{\rm
iso}(t;{J_{k}})\,.$ (A.8)

  * (iii)

The standard deviation size satisfies the doubling inequality

${\mathfrak{s}}_{2k}^{\rm iso}(t;J_{k}\cup
J_{k})\lesssim\,\ell_{t}^{1/2}\,\left({\mathfrak{s}}_{k}^{\rm
iso}(t;{J_{k}})\right)^{2}\,,$ (A.9)

where $J_{k}\cup J_{k}$ denotes the union of the multiset $J_{k}$ with itself.

  * (iv)

Monotonicity in time: for any $0\leq s\leq t$ and $\alpha\in[0,1/2]$, we have

$\mathfrak{m}_{k}^{\rm iso}(s;J_{k})\lesssim\mathfrak{m}_{k}^{\rm
iso}(t;J_{k})\qquad\text{and}\qquad{\mathfrak{s}}_{k}^{\rm
iso}(s;J_{k})\,\ell_{s}^{\alpha}\lesssim{\mathfrak{s}}_{k}^{\rm
iso}(t;J_{k})\,\ell_{t}^{\alpha}\,.$ (A.10)

Moreover, we have the following relations among ${\mathfrak{m}}^{\rm
iso}/{\mathfrak{s}}^{\rm iso}$ and ${\mathfrak{m}}/{\mathfrak{s}}$ from (5.5),
whose proofs are again immediate from (3.6) and hence omitted.

######  Lemma A.3 (${\mathfrak{m}}^{\rm iso}/{\mathfrak{s}}^{\rm