Chaos and random matrices in supersymmetric SYK

We use random matrix theory to explore late-time chaos in supersymmetric quantum mechanical systems. Motivated by the recent study of supersymmetric SYK models and their random matrix classification, we consider the Wishart-Laguerre unitary ensemble and compute the spectral form factors and frame potentials to quantify chaos and randomness. Compared to the Gaussian ensembles, we observe the absence of a dip regime in the form factor and a slower approach to Haar-random dynamics. We find agreement between our random matrix analysis and predictions from the supersymmetric SYK model, and discuss the implications for supersymmetric chaotic systems.


Introduction
A recent surge of interest in quantum chaos has revolved around a strongly-interacting quantum system called the Sachdev-Ye-Kitaev (SYK) model [1,2]. This model of N all-to-all randomly interacting Majorana fermions is solvable at strong-coupling and appears to be in the same universality class as black holes, exhibiting an emergent reparametrization invariance and an extensive ground-state entropy. More compellingly, the out-of-time order correlation function (OTOC) of the theory [1,3] saturates a universal bound on chaotic growth [4], a seemingly unique feature of gravity [5,6] and conformal field theories with a holographic dual [7]. The low-energy description of the theory in terms of a Schwarzian effective action also encapsulates dilaton gravity in AdS 2 [8,9]. This model should be seen as a valuable resource for understanding both black holes and quantum chaos.
There have already been a myriad of generalizations of the SYK model, including an extension by Fu, Gaiotto, Maldacena, and Sachdev, to a supersymmetric model of strongly interacting Majoranas [10], which has been further explored in [11][12][13][14][15]. The supersymmetric version of the model also displays many of the same holographic properties. Notably, at strong-coupling the theory has an emergent superconformal symmetry which renders it solvable and allows one to compute correlation functions. At low-energies the symmetry is broken, giving a Schwarzian-like effective action which mimics supergravity in AdS 2 [16]. Like its non-supersymmetric counterpart, the model has random matrix universality in its spectral statistics [17,18] and appears to exhibit thermalization in its eigenstates [19,20], both hinting at underlying chaotic dynamics.
Although we lack a precise definition of quantum chaos, there are still universal features one expects of quantum chaotic systems: most notably, having the spectral statistics of a random matrix [21]. Information scrambling [22,23] and chaotic correlation functions [5] have also been extolled as symptoms of chaos. Ideas from quantum information have helped make these notions more precise, quantifying how scrambling [24] and randomness [25] are encoded in OTOCs. Similarly, [26] explored the connection to random matrix dynamics, quantifying randomness and scrambling in the time evolution by random matrix Hamiltonians and computing a quantity called the frame potential. The onset of random matrix behavior can also be seen in the spectral form factor, which has been studied in the SYK model [27].
Motivated by this, we may ask the question: what are the universal features of supersymmetric SYK models, or more generally, of all supersymmetric quantum chaotic systems? And how do we quantify them from an information-theoretic standpoint?
To address this, we consider the Wishart-Laguerre ensembles, also termed random covariance matrices [28], which appeared in the random matrix classification of the supersymmetric SYK models [18]. Recall that the Hamiltonian in supersymmetric quantum mechanics is constructed as the square of a supercharge. Loosely speaking, the intution is that this random matrix ensemble arises from squaring the Gaussian random matrices, just as we might think of a chaotic supersymmetric system defined by a disordered supercharge. In this paper we consider the simplest Wishart-Laguerre ensemble, 1 the Wishart-Laguerre unitary ensemble (LUE), corresponding to supersymmetric quantum systems without additional discrete symmetries. In the following, we will quantitatively derive predictions for the spectral form factors, frame potential, and the out-of-time-ordered correlators, where a central distinction from the non-supersymmetric models arises in the spectral 1-point functions, which modifies the the early time decay of the spectral form factor. A slower decay in the LUE frame potential indicates less efficient information scrambling and the failure of the ensemble to become Haar-random. Our predictions for the LUE match those from the 1-loop partition function of the supersymmetric SYK model. The paper is organized as follows: In Section 2, we review the supersymmetric model and spectral form factor, discussing its universal features and behavior in SYK models. In Section 3, we review the basic tools in random matrix theory and then compute spectral form factors for the Wishart-Laguerre ensemble. In Section 4, we explore chaos this random matrix ensemble by computing the frame potentials and correlation functions, and comment on its complexity growth. In Section 5, we discuss chaos in supersymmetric SYK and compare with the random matrix predictions, concluding in Section 6. In Appendix A we present some numerical checks of our expressions.

