# Creation of chaos for interacting Brownian particles

Armand Bernou, Mitia Duerinckx, Matthieu Ménard

► To cite this version:

Armand Bernou, Mitia Duerinckx, Matthieu Ménard. Creation of chaos for interacting Brownian particles. Stochastic Processes and their Applications, 2025, 193, pp.104849. 10.1016/j.spa.2025.104849. hal-05417123

HAL Id: hal-05417123

https://hal.science/hal-05417123v1

Submitted on 15 Dec 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# CREATION OF CHAOS FOR INTERACTING BROWNIAN PARTICLES

ARMAND BERNOU, MITIA DUERINCKX, AND MATTHIEU MÉNARD

**ABSTRACT.** We consider a system of $N$ Brownian particles, with or without inertia, interacting in the mean-field regime via a weak, smooth, long-range potential, and starting initially from an arbitrary exchangeable $N$-particle distribution. In this model framework, we establish a fine version of the so-called creation-of-chaos phenomenon: in weak norms, the mean-field approximation for a typical particle is shown to hold with an accuracy $O(N^{-1})$ up to an error due solely to initial pair correlations, which is damped exponentially over time. Corresponding higher-order results are also derived in the form of higher-order correlation estimates. The approach is new and easily adaptable: we start from suboptimal correlation estimates obtained from an elementary use of Itô's calculus on moments of the empirical measure, together with ergodic properties of the mean-field dynamics, and these bounds are then made optimal after combination with PDE estimates on the BBGKY hierarchy.

MSC-class: Primary 60K35; Secondary 35Q70, 82C22, 82C31, 35Q84.

## CONTENTS

<table><tr><td>1. Introduction</td><td>1</td></tr><tr><td>2. Preliminary results</td><td>9</td></tr><tr><td>3. Moment estimates on empirical fluctuations</td><td>12</td></tr><tr><td>4. Creation of chaos</td><td>15</td></tr><tr><td>5. Correlation estimates</td><td>20</td></tr><tr><td>Acknowledgements</td><td>27</td></tr><tr><td>References</td><td>27</td></tr></table>

## 1. INTRODUCTION

1.1. **General overview.** We consider the underdamped Langevin dynamics for a system of $N$ Brownian particles with mean-field interactions, moving in a confining potential in $\mathbb{R}^d$, $d \ge 1$, as described by the following system of coupled SDEs: for $1 \le i \le N$,

$$ egin{cases} 
dX_t^{i,N} = V_t^{i,N} dt, \\ 
dV_t^{i,N} = -\frac{\kappa}{N} \sum_{j=1}^N \nabla W(X_t^{i,N} - X_t^{j,N}) dt - \frac{\beta}{2} V_t^{i,N} dt - \nabla A(X_t^{i,N}) dt + dB_t^i, & t \ge 0, 
\end{cases} 
\quad (1.1) $$

where $\{Z^{i,N} := (X^{i,N}, V^{i,N})\}_{1\le i\le N}$ is the set of positions and velocities of the particles in the phase space $\mathbf{X} := \mathbb{R}^d \times \mathbb{R}^d$, $W: \mathbb{R}^d \to \mathbb{R}$ is a long-range interaction potential with the action-reaction condition $W(x) = W(-x)$, $A$ is a uniformly convex confining potential, where $\{B^i\}_i$ are i.i.d. $d$-dimensional Brownian motions, and $\kappa, \beta > 0$ are given constants. For simplicity, we choose the confining potential $A$ to be quadratic,

$$ A(x) := \frac{1}{2} a |x|^2, \quad \text{for some } a > 0, \quad (1.2) $$

although perturbations can be considered as well (see Section 1.3). Turning to a statistical description of the system, in terms of a probability density $F^N$ on the $N$-particle phase space $\mathbf{X}^N$, the Langevin dynamics (1.1) is equivalent to the corresponding Liouville equation
$$ \partial_t F^N + \sum_{i=1}^{N} v_i \cdot \nabla_{x_i} F^N = \frac{1}{2} \sum_{i=1}^{N} \mathrm{div}_{v_i}((\nabla_{v_i} + \beta v_i)F^N) \\ + \frac{\kappa}{N} \sum_{i,j=1}^{N} \nabla W(x_i - x_j) \cdot \nabla_{v_i} F^N + \sum_{i=1}^{N} \nabla A(x_i) \cdot \nabla_{v_i} F^N, \quad t \ge 0. \quad (1.3) $$

In this contribution, we investigate the behavior of the system in the regime of a large number $N \gg 1$ of particles in the general setting of initial data $\{Z^{i,N}|_{t=0}\}_{1\le i\le N}$ that are only assumed to be exchangeable (that is, their law $F^N|_{t=0}$ is invariant under permutation of the particles), but that may be strongly correlated. Let us recall two standard types of results for the system:

### — Propagation of chaos:

If particles have weak correlations initially, then propagation of chaos is expected to hold, which means that only little correlations can build over time. As a consequence, the distribution of a typical particle is approximated by the solution $\mu$ of the mean-field Vlasov-Fokker-Planck equation

