Unequal time correlators of stochastic scalar fields in de Sitter space

The quantum fluctuations of a test scalar field on superhorizon scale in de Sitter spacetime can be described by an effective one-dimensional stochastic theory corresponding to a particular class of nonequilibrium dynamical systems known as the model A. Using the formulation of the latter in terms of a supersymmetric field theory, we compute various unequal-time correlators at large (superhorizon) time separations and compare with existing quantum field theory computation. This includes perturbative calculations, pushed here up to three-loop order, and a nonperturbative 1/N expansion at next-to-leading order. Exploiting the supersymmetry of the stochastic theory, we also derive a spectral representation of the field correlators and a fluctuation-dissipation relation for the infrared modes of the scalar field in de Sitter spacetime.

The quantum fluctuations of a test scalar field on superhorizon scale in de Sitter spacetime can be described by an effective one-dimensional stochastic theory corresponding to a particular class of nonequilibrium dynamical systems known as the model A. Using the formulation of the latter in terms of a supersymmetric field theory, we compute various unequal-time correlators at large (superhorizon) time separations and compare with existing quantum field theory computation. This includes perturbative calculations, pushed here up to three-loop order, and a nonperturbative 1/N expansion at next-to-leading order. Exploiting the supersymmetry of the stochastic theory, we also derive a spectral representation of the field correlators and a fluctuation-dissipation relation for the infrared modes of the scalar field in de Sitter spacetime.