Supersymmetric SYK model
We first briefly review the supersymmetric extension of the SYK model. For an in-depth discussion of the original model, see [3]. Consider N all-to-all interacting Majorana fermions ψ i with random couplings, which anticommute as {ψ i , ψ j } = δ ij . The (2q − 2)-point N = 1 supersymmetric model is constructed from the supercharge Q, a q-body Majorana interaction with odd q. The Hamiltonian is then given by the square of the supercharge as with Gaussian random couplings C i 1 ...iq of mean and variance and where J is a positive constant. We also define J as J 2 = 2 q−1 J 2 /q, with a slightly more convenient scaling in q.
In the large N limit, this model shares many of the same appealing holographic features as the SYK model, such as chaotic correlation functions, a zero-temperature entropy, and an emergent superconformal symmetry which is broken at low-energies, admitting a Schwarzianlike desciption [10]. We can compute the free energy at large N by evaluating at the saddle point, and at low temperatures find where s 0 is the zero-temperature entropy density and c is the specific heat. In the supersymmetric theory we have c = απ 2 /J with a constant α, which becomes c = π 2 /4q 2 J in the large q limit. The ground-state entropy density is computed to be s 0 = 1 2 log(2 cos π 2q ) and the ground state energy E 0 can be subtracted off.
The SYK model with N Majoranas enjoys a random matrix classification, where the symmetry class of the theory is dictated by a particle-hole symmetry [17,27]. Depending on N , the spectrum will display level statistics of one of the three Gaussian ensembles: GUE, GOE, or GSE. For the supersymmetric extension of SYK, we can similarly classify the random matrix behavior for a given number of Majoranas N , going beyond Dyson's classification to the extended 10-fold symmetry classification of Altland-Zirnbauer [30]. Understanding how anti-unitary symmetries act on the supercharge Q, we can identify the appropriate symmetry class [18]. The Hamiltonian, given as the square of the supercharge, then has random matrix description in terms of the Wishart-Laguerre ensembles. The level statistics are still those of the Gaussian ensembles, but the spectral correlations are different. Roughly, we can think of the supersymmetric SYK behaving like the square of Gaussian random matrices, which are the Wishart ensembles. For more details, see [18] as well as an extension of the classification to the N = 2 supersymmetric models [15].
Speaking generally, there a number of reasons one might wish to consider supersymmetric generalizations of SYK. For instance, much is understood about the low-energy physics in nearly AdS 2 spacetimes purportedly dual to the low-energy dynamics in SYK, but the exact holographic dual of the theory is not know. As many of the best understood examples of AdS/CFT are supersymmetric, one might hope that this particular construction might provide guidance on the correct UV completion of the SYK model. Less ambitiously, considering the supersymmetric models might be useful in contructing higher dimension analogs [12].

Spectral form factor
Quantum chaotic systems are often defined to have the spectral statistics of a random matrix. An object familiar in random matrix theory which exhibits these universal properties is the spectral form factor. We will introduce this object more precisely in our review of random matrix theory in Sec. 3.1, but the 2-point spectral form factor R 2 (t, β) can be given simply in terms of the analytically continued partition function and where the average · is taken over an ensemble of Hamiltonians (e.g. SYK, or some disordered spin system, or a random matrix ensemble). This object was discussed more recently in [27], where they studied the form factor in SYK and found that the theory revealed random matrix behavior at late times. From the bulk point of view, one motivation for studying this object was a simple version of black hole information loss [31]: 2-point functions appear to decay exponentially in terms of local bulk variables, whereas a discrete spectrum implies a finite late-time value. The same inconsistency is apparent in the spectral form factor. Some characteristic features of the time-evolved form factor R 2 (t), exhibited in both the SYK model and in random matrix theories, are: an early time decay from an initial value called the slope, a crossover at intermediate times called the dip, a steady linear rise called the ramp, and a late-time floor called the plateau. In Fig. 1 we observe these features in SYK. While the early time decay depends on the specific system, the ramp and plateau should be universal features of quantum chaotic systems. The ramp is characteristic of spectral rigidity: the long-range logarithmic repulsion of eigenvalues. The anticorrleation of eigenvalues causes the linear increase in the form factor. At late times, or at energy scales smaller than the mean spacing, the form factor reaches a plateau as degeneracies are rare and neighboring eigenvalues repel in chaotic systems.  Majoranas at inverse temperature β = 1, computed for 800 realizations of disorder. We observe the slope, dip, ramp, and plateau behaviors.