$$ \left\{ \begin{array}{l} \partial_t \mu + v \cdot \nabla_x \mu = \frac{1}{2} \mathrm{div}_v((\nabla_v + \beta v)\mu) + (\nabla A + \kappa \nabla W * \mu) \cdot \nabla_v \mu, \quad t \ge 0, \\ \mu|_{t=0} = \mu_0, \end{array} \right. \quad (1.4) $$

where $\mu_0$ stands for the initial single-particle distribution. Such a result has been established in various settings and we refer e.g. to the recent work [6, 5] in case of singular interactions. More precisely, for all $1 \le m \le N$, denoting by $F^{N,m}$ the distribution of $m$ typical particles in the system, as given by the $m$-th marginal of $F^N$,

$$ F^{N,m}(z_1, \dots, z_m) = \int_{\mathbb{X}^{N-m}} F^N(z_1, \dots, z_N) dz_{m+1} \dots dz_N, \quad 1 \le m \le N, $$

we find that $F^{N,m}$ converges to $\mu^{\otimes m}$ as $N \uparrow \infty$ for any fixed $m$. Beyond this qualitative convergence, an optimal $O(N^{-1})$ estimate for the convergence of marginals was first obtained in [11] for systems with smooth interactions (and without Brownian noise, although the argument easily extends to include diffusion). This first result was later complemented in [17], which established an analogous estimate for bounded, possibly non-smooth interactions, and in a better metric: if correlations vanish initially, then

$$ \|F_t^{N,m} - \mu^{\otimes m}\|_{TV} \le C e^{Ct} m N^{-1}, \quad \text{for all } 1 \le m \le N \text{ and } t \ge 0. \quad (1.5) $$

Furthermore, if the interaction intensity $\kappa$ is sufficiently small, ensuring uniqueness of the mean-field equilibrium, then the estimate (1.5) holds uniformly in time, without the exponential factor [18, 1]. (See also previous work relying on coupling methods, which yield a suboptimal convergence rate [4, 15].) For a detailed introduction to propagation of chaos, we refer to [7].

### — Relaxation to equilibrium:

Even if the particles are strongly correlated initially, their distribution is known to relax to Gibbs equilibrium on long times. More precisely, if the interaction intensity $\kappa$ is small enough and if the interaction kernel $K$ is Lipschitz, the following uniform-in-$N$ exponential estimate was obtained in [4],

$$ W_2(F_t^{N,m}, M^{N,m}) \le C_m e^{-\frac{1}{C}t}, \quad \text{for all } 1 \le m \le N \text{ and } t \ge 0, $$

where $M^{N,m}$ stands for the $m$-th marginal of the $N$-particle Gibbs equilibrium $M^N$. As the latter exhibits weak correlations, this entails in particular that initial correlations between the particles must be exponentially damped over time.

Both uniform-in-time propagation of chaos and convergence to Gibbs equilibrium rely essentially on the same phenomenon: due to the diffusion, correlations cannot accumulate too much over time and must somehow be damped. This leads to the natural related question of creation of chaos as first raised by Del Moral and Tugaut [9] (see also “generation of chaos” in [19, 21]): for general exchangeable initial data, is it possible to show that the mean-field approximation holds to order $O(N^{-1})$ up to an error coming from initial correlations, which is damped exponentially over time?
Such a result can in fact be deduced immediately from the coupling argument of [4], but only with a suboptimal convergence rate $O(N^{-1/2})$: more precisely, for weak Lipschitz forces, it follows from [4, Theorem 2 and Remark 12] that

$$W_2(F_t^{N,m}, \mu_t^{\otimes m}) \le C_m N^{-\frac{1}{2}} + C_m e^{-\frac{1}{C}t} W_2(F_0^{N,m}, \mu_0^{\otimes m}), \quad \text{for all } 1 \le m \le N \text{ and } t \ge 0. \quad (1.6)$$

Our first objective is to achieve the optimal rate $O(N^{-1})$ in this creation-of-chaos result, which to our knowledge is not possible via the same coupling methods. More precisely, we look for bounds of the following form, in some suitable metric,

$$\|F_t^{N,m} - \mu_t^{\otimes m}\| \le C_m N^{-1} + C_m e^{-\frac{1}{C}t} \|F_0^{N,m} - \mu_0^{\otimes m}\|, \quad \text{for all } 1 \le m \le N \text{ and } t \ge 0. \quad (1.7)$$

In the sequel, we establish estimates of this kind for weak smooth interactions, in negative Sobolev norms instead of Wasserstein distance; see Theorem 1.1 below.

While virtually any proof of uniform-in-time propagation of chaos comes with some associated creation-of-chaos result showing the damping of initial correlations, we emphasize that clean estimates of the form (1.7) are more difficult to establish. In particular, there has been a lot of recent progress based on hierarchical approaches to uniform-in-time propagation of chaos [17, 6, 18, 1, 5], starting from the BBGKY hierarchy of equations for marginals, but such methods require to deal with information on the whole hierarchy $\{F^{N,m} - \mu^{\otimes m}\}_{1 \le m \le N}$ at once and can thus a priori not be used to obtain estimates of the form (1.7) for any fixed $m$. For overdamped Langevin dynamics (see Section 1.3.1 below), modulated energy methods can also be used to establish creation of chaos, see [21], but they require much stronger initial information in the form of the relative entropy of the $N$-particle distribution with respect to a tensorized state, which again cannot lead to clean estimates of the form (1.7) for marginals.

Beyond this question on the creation and propagation of chaos, we refine our analysis by deriving corresponding estimates on correlation functions. These functions capture finer statistical information on the propagation of chaos in the system, and we shall see that the creation of chaos extends at this level as well. The two-particle correlation function is naturally defined as

$$G^{N,2} := F^{N,2} - (F^{N,1})^{\otimes 2},$$

which captures the defect to propagation of chaos (1.5) at the level of two-particle statistics. This quantity alone does not allow to reconstruct the full particle density $F^N$: we define higher-order correlation functions $\{G^{N,m}\}_{2 \le m \le N}$ as suitable polynomial combinations of marginals of $F^N$ in such a way that the full particle distribution $F^N$ is recovered in form of a cluster expansion,

$$F^N(z_1, \dots, z_N) = \sum_{\pi \vdash [N]} \prod_{B \in \pi} G^{N, \#B}(z_B), \quad (1.8)$$

where $G^{N,1} := F^{N,1}$, where $\pi$ runs through the list of all partitions of the index set $[N] := \{1, \dots, N\}$, $B$ runs through the list of blocks of the partition $\pi$, $\#B$ is the cardinality of $B$, and where for $B = \{i_1, \dots, i_l\} \subset [N]$ we write $z_B = (z_{i_1}, \dots, z_{i_l})$. As is easily checked, correlation functions are fully determined by prescribing (1.8) together with the requirement $\int_X G^{N,m}(z_1, \dots, z_m) dz_l = 0$ for all $2 \le m \le N$ and $1 \le l \le m$. If initial correlations vanish, standard propagation of chaos leads to $G_t^{N,2} = O(N^{-1})$, and a formal analysis of the BBGKY hierarchy further leads to expect

$$G_t^{N,m} = O(N^{1-m}), \quad 2 \le m \le N. \quad (1.9)$$

This is known as higher-order propagation of chaos and serves as a key tool to describe deviations from mean-field theory; see e.g. [13, 11, 16, 12], and in particular [3, Corollary 1.3] for the derivation of the Bogolyubov correction to mean field. Such estimates have been established in several settings, typically with exponential growth in time: particle systems with smooth interactions were first treated in [11], and the result was later extended to merely bounded, possibly non-smooth interactions in [16]; see also [12] for weaker bounds in case of singular forces. When the interaction intensity $\kappa$ is sufficiently small, our previous work [3] proved that (1.9) actually holds uniformly in time, a result recently
extended to bounded interactions in [23, 2]. All these results rely on the assumption that initial correlations vanish. Yet, uniform-in-time versions of (1.9) hint toward a corresponding ‘higher-order’ creation-of-chaos phenomenon at the level of correlation functions: for general exchangeable initial data, is it possible to show that the m-particle correlation $G^{N,m}$ is of order $O(N^{1-m})$ up to an error coming from initial m-particle correlations (or at least, from a finite number of correlations functions), which is damped exponentially over time? For two-particle correlations, it already follows from (1.7) that

$$ \|G_t^{N,2}\| \le CN^{-1} + Ce^{-\frac{1}{C}t}\|G_o^{N,2}\|. \quad (1.10) $$

In the sequel, we obtain similar estimates for all higher-order correlation functions; see Theorem 1.2 below. The only other result of this kind that we are aware of is the recent work of Lukkarinen and Vuoksenmaa [19] regarding energy cumulants for the stochastic Kac model. In that work, however, the authors appeal to hierarchical techniques, which require to deal with the whole hierarchy {$G^{N,m}$}$_{1\le m\le N}$ at once and do not yield simpler estimates of the form (1.10) for correlations of a fixed order.

## Notation.

- We denote by $C \ge 1$ any constant that only depends on the space dimension $d$ and on the parameters $\beta, a$, as well as possibly on other controlled quantities when specified. We use the notation $\lesssim$ for $\le C \times$ up to such a multiplicative constant $C$. We write $x \simeq y$ when both $x \lesssim y$ and $y \lesssim x$ hold. We also use the notation $\ll$ (resp. $\gg$) for $\le \frac{1}{C} \times$ (resp. $\ge C \times$) for some large enough constant $C$. We occasionally add subscripts to $C, \lesssim, \simeq, \ll, \gg$ to indicate dependence on other parameters.

- For any two integers $m \ge n \ge 0$, we use the short-hand notation $[[m, n]] := \{m, m+1, \dots, n\}$, and we set $[[m]] := [[1, m]]$.

- For $1 \le q \le \infty$, we denote by $q' := \frac{q}{q-1}$ the conjugate exponent of $q$.

- For all $z \in \mathbb{X}$, we use the notation $\langle z \rangle := (1 + |z|^2)^{\frac{1}{2}}$.

- For any measurable space $(E, \mathcal{E})$, we write $\mathcal{P}(E)$ for the space of probability measures on $E$.

- For a measure $\mu \in \mathcal{P}(\mathbb{X})$, we use the following short-hand notation for its moments,

$$ Q_r(\mu) := \int_{\mathbb{X}} \langle z \rangle^r \mu(dz), \quad r \ge 0. \quad (1.11) $$

## 1.2. Main results.
We start by defining the appropriate functional framework for our results: for $r \ge 0$, $1 < q \le 2$, and $0 < p \le 1$, we consider the following weighted negative Sobolev norms on functions $h: \mathbb{X} \to \mathbb{R}$,

$$ \|h\|_{W_p^{-r,q}(\mathbb{X})} := \sup \left\{ \left| \int_{\mathbb{X}} \varphi h \right| : \varphi \in C_c^{\infty}(\mathbb{X}), \quad \|\langle \cdot \rangle^{-p} \varphi\|_{W^{r,q'}(\mathbb{X})} \le 1 \right\} \quad (1.12) $$

and for $m \ge 1$ and functions $h: \mathbb{X}^m \to \mathbb{R}$,

$$ \|h\|_{W_p^{-r,q}(\mathbb{X})^{\otimes m}} := \sup \left\{ \left| \int_{\mathbb{X}^m} \varphi^{\otimes m} h \right| : \varphi \in C_c^{\infty}(\mathbb{X}), \quad \|\langle \cdot \rangle^{-p} \varphi\|_{W^{r,q'}(\mathbb{X})} \le 1 \right\}. \quad (1.13) $$

The choice of these spaces is dictated by the available ergodic estimates on the linearized mean-field evolution, which are not known to hold on the simpler unweighted spaces $W^{-k,1}(\mathbb{X})$; see e.g. [3] and Lemma 2.1 below. As explained in Section 1.3.1, this is considerably simplified in the overdamped setting.

Our first main result justifies the desired estimates (1.7) on the creation and propagation of chaos in the above norms.

These estimates appear to be optimal in $N$, as they match the scaling (1.5) predicted by the BBGKY hierarchy in the limit $t \uparrow \infty$, and their dependence on deviations from initial chaos seems minimal.

<footnote n="1">¹This is the norm of the injective tensor product $W_p^{-r,q}(\mathbb{X})^{\otimes m}$ on the subspace of symmetric functions, which of course differs from the norm of $W_p^{-r,q}(\mathbb{X}^m)$.</footnote>
Possible improvements may concern their dependence on $m$, the regularity assumptions, or the loss of regularity in the norms, but such refinements are currently beyond the reach of our approach.

**Theorem 1.1** (Creation of chaos). Let the interaction potential $W$ be even, smooth and decaying in the sense of $W \in W^{d+3,\infty} \cap H^s(\mathbb{R}^d)$ for some $s > \frac{d}{2} + 5$, and let the confining potential $A$ be quadratic (1.2). Given an exchangeable $N$-particle distribution $F_\circ^N \in \mathcal{P}(\mathbb{X}^N)$ with all bounded moments, consider the weak solution $F^N \in C(\mathbb{R}^+; \mathcal{P}(\mathbb{X}^N))$ of the Liouville equation (1.3) with initial data $F_\circ^N$. Given $\mu_\circ \in \mathcal{P}(\mathbb{X})$ with bounded second moments, also consider the weak solution $\mu \in C(\mathbb{R}^+; \mathcal{P}(\mathbb{X}))$ of the mean-field equation (1.4) with initial data $\mu_\circ$. Then there is $\lambda_0 > 0$ such that the following holds: given $1 < q \le 2$ and $0 < p \le \frac{1}{6}$ with $pq' \gg 1$ large enough, and provided that $0 < \kappa \ll 1$ is small enough, we have for all $1 \le m \le N$ and $t \ge 0$,

$$ \|F_t^{N,m} - \mu_t^{\otimes m}\|_{W_p^{-3,q}(\mathbb{X})^{\otimes m}} \lesssim N^{-1} + \sum_{j=1}^{\tilde{m}} e^{-j\lambda_0 pt} \|F_\circ^{N,j} - \mu_\circ^{\otimes j}\|_{W_{3p}^{-2,q}(\mathbb{X})^{\otimes j}}, \quad (1.14) $$

where we have set for abbreviation $\tilde{m} = m$ for $m$ even and $\tilde{m} = m + 1$ for $m$ odd. The exponent $\lambda_0$ only depends on $d, \beta, a$, $\|W\|_{W^{1,\infty}(\mathbb{R}^d)}$; the condition $pq' \gg 1$ only depends on $d, \beta, a$; the smallness condition $\kappa \ll 1$ further depends on $p, q, W, Q_2(\mu_\circ), Q_2(F_\circ^{N,1})$; and the multiplicative constant in (1.14) further depends on $m$ and $Q_{mp}(F_\circ^{N,1})$.

We now turn to the creation of chaos at the level of correlation estimates. The precise combinations of initial correlations that must be tracked become increasingly intricate for higher-order correlations. Although their explicit form can be readily inferred from the proof, we present here, for clarity, a simplified statement under the assumption that the initial correlations satisfy a definite scaling. A key novelty of this result is that, to estimate $m$-particle correlations, it suffices to control the initial $j$-particle correlations only for $j \le 2m - 1$. This stands in contrast with the much stronger assumptions usually required for propagation or creation of chaos: classical hierarchical methods, for instance, demand knowledge of all initial correlation functions. We note, however, that the smallness assumption on the interaction intensity $\kappa$ below is particularly restrictive, as it depends on the order $m$ of the correlation function under consideration. While this condition is necessary in our approach, it can be relaxed in the simpler setting of overdamped Langevin dynamics; see Section 1.3.1.

**Theorem 1.2** (Creation of chaos in correlation estimates). Let $W, A, F^N, \lambda_0, p, q, \kappa$ be as in the statement of Theorem 1.1, and further assume $W \in C_b^\infty \cap H^\infty(\mathbb{R}^d)$. Then the two-particle correlation function satisfies for all $t \ge 0$,

$$ \|G_t^{N,2}\|_{W_p^{-3,q}(\mathbb{X})^{\otimes 2}} \lesssim N^{-1} + e^{-\lambda_0 pt} \|G_\circ^{N,2}\|_{W_{3p}^{-2,q}(\mathbb{X})^{\otimes 2}}. \quad (1.15) $$

In addition, given $2 \le m \le N$ and $\alpha \in [0, 1]$, if we have initially

$$ \|G_\circ^{N,j}\|_{W_{3p}^{-2,q}(\mathbb{X})^{\otimes j}} \le C_j^\circ N^{\alpha(1-j)}, \quad \text{for all } 2 \le j \le 2m-1, \quad (1.16) $$

and if $\kappa$ is small enough (further depending on $m$), then we have for all $t \ge 0$,

$$ \|G_t^{N,m}\|_{W_{pm}^{-r_m,q_m}(\mathbb{X})^{\otimes m}} \lesssim N^{1-m} + e^{-\lambda_0 p_m t} N^{\alpha(1-m)}, \quad (1.17) $$

for $r_m = m + 1$, some $q'_m \simeq_m q'$, and $p_m \simeq_m p$ with $p_m q'_m = pq'$.

The multiplicative constant in (1.15) only depends on $d, \beta, a, W, Q_2(F_\circ^{N,1})$; the additional smallness condition on $\kappa$ only depends on $p, q, m, Q_2(F_\circ^{N,1})$; and the multiplicative constant in (1.17) further depends on $W, Q_{mp}(F_\circ^{N,1})$ and on the constants $(C_j^\circ)_j$ in (1.16).

1.3. Extensions. We briefly discuss some straightforward extensions of the above results.
1.3.1. *Overdamped Langevin dynamics.* Our results also apply to the overdamped Langevin dynamics, as described by the following system of coupled SDEs: for $1 \le i \le N$,

$$dY_t^{i,N} = \frac{\kappa}{N} \sum_{j=1}^{N} K(Y_t^{i,N} - Y_t^{j,N}) dt - \nabla A(Y_t^{i,N}) + dB_t^i, \quad t \ge 0, \quad (1.18)$$

where $K: \mathbb{R}^d \to \mathbb{R}^d$ is a smooth force kernel with the action-reaction condition $K(x) = -K(-x)$ (we may choose $K$ to derive from a potential, $K = -\nabla W$, but this is not necessary). In terms of a probability density $\tilde{F}^N$ on the $N$-particle phase space $(\mathbb{R}^d)^N$, this particle dynamics is equivalent to the following Liouville equation,

$$\partial_t \tilde{F}^N = \frac{1}{2} \sum_{i=1}^{N} \Delta_{x_i} \tilde{F}_t^N - \frac{\kappa}{N} \sum_{i,j=1}^{N} \mathrm{div}_{x_i}(K(x_i - x_j)\tilde{F}^N) + \sum_{i=1}^{N} \mathrm{div}_{x_i}(\nabla A(x_i)\tilde{F}^N). \quad (1.19)$$

The associated mean-field equation is now the McKean-Vlasov equation

$$\partial_t \tilde{\mu} = \frac{1}{2} \Delta \tilde{\mu} - \kappa \mathrm{div}(\tilde{\mu}(K * \tilde{\mu})) + \mathrm{div}(\tilde{\mu} \nabla A). \quad (1.20)$$

In this setting, standard parabolic theory ensures that ergodic estimates for the linearized mean-field equation hold in the *unweighted* spaces $W^{-k,1}(\mathbb{R}^d)$ for $k \ge 0$; see [10, 3]. This significantly simplifies our proof: unlike in the underdamped setting, where working in *weighted* Sobolev spaces requires delicate concentration estimates for moments of the particle dynamics, such estimates are no longer needed here. Moreover, when adapting the proof to this unweighted framework, the smallness condition on the interaction intensity $\kappa$ in correlation estimates no longer needs to depend on the order of the correlation. Indeed, in the underdamped setting, each iteration step of the correlation hierarchy incurs a slight loss of weight, necessitating the use of order-dependent weights and, consequently, an order-dependent smallness condition on $\kappa$; in the overdamped case, working in unweighted spaces removes this constraint. More precisely, we obtain the following.

**Theorem 1.3 (Creation of chaos — overdamped setting).** Let the interaction force $K$ be odd, smooth and decaying in the sense of $K \in C_b^\infty \cap H^\infty(\mathbb{R}^d)$ (say), and let the confining potential $A$ be quadratic (1.2). Given an exchangeable $N$-particle distribution $\tilde{F}_t^N \in \mathcal{P}((\mathbb{R}^d)^N)$ with bounded second moments, we consider the solution $\tilde{F}^N \in C(\mathbb{R}^+; \mathcal{P}((\mathbb{R}^d)^N))$ of (1.19) with initial data $\tilde{F}_t^N$ and denote by $\{\tilde{G}^{N,m}\}_{2 \le m \le N}$ the corresponding correlation functions. Given $\tilde{\mu}_0 \in \mathcal{P}(\mathbb{R}^d)$ with bounded second moments, also consider the weak solution $\tilde{\mu} \in C(\mathbb{R}^+; \mathcal{P}(\mathbb{R}^d))$ of the McKean-Vlasov equation (1.20) with initial data $\tilde{\mu}_0$. Then there is $\lambda_0 > 0$ such that the following hold for $0 < \kappa \ll 1$ small enough.
— For all $1 \le m \le N$ and $t \ge 0$ we have

$$\| \tilde{F}_t^{N,m} - \tilde{\mu}_t^{\otimes m} \|_{W^{-3,1}(\mathbb{R}^d)^{\otimes m}} \lesssim N^{-1} + \sum_{j=1}^{\tilde{m}} e^{-j\lambda_0 t} \| \tilde{F}_t^{N,j} - \tilde{\mu}_0^{\otimes j} \|_{W^{-2,1}(\mathbb{R}^d)^{\otimes j}}, \quad (1.21)$$

where we recall the short-hand notation $\tilde{m} = m + \mathbb{1}_{m \text{ odd}}$.

— Given $2 \le m \le N$ and $\alpha \in [0, 1]$, if we have initially

$$\| \tilde{G}_t^{N,j} \|_{W^{-2,1}(\mathbb{R}^d)^{\otimes j}} \le C_j^\alpha N^{\alpha(1-j)}, \quad \text{for all } 1 \le j \le 2m-1, \quad (1.22)$$

then we have for all $t \ge 0$,

$$\| \tilde{G}_t^{N,m} \|_{W^{-m-1,1}(\mathbb{R}^d)^{\otimes m}} \lesssim N^{1-m} + e^{-\lambda_0 t} N^{\alpha(1-m)}. \quad (1.23)$$

In addition, if the initial data are chaotic, that is, $\tilde{F}_t^N = \tilde{\mu}_0^N$, then for all $2 \le m \le N$ and $t \ge 0$ we have

$$\| \tilde{G}_t^{N,m} \|_{W^{-m-1,1}(\mathbb{R}^d)^{\otimes m}} \lesssim N^{1-m}.$$

The exponent $\lambda_0$ only depends on $d, \beta, a$, $\|W\|_{W^{1,\infty}(\mathbb{R}^d)}$; the smallness condition on $\kappa$ further depends on $W$, $Q_2(\mu_0)$, $Q_2(F_0^{N,1})$; the multiplicative constant in (1.21) further depends on $m$; the multiplicative constant in (1.23) further depends on the constants $(C_j^\alpha)_j$ in (1.22).
1.3.2. Periodic setting. Our analysis is easily adapted to particle systems on the torus $T^d$, for instance in the spatially homogeneous setting $A \equiv 0$. However, in that case, for the Langevin dynamics (1.1), ergodic estimates for the linearized mean-field equation only yield a super-polynomial decay $t^{-\infty}$ instead of exponential; see indeed [22, Theorem 56] and [3, Remark 2.14]. For the overdamped dynamics (1.18), on the other hand, the exponential decay still holds.

1.3.3. Non-quadratic confinement. Both for the underdamped and overdamped Langevin dynamics (1.1) and (1.18), our results would hold in the same form if the quadratic confinement (1.2) were replaced by $A(x) = a|x|^2 + A'(x)$ for some $a > 0$ and some smooth potential $A' \in C_c^\infty(\mathbb{R}^d)$, provided that $\|\nabla^2 A'\|_{L^\infty(\mathbb{R}^d)}$ is small enough (only depending on $\beta, W, a$). In that case, we can still appeal to [3] to ensure the validity of the ergodic estimates of Lemma 2.1 below, while the rest of our approach can be adapted directly without major difficulties.

1.4. Overview of the argument. As explained above, to obtain clean estimates on the creation of chaos at the level of marginals, we must avoid hierarchical methods and therefore develop an alternative way to capture information on correlations. A first non-hierarchical method was introduced in [11], based on Glauber and Malliavin stochastic calculus along trajectories. This approach reduced correlation estimates to a sensibility analysis of particle trajectories with respect to their initial data and driving Brownian motions. However, it fails to yield uniform-in-time estimates.

A different, non-hierarchical method was later developed in [10] and completed in our previous work [3], relying on the master equation formalism and on Lions' differential calculus. This framework allows one to directly exploit the ergodic properties of the mean-field dynamics, thereby obtaining uniform-in-time correlation estimates. More precisely, define the empirical measure $\mu^N$ of the particle system by

$$ \mu_t^N := \frac{1}{N} \sum_{i=1}^{N} \delta_{Z_t^{i,N}}, \quad t \ge 0. \qquad (1.24) $$

The starting point in [10, 3] is the following dynamical decomposition of the empirical measure $\mu^N$

$$ \int_X \varphi \mu_t^N = \mathcal{U}_t(\mu_\circ^N) + \frac{1}{N} \sum_{i=1}^N \int_0^t (\partial_\mu \mathcal{U}_{t-s})(\mu_s^N; Z_s^{i,N}) \, dB_s^i \\ + \frac{1}{2N^2} \sum_{i=1}^N \int_0^t (\partial_\mu^2 \mathcal{U}_{t-s})(\mu_s^N; Z_s^{i,N}, Z_s^{i,N}) \, ds, \quad (1.25) $$

in terms of the mean-field semigroup $\mathcal{U}_t(\mu) := \int_X \varphi \mu^t[\mu]$, with $\mu^t[\mu]$ denoting the solution of the mean-field equation (1.4) starting from $\mu$, and where $\partial_\mu$ denotes the L-derivative (see [10, 3] for details). The Brownian term in (1.25) is of order $O(N^{-1/2})$, while the last term is $O(N^{-1})$. As shown in [3], using graphical representation, this decomposition can be used to expand cumulants of the empirical measure and thereby estimate correlation functions. However, to effectively control the contribution of initial data (e.g. the term $\mathcal{U}_t(\mu_\circ^N)$ in (1.25)), one must combine this with Glauber calculus. Consequently, the method relies on nonlinear concentration properties of the initial data and thus essentially requires the system to be initially chaotic.

In the present work, in order to overcome this limitation, we develop a new, non-hierarchical approach that yields uniform-in-time correlation estimates even for non-chaotic initial data. Our method can be viewed as a linearized version of the approach in [10, 3]. Instead of decomposing the dynamics of the empirical measure through the nonlinear mean-field semigroup as in (1.25), we write

$$ \int_X \varphi \mu_t^N = \int_X (U_{0,t}\varphi) \mu_\circ^N + \frac{1}{N} \sum_{i=1}^N \int_0^t \nabla_v (U_{s,t}\varphi)(Z_s^{i,N}) \cdot dB_s^i \\ - \kappa \int_0^t \int_X \mu_s^N \nabla_v (U_{s,t}\varphi) \cdot \nabla W * (\mu_s^N - \mu_s) \, ds, $$
where {$U_{s,t}$}$_{0\le s\le t}$ stands for the backward dual linearized mean-field evolution (see (2.1) below for details). Comparing with the mean-field equation, this yields for the empirical fluctuation $\mu_t^N - \mu$,

$$ \int_X \varphi(\mu_t^N - \mu_t) = \int_X (U_{0,t}\varphi)(\mu_\circ^N - \mu_\circ) + \frac{1}{N} \sum_{i=1}^N \int_0^t \nabla_v(U_{s,t}\varphi)(Z_s^{i,N}) \cdot dB_s^i \\ - \kappa \int_0^t \int_X \mu_s^N \nabla_v(U_{s,t}\varphi) \cdot \nabla W * (\mu_s^N - \mu_s) \, ds. \quad (1.26) $$

This identity, already used for instance in [8], can be exploited as follows: the Brownian term is of order $O(N^{-1/2})$, while the last term is linear in $\mu_t^N - \mu_t$ and can therefore be treated as a perturbation of the left-hand side when $\kappa \ll 1$.

More precisely, starting from (1.26), we derive moment estimates on $\mu_t^N - \mu$. Specifically, we consider the following quantities,

$$ A^{N,m}(t) := \sup_{\varphi} \mathbb{E} \left[ \left| \int_X \varphi(\mu_t^N - \mu_t) \right|^m \right], $$

where the supremum is taken over a suitable class of smooth test functions $\varphi$. Using (1.26) together with Itô calculus and ergodic estimates for the linearized mean-field evolution {$U_{s,t}$}$_{0\le s\le t}$, and assuming that $\kappa$ is small enough for the last term in (1.26) to be treated as a perturbation, we formally obtain for all $t \ge 0$,

$$ A^{N,m}(t) \le \left(\frac{Cm}{N}\right)^{\frac{m}{2}} + C^m e^{-\lambda_0 mt} A^{N,m}(0), \quad (1.27) $$

see Proposition 3.1. Note that the scaling $N^{-m/2}$ is optimal since fluctuations of the empirical measure are known to be of order $O(N^{-1/2})$. Making this estimate rigorous requires some care in the underdamped Langevin setting. In this case, ergodic estimates for the linearized mean-field evolution hold only in weighted Sobolev spaces. As a result, one must work within such weighted spaces, which inevitably entails a loss of moments of the empirical measure in the estimates. These losses are controlled by suitable concentration bounds on the particle system. In turn, this is precisely what leads to the order-dependent smallness condition on $\kappa$ in Theorem 1.2.

Estimates of the form (1.27) have the advantage of cleanly isolating the contribution of the initial fluctuations and showing that it is damped exponentially in time. However, moments of the empirical fluctuation $\mu_t^N - \mu_t$ do not directly translate into estimates on correlation functions. In particular, as $t \uparrow \infty$, the bound (1.27) only provides an $O(N^{-1/2})$ error for the propagation of chaos, instead of the optimal $O(N^{-1})$, and an $O(N^{-m/2})$ bound for the $m$-particle correlation function, rather than the optimal $O(N^{1-m})$.

To recover the optimal scalings, we return to a hierarchical approach. Starting from the coarse correlation estimates obtained from (1.27), we can now rigorously truncate the BBGKY hierarchy and perform hierarchical analysis without requiring information on all initial correlations. We then show that the optimal rates can be rigorously recovered iteratively from the nonlinear BBGKY hierarchy satisfied by correlations: inserting suboptimal estimates into the hierarchy leads to self-improving bounds. To our knowledge, this use of the BBGKY hierarchy is new. Rather than carrying out direct estimates on the hierarchy — which in this context would require overly strong assumptions on the full family of initial correlations — we first derive coarse estimates on correlations by a different, stochastic-calculus-based method, and then use the hierarchy to bootstrap these into optimal asymptotic scalings. A similar strategy has since been employed in [2, 12] for different purposes.

Finally, compared with our previous work [3], the present analysis yields significantly stronger conclusions on the creation of chaos — with a considerably simpler and more direct proof. However, the method developed here is inherently perturbative as it requires the interaction intensity $\kappa$ to be sufficiently small. This limitation stems from the need to treat the last term in the decomposition (1.26) as a perturbation. By contrast, the framework in [3] only relies on ergodic estimates for the linearized mean-field evolution, which in some settings remain valid even for large $\kappa$.
**Plan of the paper.** Section 2 collects preliminary material, including ergodic estimates for the linearized mean-field equation and moment bounds on the empirical measure of the particle system. Building on these ingredients and using Itô calculus, Section 3 establishes a rigorous version of the key estimate (1.27) (see Proposition 3.1), providing moment bounds on empirical fluctuations.

In Section 4, we complete the proof of Theorem 1.1. We first show in Section 4.1 how the moment bounds from Proposition 3.1 can be reformulated as suboptimal mean-field estimates on the first marginal $F^{N,1}$ of the particle system. Next, in Section 4.2, we insert these preliminary estimates into the BBGKY hierarchy to recover optimal mean-field bounds on $F^{N,1}$. Finally, Section 4.3 combines these optimal bounds with the moment estimates of Proposition 3.1 to establish sharp mean-field estimates on all marginals.

Section 5 is devoted to the proof of Theorem 1.2. In Section 5.1, we first derive suboptimal bounds on correlation functions from the moment estimates of Proposition 3.1. Then, using the BBGKY hierarchy, Section 5.2 develops a nonlinear hierarchy of a priori estimates on correlations. The proof of Theorem 1.2 is completed in Section 5.3 through an iterative argument combining these ingredients.

## 2. PRELIMINARY RESULTS

In this section, we recall the available ergodic estimates on the linearized mean-field equation in weighted Sobolev spaces and we prove some concentration bounds for the empirical measure of the particle system. These ingredients are both crucial in the sequel.

### 2.1. Ergodic estimates for linearized mean-field.
As correlations are driven by the linearized mean-field operator, a key ingredient to understand the damping of initial correlations is naturally given by ergodic estimates for linearized mean field. Recall that $\mu$ stands for the solution of the mean-field equation (1.4) with initial data $\mu_0$. This equation admits a unique steady state $M$ for $0 < \kappa \ll 1$ small enough, which can be characterized as the unique solution to the fixed-point Gibbs equation

$$M(x, v) = c_M \exp \left[ -\beta \left( \frac{1}{2} |v|^2 + A(x) + \kappa W * M(x) \right) \right], \quad (x, v) \in \mathbb{X},$$

where $c_M$ is the normalizing constant such that $\int_{\mathbb{X}} M = 1$. This fixed-point equation has indeed a unique solution whenever $\kappa\beta\|W\|_{L^\infty(\mathbb{R}^d)} < 1$. For $h \in C_c^\infty(\mathbb{X})$ and $t \ge 0$, we consider the (non-autonomous) backward dual linearized mean-field evolution,

$$\left\{ \begin{array}{l} \partial_s U_{s,t}h + v \cdot \nabla_x U_{s,t}h = \left( -\frac{1}{2}\Delta_v + \frac{1}{2}\beta v \cdot \nabla_v + (\nabla A + \kappa \nabla W * \mu_s) \cdot \nabla_v \right) U_{s,t}h, \quad s \le t, \\ U_{s,t}h|_{s=t} = h. \end{array} \right. \quad (2.1)$$

In [3, Proposition 7.1], we obtained the following ergodic estimates in weighted Sobolev spaces. This is a consequence of hypocoercivity techniques together with a suitable use of the enlargement theory of Mischler and Mouhot [20]. (In [3], an estimate is obtained instead for the forward primal linearized evolution, but it is equivalent to this one by duality.)

**Lemma 2.1 (Dual ergodic estimates, [3]).** There exist $\kappa_0, \lambda_0 > 0$ (depending on $d, \beta, a, \|W\|_{W^{1,\infty}(\mathbb{R}^d)}$) such that given $\kappa \in [0, \kappa_0]$ we have for all $r \ge 0$, $1 < q \le 2$, $0 < p \le 1$ with $pq' \gg 1$ large enough (only depending on $d, \beta, a$), $t \ge s \ge 0$, and $h \in C_c^\infty(\mathbb{X})$,

$$\left\| \langle z \rangle^{-p} \left( U_{s,t}h - \int_{\mathbb{X}} hM \right) \right\|_{W^{r,q'}(\mathbb{X})} \lesssim e^{-\lambda_0 p(t-s)} \left\| \langle z \rangle^{-p} h \right\|_{W^{r,q'}(\mathbb{X})},$$

where the multiplicative constant only depends on $d, \beta, p, q, r, a, \|W\|_{W^{r+d+1,\infty}(\mathbb{R}^d)}$, and $Q_2(\mu_0)$.

We emphasize that the above ergodic estimates are not known to hold in $W^{r,\infty}(\mathbb{X})$, but only in weighted spaces. As this forces us to work with polynomial weights, we also need to control moments of the particle dynamics, which is the purpose of the next section.
2.2. **Concentration bounds for the empirical measure.** Working within weighted spaces, we will inevitably lose moments of the empirical measure in the estimates (see proof of Proposition 3.1 below). To control these losses, we need the following concentration bounds. The proof is based on Itô's formula and is reminiscent of similar computations used in the coupling argument of [4].

**Lemma 2.2 (Concentration bounds).** Let $\kappa_0 > 0$ be as in Lemma 2.1 and consider the empirical measure $\mu^N$ as defined in (1.24). There exist $\lambda_0 > 0$ (only depending on $\beta, a$) and $C_0 > 0$ (further depending on $d, \|W\|_{W^{1,\infty}(\mathbb{R}^d)}$), such that for all $2 \le m \le N$, $0 < p \le 1$, $L \ge C_0 Q_{2p}(\mu_\circ)$, and $t \ge 0$ we have

$$ \mathbb{E}[(Q_{2p}(\mu_t^N) - L)_+^m] \le \left(\frac{C_0 L m}{N}\right)^{\frac{m}{2}} + C_0^m e^{-\lambda_0 pm} \mathbb{E}\left[\left|\int_{\mathbb{X}} \langle z \rangle^{2p} (\mu_0^N - \mu_\circ)\right|^m\right]. \quad (2.2) $$

In addition, for all $r \ge 0$,

$$ \mathbb{E}[Q_r(\mu_t^N)] = Q_r(F_t^{N,1}) \le C_0(C_0 r)^{\frac{r}{2}} + C_0^r e^{-\lambda_0 r t} Q_r(F_\circ^{N,1}). \quad (2.3) $$

*Proof.* We divide the proof into two steps, starting with the concentration bound (2.2), before turning to the more classical proof of the moment estimate (2.3).

**Step 1.** Proof of (2.2).

Let $0 < p \le 1$ be fixed. In the spirit of [4], we consider the random process

$$ G_t^N := \frac{1}{N} \sum_{i=1}^{N} (R_t^{i,N})^p, \quad R_t^{i,N} := 1 + a |X_t^{i,N}|^2 + |V_t^{i,N}|^2 + \eta X_t^{i,N} \cdot V_t^{i,N}, $$

for some $\eta \in (0, 2\sqrt{a})$ to be chosen later on. Note that this range of $\eta$ ensures

$$ R_t^{i,N} \approx_{a,\eta} \langle Z_t^{i,N} \rangle^2, \quad G_t^N \approx_{a,\eta} Q_{2p}(\mu_t^N). \quad (2.4) $$

The particle dynamics (1.1) and Itô's formula yield

$$ \begin{aligned} dR_t^{i,N} = & - (a\eta |X_t^{i,N}|^2 + (\beta - \eta) |V_t^{i,N}|^2 + \frac{\eta\beta}{2} X_t^{i,N} \cdot V_t^{i,N} - 1) dt \\ & - \frac{\kappa}{N} \sum_{j=1}^{N} (2V_t^{i,N} + \eta X_t^{i,N}) \cdot \nabla W (X_t^{i,N} - X_t^{j,N}) dt + (2V_t^{i,N} + \eta X_t^{i,N}) \cdot dB_t^i. \end{aligned} $$

Provided that $0 < \eta \ll 1$ is chosen small enough (only depending on $\beta, a$), we can find $\lambda_0 > 0$ only depending on $\beta, a$ such that

$$ dR_t^{i,N} \le -2\lambda_0 R_t^{i,N} dt + C_0 dt + (2V_t^{i,N} + \eta X_t^{i,N}) \cdot dB_t^i, \quad (2.5) $$

where henceforth $C_0$ stands for a positive constant only depending on $d, \beta, a, \|W\|_{W^{1,\infty}(\mathbb{R}^d)}$, the value of which may change from line to line. By Itô's formula, using $p(p-1) \le 0$ and $(R_t^{i,N})^{p-1} \le 1$, we deduce

$$ \begin{aligned} d(R_t^{i,N})^p & \le -2\lambda_0 p (R_t^{i,N})^p dt + C_0 p (R_t^{i,N})^{p-1} dt + \frac{1}{2} p (p-1) (R_t^{i,N})^{p-2} |2V_t^{i,N} + \eta X_t^{i,N}|^2 dt \\ & \quad + p (R_t^{i,N})^{p-1} (2V_t^{i,N} + \eta X_t^{i,N}) \cdot dB_t^i \\ & \le -2\lambda_0 p (R_t^{i,N})^p dt + C_0 p dt + p (R_t^{i,N})^{p-1} (2V_t^{i,N} + \eta X_t^{i,N}) \cdot dB_t^i, \end{aligned} $$

and thus, summing over $i$ and recalling the definition of $G_t^N$,

$$ dG_t^N \le -2\lambda_0 p G_t^N dt + C_0 p dt + \frac{p}{N} \sum_{i=1}^{N} (R_t^{i,N})^{p-1} (2V_t^{i,N} + \eta X_t^{i,N}) \cdot dB_t^i. $$

Further using Itô's generalized formula for convex functions, bounding classically the expectation of the local time using the quadratic variation, we get for all $L \ge 0$ and $2 \le m \le N$,

$$ \frac{d}{dt} \mathbb{E}[(G_t^N - L)_+^m] \le \mathbb{E}\left[-2\lambda_0 pm(G_t^N - L)_+^{m-1} G_t^N + C_0 pm(G_t^N - L)_+^{m-1}\right] $$
$$ + \frac{p^2m(m-1)}{2N^2}(G_t^N - L)_+^{m-2} \sum_{i=1}^{N} (R_t^{i,N})^{2p-2} |2V_t^{i,N} + \eta X_t^{i,N}|^2 \\ \leq \mathbb{E} \left[ -2\lambda_0 pm (G_t^N - L)_+^{m-1} G_t^N + C_0 pm (G_t^N - L)_+^{m-1} + \frac{C_0 pm^2}{N} (G_t^N - L)_+^{m-2} G_t^N \right], $$

hence, decomposing $G_t^N = G_t^N - L + L$, we get for $L \geq C_0/\lambda_0$,

$$ \frac{d}{dt} \mathbb{E}[(G_t^N - L)_+^m] \leq \mathbb{E} \left[ -2\lambda_0 pm (G_t^N - L)_+^m + \frac{C_0 Lpm^2}{N} (G_t^N - L)_+^{m-2} \right]. $$

By Young's inequality, this entails

$$ \frac{d}{dt} \mathbb{E}[(G_t^N - L)_+^m] \leq -\lambda_0 pm \mathbb{E}[(G_t^N - L)_+^m] + 2\lambda_0 p \left( \frac{C_0 Lm}{\lambda_0 N} \right)^{\frac{m}{2}}, $$

and thus, by Grönwall's inequality,

$$ \mathbb{E}[(G_t^N - L)_+^m] \leq \left( \frac{C_0 Lm}{\lambda_0 N} \right)^{\frac{m}{2}} + e^{-\lambda_0 pm t} \mathbb{E}[(G_0^N - L)_+^m]. $$

Noting that

$$ \mathbb{E}[(G_0^N - L)_+^m] \leq \mathbb{E} \left[ \left( C_0 \int_X \langle z \rangle^{2p} \mu_0^N - L \right)_+^m \right] = C_0^m \mathbb{E} \left[ \left( \int_X \langle z \rangle^{2p} (\mu_0^N - \mu_\circ) + Q_{2p}(\mu_\circ) - C_0^{-1}L \right)_+^m \right], $$

and further imposing $L \geq C_0 Q_{2p}(\mu_\circ)$, we conclude

$$ \mathbb{E}[(G_t^N - L)_+^m] \leq \left( \frac{C_0 Lm}{\lambda_0 N} \right)^{\frac{m}{2}} + C_0^m e^{-\lambda_0 pm t} \mathbb{E} \left[ \left| \int_X \langle z \rangle^{2p} (\mu_0^N - \mu_\circ) \right|^m \right], $$

and the conclusion (2.2) follows.

**Step 2.** Proof of (2.3).

Using the same notation as in Step 1, and starting from (2.5), we find by Itô's formula, for any $r > 0$,

$$ \frac{d}{dt} \mathbb{E}[(R_t^{1,N})^r] \leq -2\lambda_0 r \mathbb{E}[(R_t^{1,N})^r] + C_0 r \mathbb{E}[(R_t^{1,N})^{r-1}] + \frac{1}{2} r (r-1) \mathbb{E}[(R_t^{1,N})^{r-2}] |2V_t^{1,N} + \eta X_t^{1,N}|^2. \quad (2.6) $$

In case $r \geq 1$, using $|2V_t^{1,N} + \eta X_t^{1,N}|^2 \lesssim R_t^{1,N}$ and Young's inequality, we find

$$ \frac{d}{dt} \mathbb{E}[(R_t^{1,N})^r] \leq -2\lambda_0 r \mathbb{E}[(R_t^{1,N})^r] + C_0 r^2 \mathbb{E}[(R_t^{1,N})^{r-1}] \leq -\lambda_0 r \mathbb{E}[(R_t^{1,N})^r] + \lambda_0 \left(\frac{C_0 r}{\lambda_0}\right)^r. $$

Hence, by Grönwall's inequality, recalling $r \geq 1$

$$ \mathbb{E}[(R_t^{1,N})^r] \leq \left(\frac{C_0 r}{\lambda_0}\right)^r + e^{-\lambda_0 r t} \mathbb{E}[(R_0^{1,N})^r], $$

and the conclusion (2.3) then follows (up to renaming the constants), using

$$ \mathbb{E}[(R_t^{1,N})^r] \approx r \mathbb{E}[Q_{2r}(\mu_t^N)] = \mathbb{E} \left[ \int_X \langle z \rangle^{2r} \mu_t^N \right] = Q_{2r}(F_t^{N,1}). $$

In case $0 < r \leq 1$, rather using $r(r-1) \leq 0$ and $(R_t^{1,N})^{r-1} \leq 1$ in (2.6), we find

$$ \frac{d}{dt} \mathbb{E}[(R_t^{1,N})^r] \leq -2\lambda_0 r \mathbb{E}[(R_t^{1,N})^r] + C_0 r, $$

and the conclusion follows similarly. □
### 3. MOMENT ESTIMATES ON EMPIRICAL FLUCTUATIONS

This section is devoted to the following rigorous formulation of the key estimate (1.27), providing moment bounds on the empirical fluctuations. The proof is based on a duality approach, combining Itô's formula with a Grönwall argument. The overall strategy is inspired by [8], although we work with a different dual backward evolution on test functions; in particular, we rely on the same Fourier-based reconstruction of the moment structure to close the argument.

**Proposition 3.1.** Consider the following quantities measuring the mean-field approximation error, for any $1 < q \le 2$, $0 < p \le 1$, $m \ge 1$, and $t \ge 0$,

$$A_{p,q}^{N,m}(t) := \sup \left\{ \mathbb{E} \left[ \left| \int_{\mathbb{X}} \varphi(\mu_t^N - \mu_t) \right|^m \right] : \varphi \in C_c^\infty(\mathbb{X}), \|\langle z \rangle^{-p} \varphi\|_{W^{2,q'}(\mathbb{X})} \le 1 \right\}. \quad (3.1)$$

There exists $\lambda_0 > 0$ (only depending on $d, \beta, a$, $\|W\|_{W^{1,\infty}(\mathbb{R}^d)}$) such that the following holds: given $1 < q \le 2$ and $0 < p \le 1$ with $pq' \gg 1$ large enough (only depending on $d, \beta, a$), there exist $\kappa_0, C_0 > 0$ (further depending on $p, q, Q_2(\mu_\circ)$, and $\|W\|_{W^{d+3,\infty} \cap H^s(\mathbb{R}^d)}$ for $s > \frac{d}{2} + 3$) such that given $\kappa \in [0, \kappa_0]$ we have for all $2 \le m \le N$ and $t \ge 0$,

$$A_{p,q}^{N,m}(t) \le \left(\frac{C_0 m}{N}\right)^{\frac{m}{2}} + C_0^m e^{-\lambda_0 pmt} A_{3p,q}^{N,m}(0).$$

*Proof.* Let $\lambda_0$ stand for the minimum between the corresponding exponents in Lemmas 2.1 and 2.2. Let $1 < q \le 2$ and $0 < p \le 1$ be fixed with $pq' \gg \beta, a$. 1. Instead of controlling directly the evolution of $A_{p,q}^{N,m}(t)$, we appeal to a duality method similarly as e.g. in [8]: let us consider

$$D_{p,q}^{N,m}(t) := \sup \left\{ e^{\gamma m(t-s)} \mathbb{E} \left[ \left| \int_{\mathbb{X}} (U_{s,t}\varphi)(\mu_s^N - \mu_s) \right|^m \right] : 0 \le s \le t, \varphi \in C_c^\infty(\mathbb{X}), \|\langle z \rangle^{-p} \varphi\|_{W^{2,q'}(\mathbb{X})} \le 1 \right\},$$

for some exponent $\gamma > 0$ to be chosen later, where $\{U_{s,t}\}_{0 \le s \le t}$ denotes the backward dual linearized mean-field evolution defined in (2.1). The exponential factor in this definition allows for an efficient application of the ergodic estimates for the linearized mean-field evolution. Taking $s = t$ in the supremum, we immediately see that

$$A_{p,q}^{N,m}(t) \le D_{p,q}^{N,m}(t) \quad \text{and} \quad A_{p,q}^{N,m}(0) = D_{p,q}^{N,m}(0). \quad (3.2)$$

Given a test function $\varphi \in C_c^\infty(\mathbb{X})$ with $\|\langle z \rangle^{-p} \varphi\|_{W^{2,q'}(\mathbb{X})} \le 1$, using the cancellations between the evolution of empirical fluctuations $\mu^N - \mu$ and the dual linearized mean-field evolution for the test function, we find by Itô's formula, for all $m \ge 2$ and $0 \le s \le t$,

$$\begin{aligned} e^{\gamma m(t-s)} \mathbb{E} \left[ \left| \int_{\mathbb{X}} (U_{s,t}\varphi)(\mu_s^N - \mu_s) \right|^m \right] &\le e^{\gamma m t} \mathbb{E} \left[ \left| \int_{\mathbb{X}} (U_{0,t}\varphi)(\mu_0^N - \mu_\circ) \right|^m \right] \\ &\quad - \gamma m \int_0^s e^{\gamma m(t-\tau)} \mathbb{E} \left[ \left| \int_{\mathbb{X}} (U_{\tau,t}\varphi)(\mu_\tau^N - \mu_\tau) \right|^m \right] d\tau \\ &\quad + \kappa m \int_0^s e^{\gamma m(t-\tau)} \mathbb{E} \left[ \left| \int_{\mathbb{X}} (U_{\tau,t}\varphi)(\mu_\tau^N - \mu_\tau) \right|^{m-1} \left| \int_{\mathbb{X}} \mu_\tau^N \nabla_v (U_{\tau,t}\varphi) \cdot \nabla W * (\mu_\tau^N - \mu_\tau) \right| \right] d\tau \\ &\quad + \frac{m(m-1)}{N} \int_0^s e^{\gamma m(t-\tau)} \mathbb{E} \left[ \left| \int_{\mathbb{X}} (U_{\tau,t}\varphi)(\mu_\tau^N - \mu_\tau) \right|^{m-2} \left| \int_{\mathbb{X}} |\nabla_v (U_{\tau,t}\varphi)|^2 \mu_\tau^N \right| \right] d\tau. \end{aligned} \quad (3.3)$$

From here, we split the proof into four steps and start by examining some of the terms appearing on the right-hand side of this inequality.

**Step 1.** Estimate on the penultimate term in (3.3): for any $L \ge 1$,

$$\int_0^s e^{\gamma m(t-\tau)} \mathbb{E} \left[ \left| \int_{\mathbb{X}} (U_{\tau,t}\varphi)(\mu_\tau^N - \mu_\tau) \right|^{m-1} \left| \int_{\mathbb{X}} \mu_\tau^N \nabla_v (U_{\tau,t}\varphi) \cdot \nabla W * (\mu_\tau^N - \mu_\tau) \right| \right] d\tau$$
$$ \lesssim L \int_0^s e^{\gamma m(t-\tau)} \mathbb{E} \left[ \left| \int_{\mathbb{X}} (U_{\tau,t}\varphi)(\mu_{\tau}^N - \mu_{\tau}) \right|^m \right] d\tau \\ + \frac{L}{m} \int_0^s e^{(\gamma-\lambda_0 p)m(t-\tau)} D_{p,q}^{N,m}(\tau) d\tau + \frac{1}{m} \int_0^s e^{(\gamma-\lambda_0 p)m(t-\tau)} \mathbb{E} \left[ (Q_p(\mu_{\tau}^N) - L)_+^m \right] d\tau. \quad (3.4) $$

In order to write this term as a superposition of moments, we proceed as in [8, Proposition 5.9] and decompose $\nabla W * (\mu_{\tau}^N - \mu_{\tau})$ as a superposition of Fourier modes,

$$ \int_{\mathbb{X}} \mu_{\tau}^{N} \nabla_{v}(U_{\tau,t}\varphi) \cdot \nabla W * (\mu_{\tau}^{N} - \mu_{\tau}) = C \int_{\mathbb{R}^d} i\xi \widehat{W}(\xi) \cdot \left( \int_{\mathbb{X}} \nabla_{v}(U_{\tau,t}\varphi) e_{\xi} \mu_{\tau}^{N} \right) \left( \int_{\mathbb{X}} e_{-\xi}(\mu_{\tau}^{N} - \mu_{\tau}) \right) d\xi, $$

where we have set for abbreviation $e_{\xi}(x) := e^{ix \cdot \xi}$. By the ergodic estimates of Lemma 2.1, together with the Sobolev embedding with $\ell > 2d/q'$, we find

$$ \|\langle z \rangle^{-p} \nabla_v (U_{\tau,t} \varphi) e_{\xi} \|_{L^{\infty}(\mathbb{X})} \lesssim \|\langle z \rangle^{-p} \nabla_v (U_{\tau,t} \varphi) \|_{W^{\ell,q'}(\mathbb{X})} \lesssim e^{-\lambda_0 p(t-\tau)} \|\langle z \rangle^{-p} \varphi \|_{W^{\ell+1,q'}(\mathbb{X})}. \quad (3.5) $$

Choosing $q' \ge pq' \gg 1$ large enough, we can choose $\ell > 2d/q'$ with $\ell + 1 \le 2$. Hence, recalling the choice of $\varphi$ with $\|\langle z \rangle^{-p} \varphi\|_{W^{2,q'}(\mathbb{X})} \le 1$, we get

$$ \left| \int_{\mathbb{X}} \mu_{\tau}^{N} \nabla_{v}(U_{\tau,t}\varphi) \cdot \nabla W * (\mu_{\tau}^{N} - \mu_{\tau}) \right| \lesssim e^{-\lambda_0 p(t-\tau)} Q_p(\mu_{\tau}^{N}) \int_{\mathbb{R}^d} |\xi \widehat{W}(\xi)| \left| \int_{\mathbb{X}} e_{-\xi}(\mu_{\tau}^{N} - \mu_{\tau}) \right| d\xi, $$

where the multiplicative constant only depends on $d, \beta, p, q, a$, $Q_2(\mu_\circ)$, and $\|W\|_{W^{d+3,\infty}(\mathbb{R}^d)}$. For any $L \ge 1$, decomposing $Q_p(\mu_{\tau}^N) \le L + (Q_p(\mu_{\tau}^N) - L)_+$, using $|\int_{\mathbb{X}} e_{-\xi}(\mu_{\tau}^N - \mu_{\tau})| \le 2$, and noting that the assumptions on $W$ ensure $\int_{\mathbb{R}^d} |\xi \widehat{W}(\xi)| d\xi \lesssim 1$, we can bound

$$ \left| \int_{\mathbb{X}} \mu_{\tau}^{N} \nabla_{v}(U_{\tau,t}\varphi) \cdot \nabla W * (\mu_{\tau}^{N} - \mu_{\tau}) \right| \lesssim L e^{-\lambda_0 p(t-\tau)} \int_{\mathbb{R}^d} |\xi \widehat{W}(\xi)| \left| \int_{\mathbb{X}} e_{-\xi}(\mu_{\tau}^{N} - \mu_{\tau}) \right| d\xi \\ + e^{-\lambda_0 p(t-\tau)} (Q_p(\mu_{\tau}^{N}) - L)_+ $$

Using this to estimate the penultimate term in (3.3), we find

$$ \int_0^s e^{\gamma m(t-\tau)} \mathbb{E} \left[ \left| \int_{\mathbb{X}} (U_{\tau,t}\varphi)(\mu_{\tau}^N - \mu_{\tau}) \right|^{m-1} \left| \int_{\mathbb{X}} \mu_{\tau}^N \nabla_v (U_{\tau,t}\varphi) \cdot \nabla W * (\mu_{\tau}^N - \mu_{\tau}) \right| \right] d\tau \\ \lesssim L \int_{\mathbb{R}^d} |\xi \widehat{W}(\xi)| \int_0^s e^{\gamma m(t-\tau)} e^{-\lambda_0 p(t-\tau)} \mathbb{E} \left[ \left| \int_{\mathbb{X}} e_{-\xi}(\mu_{\tau}^N - \mu_{\tau}) \right| \left| \int_{\mathbb{X}} (U_{\tau,t}\varphi)(\mu_{\tau}^N - \mu_{\tau}) \right|^{m-1} \right] d\tau d\xi \\ + \int_0^s e^{\gamma m(t-\tau)} e^{-\lambda_0 p(t-\tau)} \mathbb{E} \left[ (Q_p(\mu_{\tau}^N) - L)_+ \left| \int_{\mathbb{X}} (U_{\tau,t}\varphi)(\mu_{\tau}^N - \mu_{\tau}) \right|^{m-1} \right] d\tau. \quad (3.6) $$

Let us start by examining the first right-hand side term in this estimate. Noting that the definition of $D_{p,q}^{N,m}(t)$ entails in particular

$$ \mathbb{E} \left[ \left| \int_{\mathbb{X}} e_{-\xi}(\mu_{\tau}^N - \mu_{\tau}) \right|^{\frac{1}{m}} \right] \le \|\langle z \rangle^{-p} e_{-\xi} \|_{W^{2,q'}(\mathbb{X})} D_{p,q}^{N,m}(\tau)^{\frac{1}{m}} \lesssim \langle \xi \rangle^2 D_{p,q}^{N,m}(\tau)^{\frac{1}{m}}, $$

and noting that the assumptions on $W$ allow to estimate as follows the remaining integral with respect to $\xi$, for some $\delta > 0$,

$$ \int_{\mathbb{R}^d} \langle \xi \rangle^2 |\xi \widehat{W}(\xi)| d\xi \lesssim_{\delta} \left( \int_{\mathbb{R}^d} \langle \xi \rangle^{2\left(\frac{d+\delta}{2}+3\right)} |\widehat{W}(\xi)|^2 d\xi \right)^{\frac{1}{2}} \lesssim \|W\|_{H^{\frac{d+\delta}{2}+3}(\mathbb{R}^d)} < \infty, $$

we can use Hölder's inequality to estimate

$$ L \int_{\mathbb{R}^d} |\xi \widehat{W}(\xi)| \int_0^s e^{\gamma m(t-\tau)} e^{-\lambda_0 p(t-\tau)} \mathbb{E} \left[ \left| \int_{\mathbb{X}} e_{-\xi}(\mu_{\tau}^N - \mu_{\tau}) \right| \left| \int_{\mathbb{X}} (U_{\tau,t}\varphi)(\mu_{\tau}^N - \mu_{\tau}) \right|^{m-1} \right] d\tau d\xi \\ \lesssim L \int_0^s e^{\gamma m(t-\tau)} e^{-\lambda_0 p(t-\tau)} \mathbb{E} \left[ \left| \int_{\mathbb{X}} (U_{\tau,t}\varphi)(\mu_{\tau}^N - \mu_{\tau}) \right|^{m} \right] D_{p,q}^{N,m}(\tau)^{\frac{1}{m}} d\tau. $$
Inserting this into (3.6) and further using Young's inequality, we obtain (3.4).

**Step 2.** Estimate on the last term in (3.3):

$$ \frac{m(m-1)}{N} \int_0^s e^{\gamma m(t-\tau)} \mathbb{E} \left[ \left| \int_{\mathbb{X}} (U_{\tau,t}\varphi)(\mu_{\tau}^N - \mu_{\tau}) \right|^{m-2} \int_{\mathbb{X}} |\nabla_v (U_{\tau,t}\varphi)|^2 \mu_{\tau}^N \right] d\tau \\ \le \frac{\gamma m}{2} \int_0^s e^{\gamma m(t-\tau)} \mathbb{E} \left[ \left| \int_{\mathbb{X}} U_{\tau,t}\varphi(\mu_{\tau}^N - \mu_{\tau}) \right|^m \right] d\tau \\ + \gamma \left( \frac{Cm}{\gamma N} \right)^{\frac{m}{2}} \int_0^s e^{(\gamma-\lambda_0 p)m(t-\tau)} \mathbb{E} [Q_{2p}(\mu_{\tau}^N)^{\frac{m}{2}}] d\tau. \quad (3.7) $$

By Young's inequality, we can bound

$$ \frac{m(m-1)}{N} \int_0^s e^{\gamma m(t-\tau)} \mathbb{E} \left[ \left| \int_{\mathbb{X}} (U_{\tau,t}\varphi)(\mu_{\tau}^N - \mu_{\tau}) \right|^{m-2} \int_{\mathbb{X}} |\nabla_v (U_{\tau,t}\varphi)|^2 \mu_{\tau}^N \right] d\tau \\ \le \frac{\gamma m}{2} \int_0^s e^{\gamma m(t-\tau)} \mathbb{E} \left[ \left| \int_{\mathbb{X}} U_{\tau,t}\varphi(\mu_{\tau}^N - \mu_{\tau}) \right|^m \right] d\tau \\ + \gamma \left( \frac{2m}{\gamma N} \right)^{\frac{m}{2}} \int_0^s e^{\gamma m(t-\tau)} \mathbb{E} \left[ \left( \int_{\mathbb{X}} |\nabla_v (U_{\tau,t}\varphi)|^2 \mu_{\tau}^N \right)^{\frac{m}{2}} \right] d\tau, $$

and the claim (3.7) follows again from the ergodic estimates of Lemma 2.1 as in (3.5).

**Step 3.** Estimate on the first term in (3.3):

$$ \mathbb{E} \left[ \left| \int_{\mathbb{X}} (U_{0,t}\varphi)(\mu_0^N - \mu_\circ) \right|^m \right] \le C^m e^{-\lambda_0 pmt} D_{p,q}^{N,m}(0). \quad (3.8) $$

By definition of $D_{p,q}^{N,m}(0)$, we find

$$ \mathbb{E} \left[ \left| \int_{\mathbb{X}} (U_{0,t}\varphi)(\mu_0^N - \mu_\circ) \right|^m \right] \le \left\| \langle z \rangle^{-p} \left( U_{0,t}\varphi - \int_{\mathbb{X}} \varphi M \right) \right\|_{W^{2,q'}(\mathbb{X})}^m D_{p,q}^{N,m}(0), $$

and the claim (3.8) then follows from a direct use of the ergodic estimates of Lemma 2.1.

**Step 4.** Conclusion.

Inserting (3.4), (3.7), and (3.8) into (3.3), we are led to

$$ e^{\gamma m(t-s)} \mathbb{E} \left[ \left| \int_{\mathbb{X}} (U_{s,t}\varphi)(\mu_s^N - \mu_s) \right|^m \right] \le C^m e^{(\gamma-\lambda_0 p)m t} D_{p,q}^{N,m}(0) \\ + (C\kappa L - \frac{1}{2}\gamma)m \int_0^s e^{\gamma m(t-\tau)} \mathbb{E} \left[ \left| \int_{\mathbb{X}} (U_{\tau,t}\varphi)(\mu_{\tau}^N - \mu_{\tau}) \right|^m \right] d\tau \\ + C\kappa L \int_0^s e^{(\gamma-\lambda_0 p)m(t-\tau)} D_{p,q}^{N,m}(\tau) d\tau + C\kappa \int_0^s e^{(\gamma-\lambda_0 p)m(t-\tau)} \mathbb{E} \left[ (Q_p(\mu_{\tau}^N) - L)_+^m \right] d\tau \\ + \gamma \left( \frac{Cm}{\gamma N} \right)^{\frac{m}{2}} \int_0^s e^{(\gamma-\lambda_0 p)m(t-\tau)} \mathbb{E} \left[ Q_{2p}(\mu_{\tau}^N)^{\frac{m}{2}} \right] d\tau. $$

Provided that $\gamma \ge 2C\kappa L$, taking the supremum over $\varphi$ and $0 \le s \le t$, we deduce

$$ D_{p,q}^{N,m}(t) \le C^m e^{(\gamma-\lambda_0 p)m t} D_{p,q}^{N,m}(0) \\ + C\kappa L \int_0^t e^{(\gamma-\lambda_0 p)m(t-\tau)} D_{p,q}^{N,m}(\tau) d\tau + C\kappa \int_0^t e^{(\gamma-\lambda_0 p)m(t-\tau)} \mathbb{E} \left[ (Q_{2p}(\mu_{\tau}^N) - L)_+^m \right] d\tau \\ + \gamma \left( \frac{Cm}{\gamma N} \right)^{\frac{m}{2}} \int_0^t e^{(\gamma-\lambda_0 p)m(t-\tau)} \mathbb{E} \left[ Q_{2p}(\mu_{\tau}^N)^{\frac{m}{2}} \right] d\tau. \quad (3.9) $$
Noting that for $pq' \gg 1$ we have

$$ \mathbb{E}\left[\left|\int_X \langle z \rangle^{2p} (\mu_0^N - \mu_\circ)\right|^m\right] \le C^m D_{3p,q}^{N,m}(0), $$

and choosing $L \gg 1$ large enough (only depending on $d, \beta, a, Q_2(\mu_\circ)$, $\|W\|_{W^{1,\infty}(\mathbb{R}^d)}$), the concentration bound of Lemma 2.2 yields, for $2 \le m \le N$,

$$ \int_0^t e^{(\gamma-\lambda_0 p)m(t-\tau)} \mathbb{E}[(Q_{2p}(\mu_\tau^N) - L)_+^m] d\tau \\ \le \left(\frac{CLm}{N}\right)^{\frac{m}{2}} \int_0^t e^{(\gamma-\lambda_0 p)m(t-\tau)} d\tau + C^m D_{3p,q}^{N,m}(0) \int_0^t e^{(\gamma-\lambda_0 p)m(t-\tau)} e^{-\lambda_0 pm\tau} d\tau, $$

and thus, provided that $\gamma \le \frac{1}{2}\lambda_0 p$,

$$ \int_0^t e^{(\gamma-\lambda_0 p)m(t-\tau)} \mathbb{E}[(Q_{2p}(\mu_\tau^N) - L)_+^m] d\tau \le \left(\frac{CLm}{N}\right)^{\frac{m}{2}} + C^m e^{(\gamma-\lambda_0 p)mt} D_{3p,q}^{N,m}(0). $$

Similarly, we can also bound

$$ \left(\frac{Cm}{N}\right)^{\frac{m}{2}} \int_0^t e^{(\gamma-\lambda_0 p)m(t-\tau)} \mathbb{E}[Q_{2p}(\mu_\tau^N)^{\frac{m}{2}}] d\tau \\ \le \left(\frac{CLm}{N}\right)^{\frac{m}{2}} \int_0^t e^{(\gamma-\lambda_0 p)m(t-\tau)} d\tau + C^m \int_0^t e^{(\gamma-\lambda_0 p)m(t-\tau)} \mathbb{E}[(Q_{2p}(\mu_\tau^N) - L)_+^m] d\tau \\ \le \left(\frac{CLm}{N}\right)^{\frac{m}{2}} + C^m e^{(\gamma-\lambda_0 p)mt} D_{3p,q}^{N,m}(0). $$

Inserting these bounds into (3.9), choosing $\gamma := \frac{1}{2}\lambda_0 p$, and assuming that $\kappa$ is small enough in the sense of $2C\kappa L \le \frac{1}{2}\lambda_0 p$, we find

$$ D_{p,q}^{N,m}(t) \le \left(\frac{CLm}{N}\right)^{\frac{m}{2}} + C^m e^{-\frac{1}{2}\lambda_0 pm t} D_{3p,q}^{N,m}(0) + C\kappa L \int_0^t e^{-\frac{1}{2}\lambda_0 pm(t-\tau)} D_{p,q}^{N,m}(\tau) d\tau, $$

and thus, by Grönwall's inequality,

$$ D_{p,q}^{N,m}(t) \le \left(\frac{CLm}{N}\right)^{\frac{m}{2}} + C^m e^{-\frac{1}{2}\lambda_0 pm t} D_{3p,q}^{N,m}(0) \\ + C\kappa L e^{-\frac{1}{2}\lambda_0 pm t} \int_0^t e^{C\kappa L(t-\tau)} \left( e^{\frac{1}{2}\lambda_0 pm\tau} \left(\frac{CLm}{N}\right)^{\frac{m}{2}} + C^m D_{3p,q}^{N,m}(0) \right) d\tau. $$

As the smallness requirement on $\kappa$ precisely ensures $C\kappa L \le \frac{1}{4}\lambda_0 p$, we get after computing the remaining time integrals,

$$ D_{p,q}^{N,m}(t) \le \left(\frac{CLm}{N}\right)^{\frac{m}{2}} + C^m e^{-\frac{1}{4}\lambda_0 pm t} D_{3p,q}^{N,m}(0). $$

Recalling (3.2), the conclusion follows (up to renaming $\lambda_0$). $\square$

## 4. CREATION OF CHAOS

This section is devoted to the proof of Theorem 1.1 using the above preliminary results. Our argument shows how the estimates of Proposition 3.1 can be turned into suboptimal mean-field estimates for marginals, which in turn can be made optimal in $N$ after combination with direct estimates on the BBGKY hierarchy. A similar idea will be used again for the proof of Theorem 1.2 in Section 5.

We divide the proof into the next three subsections. Let $\lambda_0, p, q, \kappa$ be as in Proposition 3.1, and assume $p \le \frac{1}{6}$ (say). In the sequel, all multiplicative constants are allowed to depend on $d, \beta, a, p, q, Q_2(\mu_\circ), Q_2(F_\circ^{N,1})$, and $\|W\|_{W^{d+3,\infty}\cap H^s(\mathbb{R}^d)}$ for some $s > \frac{d}{2} + 5$. A subscript 'm' to constants is used to indicate further dependence on $m$ and $Q_{mp}(F_\circ^{N,1})$.
4.1. Reformulation of Proposition 3.1. In this section, we show that Proposition 3.1 implies the following suboptimal mean-field estimate, with a clean dependence on initial correlations:

$$ \|F_t^{N,1} - \mu_t\|_{W_p^{-2,q}(\mathbb{X})}^2 + \|G_t^{N,2}\|_{W_p^{-2,q}(\mathbb{X})\otimes 2} \lesssim N^{-1} + e^{-2\lambda_0 pt} \left( \|F_\circ^{N,1} - \mu_\circ\|_{W_{3p}^{-2,q}(\mathbb{X})}^2 + \|G_\circ^{N,2}\|_{W_{3p}^{-2,q}(\mathbb{X})\otimes 2} \right). \quad (4.1) $$

Noting that

$$ \mathbb{E}\left[\left(\int_{\mathbb{X}} \varphi(\mu_t^N - \mu_t)\right)\right] = \int_{\mathbb{X}} \varphi(F_t^{N,1} - \mu_t), $$

and that

$$ \begin{aligned} A_{p,q}^{N,2}(t) &\ge \sup \left\{ \left( \int_{\mathbb{X}} \varphi(F_t^{N,1} - \mu_t) \right)^2 : \varphi \in C_c^\infty(\mathbb{X}), \|\langle \cdot \rangle^{-p} \varphi\|_{W^{2,q'}(\mathbb{X})} \le 1 \right\} \\ &= \|F_t^{N,1} - \mu_t\|_{W_p^{-2,q}(\mathbb{X})}^2, \end{aligned} $$

Proposition 3.1 yields for $m=2$,

$$ \|F_t^{N,1} - \mu_t\|_{W_p^{-2,q}(\mathbb{X})}^2 \lesssim N^{-1} + e^{-2\lambda_0 pt} A_{3p,q}^{N,2}(0). \quad (4.2) $$

Also note that

$$ \mathbb{E}\left[\left(\int_{\mathbb{X}} \varphi(\mu_t^N - \mu_t)\right)^2\right] = \left(\int_{\mathbb{X}} \varphi(F_t^{N,1} - \mu_t)\right)^2 + \int_{\mathbb{X}^2} \varphi^{\otimes 2} G_t^{N,2} + \frac{1}{N} \int_{\mathbb{X}} \varphi^2 F_t^{N,1} - \frac{1}{N} \int_{\mathbb{X}^2} \varphi^{\otimes 2} F_t^{N,2}. $$

Taking the supremum over $\varphi$ in this identity, and using the Sobolev embedding together with the moment bounds of Lemma 2.2, which yields in particular

$$ \|F_t^{N,2}\|_{W_p^{-2,q}(\mathbb{X})\otimes 2} \lesssim \int_{\mathbb{X}^2} \langle z \rangle^p \langle z_* \rangle^q F_t^{N,2}(z, z_*) \, dz \, dz_* \le \int_{\mathbb{X}} \langle z \rangle^{2p} F_t^{N,1} = Q_{2p}(F_t^{N,1}), \quad (4.3) $$

we find

$$ \left| A_{p,q}^{N,2}(t) - \|G_t^{N,2}\|_{W_p^{-2,q}(\mathbb{X})\otimes 2} \right| \le \|F_t^{N,1} - \mu_t\|_{W_p^{-2,q}(\mathbb{X})}^2 + CN^{-1}Q_{2p}(F_t^{N,1}). $$

Appealing to Lemma 2.2 in form of $Q_{2p}(F_t^{N,1}) \lesssim 1 + Q_{2p}(F_\circ^{N,1}) \lesssim 1$, this means

$$ \left| A_{p,q}^{N,2}(t) - \|G_t^{N,2}\|_{W_p^{-2,q}(\mathbb{X})\otimes 2} \right| \lesssim \|F_t^{N,1} - \mu_t\|_{W_p^{-2,q}(\mathbb{X})}^2 + CN^{-1}. \quad (4.4) $$

Combining this with the result of Proposition 3.1 to get a control on $G_t^{N,2}$, we obtain

$$ \|G_t^{N,2}\|_{W_p^{-2,q}(\mathbb{X})\otimes 2} \lesssim \|F_t^{N,1} - \mu_t\|_{W_p^{-2,q}(\mathbb{X})}^2 + N^{-1} + e^{-2\lambda_0 pt} A_{3p,q}^{N,2}(0), $$

which thus yields, together with (4.2),

$$ \|F_t^{N,1} - \mu_t\|_{W_p^{-2,q}(\mathbb{X})}^2 + \|G_t^{N,2}\|_{W_p^{-2,q}(\mathbb{X})\otimes 2} \lesssim N^{-1} + e^{-2\lambda_0 pt} A_{3p,q}^{N,2}(0). \quad (4.5) $$

Now applying (4.4) at $t=0$ with $p$ replaced by $3p$, we find

$$ A_{3p,q}^{N,2}(0) \le \|F_\circ^{N,1} - \mu_\circ\|_{W_{3p}^{-2,q}(\mathbb{X})}^2 + \|G_\circ^{N,2}\|_{W_{3p}^{-2,q}(\mathbb{X})\otimes 2} + CN^{-1}, $$

and the claim (4.1) follows.
4.2. Improved bounds via BBGKY hierarchy. In this section, we use the BBGKY hierarchy to upgrade (4.1) into the following optimal-in-N mean-field bound on the first marginal:

$$ \begin{aligned} & \|F_t^{N,1} - \mu_t\|_{W_p^{-3,q}(\mathbb{X})} \\ & \lesssim N^{-1} + e^{-\frac{1}{2}\lambda_0 pt} \left( \|F_\circ^{N,1} - \mu_\circ\|_{W_p^{-2,q}(\mathbb{X})} + \|F_\circ^{N,1} - \mu_\circ\|_{W_{3p}^{-2,q}(\mathbb{X})}^2 + \|G_\circ^{N,2}\|_{W_{3p}^{-2,q}(\mathbb{X})}^2 \right) \\ & \lesssim N^{-1} + e^{-\frac{1}{2}\lambda_0 pt} \left( \|F_\circ^{N,1} - \mu_\circ\|_{W_{3p}^{-2,q}(\mathbb{X})} + \|G_\circ^{N,2}\|_{W_{3p}^{-2,q}(\mathbb{X})}^2 \right). \end{aligned} \quad (4.6) $$

Integrating the Liouville equation (1.3) with respect to its last $N-1$ coordinates, and decomposing the second marginal as $F^{N,2} = F^{N,1} \otimes F^{N,1} + G^{N,2}$, we find that the first marginal $F^{N,1}$ satisfies the following BBGKY equation,

$$ \partial_t(F^{N,1} - \mu) = R_\mu(F^{N,1} - \mu) + \kappa \nabla W * (F^{N,1} - \mu) \cdot \nabla_v F^{N,1} \\ + \kappa \int_{\mathbb{X}} \nabla W(\cdot - x_*) \cdot \nabla_v (G^{N,2} - \frac{1}{N}F^{N,2})(\cdot, z_*) dz_*, $$

where $R_\mu$ stands for the following (truncated) linearized mean-field operator (linearized at the mean-field solution $\mu$),

$$ R_{\mu}h := \frac{1}{2}\mathrm{div}_{v}((\nabla_{v} + \beta v)h) - v \cdot \nabla_{x}h + (\nabla A + \kappa \nabla W * \mu) \cdot \nabla_{v}h. $$

Given $\varphi \in C_c^\infty(\mathbb{X})$, recalling the definition (2.1) of the backward linearized evolution $\{U_{s,t}\}_{0 \le s \le t}$, we may then compute

$$ \begin{aligned} \int_{\mathbb{X}} \varphi(F_t^{N,1} - \mu_t) = & \int_{\mathbb{X}} (U_{0,t}\varphi)(F_\circ^{N,1} - \mu_\circ) - \kappa \int_0^t \left( \int_{\mathbb{X}} \nabla_v (U_{s,t}\varphi) \cdot \nabla W * (F_s^{N,1} - \mu_s) F_s^{N,1} \right) ds \\ & - \kappa \int_0^t \left( \int_{\mathbb{X}^2} \nabla_v (U_{s,t}\varphi)(z) \cdot \nabla W(x - x_*) (G_s^{N,2} - \frac{1}{N}F_s^{N,2})(z, z_*) dz dz_* \right) ds. \end{aligned} \quad (4.7) $$

Let us examine the three right-hand side terms separately. Using the ergodic estimates of Lemma 2.1, the first one is bounded by

$$ \left| \int_{\mathbb{X}} (U_{0,t}\varphi)(F_\circ^{N,1} - \mu_\circ) \right| \lesssim e^{-\lambda_0 pt} \|\langle z \rangle^{-p} \varphi\|_{W^{2,q'}(\mathbb{X})} \|F_\circ^{N,1} - \mu_\circ\|_{W_p^{-2,q}(\mathbb{X})}. $$

For the second term, using again the ergodic estimates of Lemma 2.1 along with the Sobolev embedding, we can bound

$$ \begin{aligned} & \left| \int_{\mathbb{X}} \nabla_v (U_{s,t}\varphi) \cdot \nabla W * (F_s^{N,1} - \mu_s) F_s^{N,1} \right| \\ & \lesssim Q_p(F_s^{N,1}) \|\langle z \rangle^{-p} \nabla_v (U_{s,t}\varphi)\|_{L^\infty(\mathbb{X})} \sup_{x \in \mathbb{R}^d} \left| \int_{\mathbb{X}} \nabla W(x - \cdot) (F_s^{N,1} - \mu_s) \right| \\ & \lesssim Q_p(F_s^{N,1}) e^{-\lambda_0 p(t-s)} \|\langle z \rangle^{-p} \varphi\|_{W^{2,q'}(\mathbb{X})} \|F_s^{N,1} - \mu_s\|_{W_p^{-3,q}(\mathbb{X})}, \end{aligned} $$

and thus, further using the moment bounds of Lemma 2.2,

$$ \left| \int_{\mathbb{X}} \nabla_v (U_{s,t}\varphi) \cdot \nabla W * (F_s^{N,1} - \mu_s) F_s^{N,1} \right| \lesssim e^{-\lambda_0 p(t-s)} \|\langle z \rangle^{-p} \varphi\|_{W^{2,q'}(\mathbb{X})} \|F_s^{N,1} - \mu_s\|_{W_p^{-3,q}(\mathbb{X})}. $$

For the last term in (4.7), we aim to appeal to the bound (4.1) of Step 1 for $G^{N,2}$. Note, however, that the map $(z, z_*) \mapsto \nabla_v(U_{s,t}\varphi)(z) \cdot \nabla W(x - x_*)$ is not an admissible test function for such estimates in $W_p^{-2,q}(\mathbb{X})^{\otimes 2}$, cf. (1.13). To remedy this, we decompose $W$ as a superposition of Fourier modes, which indeed allows to reduce to a tensorized test function. More precisely, by polarization, using that
$G_s^{N,2} - \frac{1}{N}F_s^{N,2}$ is symmetric in its two variables, and recalling the short-hand notation $e_\xi(x) = e^{ix \cdot \xi}$, we can bound

$$
\begin{align}
& \left| \int_{\mathbb{X}^2} \nabla_v (U_{s,t}\varphi)(z) \cdot \nabla W(x-x_*) (G_s^{N,2} - \frac{1}{N}F_s^{N,2})(z, z_*) \, dz \, dz_* \right| \nonumber \\
& \le \int_{\mathbb{R}^d} |\xi \widehat{W}(\xi)| \left| \int_{\mathbb{X}^2} e_\xi(x_*) e_{-\xi}(x) \nabla_v (U_{s,t}\varphi)(z) (G_s^{N,2} - \frac{1}{N}F_s^{N,2})(z, z_*) \, dz \, dz_* \right| d\xi \tag{4.8} \\
& \le 4 \|G_s^{N,2} - \frac{1}{N}F_s^{N,2}\|_{W_p^{-2,q}(\mathbb{X}) \otimes 2} \int_{\mathbb{R}^d} |\xi \widehat{W}(\xi)| \|\langle z \rangle^{-p} e_\xi\|_{W^{2,q'}(\mathbb{X})} \|\langle z \rangle^{-p} e_\xi \nabla_v (U_{s,t}\varphi)\|_{W^{2,q'}(\mathbb{X})} \, d\xi. \nonumber
\end{align}
$$

Noting that

$$
\|\langle z \rangle^{-p} e_{\xi}\|_{W^{2,q'}(\mathbb{X})} \lesssim \langle \xi \rangle^2, \quad \|\langle z \rangle^{-p} e_{\xi} \nabla_v (U_{s,t}\varphi)\|_{W^{2,q'}(\mathbb{X})} \lesssim \langle \xi \rangle^2 e^{-\lambda_0 p(t-s)} \|\langle z \rangle^{-p} \varphi\|_{W^{3,q'}(\mathbb{X})},
$$

we deduce

$$
\begin{aligned}
& \left| \int_{\mathbb{X}^2} \nabla_v (U_{s,t}\varphi)(z) \cdot \nabla W(x-x_*) (G_s^{N,2} - \frac{1}{N}F_s^{N,2})(z, z_*) dz dz_* \right| \\
& \lesssim e^{-\lambda_0 p(t-s)} \|\langle z \rangle^{-p} \varphi\|_{W^{3,q'}(\mathbb{X})} \|G_s^{N,2} - \frac{1}{N}F_s^{N,2}\|_{W_p^{-2,q}(\mathbb{X}) \otimes 2} \int_{\mathbb{R}^d} \langle \xi \rangle^4 |\xi \widehat{W}(\xi)| \, d\xi,
\end{aligned}
$$

where the last factor can be controlled by $\|W\|_{H^s(\mathbb{R}^d)}$ for $s > \frac{d}{2} + 5$. Combining these estimates for the different terms in (4.7), and taking the supremum over $\varphi$, we get

$$
\begin{aligned}
\| F_t^{N,1} - \mu_t \|_{W_p^{-3,q}(\mathbb{X})} & \lesssim e^{-\lambda_0 pt} \| F_o^{N,1} - \mu_o \|_{W_p^{-2,q}(\mathbb{X})} \\
& \quad + \kappa \int_0^t e^{-\lambda_0 p(t-s)} \| F_s^{N,1} - \mu_s \|_{W_p^{-3,q}(\mathbb{X})} \, ds + \kappa \int_0^t e^{-\lambda_0 p(t-s)} \| G_s^{N,2} - \frac{1}{N} F_s^{N,2} \|_{W_p^{-2,q}(\mathbb{X}) \otimes 2} \, ds.
\end{aligned}
$$

Now inserting the bound (4.1) of Step 1 for $G^{N,2}$, using (4.3) together with $Q_{2p}(F_t^{N,1}) \lesssim 1$, we obtain

$$
\begin{aligned}
\| F_t^{N,1} - \mu_t \|_{W_p^{-3,q}(\mathbb{X})} & \le C N^{-1} + C \kappa \int_0^t e^{-\lambda_0 p(t-s)} \| F_s^{N,1} - \mu_s \|_{W_p^{-3,q}(\mathbb{X})} \, ds \\
& \quad + C e^{-\lambda_0 pt} \left( \| F_o^{N,1} - \mu_o \|_{W_p^{-2,q}(\mathbb{X})} + \| F_o^{N,1} - \mu_o \|_{W_{3p}^{-2,q}(\mathbb{X})}^2 + \| G_o^{N,2} \|_{W_{3p}^{-2,q}(\mathbb{X}) \otimes 2} \right).
\end{aligned}
$$

Provided that $\kappa$ is chosen small enough so that $C\kappa < \frac{1}{2}\lambda_0 p$, the claim (4.6) follows from Grönwall's inequality.

4.3. **Conclusion.** We now use the previous steps to establish optimal-in-$N$ mean-field bounds on all marginals. For this purpose, we first express the marginals in terms of moments of the empirical measure. Let $\varphi \in C_c^\infty(\mathbb{X})$ be momentarily fixed with $\|\langle z \rangle^{-p}\varphi\|_{W^{3,q'}(\mathbb{X})} \le 1$. For all $m \ge 1$, we have

$$
\begin{align*}
\mathbb{E} \left[ \left( \int_{\mathbb{X}} \varphi (\mu_t^N - \mu_t) \right)^m \right] &= N^{-m} \sum_{i_1, \dots, i_m = 1}^{N} \mathbb{E} \left[ \prod_{l=1}^{m} \left( \varphi(Z_t^{i_l, N}) - \int_{\mathbb{X}} \varphi \mu_t \right) \right] \\
&= N^{-m} \sum_{\pi \vdash [m]} N(N-1) \dots (N-\sharp \pi + 1) \int_{\mathbb{X}^{\sharp \pi}} \left( \bigotimes_{B \in \pi} \left( \varphi - \int_{\mathbb{X}} \varphi \mu_t \right)^{\sharp B} \right) F_t^{N, \sharp \pi},
\end{align*}
$$

where we use the same notation as in (1.8) for sums over partitions. Note that the moment bounds of Lemma 2.2 together with the Sobolev embedding and the choice of $\varphi$ yield, for all $\pi \vdash [m]$, mimicking (4.3),

$$
\left| \int_{\mathbb{X}^{\sharp \pi}} \left( \bigotimes_{B \in \pi} \left( \varphi - \int_{\mathbb{X}} \varphi \mu_t \right)^{\sharp B} \right) F_t^{N, \sharp \pi} \right| \lesssim_m \int_{\mathbb{X}} \langle z \rangle^{mp} F_t^{N,1} \le (Cm)^{\frac{mp}{2}} + C^m e^{-\lambda_0 pmt} Q_{mp}(F_o^{N,1}) \lesssim_m 1.
$$
Hence, neglecting all $O(N^{-1})$ terms in the above, we get

$$ \left| \mathbb{E} \left[ \left( \int_{\mathbb{X}} \varphi(\mu_t^N - \mu_t) \right)^m \right] - \int_{\mathbb{X}^m} \left( \varphi - \int_{\mathbb{X}} \varphi \mu_t \right)^{\otimes m} F_t^{N,m} \right| \lesssim_m N^{-1}. $$

Noting that

$$ \begin{aligned} \int_{\mathbb{X}^m} \left( \varphi - \int_{\mathbb{X}} \varphi \mu_t \right)^{\otimes m} F_t^{N,m} &= \sum_{j=0}^{m} \binom{m}{j} (-1)^{m-j} \left( \int_{\mathbb{X}} \varphi \mu_t \right)^{m-j} \int_{\mathbb{X}^j} \varphi^{\otimes j} F_t^{N,j} \\ &= \sum_{j=1}^{m} \binom{m}{j} (-1)^{m-j} \left( \int_{\mathbb{X}} \varphi \mu_t \right)^{m-j} \int_{\mathbb{X}^j} \varphi^{\otimes j} (F_t^{N,j} - \mu_t^{\otimes j}), \end{aligned} $$

and focusing on the term corresponding to $j = m$ in the sum, we deduce

$$ \left| \mathbb{E} \left[ \left( \int_{\mathbb{X}} \varphi(\mu_t^N - \mu_t) \right)^m \right] - \int_{\mathbb{X}^m} \varphi^{\otimes m} (F_t^{N,m} - \mu_t^{\otimes m}) \right| \lesssim_m N^{-1} + \sum_{j=1}^{m-1} \left| \int_{\mathbb{X}^j} \varphi^{\otimes j} (F_t^{N,j} - \mu_t^{\otimes j}) \right|. \quad (4.9) $$

Taking the supremum over $\varphi \in C_c^\infty(\mathbb{X})$ with $\|\langle z \rangle^{-p} \varphi\|_{W^{3,q'}(\mathbb{X})} \le 1$, we deduce in particular for all $m \ge 1$,

$$ \|F_t^{N,m} - \mu_t^{\otimes m}\|_{W_p^{-3,q}(\mathbb{X})^{\otimes m}} \lesssim_m N^{-1} + A_{p,q}^{N,m}(t) + \sum_{j=1}^{m-1} \|F_t^{N,j} - \mu_t^{\otimes j}\|_{W_p^{-3,q}(\mathbb{X})^{\otimes j}}, \quad (4.10) $$

and thus, after a direct iteration,

$$ \|F_t^{N,m} - \mu_t^{\otimes m}\|_{W_p^{-3,q}(\mathbb{X})^{\otimes m}} \lesssim_m N^{-1} + \sum_{j=2}^{m} A_{p,q}^{N,j}(t) + \|F_t^{N,1} - \mu_t\|_{W_p^{-3,q}(\mathbb{X})}. $$

Inserting the result of Proposition 3.1 and the result (4.6) of Step 2, we then get

$$ \begin{aligned} \|F_t^{N,m} - \mu_t^{\otimes m}\|_{W_p^{-3,q}(\mathbb{X})^{\otimes m}} &\lesssim_m N^{-1} + \sum_{j=2}^{m} e^{-\lambda_0 pj t} A_{3p,q}^{N,j}(0) \\ &\quad + e^{-\frac{1}{2}\lambda_0 pt} \left( \|F_\circ^{N,1} - \mu_\circ\|_{W_{3p}^{-2,q}(\mathbb{X})} + \|G_\circ^{N,2}\|_{W_{3p}^{-2,q}(\mathbb{X})^{\otimes 2}} \right). \end{aligned} \quad (4.11) $$

From here, it remains to express the $A_{3p,q}^{N,j}(0)$'s in terms of initial correlations. Noting that

$$ \mathbb{E} \left[ \left| \int_{\mathbb{X}} \varphi(\mu_t^N - \mu_t) \right|^j \right] \le \begin{cases} \mathbb{E} \left[ \left( \int_{\mathbb{X}} \varphi(\mu_t^N - \mu_t) \right)^j \right] & : j \text{ even,} \\ \mathbb{E} \left[ \left( \int_{\mathbb{X}} \varphi(\mu_t^N - \mu_t) \right)^{j-1} \right]^{\frac{1}{2}} \mathbb{E} \left[ \left( \int_{\mathbb{X}} \varphi(\mu_t^N - \mu_t) \right)^{j+1} \right]^{\frac{1}{2}} & : j \text{ odd,} \end{cases} \quad (4.12) $$

using (4.9) at $t=0$, and taking the supremum over $\varphi \in C_c^\infty(\mathbb{X})$ with $\|\langle z \rangle^{-3p} \varphi\|_{W^{2,q'}(\mathbb{X})} \le 1$, we find for all $j \ge 2$,

$$ A_{3p,q}^{N,j}(0) \lesssim_j N^{-1} + \sum_{l=1}^{\tilde{j}} \|F_\circ^{N,l} - \mu_\circ^{\otimes l}\|_{W_{3p}^{-2,q}(\mathbb{X})^{\otimes l}}, \quad \tilde{j} := j + 1_{j \text{ odd}}. $$

Inserting this into (4.11), we get

$$ \|F_t^{N,m} - \mu_t^{\otimes m}\|_{W_p^{-3,q}(\mathbb{X})^{\otimes m}} \lesssim_m N^{-1} + e^{-\frac{1}{2}\lambda_0 pt} \|G_\circ^{N,2}\|_{W_{3p}^{-2,q}(\mathbb{X})^{\otimes 2}} + \sum_{j=1}^{\tilde{m}} e^{-\frac{1}{2}\lambda_0 pj t} \|F_\circ^{N,j} - \mu_\circ^{\otimes j}\|_{W_{3p}^{-2,q}(\mathbb{X})^{\otimes j}}. $$

Decomposing

$$ G_\circ^{N,2} = (F_\circ^{N,2} - \mu_\circ^{\otimes 2}) - (F_\circ^{N,1} - \mu_\circ) \otimes F_\circ^{N,1} - \mu_\circ \otimes (F_\circ^{N,1} - \mu_\circ), \quad (4.13) $$

the conclusion (1.14) follows for all $m \ge 1$ (up to renaming $\lambda_0$). $\square$

<figure><img src="image_2.png" /></figure>
## 5. CORRELATION ESTIMATES

This section is devoted to the proof of Theorem 1.2. Similarly as for Theorem 1.1, the argument starts from the estimates of Proposition 3.1: we show how these can be used to deduce suboptimal estimates on correlations, which in turn can be made optimal in $N$ after combination with direct estimates on the BBGKY hierarchy. Note that this time we will need to iterate BBGKY estimates multiple times to achieve optimality. We divide the proof into the next three subsections.

While in the introduction we defined correlation functions so as to satisfy the cluster expansion (1.8), we recall that they can be defined more explicitly as polynomial combinations of marginals:

$$G^{N,1} = F^{N,1},$$

$$G^{N,2} = F^{N,2} - F^{N,1} \otimes F^{N,1},$$

$$G^{N,3} = \operatorname{sym}(F^{N,3} - 3F^{N,2} \otimes F^{N,1} + 2(F^{N,1})^{\otimes 3}),$$

$$G^{N,4} = \operatorname{sym}(F^{N,4} - 4F^{N,3} \otimes F^{N,1} - 3F^{N,2} \otimes F^{N,2} + 12F^{N,2} \otimes (F^{N,1})^{\otimes 2} - 6(F^{N,1})^{\otimes 4}),$$

and so on, where the symbol ‘sym’ stands for the symmetrization of coordinates. More generally, we can write for all $2 \le m \le N$,

$$G^{N,m}(z_1, \dots, z_m) = \sum_{\pi \vdash [m]} (\sharp \pi - 1)! (-1)^{\sharp \pi - 1} \prod_{B \in \pi} F^{N, \sharp B}(z_B), \quad (5.1)$$

where we use a similar notation as in (1.8) for sums over partitions and where $\sharp \pi$ stands for the number of blocks in a partition $\pi$. We also use the following abbreviation for functional spaces: for all $r, p, q, m$,

$$W_{p,m}^{-r,q} := W_p^{-r,q}(\mathbb{X})^{\otimes m}.$$

5.1. From moment mean-field estimates to correlations. From now on, we let $\mu$ be the weak solution of the mean-field equation (1.4) with initial condition

$$\mu_0 := F_0^{N,1}.$$

In this section, we examine how the correlation functions $\{G^{N,m}\}_{2 \le m \le N}$ can be related to moments of the empirical fluctuation $\mu_t^N - \mu_t$ in terms of the $A_{p,q}^{N,m}$'s defined in (3.1), and we show that Proposition 3.1 implies the following hierarchy of estimates on correlations with suboptimal $N$-dependence ($N^{-m/2}$ instead of $N^{1-m}$, cf. (1.9)). Note the unfortunate discrepancy in the weights on the two sides of the inequality: this is natural in view of Proposition 3.1, but, since the regularity and the weight are intertwined through the condition $pq' \gg 1$, a careful choice of parameters will be required when iterating across the hierarchy of correlations.

**Proposition 5.1.** Let $\lambda_0, p, q, \kappa$ be as in Proposition 3.1 with $\mu_0 = F_0^{N,1}$. Then we have for all $2 \le m \le N, r \ge 2, 0 < p_0 \le p$, and $t \ge 0$,

$$\|G_t^{N,m}\|_{W_{p_0,m}^{-r,q}} \lesssim_m N^{-\frac{m}{2}} + \sum_{1 \le k \le l \le m-1} \sum_{j_1, \dots, j_k \ge 1} N^{l-k-m+1} \prod_{i=1}^k \|G_t^{N,j_i}\|_{W_{mp_0,j_i}^{-r,q}} \\ + e^{-\lambda_0 pmt} \times \left\{ \begin{array}{ll} \sum_{\ell=1}^{m/2} \sum_{\ell \le K \le L \le m} \sum_{j_1, \dots, j_K \ge 1} N^{L-K-m+\ell} \prod_{i=1}^K \|G_0^{N,j_i}\|_{W_{3mp,j_i}^{-2,q}} & : m \text{ even,} \\ \sum_{\ell=1}^m \sum_{\ell \le K \le L \le 2m} \sum_{j_1, \dots, j_K \le m+1} (N^{L-K-2m+\ell} \prod_{i=1}^K \|G_0^{N,j_i}\|_{W_{3mp,j_i}^{-2,q}})^{\frac{1}{2}} & : m \text{ odd,} \end{array} \right.$$

where the multiplicative constant only depends on $d, \beta, a, p_0, p, q, r, m, Q_2(F_0^{N,1})$, and $\|W\|_{W^{d+3,\infty} \cap H^s(\mathbb{R}^d)}$ for some $s > \frac{d}{2} + 3$.

*Proof.* Let us first recall that the $m$-th cumulant of a bounded random variable $X$ is defined by

$$\kappa^m(X) := \left. \frac{d}{dt} \right|_{t=0} \log \mathbb{E}[e^{tX}],$$
which can be expressed as a polynomial combination of moments of X,

$$ \kappa^m[X] = \sum_{\pi \vdash [m]} (-1)^{\sharp\pi-1}(\sharp\pi-1)! \prod_{B \in \pi} \mathbb{E}[X^{\sharp B}]. \quad (5.2) $$

We will use the standard link between cumulants of the empirical measure and correlation functions. We first used this link in [11, Section 4] and we further included a self-contained statement and proof in [3, Lemma 2.6]. It can be stated as follows: for all $\varphi \in C_c^\infty(\mathbb{X})$, $1 \le m \le N$, and $t \ge 0$,

$$ \left| \kappa^m \left[ \int_{\mathbb{X}} \varphi \mu_t^N \right] - \int_{\mathbb{X}^m} \varphi^{\otimes m} G_t^{N,m} \right| \lesssim_m \sum_{\substack{\pi \vdash [m] \\ \sharp\pi < m}} \sum_{\rho \vdash \pi} N^{\sharp\pi - \sharp\rho - m + 1} \left| \int_{\mathbb{X}^{\sharp\pi}} \left( \bigotimes_{B \in \pi} \varphi^{\sharp B} \right) \left( \bigotimes_{D \in \rho} G_t^{N, \sharp D}(z_D) \right) dz_{\pi} \right|. \quad (5.3) $$

For $m \ge 2$, as cumulants are invariant under translation, we can write $\kappa^m[\int_{\mathbb{X}} \varphi \mu_t^N] = \kappa^m[\int_{\mathbb{X}} \varphi(\mu_t^N - \mu_t)]$. From (5.2), it follows that

$$ \left| \kappa^m \left[ \int_{\mathbb{X}} \varphi \mu_t^N \right] \right| \lesssim_m \mathbb{E} \left[ \left| \int_{\mathbb{X}} \varphi (\mu_t^N - \mu_t) \right|^m \right]. $$

Hence, we can deduce for all $\varphi \in C_c^\infty(\mathbb{X})$, $2 \le m \le N$, and $t \ge 0$,

$$ \left| \int_{\mathbb{X}^m} \varphi^{\otimes m} G_t^{N,m} \right| \lesssim_m \mathbb{E} \left[ \left| \int_{\mathbb{X}} \varphi (\mu_t^N - \mu_t) \right|^m \right] + \sum_{\substack{\pi \vdash [m] \\ \sharp\pi < m}} \sum_{\rho \vdash \pi} N^{\sharp\pi - \sharp\rho - m + 1} \left| \int_{\mathbb{X}^{\sharp\pi}} \left( \bigotimes_{B \in \pi} \varphi^{\sharp B} \right) \left( \bigotimes_{D \in \rho} G_t^{N, \sharp D}(z_D) \right) dz_{\pi} \right|. $$

For $r \ge 2$, $1 < q \le 2$ with $q' > d$, and $0 < p_0 \le p \le 1$, taking the supremum over the test function $\varphi$, recognizing the definition (3.1) of $A_{p_0,q}^{N,m}(t)$, noting that for $p_0 \le p$ we have $A_{p_0,q}^{N,m} \le A_{p,q}^{N,m}$, and also noting that the Sobolev embedding implies

$$ \|\langle z \rangle^{-p_0 \ell} \varphi^{\ell}\|_{W^{r,q'}(\mathbb{X})} \lesssim \|\langle z \rangle^{-p_0} \varphi\|_{W^{r,q'}(\mathbb{X})}^{\ell}, \quad \text{for } \ell \ge 1, \quad (5.4) $$

we deduce for all $2 \le m \le N$ and $t \ge 0$,

$$ \|G_t^{N,m}\|_{W_{p_0,m}^{-r,q}} \lesssim_m A_{p,q}^{N,m}(t) + \sum_{1 \le k \le l \le m-1} \sum_{j_1, \dots, j_k \ge 1} N^{l-k+m+1} \prod_{i=1}^k \|G_t^{N,j_i}\|_{W_{mp_0,j_i}^{-r,q}}. $$

Hence, by Proposition 3.1, for suitable $\lambda_0, p, q, \kappa$,

$$ \|G_t^{N,m}\|_{W_{p_0,m}^{-r,q}} \lesssim_m N^{-\frac{m}{2}} + e^{-\lambda_0 pmt} A_{3p,q}^{N,m}(0) + \sum_{1 \le k \le l \le m-1} \sum_{j_1, \dots, j_k \ge 1} N^{l-k+m+1} \prod_{i=1}^k \|G_t^{N,j_i}\|_{W_{mp_0,j_i}^{-r,q}}. $$

It remains to express $A_{3p,q}^{N,m}(0)$ in the right-hand side in terms of initial correlations. For that purpose, we first appeal to the cluster expansion of moments in terms of cumulants: given $\varphi \in C_c^\infty(\mathbb{X})$, by definition of cumulants, we can expand for all $m \ge 1$,

$$ \mathbb{E} \left[ \left( \int_{\mathbb{X}} \varphi (\mu_0^N - \mu_\circ) \right)^m \right] = \sum_{\pi \vdash [m]} \prod_{B \in \pi} \kappa^{\sharp B} \left[ \int_{\mathbb{X}} \varphi (\mu_0^N - \mu_\circ) \right]. $$

With the choice $\mu_\circ = F_0^{N,1}$, we note that $\kappa^1[\int_{\mathbb{X}}(\mu_0^N - \mu_\circ)] = \mathbb{E}[\int_{\mathbb{X}}(\mu_0^N - \mu_\circ)] = 0$, and we may thus restrict the sum to partitions $\pi$ such that $\sharp B > 1$ for all $B \in \pi$. Appealing to (5.3) at $t = 0$ in order to control cumulants back in terms of correlation functions, and taking the supremum over $\varphi$, we deduce
$$
\sup \left\{ \mathbb{E} \left[ \left( \int_X \varphi(\mu_0^N - \mu_\circ) \right)^m \right] : \varphi \in C_c^\infty(X), \|\langle z \rangle^{-3p} \varphi\|_{W^{2,q'}(X)} \le 1 \right\} \\
\lesssim \sum_{\ell=1}^m \sum_{\substack{j_1, \dots, j_\ell \ge 2 \\ j_1 + \dots + j_\ell = m}} \prod_{i=1}^\ell \left( \sum_{1 \le k \le l \le j_i} \sum_{\substack{n_1, \dots, n_k \ge 1 \\ n_1 + \dots + n_k = l}} N^{l-k-j_i+1} \prod_{s=1}^k \|G_\circ^{N,n_s}\|_{W_{3mp,n_s}^{-2,q}} \right).
$$

In the right-hand side, we note that the sum over $j_1, \dots, j_\ell \ge 2$ with $j_1 + \dots + j_\ell = m$ requires $2\ell \le m$, so the first sum can be restricted accordingly. Reorganizing the sums, we are then led to

$$
\sup \left\{ \mathbb{E} \left[ \left( \int_X \varphi(\mu_0^N - \mu_\circ) \right)^m \right] : \varphi \in C_c^\infty(X), \|\langle z \rangle^{-3p} \varphi\|_{W^{2,q'}(X)} \le 1 \right\} \\
\lesssim \sum_{\ell=1}^{m/2} \sum_{\ell \le K \le L \le m} \sum_{\substack{j_1, \dots, j_K \ge 1 \\ j_1 + \dots + j_K = L}} N^{L-K+m+\ell} \prod_{i=1}^{K} \|G_\circ^{N,j_i}\|_{W_{3mp,j_i}^{-2,q}}.
$$

For $m$ even, this means

$$
A_{3p,q}^{N,m}(0) \lesssim \sum_{\ell=1}^{m/2} \sum_{\ell \le K \le L \le m} \sum_{\substack{j_1, \dots, j_K \ge 1 \\ j_1 + \dots + j_K = L}} N^{L-K-m+\ell} \prod_{i=1}^{K} \|G_{\circ}^{N,j_i}\|_{W_{3mp,j_i}^{-2,q}},
$$

and for $m$ odd, after combination with (4.12),

$$
A_{3p,q}^{N,m}(0) \lesssim \sum_{\ell=1}^{m} \sum_{\ell \le K \le L \le 2m} \sum_{\substack{j_1, \dots, j_K \le m+1 \\ j_1 + \dots + j_K = L}} \left( N^{L-K-2m+\ell} \prod_{i=1}^{K} \|G_{\circ}^{N,j_i}\|_{W_{3mp,j_i}^{-2,q}} \right)^{\frac{1}{2}}.
$$

This concludes the proof.

<figure><img src="image_3.png" /></figure>

5.2. **BBGKY hierarchical estimates.** Starting from the BBGKY hierarchy for marginals of $F^N$, and recalling how marginals can be expanded in terms of correlations (and vice versa), we can derive a corresponding BBGKY hierarchy of nonlinear equations for correlations. Before stating it, we start by introducing some useful short-hand notation.

**Definition 5.2.** Consider a collection ${\{h^m\}}_{1\le m\le N}$ of functions $h^m: \mathbb{X}^m \to \mathbb{R}$ such that for all $m$ the function $h^m$ is symmetric in its $m$ entries (such as ${\{F^{N,m}\}}_{1\le m\le N}$ or ${\{G^{N,m}\}}_{1\le m\le N}$).

— For $P \subset [[m]]$ with $P \neq \emptyset$, we define $h^P : \mathbb{X}^m \to \mathbb{R}$ as

$$
h^P(z_{[m]}) := h^{\sharp P}(z_P),
$$

and for $P = \emptyset$ we set $h^\emptyset := 0$.

— For $P \subset [[m]]$ with $P \neq [[N]]$, we define $h^{P \cup \{\ast\}} : \mathbb{X}^{m+1} \to \mathbb{R}$ as

$$
h^{P \cup \{\ast\}}(z_{[m]}, z_{\ast}) := h^{\sharp P+1}(z_P, z_{\ast}),
$$

and for $P = [[N]]$ we set $h^{[N] \cup \{\ast\}} = 0$.

— For $P \subset [[m]]$ and $k, \ell \in P$, we define

$$
S_{k,\ell} h^P := \nabla W(x_k - x_\ell) \cdot \nabla_{v_k} h^P.
$$

— For $P \subset [[m]]$ and $k \in P$, we define

$$
H_k h^{P \cup \{\ast\}}(z_{[m]}) := \int_{\mathbb{X}} \nabla W(x_k - x_\ast) \cdot \nabla_{v_k} h^{P \cup \{\ast\}}(z_{[m]}, z_\ast) \, dx_\ast.
$$
With this notation, we now formulate the BBGKY hierarchy of nonlinear equations satisfied by correlation functions. To the best of our knowledge, such hierarchies for correlations were first written down by Ernst and Cohen [14] in the context of Boltzmann-type systems. We refer to [16, Section 4] for the derivation in the case of the overdamped dynamics (1.18): lengthy but straightforward algebraic manipulations are needed to collect the different factors in the correlation hierarchy. For convenience, we shall write $A - B = A \setminus B$ for set difference, with for instance $A - B - C = A \setminus (B \cup C)$ and $A \cup B - C = (A \cup B) \setminus C$.

**Lemma 5.3** (Correlation hierarchy, see e.g. [14, 16, 23]). For fixed $N$, correlation functions satisfy the following hierarchy of equations: for all $1 \le m \le N$,

$$ \partial_t G^{N,m} = L_{N,m} G^{N,m} + \kappa \frac{N-m}{N} \sum_{k=1}^{m} H_k G^{N,[m]\cup\{*\}} \\ - \kappa \sum_{k=1}^{m} \sum_{A \subset [m] - \{k\}} \frac{m-1-\sharp A}{N} H_k(G^{N,A \cup \{k,*\}} G^{N,[m]-\{k\}} - A) \\ + \kappa \frac{N-m}{N} \sum_{k=1}^{m} \sum_{\substack{A \subset [m] - \{k\} \\ A \ne \emptyset}} H_k(G^{N,A \cup \{k\}} G^{N,[m]\cup\{*\}} - A - \{k\}) \\ - \kappa \sum_{k=1}^{m} \sum_{A \subset [m] - \{k\}} \sum_{B \subset [m] - \{k\} - A} \frac{m-1-\sharp A - \sharp B}{N} H_k(G^{N,A \cup \{k\}} G^{N,B \cup \{*\}} G^{N,[m]-A-B-\{k\}}) \\ + \kappa \frac{N-m}{N} \sum_{k,\ell=1}^{m} S_{k,\ell} G^{N,[m]} + \kappa \frac{N-m}{N} \sum_{k \ne \ell} \sum_{A \subset [m] - \{k,\ell\}} S_{k,\ell} (G^{N,A \cup \{k\}} G^{N,[m]-A-\{k\}}) \\ + \kappa \frac{N-m}{m} \sum_{k=1}^{m} H_k \left( (F^{N,\{*\}} - \mu(z_*)) G^{N,[m]} + (F^{N,\{k\}} - \mu(z_k)) G^{N,[m]\cup\{*\}-\{k\}} \right), $$

where $L_{N,m}$ stands for the following (non-autonomous) $m$-particle linearized mean-field operator,

$$ L_{N,m} := \sum_{k=1}^{m} \mathrm{Id}^{\otimes k-1} \otimes L_{N,m}^{\circ} \otimes \mathrm{Id}^{\otimes m-k}, \qquad (5.5) $$

$$ L_{N,m}^{\circ} h := \frac{1}{2} \mathrm{div}_v ((\nabla_v + \beta v)h) - v \cdot \nabla_x h + \left( \nabla A + \kappa \frac{N-m}{N} \nabla W * \mu \right) \cdot \nabla_v h \\ + \kappa \frac{N-m}{N} \nabla_v \mu \cdot \int_{\mathbb{X}} \nabla W(x - x_*) h(z_*) dz_* $$

Starting from this hierarchy of nonlinear equations for correlations, combined with ergodic estimates for the linearized mean-field evolution, we establish the following a priori estimates. While relating norms with similar weights, a loss of regularity, due to the derivatives in the hierarchy, can be observed here, creating additional constraints on the choice of regularity and weights when combining this result with Proposition 5.1 in the proof of Theorem 1.2.

**Proposition 5.4.** There is $\lambda_0 > 0$ (only depending on $d, \beta, a$, $\|W\|_{W^{1,\infty}(\mathbb{R}^d)}$) such that the following holds: given $1 < q \le 2$ and $0 < p \le \frac{1}{6}$ with $pq' \gg 1$ large enough (only depending on $d, \beta, a$), provided that $0 < \kappa \ll 1$ is small enough (further depending on $p, q$, $Q_2(F_\circ^{N,1})$, and on $\|W\|_{W^{d+3,\infty} \cap H^s(\mathbb{R}^d)}$ for some $s > \frac{d}{2} + 5$), we have for all $2 \le m \le N$, $r \ge 4$, and $t \ge 0$,

$$ \|G_t^{N,m}\|_{W_{p,m}^{-r,q}} \lesssim e^{-\lambda_0 pmt} \|G_{\circ}^{N,m}\|_{W_{p,m}^{-r,q}} + \int_0^t e^{-\lambda_0 pm(t-s)} \left( \|G_s^{N,m+1}\|_{W_{p,m+1}^{1-r,q}} \right. $$
$$ + \sum_{\substack{i_1+i_2=m+1 \\ 1 \le i_1, i_2 \le m-1}} \|G_s^{N,i_1}\|_{W_{p,i_1}^{1-r,q}} \|G_s^{N,i_2}\|_{W_{p,i_2}^{1-r,q}} + N^{-1} \sum_{\substack{i_1+i_2=m \\ 1 \le i_1, i_2 \le m-1}} \|G_s^{N,i_1}\|_{W_{p,i_1}^{1-r,q}} \|G_s^{N,i_2}\|_{W_{p,i_2}^{1-r,q}} \\ + N^{-1} \sum_{\substack{i_1+i_2+i_3=m+1 \\ 1 \le i_1, i_2, i_3 \le m-1}} \|G_s^{N,i_1}\|_{W_{p,i_1}^{1-r,q}} \|G_s^{N,i_2}\|_{W_{p,i_2}^{1-r,q}} \|G_s^{N,i_3}\|_{W_{p,i_3}^{1-r,q}} \\ + (N^{-1} + e^{-\lambda_0 ps} \|G_\circ^{N,2}\|_{W_{3p,2}^{-2,q}}) \|G_s^{N,m}\|_{W_{p,m}^{1-r,q}} \right) ds, \quad (5.6) $$

where the multiplicative constant further depends on $m, r$, $Q_{mp}(F_\circ^{N,1})$, and on $\|W\|_{W^{d+r+1,\infty}\cap H^s(\mathbb{R}^d)}$ for some $s > 2r + \frac{d}{2} - 1$.

**Proof.** We start by recalling the available ergodic estimates for the linearized mean-field operator $L_{N,m}$. Given $m \ge 1$ and $u_\circ^m \in C_c^\infty(X^m)$ with $\int_X u_\circ^m(z_{[m]}) dz_j = 0$ for all $1 \le j \le m$, let us consider the solution $u^m$ of

$$ \left\{ \begin{array}{ll} \partial_t u^m = L_{N,m} u^m, & t \ge 0, \\ u^m|_{t=0} = u_\circ^m, & \end{array} \right. $$

where in the definition (5.5) of $L_{N,m}$ we recall that we have chosen $\mu$ as the solution of the mean-field equation (1.4) with initial condition $\mu_\circ = F_\circ^{N,1}$. By definition of $L_{N,m}$ as a Kronecker sum, we can write

$$ u_t^m = (V_{0,t}^{N,m})^{\otimes m} u_\circ^m, $$

where $\{V_{s,t}^{N,m}\}_{0 \le s \le t}$ stands for the one-particle evolution defined by

$$ \left\{ \begin{array}{ll} \partial_t V_{s,t}^{N,m} h = L_{N,m}^\circ V_{s,t}^{N,m} h, & t \ge s, \\ V_{s,t}^{N,m} h|_{t=s} = h. & \end{array} \right. $$

In each component, we may then appeal to the ergodic estimate derived in [3, Theorem 2.13]: there exist $\kappa_0, \lambda_0 > 0$ (only depending on $d, \beta, a, \|W\|_{W^{1,\infty}(\mathbb{R}^d)}$) such that, given $\kappa \in [0, \kappa_0]$, $1 < q \le 2$, and $0 < p \le 1$ with $pq' \gg 1$ (only depending on $d, \beta, a$), we have for all $r \ge 2$ and $t \ge 0$,

$$ \|u_t^m\|_{W_{p,m}^{-r,q}} \lesssim e^{-\lambda_0 pm t} \|u_\circ^m\|_{W_{p,m}^{-r,q}}, $$

where the multiplicative constant only depends on $d, \beta, a, p, q, r, m$, $Q_2(F_\circ^{N,1})$, $\|W\|_{W^{d+r+1,\infty}(\mathbb{R}^d)}$. Testing the definition of $S_{k,l}h^P$ with a tensorized test function, decomposing $W$ as a superposition of Fourier modes and using polarization similarly as in (4.8) to rewrite the resulting expression as a superposition of integrals of $h^P$ tested with a tensorized test function (recall the definition (1.13) of functional spaces), we are led to

$$ \|S_{k,l}h^P\|_{W_{p,\#P}^{-r,q}} \lesssim \|h^P\|_{W_{p,\#P}^{1-r,q}}, $$

and similarly,

$$ \|H_k h^{P \cup \{\ast\}}\|_{W_{p,\#P}^{-r,q}} \lesssim \|h^{P \cup \{\ast\}}\|_{W_{p,\#P+1}^{1-r,q}}, $$

where the multiplicative constants only depend on $d, r, \#P$, and $\|W\|_{H^s(\mathbb{R}^d)}$ for some $s > 2r + \frac{d}{2} - 1$. Now applying this to the equation for $G^{N,m}$ in Lemma 5.3, we get for all $2 \le m \le N$, $r \ge 2$, and $t \ge 0$,

$$ \|G_t^{N,m}\|_{W_{p,m}^{-r,q}} \lesssim e^{-\lambda_0 pm t} \|G_\circ^{N,m}\|_{W_{p,m}^{-r,q}} \\ + \int_0^t e^{-\lambda_0 pm(t-s)} \left( \|G_s^{N,m+1}\|_{W_{p,m+1}^{1-r,q}} + N^{-1} \sum_{\substack{i_1+i_2=m+1 \\ 1 \le i_1, i_2 \le m}} \|G_s^{N,i_1}\|_{W_{p,i_1}^{1-r,q}} \|G_s^{N,i_2}\|_{W_{p,i_2}^{1-r,q}} \right. \\ \left. + \sum_{\substack{i_1+i_2=m+1 \\ 1 \le i_1, i_2 \le m-1}} \|G_s^{N,i_1}\|_{W_{p,i_1}^{1-r,q}} \|G_s^{N,i_2}\|_{W_{p,i_2}^{1-r,q}} + N^{-1} \sum_{\substack{i_1+i_2=m \\ 1 \le i_1, i_2 \le m-1}} \|G_s^{N,i_1}\|_{W_{p,i_1}^{1-r,q}} \|G_s^{N,i_2}\|_{W_{p,i_2}^{1-r,q}} \right) ds, $$
$$ + N^{-1} \sum_{\substack{i_1+i_2+i_3=m+1 \\ 1 \le i_1, i_2, i_3 \le m-1}} \|G_s^{N, i_1}\|_{W_{p, i_1}^{1-r, q}} \|G_s^{N, i_2}\|_{W_{p, i_2}^{1-r, q}} \|G_s^{N, i_3}\|_{W_{p, i_3}^{1-r, q}} \\ + (N^{-1} + \|F_s^{N, 1} - \mu_s\|_{W_{p, 1}^{1-r, q}}) \|G_s^{N, m}\|_{W_{p, m}^{1-r, q}} \) ds. $$

Note that in the first sum in the right-hand side the terms with $(i_1, i_2) = (m, 1)$ or $(1, m)$ can be bounded by $N^{-1}\|G_s^{N,m}\|_{W_{p,m}^{1-r,q}}$, which already appears in the last term. Hence, this first sum can be restricted to $1 \le i_1, i_2 \le m-1$, which is then bounded by the second sum. For the last right-hand side term, involving $F^{N,1}-\mu$, we appeal to Theorem 1.1 with initial condition $\mu_0 = F_0^{N,1}$: provided that $\lambda_0$ is chosen small enough (only depending on $d, \beta, a$, $\|W\|_{W^{1,\infty}(\mathbb{R}^d)}$), for all $1 < q \le 2$ and $0 < p \le \frac{1}{6}$ with $pq' \gg 1$ large enough (only depending on $d, \beta, a$), provided that $\kappa \ll 1$ is small enough (further depending on $p, q$, $Q_2(F_0^{N,1})$, $\|W\|_{W^{d+3,\infty}\cap H^s(\mathbb{R}^d)}$ for some $s > \frac{d}{2}+5$), we have

$$ \|F_t^{N,1} - \mu_t\|_{W_{p,1}^{-3,q}} \lesssim N^{-1} + e^{-2\lambda_0 pt} \|G_0^{N,2}\|_{W_{3p,2}^{-2,q}}, $$

which then yields the conclusion.

<figure><img src="image_4.png" /></figure>

5.3. **Proof of Theorem 1.2.** The claim (1.15) already follows from (4.1) in the proof of Theorem 1.1 with the choice $\mu_0 = F_0^{N,1}$. It remains to prove (1.17). Let $1 < q \le 2$ and $0 < p \le \frac{1}{6}$ with $pq' \gg 1$ large enough (only depending on $d, \beta, a$), let $2 \le m \le N$ and $\alpha \in [0, 1]$ be fixed, and assume that initially

$$ \|G_0^{N,j}\|_{W_{3p,j}^{-2,q}} \le C_j^\alpha N^{-\alpha(j-1)}, \quad \text{for all } 1 \le j \le 2m-1. \qquad (5.7) $$

We split the proof into two main steps.

**Step 1:** Suboptimal estimates.

Let $\lambda_0, \kappa$ be as in Proposition 3.1. As already emphasized, the suboptimal hierarchical estimates on correlations in Proposition 5.1 involve a loss in the weights. To accommodate this in the iterated use of the hierarchy, we choose the initial weight carefully, setting $p_0 = p/n!$ for all $2 \le n \le 2m-1$. A straightforward computation then shows that, under the initial assumption (5.7), Proposition 5.1 yields for all $r \ge 2$, and $t \ge 0$,

$$ \|G_t^{N,n}\|_{W_{p/n!,n}^{-r,q}} \lesssim N^{-\frac{n}{2}} + e^{-\frac{\lambda_0 p}{n!}t} N^{-\alpha \frac{n}{2}} + \sum_{1 \le k \le l \le n-1} \sum_{j_1, \dots, j_k \ge 1} N^{l-k-n+1} \prod_{i=1}^{k} \|G_t^{N,j_i}\|_{W_{p/(n-1)!,j_i}^{-r,q}}. $$

By a direct iteration of these estimates, inductively eliminating correlations in the right-hand side, and using $\|G_t^{N,1}\|_{W_{p,1}^{-r,q}} \lesssim 1$, we deduce for all $2 \le n \le 2m-1$, $r \ge 2$, and $t \ge 0$,

$$ \|G_t^{N,n}\|_{W_{p/n!,n}^{-r,q}} \lesssim N^{-\frac{n}{2}} + e^{-\frac{\lambda_0 p}{n!}t} N^{-\alpha \frac{n}{2}}. \qquad (5.8) $$

**Step 2:** Conclusion via BBGKY analysis.

We leverage the suboptimal correlation estimates of Step 1, using the BBGKY estimates of Proposition 5.4 and arguing by induction. As the desired estimates (1.17) are already known to hold for two-particle correlations, cf. (1.15), we can assume that there is some $3 \le n \le m$ such that the following holds: provided that $\kappa$ is small enough (also depending on $n$), we have for all $t \ge 0$,

$$ \|G_t^{N,j}\|_{W_{p_j,j}^{-r_j,q_j}} \lesssim N^{1-j} + e^{-\lambda_0 p_j t} N^{\alpha(1-j)}, \quad \text{for all } 2 \le j \le n-1, \qquad (5.9) $$

where we have set

$$ r_j := j^2, \quad p_j := \frac{p}{(2j-2)!}, \quad q_j := \frac{(2j-2)!q'}{(2j-2)!q'^{-1}}. $$

This set of parameters may appear somewhat intricate at first sight. We will explain the constraints leading to the choice of $p_j$ and $r_j$ throughout the argument. The exponent $q_j$ is chosen so that $q'_j = (2j-2)!q'$, ensuring that $p_j q'_j = pq' \gg 1$ and thereby allowing us to apply the ergodic estimates from Lemma 2.1. Under this induction assumption (5.9), we shall show that the same estimate is also
automatically valid for $j = n$ up to further restricting the range of $\kappa$, which then concludes the proof. We split the argument into two further substeps.

**Substep 2.1:** Proof that for $\kappa$ small enough (also depending on $n$) we have for all $0 \le k \le n - 2$, $r \ge r_{n-1} + 1$, and $t \ge 0$,

$$ \|G_t^{N,n+k}\|_{W_{p_n,n+k}^{-r,q_n}} \lesssim N^{1-n-k} + e^{-\lambda_0 p_n t} N^{\alpha(1-n-k)} \\ + \int_0^t e^{-\lambda_0 p_n(n+k)(t-s)} \left( \|G_s^{N,n+k+1}\|_{W_{p_n,n+k+1}^{1-r,q_n}} + (N^{-1} + e^{-\lambda_0 p_n s} N^{-\alpha}) \|G_s^{N,n+k}\|_{W_{p_n,n+k}^{1-r,q_n}} \right. \\ \left. + \sum_{i=2}^{k+1} (N^{1-i} + e^{-\lambda_0 p_n s} N^{\alpha(1-i)}) \|G_s^{N,n+k+1-i}\|_{W_{p_n,n+k+1-i}^{1-r,q_n}} \right) ds. \quad (5.10) $$

We apply Proposition 5.4 with $m = n + k$ and with exponents $p_n, q_n$ satisfying $p_n q'_n = pq' \gg 1$. Provided that $\kappa$ is small enough (depending on $n$), we get from Proposition 5.4, for all $0 \le k \le n - 2$, $r \ge r_{n-1} + 1$, and $t \ge 0$,

$$ \|G_t^{N,n+k}\|_{W_{p_n,n+k}^{-r,q_n}} \lesssim e^{-\lambda_0 p_n t} \|G_{\circ}^{N,n+k}\|_{W_{p_n,n+k}^{-r,q_n}} + \int_0^t e^{-\lambda_0 p_n(n+k)(t-s)} \left( \|G_s^{N,n+k+1}\|_{W_{p_n,n+k+1}^{1-r,q_n}} \right. \\ + \sum_{\substack{i_1+i_2=n+k+1 \\ 1 \le i_1, i_2 \le n+k-1}} \|G_s^{N,i_1}\|_{W_{p_n,i_1}^{1-r,q_n}} \|G_s^{N,i_2}\|_{W_{p_n,i_2}^{1-r,q_n}} + N^{-1} \sum_{\substack{i_1+i_2=n+k \\ 1 \le i_1, i_2 \le n+k-1}} \|G_s^{N,i_1}\|_{W_{p_n,i_1}^{1-r,q_n}} \|G_s^{N,i_2}\|_{W_{p_n,i_2}^{1-r,q_n}} \\ + N^{-1} \sum_{\substack{i_1+i_2+i_3=n+k+1 \\ 1 \le i_1, i_2, i_3 \le n+k-1}} \|G_s^{N,i_1}\|_{W_{p_n,i_1}^{1-r,q_n}} \|G_s^{N,i_2}\|_{W_{p_n,i_2}^{1-r,q_n}} \|G_s^{N,i_3}\|_{W_{p_n,i_3}^{1-r,q_n}} \\ \left. + (N^{-1} + e^{-\lambda_0 p_n s} \|G_{\circ}^{N,2}\|_{W_{3p_n,2}^{-2,q_n}}) \|G_s^{N,n+k}\|_{W_{p_n,n+k}^{1-r,q_n}} \right) ds. $$

Note that the condition $r \ge r_{n-1} + 1$ ensures $r - 1 \ge r_j$ for all $1 \le j \le n - 1$. Using $p_j q'_j = pq' > 2d$, by Jensen's inequality, we get

$$ \|G^{N,j}\|_{W_{p_n,j}^{1-r,q_n}} \lesssim \|G^{N,j}\|_{W_{p_j,j}^{-r_j,q_j}}. $$

Using the initial assumption (5.7), as well as the induction assumption (5.9) for $G^{N,j}$ with $j \le n - 1$, the claim (5.10) follows from the above estimate after straightforward simplifications.

**Substep 2.2:** Proof that for $\kappa$ small enough (depending on $n$) we have for all $0 \le \ell \le n - 2$ and $t \ge 0$,

$$ \|G_t^{N,n+k}\|_{W_{p_n,n+k}^{n-2-r_n-\ell,q_n}} \lesssim N^{-\frac{n+k+\ell}{2}} + e^{-\lambda_0 p_n t} N^{-\alpha \frac{n+k+\ell}{2}}, \quad \text{for all } 0 \le k \le n-2-\ell. \quad (5.11) $$

This will conclude the proof: indeed, choosing $\ell = n - 2$ and $k = 0$, the above estimate precisely yields (5.9) with $j = n$ as desired. We turn to the proof of (5.11) and argue by induction on $\ell$. First, for $\ell = 0$, the claim (5.11) follows from (5.8). Our choice for $p_n$ is easily understood here, as the use of (5.8) requires $p_n \le p/(n+k)!$ for all $0 \le k \le n - 2$, which is guaranteed for $p_n = p/(2n - 2)!$. Next, we can assume that there is some $1 \le L \le n - 2$ such that (5.11) holds with $\ell = L - 1$, that is,

$$ \|G_t^{N,n+k}\|_{W_{p_n,n+k}^{n-1-r_n-L,q_n}} \lesssim N^{-\frac{n+k+L-1}{2}} + e^{-\lambda_0 p_n t} N^{-\alpha \frac{n+k+L-1}{2}}, \quad \text{for all } 0 \le k \le n-1-L, \quad (5.12) $$

and it remains to show that (5.11) then automatically also holds for $\ell = L$. This follows directly from the result (5.10) in Step 2.1 combined with the induction assumption (5.12).

To apply (5.10), we require $r_n + 2 - n \ge r_{n-1} + 1$, which holds for $r_n = n^2$, thus explaining our choice for this parameter. $\square$
ACKNOWLEDGEMENTS

The authors thank Arnaud Guillin for pointing out the validity of the creation-of-chaos estimate (1.6) in Wasserstein distance as a consequence of [4]. AB acknowledges financial support from the AAP “Accueil” from Université Lyon 1 Claude Bernard. MD and MM acknowledge financial support from the European Union (ERC, PASTIS, Grant Agreement n°101075879).<fnref n="2" />

REFERENCES

[1] A. Béjar-López, A. Blaustein, P.-E. Jabin, and J. Soler. Longtime and chaotic dynamics in microscopic systems with singular interactions. *Preprint arXiv:2411.08614*, 2024.

[2] A. Bernou and M. Duerinckx. Mean-field approximation, Gibbs relaxation, and cross estimates. Preprint, arXiv:2507.23123.

[3] A. Bernou and M. Duerinckx. Uniform-in-time estimates on corrections to mean field for interacting Brownian particles. *Probability Theory and Related Fields*, May 2025.

[4] F. Bolley, A. Guillin, and F. Malrieu. Trend to equilibrium and particle approximation for a weakly selfconsistent Vlasov-Fokker-Planck equation. *M2AN Math. Model. Numer. Anal.*, 44(5):867–884, 2010.

[5] D. Bresch, M. Duerinckx, and P.-E. Jabin. A duality method for mean-field limits with singular interactions. Preprint, arXiv:2402.04695.

[6] D. Bresch, P.-E. Jabin, and J. Soler. A new approach to the mean-field limit of Vlasov-Fokker-Planck equations. *Anal. PDE*, 18(4):1037–1064, 2025.

[7] L.-P. Chaintron and A. Diez. Propagation of chaos: A review of models, methods and applications. I. Models and methods. *Kinetic and Related Models*, 15(6):895, 2022.

[8] F. Cornalba, J. Fischer, J. Ingmanns, and C. Raithel. Density fluctuations in weakly interacting particle systems via the Dean-Kawasaki equation. Preprint, arXiv:2303.00429.

[9] P. Del Moral and J. Tugaut. Uniform Propagation of Chaos and Creation of Chaos for a Class of Nonlinear Diffusions. *Stoch. Anal. Appl.*, 37(6):909–935, 2019.

[10] F. Delarue and A. Tse. Uniform in time weak propagation of chaos on the torus. *Annales de l'Institut Henri Poincaré, Probabilités et Statistiques*, 61(2):1021–1074, 2025.

[11] M. Duerinckx. On the Size of Chaos via Glauber Calculus in the Classical Mean-Field Dynamics. *Commun. Math. Phys.*, 382:613–653, 2021.

[12] M. Duerinckx and P.-E. Jabin. Correlation estimates for Brownian particles with singular interactions. Preprint, arXiv:2510.01507.

[13] Mitia Duerinckx and Laure Saint-Raymond. Lenard-Balescu correction to mean-field theory. *Probab. Math. Phys.*, 2(1):27–69, 2021.

[14] M. H. Ernst and E. G. D. Cohen. Nonequilibrium fluctuations in mu space. *J. Stat. Phys.*, 25(1):153–180, 1981.

[15] A. Guillin, P. Le Bris, and P. Monmarché. Convergence rates for the Vlasov-Fokker-Planck equation and uniform in time propagation of chaos in non convex cases. *Electronic Journal of Probability*, 27(none), January 2022.

[16] E. Hess-Childs and K. Rowan. Higher-order propagation of chaos in $L^2$ for interacting diffusions. *Probab. Math. Phys.*, 6(2):581–646, 2025.

[17] D. Lacker. Hierarchies, entropy, and quantitative propagation of chaos for mean field diffusions. *Probab. Math. Phys.*, 4(2):377–432, 2023.

[18] D. Lacker and L. Le Flem. Sharp uniform-in-time propagation of chaos. *Probab. Theory Related Fields*, 187:443–480, 2023.

[19] J. Lukkarinen and A. Vuoksenmaa. Generation of chaos in the cumulant hierarchy of the stochastic Kac model. Preprint, arXiv:2407.17068.

[20] S. Mischler and C. Mouhot. Exponential stability of slowly decaying solutions to the kinetic-Fokker-Planck equation. *Arch. Ration. Mech. Anal.*, 221(2):677–723, 2016.

[21] M. Rosenzweig and S. Serfaty. Modulated logarithmic Sobolev inequalities and generation of chaos. *Ann. Fac. Sci. Toulouse, Math. (6)*, 34(1):107–134, 2025.

[22] C. Villani. Hypocoercivity. *Mem. Amer. Math. Soc.*, 202(950):iv+141, 2009.

[23] P. Xie. Uniform-in-time estimates on the size of chaos for interacting particle systems. *SIAM J. Math. Anal.*, 57(5):5358–5381, 2025.

<footnote n="2">²Views and opinions expressed are however those of the authors only and do not necessarily reflect those of the European Union or the European Research Council Executive Agency. Neither the European Union nor the granting authority can be held responsible for them.</footnote>
(Armand Bernou) UNIVERSITÉ DE LYON, UNIVERSITÉ CLAUDE BERNARD LYON 1, LABORATOIRE DE SCIENCES ACTUARIELLE ET FINANCIÈRE, 50 AVENUE TONY GARNIER, F-69007 LYON, FRANCE

Email address: armand.bernou@univ-lyon1.fr

(Mitia Duerinckx) UNIVERSITÉ LIBRE DE BRUXELLES, DÉPARTEMENT DE MATHÉMATIQUE, 1050 BRUSSELS, BELGIUM

Email address: mitia.duerinckx@ulb.be

(Matthieu Ménard) UNIVERSITÉ LIBRE DE BRUXELLES, DÉPARTEMENT DE MATHÉMATIQUE, 1050 BRUSSELS, BELGIUM

Email address: matthieu.menard@ulb.be