I. INTRODUCTION
The dynamics of quantum fields in an expanding spacetime is a subject of primordial importance for cosmology and inflation. In this context, the usual approach is a semi-classical treatment where spacetime is treated classically and interacts with a matter content which is of quantum nature and can have a possible backreaction on the geometry [1]. De Sitter spacetime is of particular interest both physically, as it is a good approximation of the inflationary phase, and mathematically, because of its high degree of symmetry.
The computation of quantum corrections in the presence of interactions is a lot more complicated in a curved background as usual pertubative tools are not always available. One particular setup in which non trivial effects arise is the case of light scalar fields in the expanding Poincaré patch of de Sitter spacetime, particularly relevant for inflationnary cosmology. The scalar fields mode functions are significantly modified by the curvature with, in particular, a strong amplification of the infrared modes, which can be viewed as intense particle production from the gravitational field [2][3][4]. This effect is at the origin of infrared and secular divergences in loop computations that limit the use of perturbation theory [5,6].
A variety of non perturbative treatments exists to address the question of the nonlinear effects (e.g. selfinteractions), see Refs. [4,[6][7][8][9][10][11][12][13][14][15][16][17][18][19][20][21][22][23] for various examples. The most prominent one is certainly the stochastic approach, developed in Ref. [6]. It gives an effective description of the dynamics of the infrared, long wavelength, modes in terms of an effective Langevin equation. The infrared modes of the scalar fields behave classically as a result of the aforementioned gravitational amplification and experience a random noise which encodes the effect of the ultraviolet modes crossing the horizon during expansion. The Langevin dynamics can be treated through the equivalent Fokker-Planck equation. This gives access, for example, to the late-time, equilibrium probability dis-tribution for the fields, from which one can compute various equal-time correlators, often analytically for simple enough potentials. Unequal-time correlators or genuine nonequilibrium properties, which contain important information about the long time/distance properties of the theory (dynamical time scales, spectral indices, etc.), are more difficult to access analytically and even numerically in some situations. For instance, for a simple quartic potential, the cases of vanishing or of negative square mass are intrinsically nonperturbative.
The stochastic Langevin equation is a particular case of the so-called model A in the Halperin et al. classification of nonequilibrium dynamical systems [24]. In the present article, we shall use tools developed in this context to compute various unequal-time correlators at large time separation, which gives access to different autocorrelation and relaxation time scales. In the stationary state, the problem can be formulated as a supersymmetric onedimensional field theory [23,25,26], free of ultraviolet divergences and which is a lot easier to manage than the original D-dimensional quantum field theory (QFT).
This one-dimensional field theory gives analytic access to properties of the stationary state reached by the scalar fields in the late-time limit. Diagram resummations, previously performed in the complete four dimensional field theory [18,27], can be done here in a simpler way which reproduces the leading infrared behavior. We compute various correlators in two approximations schemes. First, in a perturbative expansion in the self-interaction coupling constant, which is, however, limited to not too light fields. The second approximation scheme is the 1/N expansion, where N is the number of scalar fields. The latter allows to consider the interesting case of massless fields and of a symmetry breaking potential [27,28].
Along with the path integral formulation of the model A comes some interpretation of the different correlators and specific relations which are usually formulated in a statistical physics language. Using this analogy allow us here to reformulate these results in terms of our particular model and discuss some consequences for the scalar field correlator.
After briefly reviewing the effective stochastic approach, we present the functional formulation of the Langevin equation and discuss the supersymmetry of the resulting field theory in Sec. II. Various properties of the field correlators, independent of any approximation scheme are discussed in Sec. III. Our calculations in the perturbative and the 1/N expansions are presented in Sec. IV. We conclude in Sec. V. Additional calculations and technical details are presented in the various appendices.

II. GENERAL SETUP
We briefly recall the effective stochastic theory for the superhorizon modes of light scalar fields in de Sitter spacetime and review the functional formulation of the resulting one-dimensional model A as a supersymmetric field theory. We consider an O(N ) symmetric scalar field theory on the expanding Poincaré patch of a D-dimensional de Sitter spacetime with d spatial dimensions (D = d+1). The metric reads ds 2 = −dt 2 +a(t)d x 2 , with a(t) = e Ht , where t is the cosmological time and we set the Hubble rate H = 1. The classical action reads whereφ 2 =φ aφa and x denotes the appropriate, invariant integration measure.

A. Effective stochastic approach
For light fields in units of H, the (quantum) fluctuations of long wavelength, superhorizon modes are well described by the effective Langevin equation [6] where the dot denotes a time derivative and we used the notationV ,â = ∂V /∂φ a . Here, the infrared fieldŝ ϕ a , spatially smeared over a Hubble patch, effectively behave as classical stochastic fields whose fluctuations mimic those of the long wavelength modes of the original quantum fields. Those stochastic fluctuations are driven by the random kicks from the (quantum) subhorizon modes which cross the horizon at a constant rate due to the gravitational redshift. This is represented by the noise termξ a , whose stochastic properties reflect the quantum state of the system. For the Bunch-Davies (BD) vacuum, and treating the ultraviolet modes in the linear approximation, one finds [6] ξ a (t, with Ω n = 2π n/2 /Γ(n/2) and where the function F reflects the spatial smearing: it can always be normalized as F(0) = 1 and it vanishes rapidly for spatial separations | x − x | 1. Its precise form depends on the smearing procedure. Within a single Hubble patch, F ≈ 1 and the time evolution of the infrared fields is described by an effective one-dimensional Langevin equation with a Gaussian white noise. At sufficiently late times, the system is driven towards a stationary regime where, e.g., the equilibrium distribution of field values is given by The latter describes the equal-time statistical properties of the stochastic process and reflects the quantum fluctuations of the infrared modes of the original quantum fields in the BD vacuum. It can be seen as the Boltzmann distribution for a thermal system. Introducing the Hamiltonian for the superhorizon field in the Hubble patch under consideration asĤ Hubble patch (of radius H −1 = 1), the distribution (4) reads P ∝ e −βĤ , with β = Ω D+1 /V d = 2π the inverse Gibbons-Hawking temperature [29]. It is useful to rescale the variables so as to absorb the various volume factors. Defininĝ we get, denoting the time derivative with a dot, This is a particular, (0+1)-dimensional case of the model A in the classification of Halperin et al. [24], which has been widely studied in the context of out-of-equilibrium statistical physics. It can be given an elegant functional formulation by means of the Janssen-de Dominicis (JdD) procedure [25,26], which provides an efficient starting point for implementing various field techniques [26]. Recent examples in the present context include diagrammatic methods [30,31] or renormalization group techniques [23]. We now briefly review the JdD procedure.

B. Path integral formulation
The expectation value of an operator O(ϕ) can be formally expressed as where ϕ ξ is a solution of solution of Eq. (6) with given initial conditions and is the normalized probability distribution of the noise, with t = +∞ −∞ dt. In general, one should also average over initial conditions in Eq. (8). However, the latter become irrelevant if we restrict our considerations to the stationary regime. Assuming the uniqueness of the solution of Eq. (6) for a given realization of the noise (and given initial conditions), one writes where J [ϕ] = |Det[δ ab ∂ t + V ,ab ]| is the appropriate functional Jacobian. Under the above uniqueness assumption, one can forget the absolute value on the determinant and exponentiate the latter in terms of Grassmann fields Similarly, one exponentiates the functional delta as where the so-called response fieldsφ a are purely imaginary. Integration over the Gaussian noise ξ a finally gives, up to an irrelevant constant factor N , with the following action (14) This one-dimensional statistical field theory with 4N fields describes the leading infrared behaviour of the underlying QFT in de Sitter spacetime. Alternatively, we can use a more symmetric form of the action by changing the variableφ a → F a = i(φ a −φ a ). The action rewrites as where we neglect the boundary term t 2φ a V a = tV in the stationary state. This form of the action makes clear another link, namely, it relates to a supersymmetric quantum mechanics after the Wick rotation t → iτ [32].

C. Supersymmetry
The action (14) or, equivalently, (15), possesses various symmetries, such as the time-translation and the timereversal symmetries of the stationary regime, which can be conveniently encoded in a supersymmetry that mixes the bosonic and fermionic degrees of freedom [26,32]. To exhibit the latter, it is convenient to recast the various fields into the superfield Φ a (t, θ,θ) = ϕ a (t) +θψ a (t) +ψ a (t)θ +θθF a (t), (16) living on the superspace (t, θ,θ), with Grassmann directions θ andθ. The generators of the supersymmetry can be written as Q = i∂θ + θ∂ t andQ = i∂ θ +θ∂ t , and the covariant derivatives D = i∂θ − θ∂ t ,D = i∂ θ −θ∂ t allow to write the action in the following form with 1 z = (t,θ, θ), dz = dt dθdθ and K = 1

III. GENERAL PROPERTIES OF THE CORRELATOR
The general form of the superfield correlators is constrained by various considerations, most prominently the symmetries and causality. In this Section, we detail the case of the connected 2 two-point correlator, with the notation For simplicity, we consider a single field (N = 1). The generalization to arbitrary N is trivial.

A. Supersymmetry constraints
The dependence of the inverse propagator Γ (2) and the propagator G on the Grassmann variables are strongly constrained by the supersymmetry of the action. First, the anticommutator Q,Q = 2i∂ t generates the timetranslation invariance, so that it proves more convenient to work in frequency space and similarly for the two-point vertex Γ 12 (ω). The general dependence of the latter in the Grassmann variables involves a priori six independent functions: Supersymmetry implies the Ward identities where the numerical index indicates the Grassmann variable each operator Q orQ is acting on. These yield four independent constraints which are solved as Renaming A(ω) = η(ω) and B(ω) = iγ(ω), the general structure of the two-point vertex is [26,32] where the two Grassmann structures denote, respectively, the Dirac function in Grassmann coordinates and the supersymmetric d'Alembertian op- The superfield propagator is obtained by inversion, 12 (ω)G 23 (ω) = δ 13 , with 2 = dθ 2 dθ 2 , and reads Using the decomposition (16) of the superfield, we obtain the various correlators 3 . Now, from the path integral representation (13), we see that both ϕ(t)ϕ(t ) and ϕ(t)φ(t ) are real (despiteφ being imaginary) and thus G ϕϕ (t) ∈ R and G ϕF (t) ∈ iR. Using also the permutation identity of the superfield correlator, G 12 (t) = G 21 (−t), we conclude, in frequency space, that both the functions γ(ω) and η(ω) are real and even.

B. Fluctuation-dissipation relation
The stationary, equilibrium state of the system is characterized by a fluctuation-dissipation relation which directly follows from the above constraints. This relates the statistical correlator G ϕϕ (ω) (fluctuation) to the response function 4 G ϕφ (ω) (dissipation) or, more precisely, to the stochastic spectral function ρ, which we now introduce. The response function is given by and we define the stochastic spectral function as where the second equality follows from Eqs. (31)-(33).
In real time, this reads This is the announced fluctuation-dissipation relation characteristic of a thermal state in the high temperature (classical field) regime as discussed in the Appendix A. An interesting consequence of the above relation is the exact identity which can be proven as follows. In the limit 5 t → 0 + , we have, using the relation (35) and Eq. (6) The equal-time average in the last equality can be computed with the one-point equilibrium distribution (4) with the proper rescaling (5). The result (36) follows from the identity obtained after integration by parts.

C. Causality
Further interesting information can be obtained from causality. The latter implies, in particular, that the response function vanishes identically for negative times, [26]. From the definition (34) of the spectral function and the fact that G ϕφ (t) ∈ R, we easily deduce that ρ(t) = G ϕφ (t) − G ϕφ (−t) and thus that 6 or, equivalently, in frequency space, which implies that G ϕφ is analytic in the upper half complex frequency plane. Also, using the fuctuation-dissipation relation (34), we deduce This yields an exact expression for the so-called 7 dynamical mass m dyn , which measures the amplitude of the equal-time fluctuations of the stochastic field within a Hubble patch as Using Eqs. (33) and (41), we deduce Such a relation is reminiscent of the concept of screening mass, or susceptibility in thermal (quantum/statistical) field theory, which are related to the value of the (inverse) propagator at vanishing momentum and frequency and typically measure the overall response of the system to a static perturbation. These are to be distinguished from the so-called pole masses, or correlation lengths, which are associated to the poles of the response function and describe correlations between different spacetime points. The latter have their analogs in the present stochastic model, which we now discuss.

D. Mass hierarchy
Using the Fokker-Planck formulation of the Langevin equation (6), one shows that the unequal-time (connected) correlator for a given local function A(ϕ) of the field can be written as [6,34] where the Λ n, 's are the eigenvalues of the (properly rescaled) Fokker-Planck operator and the C A n, 's are appropriate coefficients. Because of the O(N ) symmetry, the latter can be labelled in terms of the eigenvalues ∈ N of the N -dimensional angular momentum and another possible index n. In the case of a quadratic potential, the latter is a single positive integer and the possible values of are constrained such that n − is even and positive. We expect this to remain true for λ = 0.
The eigenvalues are nonnegative real numbers 8 . Of course, some C A n, may vanish, e.g., due to symmetry selection rules [34]. For instance, the case A = ϕ only involves the vector channel = 1, so that the only nonvanishing coefficients in the decomposition (44) are C ϕ 2n+1,1 . Similarly, for the composite field χ = ϕ 2 /(2N ) in the scalar ( = 0) channel, the only nonvanishing terms are C χ 2n,0 . The correlations of various quantities of interest at large time separations are thus governed by the lowest eigenvalues contributing to the decomposition (44).
Below, we shall compute the ϕϕ and χχ correlators in various approximation schemes, from which we can extract the eigenvalues Λ 2n+1,1 and Λ 2n,0 , respectively, at each approximation order. Introducing the redefinitions ,0 ), and the following notation for the tree-level correlator of a field of mass m The eigenvalues Λ n, and the coefficients c ϕ,χ n are directly obtained as the poles and residues of the relevant response function, e.g., An obvious relation is which directly follows from the definition (42). Another constraint on the coefficients c ϕ 2n+1 is the following sum rule which directly follows from Eqs. (35) and (36).

E. Effective noise correlator
Finally, we mention that the η component of the selfenergy (27) can be interpreted as the effective noise correlator dressed by the nonlinear effect of the infrared modes themselves. Indeed, as recalled in the Appendix A, the general expression of the correlator of a Langevin process with a colored noise is, in frequency space, Using the exact relations (31) and (33), we deduce that can be interpreted as an effective colored noise kernel as announced. The tree-level exppression η free (ω) = 1 corresponds to the white noise contribution (7) from the ultraviolet modes in the present effective stochatic theory. As we shall see below, nonlocal loop corrections bring a nontrivial frequency dependence which corresponds to the effective dressing of the noise kernel from the nonlinear infrared dynamics.

IV. EXPLICIT CALCULATIONS
We now turn to explicit computations of the ϕϕ and χχ correlators in two approximation schemes previously studied in the D-dimensional QFT [18,27], namely, the perturbative expansion and the 1/N expansion. We consider an O(N )-symmetric scalar theory with quartic self interaction, whose superpotential is given by There is no possibility of spontaneously broken symmetry in the present low dimensional system [17,35,36]. We thus have Φ a = 0 and G ab 12 (ω) = G 12 (ω)δ ab , including in the case m 2 < 0.
where the first two terms on the right-hand side correspond to the free field case. We denote the tree-level superpropagator for a field with mass m as We also introduce the supercorrelator of the composite field X = Φ 2 /(2N ), which, in the free theory, is simply given by the oneloop diagram of Fig. 1. This is easily computed as [see Eq. (B11)] The component at θ 1,2 =θ 1,2 = 0 is From the decompositions (46) and (47) and the free-field expressions (45) and (59), we read Λ free 1,1 = m 2 , c free 2n+1 = δ n,0 , Λ free 2,0 = 2m 2 , and c free 2n = δ n,1 /(2N m 2 ). This agrees with the known spectrum of the free case, which is just that of a O(N )-symmetric harmonic oscillator [6,34] Λ free n, = nm 2 .
A. The perturbative expansion We first compute the self-energy at two-loop order in a perturbative expansion (the three-loop order is computed in Appendix C). The relevant diagrams are shown in Fig. 2. Their explicit evaluation is straightforward and we shall only give the resulting expressions here. The details can be found in Appendix B. The one-loop contribution, diagram (b), yields which simply corresponds to a constant shift of γ(ω), that is, a mere mass renormalization. The same is true for the two-loop local 9 contribution given by diagram (c) in Fig. 2, which reads Σ (c) A nontrivial frequency dependence appears with the nonlocal contribution, diagram (d), which can be written as Altogether, we obtain, for the functions γ and η in Eq. (27), where we have introduced the dimensionless couplinḡ and the renormalized mass We immediately obtain the expression of the dynamical mass as As explained in Sec. III, the relevant mass hierarchy can be directly read off the response function. Using the expressions (64) and (65), the latter can be written as with the poles given by and the residues 9 Here, local means that both external legs are attached to the same vertex. In particular, we verify the sum rule (50) at this order. The two-pole structure (69) at the present order of approximation precisely coincides to the splitting of the propagator obtained in the QFT calculation of Ref. [18], which reads with G m 2 given in Eq. (45). The expressions of the various masses and coefficients exactly agree, with the iden- , and with the rescaling (5), that is, 10 We now come to the two-loop correction to the χχ correlator, given by the two diagrams in Fig. 3. The diagram (e) simply corresponds to the effect of the one-loop mass renormalisation of one propagator line (the same is true for the diagram (c) of Fig. 2) and can be easily computed. Equivalently, we can treat this diagram with the following trick [18,27]. We implicitly include it in the one-loop diagram (a) of Fig. 1 by using effective propagator lines with an effective mass M . We then replace the latter by its expression (67) and systematically expand at the relevant order of approximation.
Each loop in the diagram (a) of 1 and the diagram (b) of 3 is given by Eq. (58), with m 2 → M 2 and the sum reads and extracting the component at vanishing Grassmann variables, we obtain, in terms of the renormalized mass (78) 10 In particular, the quantity namedλ in Ref. [18] is the same as here. In the last equation, we have used the knowledge of the general structure (44) of the correlator to resum the twoloop correction to the propagator in the appropriate form (i.e., a correction to the corresponding self-energy). We can directly read off the expressions We note that the perturbative calculation of the propagator at orderλ 2 only gives access to the leading-order expression of the infrared subleading eigenvalue Λ 3,1 because the corresponding coefficient c ϕ 3 is, itself, of order λ 2 . It is interesting to push our perturbative calculation to three-loop order so as to obtain the first correction to Λ 3,1 and compare to the perturbative results of Ref. [34] obtained by directly solving the Fokker-Planck equation. We present this calculation in the Appendix C. The three-loop expressions of m 2 dyn , Λ 1,1 , c ϕ 1 and c ϕ 3 can be found there. Here, we simply gather the next-toleading results for the lowest eigenvalues: which reproduce (and generalize to arbitrary D and N ) the perturbative results of Ref. [34] for D = 4 and N = 1 (in that case, Λ n, = Λ n ). The present perturbative calculations are controlled by the dimensionless expansion parameterλ ∝ λ/m 4 and are thus invalid in the zero mass limit as well as in the negative square mass case. These cases require a nonperturbative treatment, such as the 1/N expansion, studied in Ref. [27] in the QFT context and that we now describe in the present stochastic framework.

B. The 1/N expansion
We closely follow Ref. [27] for the diagrammatic formulation of the 1/N expansion, which we adapt to the = + • + . . . present (supersymmetric) theory. In particular, we separate the local and nonlocal contributions 11 to the selfenergy Σ and grab the former in an effective square mass M 2 , which satisfies the following exact gap equation where σ is given by the diagram (b) of Fig. 2, but computed with the full propagator, namely, Here, we have used G 11 (ω) = G ϕϕ (ω) together with Eqs. (33) and (41). In the spirit of the 1/N expansion, we write M 2 = M 2 0 +O(1/N ). At LO, there are no nonlocal contributions to the self-energy and the propagator is simply given by a tree-level-like propagator G To compute the NLO propagator, we first compute the nonlocal contributions to the self-energy Σ at NLO in terms of the LO propagator G the topology depicted in Fig. 4(g). Each one-loop bubble, corresponding to the diagram (h), gives a contribution and summing the infinite sum of bubbles is achieved by solving the integral equation The function I resums the infinite chain of bubble diagrams, as is depicted in Fig. 5, where it is represented as a wiggly line. In terms of the latter the nonlocal contribution to the NLO self-energy is obtained as the diagram (i) of Fig. 5, which gives the one-loop expression Again, we skip the details of the calculations and refer the reader to the Appendix D for details. The calculation of the one-loop bubble follows the same lines as that of diagram (a) above. It can be written as and we get, for the infinite series of bubbles, where we defined the effective dimensionless coupling which is the large-N analog ofλ defined in Eq. (66). The nonlocal self-energy at NLO reads and has a similar structure as the two-loop nonlocal selfenergy in the previous perturbative calculation, Eq. (63). We finally get, for the functions γ(ω) and η(ω), As in the previous case, the response function and the field correlator can be decomposed as a sum of two poles, see Eq. (69). At the present order of approximation, we get and Similarly to the previous perturbative calculation, the coefficient c ϕ 3 being of order 1/N , we only obtain the LO expression for Λ 3,1 .
Let us now consider the χχ correlator which, at LO, is simply given by the infinite chain of bubbles. Indeed, one easily shows (see Appendix D) that From this, we get the (connected) correlator of the composite field χ = ϕ 2 /(2N ) and we deduce the LO expressions We finally need to solve Eq. (84) for the local contribution M 2 at NLO. To this aim, we use from which we obtain Collecting the previous results, we have, for the dynamical mass, and for the lowest eigenvalues As for the previous perturbative expressions, the above results exactly agree with those of the direct QFT calculations in Ref. [27]. In fact the agreement concerns all the intermediate quantities Π, I and Σ, using the rescalings (5) of the parameters and 2d Σ (110) for the different two-point functions. The very same results have also been recently obtained from a QFT calculation in Euclidean de Sitter in Ref. [28]. That such very different calculations agree is a nontrivial result. Such an agreement between the stochastic approach and direct QFT calculations on either Lorentzian or Euclidean de Sitter was already well-known for equal-time correlators, e.g., ϕ n , which measure the local field fluctuations [4,9,13]. Although expected on the basis of general arguments [4,30,31], the agreement mentioned here for unequal-time (nonlocal) correlators is far less trivial, in particular, for nonperturbative approximation schemes, and the present results, together with those of Refs. [27] and [28] provide an explicit nontrivial check.

C. Discussion
We now discuss the results we obtained for the eigenvalues and associated correlators in several regimes. First, we check that the expressions we have for the Λ n, and c ϕ,χ n coincide in the limit where we take both N large and a the couplingλ small. In this regime, introducingλ ∞ = lim N →∞λ = λ/(12m 4 ), we have ∞ , thus the two effective coupling coincide. For example, it is easy to check that Eq. (C14) coincides with the first member of Eq. (107) to give The 1/N expansion allows the study of the nonperturbative regime inλ, which correspond to either small or negative m 2 [27]. This is illustrated in Fig. 6, where we show the effective couplingλ as a function of m 2 for fixed coupling λ. The latter is of order one in the small mass regime m 2 = 0 and becomes large for m 2 < 0. Let us discuss the leading-order eigenvalues (107)-(109) in these regimes. The latter rewrite, in terms of the parameters m 2 and λ, and are plotted as functions of m 2 for λ = 1 in Fig. 7.
In the small mass regime, we have where we see that all eigenvalues are of the same order, so that all the correlators computed here have relatively large autocorrelation times, in particular, in the case of small coupling λ 1. This reflects the fact that the potential is very flat in that case.
Instead, in the case of a steep symmetry breaking treelevel potential, with m 2 < 0 and λ/m 4 1, we have The presence of a small (Λ 1,1 ) and a large (Λ 3,1 ) eigenvalue in the correlator of the vector field ϕ reflects the existence of a flat (Goldstone mode) and a steep (Higgs mode) direction in the tree-level potential. 12 The eigenvalue Λ 2,0 is, again the longitudinal Higgs mode, the only one which contributes to the correlator of the scalar field χ.
Interestingly, the present large-N results share similarities with similar analytical results in the case N = 1 in the limit of a steep double-well potential [6]. When the two minima are far apart, the situation can be described as a superposition of two single-well spectra with tunnel effect yielding infinitesimally split energy levels. Because the ground (equilibrium) state has Λ 0 = 0, this results in an exponentially suppressed, instanton-like value of Λ 1 ∝ exp(−a/λ), with a a constant. Higher eigenvalues are essentially those of the unperturbed separate Gaussian wells of square mass 2|m 2 |, that is, Λ 2n ≈ Λ 2n+1 ≈ 2n|m 2 |. Not surprisingly, our result (117) for Λ 2,0 and Λ 3,1 agree with the case N = 1 in the limit of a steep symmetry breaking potential since these correspond to the heavy longitudinal mode with tree-level square mass 2|m 2 |. The main difference in the case of a continuous symmetry N > 1 is the presence of flat directions (Golsdtone modes) in the potential, which result in a milder, power law suppression for Λ 1,1 . In fact, we observe that the 1/N expansion becomes singular for the latter when N → 1. In the limit of a steep symmetry breaking potentialλ 1, we have and similarly for the coefficient c ϕ 1 = 1−1/N +O 1/N 2 . In terms of correlation functions, this implies that in the case m 2 < 0 the correlation time in the vector channel (ϕ) is considerably larger than in the scalar channel (χ), which does not see the flat transverse direction. It is to be expected that such large correlation time also occur for composite fields in higher representations (tensor channels). These correlation times are related to other quantities of physical interest, such as the relaxation (or equilibration) times from an excited state to the BD vacuum, decoherence time scales [39,40], or, closest to standard phenomenological interest, to the spectral index of various observables [6,34]. Exploiting de Sitter invariance, the spectral index of a given field A can be read off the decomposition (44) as n A − 1 = 2Λ A , with Λ A the lowest eigenvalue contributing to the sum (44). For instance, the spectral index of the field ϕ is given by n ϕ − 1 = 2Λ 1,1 . Similarly, Λ 2,0 is related to the spectral index of the field χ or of other typical O(N ) scalar quantity. For instance, as discussed in Ref. [34], the density contrast δ = (V − V )/ V has a spectral index n δ − 1 = 2Λ 2,0 .
Finally, we mention an interesting role played by the eigenvalue Λ 3,1 based on the discussion in Sec. III E. According to Eqs. (65) and (95), we have both in the perturbative and in the 1/N expansions, with In real time, this gives  With Eq. (53), we see that Λ −1 3,1 is the autocorrelation time of the effective colored noise correlation in the vector channel due to the infrared modes while cΛ 1,1 controls the amplitude of the colored contribution. In the perturbative regime m 2 > 0, the autocorrelation time is small ∼ 1/m 2 with small amplitude ∼ λ 2 /m 2 . However the autocorrelation time can be either parametrically large ∼ 1/ √ λ with a small amplitude ∼ √ λ for m 2 = 0, or small ∼ 1/|m 2 | with "large" amplitude ∼ 1 for m 2 < 0.
We close this Section by comparing the expressions of Λ 1,1 at leading and next-to-leading orders in the 1/N expansion as a function of the parameters of the theory for the extreme case N = 2 in Figure 8.

V. CONCLUSIONS
To conclude, we used the JdD path integral formulation of the stochastic equation that describes the infrared dynamics of an O(N ) theory of test scalar fields to study the two-point unequal-time correlators of various operators. The resulting field theory is a one-dimensional supersymmetric theory with N scalar superfields. This supersymmetry is a mere consequence of the symmetries of the original system in stationary state, and was used to show that the correlators of the various fields are not independent. One of the obtained relation can be interpreted as a fluctuation-dissipation relation, showing the analogy of the system with a Brownian motion with a thermal noise at de Sitter temperature [41].
Having in mind the result of the Fokker-Planck formulation, which is usually solved as an eigenvalue problem [6,34], we then discussed the general structure of the unequal-time two-point correlator of composite operators of the scalar field. It can be expressed as a sum of free propagator with a hierarchy of mass scales, which corresponds to a subset of the tower of eigenvalues previously mentioned. We have computed explicitly the ϕϕ and χχ correlators in two specific limits, the perturbative case up to three loop and the 1/N expansion at NLO, and we have obtained the values of the first three eigenvalues in both cases. We have checked that our results coincide with other computations, done either in the Lorentzian [18,27] or Euclidean [28] field theory.
The result from the 1/N expansion is particularly interesting as it allows us to probe nonperturbative regimes. Such regimes corresponds to the massless and symmetry breaking case, the latter being difficult to probe numerically. In the limit of a deeply broken initial potential, we find that the lowest nonzero eigenvalue is strongly suppressed with respect to higher order ones. This has direct physical consequence, e.g., in terms of equilibration times or power spectra of fields in the different representations of the O(N ) group. For instance, the vector channel = 1 has typically long range spacetime correlations whereas the scalar channel = 0 is typically (sometimes significantly) of shorter range.
There are several directions to extend the present analysis. First, we used the computation of correlators to access the mass scale hierarchy. When combined with our expansion schemes, this only gives the first eigenvalues due to the coefficients appearing in the sums of free propagators of the correlators. Alternative formulations, directly at the level of the Fokker-Planck equation, may be able to grasp the full hierarchy directly.
On a more speculative level, this work is limited to test scalar fields, and an important question would be to extend the present considerations to a more realistic inflationary setup.
Here, we allow for an arbitrary normalization of the smeared field. The retarded Green function is defined asĜ and solves the following equation In frequency space, we havê and, for the spectral function, In real time, this gives the standard expression for the (over)damped harmonic oscillator It is useful to rewrite these expressions in terms of the roots of the inverse retarded function G −1 R (−iω ± ) = 0, that is, ω ± = d/2 ± ν. We haveĜ

Statistical correlator
The statistical correlator of the quantum field is defined asF The effective stochastic theory relies on the assumption that the smeared infrared field is essentially classical. At the level of the field correlators, this means that the spectral function (which encodes the noncommuting aspects of the quantum fields) is small compared to the statistical correlator (which encodes the typical occupation number) [42], that is,ρ F . In that case, the statistical correlator (A13) reduces to the stochastic correlator In the present simple model, the latter is obtained as follows. Equation (A1) is solved aŝ where the solution of the homogeneous equation readŝ with A ± arbitrary constants of integration. The latter describes the transient regime from generic initial conditions to the late-time equilibrium. The equilibrium stochastic correlator is obtained as, in frequency space, whereN (ω) = 2/(dΩ D+1 ) is the Fourier transform of the noise correlator (A2). We get or, equivalently, in real time,

Slow-roll limit
We are now in a position to clearly identify the necessary requirement for the slow-roll limit, which, we recall, amounts to neglecting the termφ in Eq. (A1). We see that, in that case, the only root of the homogeneous equation is ω = −im 2 /d. This corresponds to, first, taking the small mass limitm 2 1, so that ω − ≈m 2 /d ω + ≈ d and, second, keeping only the contribution from the lowest pole in Eqs. (A10), (A11), and (A18). In real time, this amounts to neglecting e −ω+|t| e −ω−|t| , which amounts to a coarse graining on time scales ∼ ω −1 − . We define the stochastic response and spectral functions 14 G ϕφ and ρ aŝ (A23) 14 From their definitions, these are insensitive to a change of normalization so that, e.g., G ϕφ = Gφφ. In other words, the response fieldφ rescales as the inverse of ϕ.
In the slow-roll limit we thus have G ϕφ (t) = θ(t)ρ(t), with It is important to notice that in this coarse graining process, we lose the property (A4), for which the presence of both roots ω ± is essential. In other words, the coarse grained, stochastic spectral function does not resolve the short-time structure close to t = 0. However, the constraint (A4) from quantum mechanics is replaced by another one in the coarse-grained theory. Indeed, one easily checks in the present simple example that form 2 → 0, the maximal value of the QFT spectral function isρ max → Z/d, that is, ρ max → 1, with occurs at t max → d −1 ln d 2 /m 2 . In units of the stochastic time scale ω −1 − = d/m 2 , we have ω − t max = (m 2 /d 2 ) ln d 2 /m 2 → 0 + . We thus have where t → 0 + is to be understood in units of the relevant time scale ω −1 − . The identity (A25) also follows directly from the expression (A24).
Finally, the stochastic correlator (A19) becomes, in the slow-roll limit, which corresponds to Eq. (45) after the proper rescaling of the field Gφφ(t) = 2 dΩ D+1 G ϕϕ (t). Here, we see that in the slow-roll limit, we automatically satisfy the classicality condition Gφφ ρ mentioned above. We illustrate the coarse-graining process in Fig. 9 by showing the stochastic spectral function and correlator for various values ofm 2 in units of the stochastic time scale ω − . We clearly see how the ultraviolet, short-time structure near t = 0 is washed out, leading to a singular behavior. We also see the condition (A25) emerging.
We end by remarking that the fluctuation-dissipation relation discussed above remains valid in the slow-roll limit and reads which is nothing but the relation (35) for the particular model considered here, again with the appropriate rescaling of the fields. Using the expression (56) of the tree-level propagator and the relations δ 2 12 = 0, (B2) δ 12 (K ω δ 12 ) = δ 12 , (B3) (K ω δ 12 )(K ω δ 12 ) = K ω+ω δ 12 , we easily perform the Grassmann algebra to get where we defined This generalizes to all the components of the superpropagator G 12 (t). Using this relation twice on the expression (B12) gives Eq. (63).
This is the only combination compatible with the decomposition (48) in simple fractions. Computing the poles and residue yields Λ 3,1 = 3M 2 (1 + 2βλ), (C7) We now need to compute the effective square mass M 2 at three-loop order. This can be done either by a direct calculation of the relevant local contributions to the self-energy, or by following the strategy adopted in Sec. IV B, that is, by solving the gap equation (84) at the appropriate order of approximation. We implement the latter here since this does not involve computing any new diagram. The exact gap equation for M 2 reads where, at the present order of approximation, This is readily solved as we finally get