SYK form factor and GUE
Recently, [27] studied the form factor in SYK and found agreement with random matrix theory, showing analytically and numerically the aspects of the dip, ramp, and plateau of SYK agree with those of the Gaussian unitary ensemble (GUE), an ensemble of L × L random Hermitian matrices. We will avoid explicitly introducing and defining the original Majorana, instead simply mentioning a few details to better frame the discussion of the model's supersymmetric extension. The emergent reparamentrization invariance of SYK at strongly-coupled is broken spontaneously and explicitly at low-energies, yielding an effective description in terms of the Schwarzian derivative [1,3]. The 1-loop partition function of the Schwarzian theory Z Sch 1-loop ∼ e cN/2β /β 3/2 , can be analytically continued to β + it to study the form factor of SYK. At early times, R 2 (t, β) is dominated by the disconnected piece which gives a 1/t 3 power law decay, normalized by its initial value we have for times greater than t √ N when the time dependence in the exponent disappears and where c is the specific heat of the theory. To isolate this contribution, [27] considered a special limit (a 'triple scaled' limit) where only the Schwarzian contributes. Moreover, [32] showed that the Schwarzian theory is 1-loop exact and recieves no higher-order corrections, indicating that the power-law decay predicted by the Schwarzian should dominate the disconnected form factor for long times. 2 This power law decay is simply the Laplace transform of the statement that the spectrum has a square-root edge 3 Knowing the free energy in the large N limit, we can also show that the form factor of SYK transitions to a ramp at a dip time t d ∼ e N s 0 /2 , growing linearly until a plateau time of t p ∼ e N s 0 +cN/2β , where s 0 is the zero-temperature entropy density. Many of these features of the SYK form factor agree with the universal predictions from GUE. The form factor for GUE has been studied extensively in the random matrix literature [37][38][39] and references therein, and revisited more recently in the context of SYK and black holes in [26,27,40]. Simply stating the results, the early-time decay of the GUE form factor transitions to a linear ramp at a dip time of t d ∼ √ L, growing linearly until the plateau time t p ∼ L. We note that around the plateau time the ramp is not quite linear as nonperturbative effects become important as we transition to the plateau [41]. The non-universal early time decay also has the same power law 1/t 3 , due to the fact the Wigner semicircle law for Gaussian random matrices ρ(λ) = 1 2π √ 4 − λ 2 , also exhibits a square-root edge.

Supersymmetric SYK form factor
From the large N partition function of the supersymmetric theory, we can also make predictions as to the behavior of the spectral form factor. We will present a more explicit treatment of this in Sec. 5. At low-energies, the fluctuations around the large N saddle point of the supersymmetric theory break superconformal symmetry; the action for these reparametrizations is a super-Schwarzian [10], where the action integrates over τ and a superspace coordinate θ and the super-Schwazian acts just like the standard Schwarzian derivative except as a super-derivative, respecting a similar chain rule. The action gives a 1-loop partition function which differs in the 1-loop determinant from the SYK model. The super-Schwarzian theory is also 1-loop exact [32], ensuring its validity away from very early times. Analytically continuing the partition function β → β + it, disconnected piece of the form factor which dominates at early times, is exhibiting a 1/t decay in the slope, slower with than the decay in SYK. This can also be understood as the contribution from the edge of the spectrum, where the Laplace transform of the 1-loop partition function gives observing a square-root growth at the edges of the spectrum. As we discuss later, computing the ramp function for supersymmetric SYK, we find the ramp and slope intersect at a dip time t d ∼ e N s 0 , which is the same time scale as the ramp's transition to the plateau t p ∼ e N s 0 . The slow decay at early times means that the slope transitions to ramp behavior at the same time-scale as the plateau time, i.e. the ramp is hidden beneath the slope. We plot the 2-point form factor for the model in Fig. 2. Subtracting the disconnected contribution reveals the ramp in the connected form factor, also plotted. The lack of a dip in the supersymmetric model will have implications for our discussion of the frame potential and randomness.  Notation A brief comment on notation. In recent work studying the spectral form factor, the normalized 2-point form factor is often denoted as g(t, β), and its connected component as g c (t, β): While in [26], we denoted the 2-point form factor as R 2 (t, β), and more generally the 2k-th form factor as R 2k (t, β). Just to be clear For us, working directly with the numerator turns out to be more convenient when discussing the frame potential and correlation functions, and avoids subtleties regarding the appropriate or tractable normalization, i.e. 'quenched' vs 'annealed'.

Basic setup in random matrix theory
In this paper, we consider the Wishart-Laguerre Unitary Ensemble (LUE), an ensemble of L × L random matrices which can be generated as H † H, where H is a complex Gaussian random matrix with normally distributed complex entries drawn with mean 0 and variance σ 2 = 1/L. This is the 'physics normalization', where the spectrum does not scale with system size. 4 The joint probability distribution of LUE eigenvalues is given by where ∆(λ) is the Vandermonde determinant and the constant factor is define such that the distribution integrates to unity. One can think of LUE matrices as square of a Gaussian random matrix. More generally, we could define L × L Wishart matrices generated by L × L Gaussian matrices, where L ≥ L, which gives a slightly more general eigenvalue distribution. But given the supersymmetric Hamiltonians we consider defined as the square of the supercharge, we just consider Wishart matrices generated by square matrices with L = L . We average over the random matrix ensemble as The spectral density is given by integrating the joint probability P (λ) over L − 1 variables, More generally, we can define the k-point spectral correlation function by integrating over all but k arguments Recall that for the Gaussian ensembles, we may take the large L limit famously recover Wigner's semicircle law for the distribution of eigenvalues. Instead in the LUE, we take the large L limit and find [42] which is referred to as the Marčenko-Pastur distribution. Just as in the GUE, the LUE is a determinantal point process, which means the k-point spectral correlators are given by a kernel K as Demonstrating the universality of Dyson's sine kernel [43], the Wishart ensemble has sine kernel statistics in the large L limit [28,44], meaning where u is an arbitrary constant valued in [0, 4]. We will fix the value of u numerically. 5 The spectral form factor, defined as the Fourier transform of the spectral correlation functions, is a standard quantity to consider in random matrix theory; see [37] for an overview. We define the 2-point spectral form factor in terms of the analytically continued partition function Z(β, t) as 6 More generally, we can consider k-point spectral form factors could be defined as In the following subsections, we will compute the LUE spectral form factors and compare analytical results with numerical observations. At large L, we compute the spectral form factors by Fourier transforming the determinant of kernels in Eq. (3.6). We integrate the products of K as [37] The analogous constant in considering the GUE would be fixed to u = 0, given the symmetry of the spectrum. However, for the LUE u = 0 it is divergent. The value of u specifies the center of the two eigenvalues λi and λj. 6 This is slightly different than the standard presentation in the RMT literature, where the form factor is usually given as the Fourier transform of a connected form factor, called the cluster function. Here we work with both connected and disconnected pieces.
where the Fourier transform of the sine kernel is g(k) ≡ dr e 2πikr sin(πr) πr = 1 for |k| < 1/2 0 for |k| > 1/2 . (3.13) and where α L ≡ Lρ(u). The integral over the sine kernel is unbounded and can be treated by imposing a cutoff. We use the box approximation [26] α L dλ e i n j=1 k j λ → α L , (3.14) fixed such that Eq. (3.12) over the truncated range with k i = 0 integrates to L. This will be helpful in computing the higher-point spectral form factors, for instance, R 4 . It will also be convenient to define the following functions which will appear in computing the LUE form factors

Two point form factor at infinite temperature
Let us start with the simplest case, the two point spectral form factor at infinite temperature β = 0. Pulling out coincident eigenvalues, we have The determinant of kernels in Eq. (3.6) gives a squared 1-point function and 2-point function contribution. Using the integration formula in Eq. (3.12), we obtain in terms of the functions defined above, and where In Fig. 3, we plot the infinite temperature LUE 2-point form factor as derived in Eq. (3.17) along side the GUE form factor (see [26]). Note that unlike in the GUE case there is no dip or ramp. The lack of a intermediate time scale at which the initial slope decay transitions at the dip to a linear growth to a plateau, is due to the slow decay of the 1-point functions which gives the slope. Subtracting off the contribution from the 1-point functions defines the connected piece of the 2-point form factor which exposes the linear growth before the plateau. The connected components are also plotted in Fig. 3. The transition point in the function of r 2 is defined as the plateau time t p = 2πα L , where α L = Lρ(u). The value of 2πα L is not straightforwardly fixed given the unbounded support when integrating over kernels. The constant also determines the linear slope of the ramp function r 2 prior to the plateau. As we discuss in App. A, the constant u is fixed by numerically fitting to the ramp. We find a plateau time of t p ∼ πL/2 for the LUE 2-point form factor.
Using the asymptotic form of the Bessel function, we conclude that the disconnected piece decays at early times (for t much smaller than L but larger than O(1)) as This O(1/t) decay of the LUE form factor is to be contrasted with the slower O(1/t 3 ) decay in both the GUE and the SYK model [26,27]. However, the connected piece, dominated by the universal sine kernel in the large L limit, still sees the steady linear rise O(t) at intermediate time scales. This fact reaffirms the expectation that the decay in the disconnected piece, the Fourier transformed one-point functions, is model dependent. However, the ramp in the connected 2-point function is a universal feature of quantum chaotic systems. In addition to a hidden dip, another difference with the GUE result is the lack of an oscillating decay in the LUE at infinite temperature. In the GUE, the Bessel function decay at β = 0 gives a true dip time O(1). The envelope of this decay was what we considered as the decay to a dip given that a finite β smoothed out the oscillations.

Two point form factor at finite temperature
Now let us consider the two point form factor at the finite temperature. For small β, one may effectively insert the one point distribution in the integration formula. We walk through the computation in some detail as it will mimic the calculation of the supersymmetric SYK form factor in Sec. 5. To be concrete, we write simply integrating the kernels as specified above. For the final integral, we make the change of variables which allows us to compute where we regulate the unbounded integral with the insertion of ρ(u 1 ). The 2-point spectral form factor at finite temperature is We plot the analytic result in Fig. 4 and observe that at finite temperature there is still no clear dip time in LUE, unlike for the GUE, and that the plateau time t p does not depend on β. For the LUE, we define h 1 (β) ≡ r 1 (2iβ), a purely real function of the inverse temperature, with the plateau value R 2 (t p , β) = h 1 (2β)L .  At small but finite β we have compared to the GUE result 1 + 2β 2 + O(β 4 ) [26], one can see that the LUE plateau value is smaller than GUE, which we also observed in numerics.

Four point form factor at infinite temperature
As an example of a higher point form factor, we compute the 4-point R 4 at infinite temperature. By definition we have To evaluate the expression we must consider all possible ways in which the eigenvalues can collide in the sum, i.e. all equal, λ i = λ j , λ k = λ , etc, and treat them separately. Making use of the 2-point form factors we derived above, and computing the 3 and 4-point function contributions by expanding the determinant and integrating products of kernels as Eq. (3.12), we obtain (3.29) In the large L limit, some of the terms above are subdominant or suppressed in L at all times, allowing us to simplify the expression as R 4 (t) ≈ L 4 |r 1 (t)| 4 + 2L 2 r 2 2 (t) − 4L 2 r 2 (t) + 2L 2 − 7Lr 2 (2t) + 4Lr 2 (3t) + 4Lr 2 (t) − L , (3.30) similar to the result we derived for the GUE [26]. At times much earlier than the plateau time, we have Again, we find a slow decay of O(1/t 2 ) and thus no visible dip at large L. The plateau time is still 2πα L , with a plateau value R 4 (t p ) = 2L 2 − L ∼ 2L 2 .

Chaos and Wishart matrices
We want to study the chaotic nature of time-evolution by LUE Hamiltonians. Consider the ensemble of unitary time-evolutions generated by LUE random matrices We want to understand how random LUE time-evolution is by asking when the ensemble forms a k-design. Computing the frame potential for the ensemble will quantifies a distance to Haar-randomness. We also compute correlation functions of operators evolved by the LUE to look at early-time chaos in the chaotic decay of 2k-point functions.

QI overview
Before discussing the frame potential and measures of chaos for the random matrix ensemble, we will briefly overview the quantum information theoretic concepts and tools we use, namely the notion of a unitary k-design and the frame potential. For a more in-depth review of these in the context of information scrambling in chaotic systems, see [25,26]. For a finite dimensional quantum mechanical system, with Hilbert space H of dimension L, the unitary group U (L) can be equipped with the Haar measure, the unique left/right invariant measure on U (L). Given some ensemble of unitary operators E, we say that the ensemble forms a unitary k-design if it reproduces the first k-moments of Haar for any operator. More intuitively, we should think of this as capturing how random the ensemble is, in that the ensemble is sufficiently spread out over the unitary group to reproduce its statistics. A precise measure of Haar-randomness is the frame potential [45], which measures the 2-norm distance between the k-th moments of an ensemble E and Haar. The k-th frame potential is defined with respect to an ensemble E as The frame potential for any ensemble E is lower bounded by the Haar value Haar , (4.4) with equality iff E forms a k-design. The k-frame potential for the Haar ensemble is simply F (k) The frame potential appeared in the context of information scrambling and black holes as the average of all out-of-time ordered correlators [25] 1 L 4k where "B(t)" = U BU † and U ∈ CE for any ensemble of unitaries E, with each A i and B i summed over all Pauli operators. This makes precise an approach to randomness, where the chaotic decay of correlators at late-times means the frame potential becomes small and the ensemble forms a k-design.

Frame potentials
First frame potential at β = 0 We start by computing the first frame potential at infinite temperature F (k) E for the ensemble of LUE time-evolutions. Following [26], we have (4.6) Using the unitary invariance of the ensemble and integrating using the second moment of the Haar ensemble, we find with the same dependence on the form factors as in the GUE case. In Fig. 5 we plot our analytic form of the first frame potential of the LUE at infinite temperature. We can see that there are significant differences between the supersymmetric and non-supersymmetric cases. The slow decay of the LUE means there the ensemble does not form a k-design at the dip. At late-times after the plateau, we find the frame potential approaches a value of 2.
First frame potential at finite β We can also generalize the frame potential to finite temperature by averaging over all thermal 2k-point functions with operators spaced equidistant on the thermal circle (i.e. inserting ρ 1/2k between operators in the 2k-OTOC). Averaging over operators, we find [25] F (k) Tr(e −βH 1 )Tr(e −βH 2 )/L 2 (4.8) with the normalization that gives the standard frame potential as β → 0. For the LUE, we compute the finite temperature frame potential just as above, Haar integrating to find where, we define a slightly more conveniently normalized form factor As it is more analytically tractable, we opt to separately average the numerator and denominator (the 'quenched' version), and checked numerically that the results are in good agreement. We see that at early times, near t = 0, we have the β-dependent value While at late times, after the plateau time, we have F Second frame potential at β = 0 The second frame potential for the LUE at infinite temperature in expressed in terms of the spectral form factors as [26] F where we have defined The 4-point form factor with two coincident eigenvalues, R 4,2 (t), is simply R 2 (2t). The 3-point form factor R 4,1 (t) for the LUE can be computed just as in Sec. 3, where we find R 4,1 (t) = L 3 Re r 1 (2t)r * 2 1 (t) − L 2 Re(r 1 (2t))r 3 (2t)r 2 (t) − 2L 2 Re(r * 1 (t))r 3 (t)r 2 (2t) + L 2 |r 1 (2t)| 2 + 2L 2 |r 1 (t)| 2 + 2Lr 2 (3t) − Lr 2 (2t) − 2Lr 2 (t) + L . (4.14) We plot the second frame potential for LUE alongside the first frame potential in Fig. 5. The second frame potential has an initial value of L 4 and late-time value of 10, just as for the GUE. But again the difference arises at intermediate time scales, where the LUE fails to form a k-design.

Correlation functions
As we discussed before, the recent interest in quantum chaos has involved extensive discussion of out-of-time order correlation functions (OTOCs). Namely, the following 4-point functions of pairs of operators in thermal states We consider OTOCs with operators evolved by LUE Hamiltonians and averaged over the random matrix ensemble. In [26], we studied 2k-OTOCs and related them to spectral quantities, both by averaging over the operators in the correlation function or over an ensemble of Hamiltonians. In that work, we averaged 2k-OTOCs over the GUE and related the correlators to spectral quantities using the unitary invariance of the measure. As the LUE is similarly invariant, the relation between correlation functions averaged over the random matrix ensemble and the form factors will be the same as thus parts of the discussion here will closely follow [26]; the differentiating aspects of LUE time-evolution thus lie in the spectral form factors themselves. First we look at the 2-point function and integrate over Hamiltonians drawn from the LUE, using the unitary invariance of the measure and Haar integrating in the eigenvalue basis where AB c denotes the connected correlator. For non-identity Paulis, the expression is nonzero for B = A † , and thus LUE average : We next compute the 4-point OTOC averaged over the LUE, using the fourth moment of Haar and looking at the leading order behavior (4.21) The characteristic time-scale for decay of LUE 2-point functions is t 2 ∼ O(1), or for systems at finite temperature O(β). The time-scale for 4-point function decay is also order 1, but faster than the decay of 2-point functions t 4 ∼ t 2 /2. Although the decay is slower than for GUE, unsurprisingly, the conclusion about the LUE's perception of early-time chaos is the same: the LUE 4-point OTOCs decay faster than the LUE 2-point functions, which means the random matrix ensemble fails to describe scrambling at early times.

Complexity
Lastly, we briefly comment on the complexity growth under time-evolution of LUE Hamiltonians. Here we simply discuss the results; details and definitions of ensemble complexity and its relation to the frame potential are given in [25,26]. The gate complexity of an ensemble E, i.e. the number of gates needed to generate E, is lower bounded by the frame potential as At early times before the dip time t t d , the dominant contribution to the k-th frame potential is F [26]. For k L, the 2k-th form factor goes as R 2k ∼ r 2k 1 , the function defined in Eq. (3.15) in terms of Bessel functions. The decay r 2 1 ∼ 1/t, gives a lower bound on the growth of the circuit complexity where the slower decay for LUE still gives the same logarithmic lower bound as GUE. Interestingly, in GUE the 1-point function contribution to the form factor at early times is an oscillating Bessel function decay J 2 0 (2t)/t 2 , which formally gives a dip time O(1). As these oscillations are not present in the LUE, we can bound the complexity up to the dip time even at infinite temperature. But for large k, we recover the quadratic growth of complexity: C ≥ t 2 / log n, hinting again at the unphysical nature of LUE evolution at early times.

Chaos in supersymmetric SYK
The supersymmetric SYK model admits a classification by Wishart-Laguerre random matrix ensembles and has a density of states which closely follows a Marčenko-Pastur distribution [18]. Having discussed the properties of LUE random matrices, we turn to the supersymmetric SYK model and check that the form factor acts similarly. From the frame potential, we then discuss the Haar-randomness of the model's time evolution.
Assuming that the spectral statistics of the theory are Gaussian, as both SYK and the Wishart matrices are, allows us to use the sine kernel to compute the spectral n-point functions. We note that if the statistics are GUE/GOE/GSE, the sine kernel is slightly modified and the ramp function differs as we approach t ∼ L, but the universal growth of the ramp is still present. Knowing that the supersymmetric SYK model has Gaussian spectral statistics [18], we can compute the finite temperature form factor for the theory just as in Eq. (3.22), and find where r 2 (t) is the ramp function from the LUE and we define E = 1 2 (λ 1 + λ 2 ). Continuing, we find the finite temperature form factor As a sanity check, the late-time value Z(2β) here matches the infinite-time average of the spectral form factor. As we discussed in Sec. 2, the 1-loop partition function from the super-Schwarzian theory is where s 0 is the ground-state entropy density and c is the specific heat. At early times, the form factor is dominated by its disconnected component, decaying as 1/t Early : for times greater than t ∼ √ N = log L/2, but shorter than t ∼ √ L. Computing the connected form factor, we find Equating the 1/t decay with the ramp gives a dip time t d ∼ e N s 0 , the same order as the plateau time t p . Even in light of the exactness of the super-Schwarzian theory, we should be cautious in extrapolating to very late times. It is possible that in the large N theory the slope is not well-described by the effective theory at late times and, in turn, decays faster at an intermediate time scale. Lastly, to get a hint at the nature of scrambling and an approach to randomness in SYK and its supersymmetric extension, we numerically plot the first frame potential for each in Fig. 6 at infinite temperature and for N = 16 Majoranas. The faster decay and dip that appears for SYK means the frame potential decays quickly, forming an approximate k-design at the dip time. Although the dip value of the SYK frame potential for N = 16 is larger than the Haar value, we checked that as we increase N the dip value decreases and expect that SYK forms an approximate k-design in the large N limit. The frame potential for the supersymmetric model exhibits a much more gradual approach to its minimal value which is larger than in SYK, indicating less effective information scrambling and a greater distance of the ensemble to forming a k-design. It would be interesting to see, either numerically or analytically, if these behaviors persist at large N . Both theories, like their random matrix counterparts, become less random and increase after the dip, deviating further from an approximate design, which suggests that k-invariance [26] might provide a better insight in how information scrambles in SYK models.

Conclusion and outlook
In this paper, we considered the Wishart-Laguerre unitary ensemble in order to understand universal features of supersymmetric quantum mechanical systems. We computed the 2-point spectral form factor for the LUE and found the one-point function contribution gives a 1/t power law decay at early times, hiding the dip and transitioning directly into the plateau. This is relatively slow compared to the ∼ 1/t 3 decay seen in both SYK and the GUE. The universal ramp behavior from the sine kernel can be seen in the connected LUE 2-point form factor. These results agrees with the prediction from the 1-loop partition function in supersymmetric SYK. This slow decay implies the onset of a random matrix description occurs at much later times. This can best be seen from the frame potential, where we find a more gradual decay to Haar-random dynamics. Moreover, the frame potential for the LUE, unlike that of the GUE, does not reach the Haar value and does not form an approximate k-design. This is also what we predict and observe numerically in the supersymmetric SYK model, where the slower decay and larger dip value imply less effective information scrambling.
The supersymmetric model, while maximally chaotic, sees a slower onset of random matrix behavior-made evident by the lack of a dip in the form factor and by the slow approach to Haar-randomness in the frame potential. The apparent distinction here between early-time chaos, in terms of chaotic correlation functions, and late-time chaos, in terms of scrambling and Haar-randomness, demands a deeper understanding.
Note added: In the preparation of this draft, [46] appeared which also considers the infinite temperature 2-point spectral form factor for Wishart matrices in a different context. Namely, they study the statistical properties of the reduced density matrix on spatial regions in quantum many-body systems. They also comment on universal features of Wishart matrices in Floquet systems. As there is a sense in which Floquet systems may be thought of as supersymmetric quantum mechanics [47], where the Floquet unitary is built from two 'supercharges', it would be interesting to explore further connections with our work. JL acknowledges support from the U.S. Department of Energy, Office of Science, Office of High Energy Physics, under Award Number DE-SC0011632.

A Numerics
In this appendix we discuss numerics to fix an analytic form of the form factors for LUE and to further provide checks on the expressions we derived for the form factors and frame potentials. As we mentioned in Sec. 3, there was a free parameter u in the expressions we derived for the k-point form factors. This dependence appears in the ramp function r 2 (t), defined in Eq. (3.15), and determines both the slope of the linear ramp in R c 2 (t) and the plateau time. Numerically computing the connected 2-point form factor for L = 500, we fix u by fitting the ramp between times ∼1 and √ L/2. We know that the early time behavior of the ramp is quadratic before t ∼ 1 and expect a loss of analytic control as we approach the plateau time. We thus linearly fit points in this intermediate regime and find u = 1.156. We hope to derive this result more rigorously in the future. We also present some numerical checks of our expressions for the LUE 2-point form factor in Fig. 7, where we find good agreement in the slope, ramp, and plateau. Our results were derived for LUE at large L and thus should capture the perturbative behavior. But in the transition to the plateau, nonperturbative effects [41] become important and our results deviate from numerics in this regime. After the plateau time, we return to contributions from the 1-point function. At very early times, before t ∼ O(1), the connected component grows as R c 2 (t) ∼ t 2 . This quadratic growth can be derived from an impressive integral representation of the connected 2-point form factor [39]. We have also checked our expressions of the finite temperature and higher point LUE spectral functions and found good agreement with numerics.