Global Stability of Minkowski Space for the Einstein–Vlasov System in the Harmonic Gauge

Minkowski space is shown to be globally stable as a solution to the massive Einstein–Vlasov system. The proof is based on a harmonic gauge in which the equations reduce to a system of quasilinear wave equations for the metric, satisfying the weak null condition, coupled to a transport equation for the Vlasov particle distribution function. Central to the proof is a collection of vector fields used to control the particle distribution function, a function of both spacetime and momentum variables. The vector fields are derived using a general procedure, are adapted to the geometry of the solution and reduce to the generators of the symmetries of Minkowski space when restricted to acting on spacetime functions. Moreover, when specialising to the case of vacuum, the proof provides a simplification of previous stability works.


The Einstein-Vlasov System
The Einstein-Vlasov system provides a statistical description of a collection of collisionless particles, interacting via gravity as described by Einstein's general theory of relativity. A fundamental problem is to understand the long time dynamics of solutions of this system. The problem is a great challenge even in the absence of particles, and global works on the vacuum Einstein equations all either involve simplifying symmetry assumptions or solutions arising from small data. In the asymptotically flat setting, small data solutions of the vacuum Einstein equations were first shown to exist globally and disperse to Minkowski space in the monumental work of Christodoulou-Klainerman [14]. An alternative proof of the stability of Minkowski space was later given by Lindblad-Rodnianski [37] in which a global harmonic coordinate system was constructed, described below. For the Einstein-Vlasov system, the global properties of small data solutions were first understood when the initial data were assumed to be spherically symmetric, an assumption under which the equations simplify dramatically, by Rein-Rendall [41] in the massive case, when all particles are assumed to have mass 1, and by Dafermos [15] in the massless case, when all particles are assumed to have mass 0. See also work of Andréasson-Kunze-Rein [4] for a global existence result in spherical symmetry for a class of large "outgoing" data. Without the simplifying assumption of spherical symmetry, small data solutions of the massless Einstein-Vlasov system were later understood by Taylor [46].
In this work the problem of the stability of Minkowski space for the Einstein-Vlasov system, without any symmetry assumptions, is addressed in the case that all particles have mass 1 (and can easily be adapted to the case that all particles have any fixed mass m > 0). The system takes the form 2) where the unknown is a 3 + 1 dimensional manifold M with Lorentzian metric g, together with a particle distribution function f : P → [0, ∞), where the mass shell P is defined by Here Ric(g) and R(g) denote the Ricci and scalar curvature of g respectively. A coordinate system (t, x) for M, with t a time function (that is the one form dt is timelike with respect to the metric g), defines a coordinate system (t, x, p 0 , p) for the tangent bundle T M of M conjugate to (t, x), where (t, x i , p 0 , p i ) denotes the point The mass shell relation in the definition of P, should be viewed as defining p 0 as a function of (t, x 1 , x 2 , x 3 , p 1 , p 2 , p 3 ). Here Greek indices run over 0, 1, 2, 3, lower case Latin indices run over 1, 2, 3, and often the notation t = x 0 is used. The vector X is the generator of the geodesic flow of M which, with respect to the coordinate system (t, x, p) for P, takes the form The volume form ( √ | det g|/ p 0 ) d p 1 d p 2 d p 3 in (1.2) is the induced volume form of the spacelike hypersurface P (t,x) ⊂ T (t,x) M when the tangent space T (t,x) M is endowed with the metric g μν (t, x) d p μ d p ν induced by g on M.
Einstein equations, (see also the recent textbook of Ringström [42]) guarantees that, for any initial data set as above, there exists a globally hyperbolic solution (M, g, f ) of the system (1.1)-(1.3) which attains the initial data, in the sense that there exists an imbedding ι : → M under the pullback of which the induced first and second fundamental form of g are g and k respectively, and the restriction of f to the mass shell over ι( ) is given by f 0 . As in [37], the proof is based on the harmonic gauge (or wave gauge), that is a system of coordinates x μ satisfying g x μ = 0, (1.6) where g denotes the geometric wave operator of the metric g.
The initial data set where e denotes the Euclidean metric on R 3 , satisfies the constraint equations and gives rise to Minkowski space, the trivial solution of the system (1.1)-(1.3). The main result of this work concerns solutions arising from initial data sufficiently close to the trivial initial data set (1.7). The initial data will be assumed to be asymptotically flat in the sense that is diffeomorphic to R 3 and there exists a global coordinate chart (x 1 , x 2 , x 3 ) of and M ≥ 0 and 0 < γ < 1 such that as r = |x| → ∞. For given γ > 0 and such an initial data set, define, for any N ≥ 0, the initial energy (1 + r ) 1/2+γ +|I | ∇∇ I h 1 L 2 ( ) + (1 + r ) 1/2+γ +|I | ∇ I k L 2 ( ) , (1.9) where ∇ = (∂ x 1 , ∂ x 2 , ∂ x 3 ) denotes the coordinate gradient, and the initial norms for f 0 : (1.10) By the Sobolev inequality there exists a constant C such that D N −4 ≤ CV N for any N ≥ 0.
The main result of this work is the following: Theorem 1.1. Let ( , g, k, f 0 ) be an initial data set for the Einstein-Vlasov system (1.1)-(1.3), asymptotically flat in the above sense, such that f 0 is compactly supported. For any 0 < γ < 1 and N ≥ 11, there exists ε 0 > 0 (depending on the size of supp( f 0 )) such that, for all ε ≤ ε 0 and all initial data satisfying there exists a unique future geodesically complete solution of (1.1)- (1.3), attaining the given data, together with a global system of harmonic coordinates (t, x 1 , x 2 , x 3 ), relative to which the solution asymptotically decays to Minkowski space.
The analogue of Theorem 1.1 for the massless Einstein-Vlasov system (which is the system (1.1)-(1.3) but with the mass shell relation (1.4) replaced with the relation g μν p μ p ν = 0, so that all particles travel through spacetime along null geodesics, as opposed to unit timelike geodesics in the massive case considered here) was resolved by Taylor [46]. It should be noted that the difficulties encountered in [46] are very different to the difficulties encountered in the present work. The work [46] involves a double null gauge (in contrast to the harmonic gauge employed here) and, as is shown as part of the bootstrap argument in the proof, f , being initially supported in a spatially compact set, is supported only in the wave zone |x| ∼ t in a region of finite retarded length. Due to the slow decay of the metric in the wave zone, the work [46] relies crucially in exploiting a null structure 1 present in the Vlasov equation. 2 In contrast, the main difficulties in the present work arise in the interior region |x| < t. The "Minkowski vector fields" (see Section 1.5.1 below), which one would use to control the decoupled Vlasov equation on a fixed Minkowski background, are insufficient by themselves for the proof of Theorem 1.1 and have to be further adapted to the geometry of the spacetimes considered. See Section 1.5 below for a further discussion of the vector fields used in the proof.
The assumption in Theorem 1.1 that f 0 is compactly supported is made for simplicity and, we believe, our method can be extended to relax this assumption.
The assumption implies that the solution f is supported, at late times, away from the wave zone |x| ∼ t. The issue of the slow decay of the metric g in the wave zone is therefore not relevant in the proof of Theorem 1.1 when controlling the solution of the Vlasov equation f , and its derivatives. Were f 0 not compactly supported, this slow decay of g would be relevant and it would be necessary to exploit a form of null structure present in the Vlasov equation, as in the massless case discussed above, together with the weak null structure of Einstein's equations in wave coordinates.
A similar stability result to Theorem 1.1 was shown independently by Fajman-Joudioux-Smulevici [19]. We remark on some of the differences between the two proofs in Section 1.5 after presenting our argument.
There have been a number of related stability works on the Einstein-Vlasov system. In addition to those discussed earlier, there has been related work by Ringström [42] on the Einstein-Vlasov system in the presence of a positive cosmological constant, where the analogue of the Minkowski solution is the de Sitter spacetime. See also [5]. A stability result for a class of cosmological spacetimes was shown in 2 + 1 dimensions with vanishing cosmological constant by Fajman [16], and later extended to 3 + 1 dimensions by Andersson-Fajman [2]. See also the recent work of Moschidis [39] on the instability of the anti-de Sitter spacetime for a related spherically symmetric model in the presence of a negative cosmological constant. A much more comprehensive overview of work on the Einstein-Vlasov system can be found in the review paper of Andréasson [3].
There has also been work on the problem of the stability of Minkowski space for the Einstein equations coupled to various other matter models [7,20,22,[28][29][30]43,48].
Small data solutions of the Vlasov-Poisson system, the non-relativistic analogue of the Einstein-Vlasov system, were studied in [6,23], and those of the Vlasov-Maxwell system in [21]. We note in particular works of Smulevici [44] and Fajman-Joudioux-Smulevici [18] (following [17]) on the asymptotic properties of small data solutions of the Vlasov-Poisson and the related Vlasov-Nordström systems respectively, where issues related to those described in Section 1.5.2 arise and are resolved using an alternative approach. Moreover, there are global existence results for general data for Vlasov-Poisson [38,40] and Vlasov-Nordström [8].
Theorem 1.1 follows as a corollary of the following theorem: Theorem 1.2. For any 0 < γ < 1, K, K > 0 and N ≥ 11, there exists ε 0 > 0 such that, for any data (g μν , ∂ t g μν , f )| t=0 for the reduced Einstein-Vlasov system (1.2), (1.3), (1.11) which satisfy the smallness condition for any ε ≤ ε 0 and the wave coordinate condition (1.12), and such that for there exists a global solution attaining the data such that (E N (t)) 1
Note that the wave coordinate condition (1.12), if satisfied initially, is propagated in time by the reduced Einstein-Vlasov system (1.2), (1.3), (1.11) (this fact is standard; see, for example, Section 4 of [36], which requires only minor modifications for the presence of matter).
One can show that, with this choice, (E N (0)) 1/2 E N , where E N given by (1.9) is the norm of geometric data, and moreover (g μν , ∂ t g μν )| t=0 satisfy the wave coordinate condition (1.12), see, for example, [36,37]. It is therefore clear that Theorem 1.1 follows from Theorem 1.2 (the future causal geodesic completeness can be shown as in [36]) and so the goal of the paper is to establish the proof of Theorem 1.2.

Estimates for the Vlasov Matter
In what follows it is convenient, instead of parameterising the mass shell P by (t, x, p), to instead parameterise it by (t, x, p ), where p i = p i / p 0 for i = 1, 2, 3. Note that, by the mass shell relation (1.4), in Minkowski space ( p 0 ) 2 = 1 + ( p 1 ) 2 + ( p 2 ) 2 + ( p 3 ) 2 and, under a mild smallness condition on g − m, | p | < 1. Abusing notation slightly, we will write f (t, x, p ) for the solution of the Vlasov equation (1.3).
Let { t } denote the level hypersurfaces of the time coordinate t, and let X (s, t, x, p ) i , P (s, t, x, p ) i denote solutions of the geodesic equations (see Section 2) (1.17) normalised so that (s, X (s, t, x, p )) ∈ s , with Here where α βγ are the Christoffel symbols of the metric g with respect to a given coordinate chart (t, The notation X (s), P (s) will sometimes be used for X (s, t, x, p ), P (s, t, x, p ) when it is clear from the context which point (t, x, p ) is meant, and the notation X (s) = (s, X (s)) will sometimes be used.
It follows that the Vlasov equation (1.3) can be rewritten as 18) for all s. The notation (y, q) will be used to denote points in the mass shell over the initial hypersurface, P| t=0 . In Theorem 1.2 it is assumed that f 0 has compact support; |y| ≤ K and |q| ≤ K for (y, q) ∈ supp( f 0 ), for some constants K , K .
Under relatively mild smallness assumptions on h = g − m, see Proposition 2.1, it follows that there exists c < 1, depending only on K , such that solutions of the Vlasov equation satisfy The main new difficulties in the proof of Theorem 1.2, arising from the coupling to the Vlasov equation, are resolved in the following theorem, which is appealed to in the proof of Theorem 1.2: Theorem 1.3. For a given t ≥ 0 and N ≥ 1, suppose that g is a Lorentzian metric such that the Christoffel symbols of g with respect to a global coordinate system for all t ∈ [0, t]. Then there exists ε 0 > 0 such that, for ε < ε 0 and any solution f of the Vlasov equation (1.3) satisfying and metric satisfying (1.20) and |g − m| ≤ ε, the components of the energy momentum tensor T μν (t, x) satisfy where k = |I | and k = k/2 + 1, and Here the constants D k depend only on C N , K, K and c, and ε 0 depends only on c.
In the proof of Theorem 1.2, the L 2 estimates of Theorem 1.3 are used to prove energy estimates for the metric g; see the discussion in Section 1.4.2. The L 1 estimates of Theorem 1.3 are used to recover the pointwise decay (1.20) for the lower order derivatives of the Christoffel symbols ; see the discussion in Section 1.4.7 below. Remark 1.4. The proof of Theorem 1.3 still applies when a ≥ 1, though the theorem is only used in the proof of Theorem 1.2 for some fixed 1 2 < a < 1. The case of a = 1 is omitted in order to avoid logarithmic factors. The proof of an appropriate result when a > 1 is much simpler, although, when Theorem 1.3 is used in the proof of Theorem 1.2, one could not hope for the assumptions (1.20) to hold with a > 1, see Section 1.5.2. Remark 1.5. In Section 5 a better L 2 estimate in terms of t behaviour, compared with the L 2 estimate of Theorem 1.3, is shown to hold for Z I T μν , which involves one extra derivative of ; see Proposition 5.8. It is important however to use the L 2 estimate which does not lose a derivative in the proof of Theorem 1.2.

Overview of the Global Existence Theorem for the Reduced Einstein Equations
It should be noted from the outset that the reduced Einstein-Vlasov system (1.2), (1.3), (1.11) is a system of quasilinear wave equations coupled to a transport equation. It is well known that the general quasilinear wave equation does not necessarily admit global solutions for all small data [24,25]. The null condition, an algebraic condition on the nonlinearity of such equations, was introduced by Klainerman [26], and used independently by Klainerman [27] and Christodoulou [13], as a sufficient condition that small data solutions exist globally in time and are asymptotically free. However, as was noticed by Lindblad [32,33] and Alinhac [1], there are quasilinear equations that do not satisfy the null condition but still admit global solutions for all sufficiently small data. In fact, the classical null condition fails to be satisfied by the vacuum Einstein equations in the harmonic gauge ((1.11) with T ≡ 0), though it was noticed by Lindblad-Rodnianski [35] that they satisfy a weak null condition, which they used to prove a small data global existence theorem [36,37].
The proof of Theorem 1.2 follows the strategy adopted in [37]. The new difficulties, of course, arise from the coupling to the Vlasov equation. A fundamental feature of the problem arises from the fact that, whilst the slowest decay of solutions to wave equations occurs in the wave zone, where t ∼ r , the slowest decay of solutions of the massive Vlasov equation occurs in the interior region t > r . The most direct way to exploit this fact is to impose that f 0 has compact support, in which case, as will be shown in Proposition 2.1, the support condition (1. 19) holds and f actually vanishes in the wave zone at late times.
Since they have been described at length elsewhere, the difficulties associated to the failure of (1.11) to satisfy the classical null condition are only briefly discussed here. Suffice it to say that there is a rich structure in the equations (1.11) which is exploited heavily (see the further discussion in the introductions to [36,37]). The main new features of this work are contained in the proof of Theorem 1.3. Indeed, for a given inhomogeneous term T in (1.11) which satisfies the support conditions and estimates of Theorem 1.3, the small data global existence theorem of [37] mostly goes through unchanged. An outline is given below, including a discussion of some observations which lead to simplifications compared with the proof in [37] (most notably the fact that the structure of the equations are better preserved under commuation with modified Lie derivatives, described in Section 1.4.4 below, but also the use of a simpler L ∞ -L ∞ estimate, described in Section 1.4.6).
The proof of Theorem 1.2 is based on a continuity argument. One assumes that the bounds hold for all t ∈ [0, T * ] for some time T * > 0 and some fixed constants C N and δ, and the main objective is to use the Einstein equations to prove that the bounds (1.21) in fact hold with better constants, provided the initial data are sufficiently small.  [45,47] that the constraint equations imply that the solution is trivial. The contribution of the mass is therefore identified explicitly using the decomposition (1.13), and the reduced Einstein equations are recast as a system of equations for h 1 : The term g h 0 μν is treated as an error term. Note that h 0 is defined so that h 0 μν , which is a good approximation to g h 0 μν , is supported away from the wave zone t ∼ r and so only contributes in the interior region, where h 1 will be shown to have fast decay.

Energy Inequality with Weights
An important ingredient in the procedure to recover the assumption on the energy (1.21) is the energy inequality with weights, which holds, for any suitably regular function φ, under mild assumptions on the metric g (see Lemma 6.29), where the weight w is as in (1.14) and (see Section 7 for details of how estimates for T μν follow from estimates for T μν ).
By the reduced Einstein equations (1.22) and the energy inequality (1.23), The first four terms on the right hand side arise already in [37] and so, combining estimates which will be shown for these terms in Section 6 (see also the discussion below) with (1.24), the bound (1.25) implies that and so the Grönwall inequality yields In the proof of Theorem 1.2, such a bound for Q N will lead to a recovery of the assumption (1.21) on E N with better constants provided C N is chosen to be sufficiently large and ε is sufficiently small.

The Structure of the Nonlinear Terms
As discussed above, whether a given quasilinear wave equation admits global solutions for small data or not depends on the structure of the nonlinear terms (moreover, the main analysis of the nonlinear terms is relevant in the wave zone, where T μν vanishes at late times and so plays no role in this discussion). A closer inspection of the nonlinearity in (1.22) reveals (see [11,37]) that where Q μν (∂h, ∂h) is a linear combination of null forms (satisfing |Q μν (∂h, ∂h)| |∂h||∂h|), |G μν (h)(∂h, ∂h)| |h||∂h| 2 denote cubic terms and Clearly the failure of the semilinear terms of the system (1.22) to satisfy the classical null condition arises in the P(∂ μ h, ∂ ν h) terms. In [35] it was observed that the semilinear terms of (1.22), after being decomposed with respect to a null frame N = {L, L , S 1 , S 2 }, where possess a weak null structure. It is well known that, for solutions of wave equations, derivatives tangential to the outgoing light cones ∂ ∈ T = {L , S 1 , S 2 } decay faster and so, neglecting such ∂h derivatives of h, (1.29) and, neglecting cubic terms and quadratic terms involving at least one tangential derivative, With respect to the null frame (1.28), the Einstein equations (1.22) become Decomposing with respect to the null frame (1.28), (see [36]). Except for the ∂ q h L L ∂ q h L L term, P N (∂ q h, ∂ q h) only involves the non L L components of h. The wave coordinate condition (1.12) with respect to the null frame becomes neglecting tangential derivatives and quadratic terms, see [36]. In particular, the ∂ q h L L ∂ q h L L term in P N (∂ q h, ∂ q h) can be neglected. The asymptotic identity (1.33) moreover implies (see equation (6.47)) that the leading order behaviour A decoupling therefore occurs in the semilinear terms of (1.22), modulo terms which are cubic or involve at least one "good" ∂ derivative, and the right hand side of the second identity in (1.31) only depends on components we have better control on by the first identity in (1.31). A further failure of (1.22) to satisfy the classical null condition arises in the quasilinear terms. Expressing the inverse of the metric g μν as g μν = m μν + H μν , the reduced wave operator takes the form which differs from the Minkowski wave operator only by the term H L L ∂ 2 q , plus terms which involve at least one tangential ∂ derivative. This main quasilinear term is controlled by first rewriting the wave coordinate condition (1.12) as (1.35) where |W ν (h, ∂h)| |h||∂h| is quadratic, and using the formula for any F μν , to rewrite ∂ μ H μν in terms of the null frame.
that is ∂ q H L L is equal to quadratic terms plus terms involving only tangential ∂ derivatives. Integrating ∂ q H L L from initial data {t = 0} then gives that H L L is approximately equal to the main contribution of its corresponding initial value, 2M/r .

Commutation
In order to apply the energy inequality (1.23) to improve the higher order energy bounds (1.21), it is necessary to commute the system (1.22) with the vector fields Z . Instead of commuting with the vector fields Z directly as in [37], notice that, for any function φ, where the fact that ∂ μ ∂ ν Z λ = 0 for each Z and μ, ν, λ = 0, 1, 2, 3 has been used.
Here L Z denotes the Lie derivative along the vector field Z (see Section 6.3 for a coordinate definition). The procedure of commuting the system (1.22) therefore becomes computationally much simpler if it is instead commuted with the Lie derivatives along the vector fields, L Z . In fact, the procedure simplifies further by commuting with a modified Lie derivative L, defined in the (t, x) coordinates by the formula The modified Lie derivative has the property that L Z m = 0 for each of the vector fields Z and moreover a computation shows that, in the case that φ μν is a (0, 2) tensor, the commutation property (1.36) holds for each of the vector fields Z . The commutation error in (1.36) can be controlled by ( L Z H L L )∂ 2 φ μν , plus terms which involve at least one tangential derivative of φ μν The Lie derivative along any of the vector fields Z commutes with partial derivatives ∂ (see Proposition 6.24). This fact leads to the commutation formula involving the modified Lie derivative. The term L Z H L L is then controlled easily by using the formula (1.35) and repeating the argument, described in Section 1.4.3, used to control H L L itself. When applying the commutation formula (1.36) to the reduced Einstein equations (1.22), it in particular becomes necessary to estimate the Lie derivative of the nonlinear terms, L I Z F μν (h)(∂h, ∂h) . Recall the nonlinear terms take the form (1.26). The modified Lie derivative also simplifies the process of understanding derivatives of the nonlinear terms due to the following product rule. Let h αβ and k αβ be (0, 2) tensors and let S μν (∂h, ∂k) be a (0, 2) tensor which is a quadratic form in the (0, 3) tensors ∂h and ∂k with two contractions with the Minkowski metric (in particular P(∂ μ h, ∂ ν k) or Q μν (∂h, ∂k)). Then and so the desirable structure of the nonlinear terms P(∂ μ h, ∂ ν h) and Q μν (∂h, ∂k) described in Section 1.4.3 is preserved after applying Lie derivatives.
The Lie derivatives of the energy momentum tensor are controlled by Theorem 1.3 since, for any function φ, the quantities |I |≤N |Z I φ| and |I |≤N | L I Z φ| are comparable.

The Klainerman-Sobolev Inequality with Weights
In order to control the derivatives of the nonlinear terms and the error terms arising from commuting the system (1.22), described in Section 1.4.4, when using the energy inequality (1.23), pointwise estimates for lower order derivatives of the solution are first shown to hold. The Klainerman-Sobolev Inequality can be used to derive non-sharp bounds for |∂ Z I h 1 | for |I | ≤ N − 3 (see equation (6.10)) directly from the bound on the energy (1.21). These pointwise bounds can be integrated from {t = 0} to also give pointwise bounds for |Z I h 1 | for |I | ≤ N − 3 which, using the fact that |∂φ| ≤ C 1+t+r |I |=1 |Z I φ| for any function φ, lead to strong pointwise estimates for all components of ∂ Z I h 1 for |I | ≤ N − 4. See Section 6.1 for more details.
Since, without restricting things to tangential derivatives, it is only true that |∂φ| ≤ C 1+|t−r | |I |=1 |Z I φ| for any function φ, the Klainerman-Sobolev Inequality does not directly lead to good pointwise estimates for all derivatives of h 1 . Some further improvement is necessary to control the terms in the energy estimate and recover the inequality (1.21).

L ∞ -L ∞ Estimate for the Wave Equation
The pointwise decay obtained for the transverse derivative of certain components of h 1 "for free" from the wave coordinate condition is not sufficient in the wave zone t ∼ r to close the energy estimate. The decay in this region is further improved by an L ∞ -L ∞ estimate, obtained by integrating the equations along the outgoing characteristics of the wave equation. In fact, instead of using the estimate for the full wave operator g , as in [37], it suffices to use the estimate for the operator 0 = m αβ − M r χ r 1+t δ αβ ∂ α ∂ β . Moreover, using the pointwise decay obtained from the Klainerman-Sobolev inequality and the wave coordinate condition, it can be seen that, for the purposes of this estimate, the essential contribution of the failure of (1.11) to satisfy the classical null condition is present in the P S (∂ q h, ∂ q h) terms, defined by (1.34). See Proposition 6.20 for a precise statement of this, and Lemma 6.21 for a proof of the L ∞ -L ∞ inequality. The fact that f is supported away from the wave zone can be shown using only the decay obtained from the Klainerman-Sobolev Inequality, and so the T term in (1.22) plays no role in Lemma 6.21. The pointwise decay of higher order Lie derivatives of h 1 is similarly improved in Section 6.4.

The
Hörmander L 1 -L ∞ Inequality Whilst the pointwise decay for lower order derivatives of h described above is sufficient to recover the assumptions (1.21) in the vacuum (when T μν ≡ 0), the interior decay is not sufficient to satisfy the assumptions of Theorem 1.3. In Proposition 6.11 the Hörmander L 1 -L ∞ inequality, Lemma 6.7, is used, together with the assumptions (1.21) on the L 1 norms of Z I T μν (t, ·) and on the energy E N (t) of h 1 , in order to improve the interior decay of h 1 and lower order derivatives. This improved decay ensures that the assumptions of Theorem 1.3 are satisfied and hence the theorem can be appealed to in order to recover the assumptions (1.21) on the L 1 norms of Z I T μν , and to control the L 2 norms of Z I T μν arising when the energy inequality (1.23) is used to improve the assumptions (1.21) on the energy E N . See Section 7 for further details on the completion of the proof of Theorem 1.2.

Vector Fields for the Vlasov Equation
The remaining difficulty in the proof of Theorem 1.2 is in establishing the L 1 and L 2 estimates of the vector fields applied to components of the energy momentum tensor, Z I T μν (t, x), of Theorem 1.3. For simplicity, we outline here how bounds are obtained for Zρ(t, x), for Z = i j , B i , S, where ρ(t, x) is the momentum average of f , defined by The bounds for Zρ(t, x) will follow from bounds of the form for a suitable collection of vector fields Z , acting on functions of (t, x, p ), which reduce to the Z = i j , B i , S vector fields when acting on spacetime functions, that is functions of (t, x) only. Throughout this section, and in Sections 4 and 5, it is convenient, instead of considering initial data to be given at t = 0, to consider initial data for the Vlasov equation to be given at t = t 0 for some t 0 ≥ 1. 4 It follows from the form of the Vlasov equation (1.18) that (1.38) (where X (t 0 , t, x, p ), P (t 0 , t, x, p ) are abbreviated to X (t 0 ), P (t 0 ) respectively) for any vector Z . Since derivatives of f | t=t 0 are explicitly determined by initial data and behave like f , an estimate for |Z f (t, x, p )| will follow from appropriate estimates for |Z X (t 0 , t, x, p ) i | and |Z P (t 0 , t, x, p ) i |.

General Procedure and Vector Fields in Minkowski Space
A natural way to extend a given vector field Z on M to a vector field on P, which by construction will have the property that Z (X (t 0 ) i ) satisfy good bounds, is as follows. For a given vector field Z on M, let Z λ : M → M denote the associated one parameter family of diffeomorphisms, so that Under a mild assumption on g, for fixed τ any point (t, x, p ) ∈ P with t > τ can be uniquely described by a pair of points is the point where the geodesic emanating from (t, x) with velocity p intersects the hypersurface τ (recall that { t } denotes the level hypersurfaces of the function t), that is (t, x, p ) ∈ P can be parameterised by {(t, x), (τ, y)} to get (t, x, p ) = (t, x, p X (t, x, τ, y)). 5 The subscript X is used in p X (t, x, τ, y) to emphasise that p is parametrised by y using the geodesics X (in contrast to suitable approximations to the geodesics, as will be considered later). Now the action of Z λ on (t, x) and (τ, y) induces an action on P at time t, given by For fixed t 0 we define the vector field Z by (1.40) which results in a good bound for |Z (X (t 0 ) i )|. In particular, if X i (t 0 , t, x, p ) and To see that (1.40) indeed holds, first note that the left hand side is the derivative of X t 0 , i with respect to λ at λ = 0. Also the first term on the right hand side is the derivative of Z λ (t 0 , y) i at λ = 0. The equality (1.40) follows from taking the derivative with respect to λ, and setting λ = 0, of both sides of the identity In Minkowski space, y has the explicit form y = x − (t − τ ) p and, when Z is chosen to be i j , B i , S, a straightforward computation, see Section 3, shows that the resulting vectors Z M take the form, Let q X (t, x, t 0 , y) = P(t 0 , t, x, p) be the initial momentum of the geodesic going through (t, x, p). By differentiating the equality with respect to λ, one can similarly obtain an estimate for Z ( P i (t 0 )). In the simple case of Minkowski space, q X has the explicit form The rotation vector fields M i j and a form of the scaling vector field S M were used in [46] for small data solutions of the massless Einstein-Vlasov system (note though that the above procedure of using Z to define Z breaks down when the mass shell P becomes the set of null vectors, as is the case for the massless Einstein-Vlasov system). 6 The vector fields where 0 < c < 1 and K ≥ 0. 7 On a spacetime whose Christoffel symbols decay as such, it can be shown that the Minkowski vector fields Z M of the previous section only satisfy 6 The proof in [46] is based on a double null foliation, and an associated double null coordinate system (u, v, θ 1 , θ 2 ), of the spacetimes which are constructed, and so the language used there is slightly different. In the coordinate system (u, v, θ 1 , θ 2 , p θ 1 , p θ 2 , p v ) conjugate to the double null coordinate system for M, the vector fields ∂ θ A , for A = 1, 2, are used. Defining appropriate Cartesian coordinates, one can show that ∂ θ A take the form of M i j . The proof in [46] in fact reduces to a semi global problem since the matter is shown, as part of the bootstrap argument in the proof, to be supported in a strip of finite retarded u length. The vector (v − u)∂ v is also used which, since u remains of size 1 in the support of the matter, agrees to leading order with the vector field u∂ u + v∂ v which, when written with respect to an appropriate Cartesian coordinate system, is seen to be equal to S M . 7 It should be noted that there are two contributions to this slow interior decay. The first arises from the failure of the Einstein equations in the harmonic gauge to satisfy the classical null condition of [26]. Indeed, it was recently shown by Lindblad [34] that small data solutions of the vacuum Einstein equations in the harmonic gauge satisfy this decay rate (compare with [14] where the Ricci coefficients associated to the maximal-null foliation decay in the interior at a faster rate). The second contribution arises from the presence of the Vlasov matter, in the form of the energy momentum tensor as a source term in the Einstein equations. This fact can be more easily seen in a simplified setting. Indeed, if T (t, x) denotes a function which decays at rate t −3 for |x| ≤ ct + K and vanishes for |x| ≥ ct + K -the sharp behaviour of the components of the energy momentum tensor associated to solutions of the Vlasov equation on Minkowski space-the sharp interior behaviour of solutions of for solutions f of the Vlasov equation. This logarithmic loss compounds at higher orders, and cannot be used to recover the sharp bounds (1.42) in the context of Theorem 1.2. The proof of Theorem 1.2 is therefore based on a different collection of vector fields, Z , which are adapted to the geometry of the background spacetime and again reduce to Z = i j , B i , S when acting on spacetime functions, and satisfy a good bound of the form (1.37) when applied to solutions f of the Vlasov equation. The vector fields can be derived using the procedure described in Section 1.5.1, which in fact did not rely on the background spacetime being Minkowski space. Instead of the expression (1.39), the components y i are defined using approximations X 2 (s, t, x, p ) to the true geodesics X (s, t, x, p ) of the spacetimes to be constructed 8 . One could also use the geodesics themselves, however doing so involves estimating the components T μν at the top order in a slightly different way to how they are estimated at lower orders. We choose to use the approximations to the geodesics so that T μν can be estimated in the same way at all orders. A derivation of the vector fields obtained using this procedure is given in Section 3, but here the failure of the Minkowski vector fields Z M are identified explicitly and shown how to be appropriately corrected. The two procedures agree up to lower order terms. Instead of the sharp interior bounds (1.42), the proof of Theorem 1.3 requires only the weaker bounds for |I | ≤ N /2 + 2, where 1 2 < a < 1. Consider first the rotation vector fields, and recall the expression (1.38). The rotation vector fields i j are defined using approximations to the geodesics. The geodesics take the form Using the fact that where each of the first approximations arise by replacing P (s, t, x, p ) k and X (s, t, x, p ) k by their respective values in Minkowski space, and the second arise from the fact that x i t ∼ p i , which holds asymptotically along each given geodesic (see Proposition 2.2 below), the approximations to the geodesics are defined as for t 0 ≤ s ≤ t. It will be shown in Section 2 that X 2 (s, t, x, p ) k are good approximations to the geodesics X (s, t, x, p ) k in the sense that for all t 0 ≤ s ≤ t and k = 1, 2, 3. The idea is now to construct vector fields so that the vector fields applied to X 2 are bounded. Then one can show that (1.45) is true with X 2 − X replaced by Z (X 2 − X ). See Section 4 for more details. The approximations X 2 have the desirable property, which will be exploited below, that vanishes (and in particular does not involve derivatives of ).
Applying the Minkowski rotation vector fields to the approximations X 2 gives In the final equality, for (t, x, p ) ∈ supp( f ), the first two terms are bounded (as can be seen from (1.45) and the fact that f | t=t 0 has compact support). However on a spacetime only satisfying the bounds (1.43), the terms on the last line in general grow in t. The vector fields i j are defined so that these terms are removed: where if the functions˚ l i j are defined as it follows from (1.46) and (1.47) that and so for (t, x, p ) ∈ supp( f ) and k = 1, 2, 3. It can similarly be shown (see Section 4.2) that which, by (1.38), leads to a good bound for i j f (t, x, p ). We remark that the identity (1.48) can be expressed as which is exactly (1.40) when Z = i j and X is replaced by X 2 . The rotation vector fields i j therefore arise by following the procedure in Section 1.5.1 with X replaced by X 2 ; see Section 3.
Since the functions˚ l i j do not depend on p , and so the good bounds for i j f (t, x, p ) lead to good bounds for i j ρ(t, x). Had the true geodesics been used to define the vector fields i j , the functions˚ l i j would have involved a term of the form P (s ) ds and so ∂ p k˚ l i j would involve second order derivatives of . The average ˚ l i j ∂ p l f (t, x, p ) d p has better decay than the term˚ l i j ∂ p l f (t, x, p ) does pointwise and so, to exploit this fact when obtaining a bound for i j ρ(t, x), it is necessary to integrate this term by parts If the true geodesics are used to define i j then, in the setting of Theorem 1.2, the fact that ∂ p k˚ l i j involves second order derivatives of would mean that this integration by parts cannot be used when estimating derivatives of T μν at the top order, and so T μν would be estimated in a slightly different way at the top order.
The approximations to the geodesics X 2 (s, t, x, p ) are used in a similar manner in Section 4.1 to define vector fields, B i and S. 9 In Hwang-Rendall-Velázquez [23] derivatives of the average f (t, x, p) d p for solutions of the Vlasov-Poisson system are controlled and L ∞ analogues of the estimates of Theorem 1.3 are obtained. The approach is different to that taken in the present work, though the estimates in [23] are obtained by similarly first controlling derivatives of the analogues of the maps X , P. The re-parameterisation (t, x, y) of P is also used, though whilst here it is only formally used to motivate the definition of the Z vector fields, in [23] ∂ x derivatives, in the (t, x, y) coordinate system, of the analogue of the maps X , P are controlled.
It is also interesting to compare the present work with the independent work of Fajman-Joudioux-Smulevici [19]. The proof in [19] is also based on the vector field method (following [29,37] in the vacuum) and, accordingly, the new elements of the proof also involve controlling vector fields applied to the components of the energy momentum tensor T μν . The issue of the failure of the Minkowski vector fields applied to f to be bounded arises. A key step in [19] is therefore also to introduce a new collection of vector fields further adapted to the geometry of the spacetimes under consideration. The vector fields introduced are different to those introduced here. The construction in [19] proceeds roughly by, in the (t, x, p) coordinate system for P, considering the vectors

Outline of the Paper
In Section 2 the system (1.17) is used to prove the property (1.19) regarding the support of the matter, and the maps (1.44) are shown to be good approximations to the geodesics. In Section 3, which is not required for the proof of Theorem 1.2 or Theorem 1.3, the discussion of the vector fields in Section 1.5 is expanded on and a derivation of the vector fields used in the proof of Theorem 1.3 is given. In Section 4 combinations of the vector fields applied to the geodesics are estimated. In Section 5 the estimates of Section 4 are used to obtain estimates for spacetime vector fields applied to the components of the energy momentum tensor and hence prove Theorem 1.3. In Section 6 the solution of the reduced Einstein equations is estimated in terms of the components of the energy momentum tensor. The results of the previous sections are combined in Section 7 to give the proof of Theorem 1.2.

The Support of the Matter and Approximations to the Geodesics
In this section the Vlasov equation on a fixed spacetime is considered. It is shown that under some assumptions on the metric the solution is supported, for large times, away from the wave zone x ∼ t provided f 0 is compactly supported. Curves X 2 , which approximate the timelike geodesics X , are introduced, which are later used to define the Z vector fields.
Note that the characteristics of the operator 1 p 0 X, where X is given by (1.5), and so d ds

Properties of the Support of the Matter
The results of this section will rely on the Christoffel symbols satisfying the assumptions for various small N . Then if (2.7) below holds in supp( f ), the assumption (2.5) implies that (2.6) in supp( f ). The notation (y, q) will be used for points in the mass shell over the initial hypersurface {t = 0}, so y ∈ 0 , q ∈ P (0,y) . The next result (Proposition 2.8) guarantees that for (t, x, p ) ∈ supp( f ), for some constants K and c.
Proposition 2.1. Suppose that |y| ≤ K and |q| ≤ K , for some K, K ≥ 1, and that for some fixed a, δ > 0 Proof. Let s 1 be the largest number such that |P(s, 0, y, q)| ≤ 2K for 0 ≤ s ≤ s 1 . We will show that then it follows that |P(s, 0, y, q)| < K + 1, for 0 ≤ s ≤ s 1 , contradicting the maximality of s 1 . Let p = P(s, 0, y, q). Since g αβ p α p β = −1 it follows that 1 which proves the second and hence third part of (2.9), assuming the weaker bound of the first. It follows that by assumption. Hence, using the weak inductive assumption, It follows that |P(s, 0, y, q)| ≤ K + 1/4 in the support of f , and (2.9) follows.
for all s ≥ 0, using the bound (2.8). The proof follows by integrating forwards from s = 0 and using the fact that |y| + |q| ≤ C, and dividing by s.

Translated Time Coordinate
It is convenient to use this translated time coordinatet in what follows. In particular, the vector fields of Section 4 will be defined using this variable. The main advantage is that the spacetime vector fieldsZ for some smooth functions I J . Estimates for ∂ IZ J T μν will then follow directly from estimates forZ I T μν , which are less cumbersome to obtain; see Section 5.4. For simplicity the˜will always be omitted, and statements just made for t ≥ t 0 .

Approximations to Geodesics
Define, for (t, x, p ) ∈ supp( f ) and t 0 ≤ s ≤ t, and set The geodesic equations (2.2) can be used to derive the following equations for X and P: It follows from the next proposition that the curves s → X 2 (s, t, x, p ) are good approximations to the geodesics s → X (s, t, x, p ). Recall, from Section 1.2.2, the notation .
Proof. First note that equation (2.2) and the bounds assumed on imply where y = X (0, t, x, p ), q = P (0, t, x, p ). Proposition 2.2, and the integration of (2.17) from s to t and (2.11) then gives Moreover, Integrating equation (2.15) backwards from s = t, and using the fact that P(t, t, x, and integrating the equation (2.14) backwards from s = t and using X (t, t, x, p ) i = 0 gives Proof. The proof is an immediate consequence of the first bound of (2.7), and Proposition 2.3.

The Vector Fields
The general procedure for using approximate geodesics to lift a vector field Z on M to a vector field Z on P, outlined in Section 1.5 for geodesics, is described in more detail in Section 3.1. Here a formula for the action of the lifted vector fields on initial conditions for the approximate geodesics is derived. The essential property of the lifted vector fields is that they are bounded when applied to these initial conditions, which one can deduce from this formula. It is however computationally more convenient to define the vector fields from their action on initial conditions for the approximate geodesics, which we do in Section 3.2. Section 3.2 is independent of Section 3.1, but we include Section 3.1 for the purpose of conceptual justification of the vector fields and the formula mentioned above.

Parametrization of Momentum Space with Physical Initial and Final
Coordinates for the Approximate Geodesics Given first approximations to the geodesics X 1 (s, t, x, p ) and P 1 (s, t, x, p ) we define the second approximations X 2 (s, t, x, p ) and P 2 (s, t, x, p ) to the geodesics through (t, x, p ), to be the solutions of the system d ds Under some mild assumptions on the metric g, 10 for fixed τ any point (t, x, p ) ∈ P with t > τ can be described uniquely by the pair of points is the point where the approximate geodesic X 2 emanating from (t, x) with velocity p intersects the hypersurface τ , that is (t, x, p ) ∈ P can be parameterised by where the subscript X 2 is used in p X 2 (t, x, τ, y) to emphasise that p is now parametrised by y using the approximations X 2 to the geodesics. Integrating (3.1) backwards from t to t 0 gives that p = p X 2 (t, x, t 0 , y) is implicitly given by Here we used Taylor's formula with integral remainder The first approximate geodesics in the previous section are independent of p so then is independent of p and the above relation in fact gives p X 2 explicitly.

Lifting of Physical Vector Fields to Vector Fields on Momentum Space
Adapted to the Approximate Geodesics For a given vector field Z on M, let Z λ : M → M denote the associated one parameter family of diffeomorphisms, Now the action of Z λ on (t, x) and (τ, y) in M induces an action on (t, x, p ) in P, given by For fixed t 0 we define the vector field Z by We have where

The Lifted Vector Fields Applied to the Initial Conditions of the Approximate Geodesics Parameterized by the Final Momentum Space A computation shows that
To see that (3.7) indeed holds, first note that the left hand side is the derivative of i with respect to λ at λ = 0. Also the first term on the right hand side is the derivative of Z λ (t 0 , y) i at λ = 0. The equality (3.7) follows from taking the derivative, with respect to λ at λ = 0 of both sides of the identity The identity (3.8) is just two different ways of expressing that we apply the transformation λ to the initial point of the path (t 0 , y). The left hand side is that we apply the transformation to the approximate geodesic which will lead to the transformation of the initial point.

Computation of the Lifted Vector Fields from Their Action on Initial Conditions for the Approximate Geodesics
In this section we will compute the vector fields from how they act of the initial conditions. We will use a slight modification of the formula for how they act on the initial conditions derived in the previous section and we will hence get a slightly modification of vector fields in the previous section. They are computationally slightly simpler to use than the vector fields of the previous section, though they are mildly singular for t close to t 0 . It is therefore assumed throughout this section that t ≥ t 0 + 1, and the vector fields will only be used in the proof of Theorem 1.3 under this assumption.
The vector fields can be derived by imposing that they have the form where Z i are to be determined, and insisting that the relation holds, instead of the relation (3.7) which involves P 2 (t 0 , t, x, p ) i instead of p i and is satisfied by the vector fields of the previous section. Equation (3.10) is in particular true for the Minkowski vector fields Z M .
We now make the further assumption that the first approximate geodesics X 1 and P 1 in (3.1) are independent of p in which case by (3.3) the second approximation to the geodesics are given by with given by (3.4) now independent of p . Applying the expression (3.9) to (3.11) gives Substituting the definition (3.10) into (3.12) and solving for Z i we obtain The coefficients of the vector fields we are considering are linear functions where p 0 = 1 and we defined 0 = 0. Hence Proof. Changing variables gives Integrating by parts we have and adding things up we get the lemma.

Remark 3.2.
We note Z is of exactly the same form of and hence can be estimated in the same way, and moreover a vector field Z applied to Z will produce a ( X ) Z of the same form.
One can further obtain bounds for Z applied to the initial conditions for P as follows: we have which again is bounded, provided the components of Z grow at most like t.

The Rotation Vector Fields
Since the rotations satisfy Z 0 = 0 it easily follows from the above lemma that where is given by (3.17) and

The Scaling Vector Field
By the above lemma we have

Vector Fields Applied to Geodesics
Recall from Section 2 the definitions (2.13) of X 2 (s, t, x, p ) i and P 2 (s, t, x, p ) i , and In this section estimates for combinations of appropriate vector fields, Z I (see Sections 3 and 4.1 below), applied to the geodesics X (t 0 , t, x, p ) and P (t 0 , t, x, p ) of a fixed spacetime satisfying the pointwise bounds The vector fields Z are defined so that such good estimates indeed hold (see Proposition 4.1 for the case |I | = 1, and Proposition 4.9 for general I ).
In Section 4.1 the vector fields are recalled from Section 3. In Section 4.2 the basic idea of the remainder of the section is illustrated by considering only one rotation vector field i j applied to the geodesics. In Section 4.3 schematic ex- In Sections 4.1-4.9 it will be assumed that t ≥ t 0 + 1. Section 4.10 is concerned with the case t 0 ≤ t ≤ t 0 +1. It will be assumed throughout this section that |x| ≤ ct, where 0 < c < 1 (recall the discussion in Section 2.2).

Vector Fields
The proof of the main result uses the following collection of vector fields, introduced in Section 3, which are schematically denoted Z for i, j = 1, 2, 3, i < j, and reduce to the standard rotations, boosts and scaling vector field of Minkowski space when acting on spacetime functions. Recall the notation, for i = 1, 2, 3, Recall that the vector fields i j , B i , S take the following form: first, for i = 1, 2, 3. Finally, By the bound (2.16), the maps X 2 (s, t, x, p ) i are good approximations to the true geodesics, X (s, t, x, p ) i . Recall that the vector fields i j , B i , S are defined so that, when applied to X 2 (s, t, x, p ) i the following hold at s = t 0 : for i, j, k = 1, 2, 3.
Proof. The proof is a straightforward computation. See also Section 3.2, where the above form of the vector fields i j , B i , S are derived by insisting that the result of this proposition is true. (Note that (4.7)-(4.9) are nothing other than equation Since, in Minkowski space, the maps X (s, t, x, p ) i are simply given by it is easy to see that Proposition 4.1 in fact holds in Minkowski space with t 0 replaced by any s:

Estimates for One Rotation Vector Field Applied to the Geodesics
Recall that the main goal of Section 4 is to control combinations of vector fields applied to the components of X (t 0 , t, x, p ) and P (t 0 , t, x, p ), that is to prove Proposition 4.24 and Corollary 4.27. Proposition 4.24 and Corollary 4.27 require introducing schematic notation and first obtaining preliminary results (see Sections 4.3-4.7). The main idea, however, is straightforward and can already be understood. In order to illustrate the idea, in this section it is shown how to estimate one rotation vector applied to the components of X (t 0 , t, x, p ) and P (t 0 , t, x, p ). This section is included for the purpose of exposition and the results are not directly used in the proof of Proposition 4.24 and Corollary 4.27 or elsewhere.
By Proposition 4.1, in order to control i j applied to the components of X (t 0 , t, x, p ) and P (t 0 , t, x, p ) it suffices to control i j applied to the components of X (t 0 , t, x, p ) and P(t 0 , t, x, p ) or, more generally, to show the following: The proof of Proposition 4.3 relies on two facts. The first is the fact that, for which follows from the fact that X (t, t, x, p ) = P(t, t, x, p ) = 0, and that i j does not involve ∂ t derivatives. The analogue of (4.10) is not true, due to the presence of the ∂ t derivatives, when i j is replaced by B i or S. These "final conditions" for B i and S, and higher order combinations of all of the Z vector fields, applied to X and P are estimated in Section 4.7.
The second fact required for the proof of Proposition 4.3 is the following estimate for i j applied to the right hand side of equation (2.15).
for all t 0 ≤ s ≤ t.
Proof. Recall that The bound for the first term follows from writing X (s) l −s x l t = X (s) l +X 2 (s) l −s x l t and, using the definition (2.13) for X 2 and the definition of˚ l i j (see (4.3)), Using the assumptions (4.1) on , the fact that For the second term, note that and, as above, Hence, Similarly, using the fact that it follows that and similarly, In a similar manner it follows that from which (4.11) then follows.
A higher order analogue of Proposition 4.4, which includes also the boosts and scaling and is used in the proof of Proposition 4.24 and Corollary 4.27, is obtained in Section 4.5.
Proof of Proposition 4.3. The proof proceeds by applying i j to the system (2.14)-(2.15). Applying i j to the equation (2.14) and integrating backwards from s = t, using (4.10), that and, summing over k = 1, 2, 3, the Grönwall inequality (see Lemma 4.25) gives Inserting this bound into the equation (2.14) for X , after applying i j , integrating backwards from s = t and using (4.10) gives The proof follows after inserting this bound back into (4.12).

Repeated Vector Fields Applied to the Initial Conditions for Approximations to Geodesics
Recall the discussion at the beginning of Section 4. In order to motivate the results of this section and the next section note that, after applying Z I to the equation (2.15), the term (amongst others) appears. Consider the first of these summands. Rewriting the first term Z I (X k ) is controlled, in the proofs of Propositions 4.23 and 4.26, by the Grönwall inequality. In Section 4.4, the term Z I X k 2 (s) − sx k t is controlled by first rewriting t, x, p )).
The second term is computed schematically in Proposition 4.5 below, and is controlled in Section 4.4. For the first term, the fact that is used, along with the schematic expressions of Propositions 4.6 and 4.7. The expression 14) similarly appears after applying Z I to equation (2.15), where and is therefore also estimated in Section 4.4.
The following generalises Proposition 4.1 to higher orders. Recall it is assumed that t ≥ t 0 + 1. Section 4.10 is concerned with the case t 0 ≤ t ≤ t 0 + 1. Proof. The result is clearly true for |I | = 1 by Proposition 4.1. The result for |I | ≥ 2 then follows from a straightforward induction argument after noting that, for any multi index J j ,

Proposition 4.5. For any multi index I , there exist smooth functions i I , i I, j such that
where |L| = |J j | + 1, and also noting that where the equalities (4.4), (4.5), (4.6) have been used. Hence, for any smooth function for some smooth functions i,J 1 ,...,J k I,i 1 ,...,i k .
The result for |I | ≥ 2 follows from a straightforward induction, as in the proof of Proposition 4.5, now also using the fact that i j

Preliminary Estimates for Repeated Vector Fields Applied to Approximations to Geodesics
Proposition 4.8. Suppose t ≥ t 0 + 1, |x| ≤ ct and the bounds (4.1) hold. Then, for |I | ≤ N , where Z I is a product of |I | of the vector fields i j , B i , S. Moreover, if |I | ≤ N 2 + 2, then Proof. Recall the definition (4.2) of . Note that and that i j μ αβ s, s and S μ αβ s, s Using the equalities (4.18), the result for |I | = 1 clearly follows from the L ∞ bounds (4.1) for μ αβ (recall that, when restricted to spacetime functions, the Z vector fields are equal to the Z vector fields). The proof for |I | ≥ 2 is a straightforward induction. Proposition 4.9. Supposing that t ≥ t 0 + 1, |x| ≤ ct, (t, x, p ) ∈ supp( f ) and the bounds (4.1) hold, then, for |I | ≤ N , for i = 1, 2, 3 and for all t 0 ≤ s ≤ t. In particular, for |I | ≤ N 2 + 2, Proof. Consider first the bound (4.19). Proposition 4.5 implies that where the fact that |X 2 (t 0 ) i | 1 has been used (see Corollary 2.4). If |J l | ≤   for all t 0 ≤ s ≤ t, i = 1, 2, 3. In particular, for |I | ≤ N 2 + 2, Moreover, In particular, for |I | ≤ N 2 + 2, Proof. After writing Proposition 4.9 implies that The proof of the bound (4.24) follows after dividing by s and using the fact that if s ≤ s (recall also that t ≥ t 0 + 1). The proof of (4.25) follows similarly by writing and using the bound (4.20). The bound (4.27) similarly follows from Proposition 4.9 after writing The lower order estimates (4.26) and (4.28) follow from rewriting and using the pointwise bounds (4.1) for lower order derivatives of .

Schematic Notation and Repeated Vector Fields Applied to Differences
To make long expressions more concise, the s dependence of many quantities is suppressed throughout this section.
The proofs of Propositions 4.23 and 4.26 below follow from applying vector fields to the system (2.14)-(2.15). It is therefore necessary to estimate vector fields applied to the difference which appears on the right hand side of equation (2.15), where X 1 (s) = sx/t and P 1 = dX 1 ds = x/t, and is defined in (2.4). Apart from controlling vector fields applied to X 1 and P 1 we also control vector fields applied to X 2 and P 2 , and moreover the differences X 1 − X 2 and P 1 − P 2 decay, see (4.13) and (4.14) and the propositions to follow. Furthermore, because of the definition of X 2 and P 2 in terms of X 1 and P 1 the differences X = X − X 2 and P = P − P 2 are small, see Proposition 2.3. It may be tempting to write this as differences with evaluated at (X 2 , P 2 ). However we do not want to involve an estimate of applied to X 2 so instead we will first differentiate (4.29) and use the decomposition X 1 − X = X 1 − X 2 − X to the factors that come out when we differentiated. The following result is straightforward to show: Proof. First one differentiates to get a sum of terms of the form and then one makes a different decomposition depending on the size of k. Then one proceeds by organizing them in order so |I 1 | is smallest. If k < |J |/2, one writes and replaces them one by one: This produces the first two sums with k < |J |/2. For k ≥ |J |/2 one simply does the same thing but with X 1 and X interchanged.
Note that, for each term in the equality (4.30), |J i | ≤ |J |/2 if i = k and therefore, in the applications of Lemma 4.11 below, the factors Z J i X can be estimated by the their L ∞ norms using induction by previous estimates. Given that X 1 is known the expression (4.30) can be thought of as linear in the unknown X − X 1 . However when we prove L 2 estimates for high derivatives we also have to take into account how Z J k X 1 depends on high derivatives of . Either k ≤ |J |/2 is small, in which case we can use L ∞ estimates for ∂ k F and ∂ k+1 F, or k ≥ |J |/2 is large, and as a result |J i | ≤ |J |/2 are small for all i, and we can use L ∞ estimates for all the other factors. Applying Lemma 4.11 we get the following estimates: Proof. The bound (4.31) follows from Lemma 4.11 applied to ( P ), after noting that |Z I ( P 1 )| 1 for any I , using the form of the vector fields Z and the fact that since is smooth. The bound (4.32) is even simpler.
Proof. The Lemma is again a straightforward application of Lemma 4.11, noting again that Z I X 1 s 1 for any I since X 1 (s, t, x) = s x t . Note that for some homogeneous functions A that are smooth when |x|/t ≤ c < 1, and hence, In the application we will estimate low derivatives of with its L ∞ norm so that the result only depends on the differences of functions evaluated at X and X 1 and functions evaluated at X 1 but not at X . Lemma 4.14. Suppose t ≥ t 0 + 1, |x| ≤ ct, (t, x, p ) ∈ supp( f ) and suppose, for some t 0 ≤ s ≤ t, that |Z K X | ≤ Cs for |K | ≤ |J |/2. Then, at (s, t, x, p ), Instead of applying vector fields to the decomposition (4.29) we first apply vector fields and then apply this decomposition to the terms with more derivatives falling on , but if more fall on we apply the decomposition with ( P) interchanged with ( P 1 ) and (X ) with (X 1 ): Using this decomposition and the previous lemmas we obtain Proof. The above decomposition and Lemma 4.12 give The first bound then follows from Lemma 4.14. The proof of the second bound is straightforward.
The following corollary of Proposition 4.15 is used at lower orders: Proof. The proof follows from Proposition 4.15 if we note that in the support of X and X 1 ,

Parameter Derivatives of the Equations and Vector Fields Applied to Their Differences
The vector fields applied to the system (2.14)-(2.15) will be estimated by integrating from the final time t and in order to do this we need to control the final conditions for Z I X and Z I P at time t. Note that, for i, j = 1, 2, 3, t, x, p )) (similarly for ∂ p i in place of ∂ x i and for X 2 , P , P 2 respectively in place of X ) and so, if Z is a vector which does not involve ∂ t derivatives, these final conditions vanish. Some of our vector fields, however, also involve ∂ t derivatives. Therefore we need to estimate higher d ds derivatives of the system, which can be recast as spacetime derivatives. Recall that X (s) = (s, X (s)). We The structure of higher order d ds derivatives is very simple. Either the derivative falls on ( P ), in which case we can substitute the second equation for d P/ds and get another factor of ( X ), or the derivative falls ( X ), which produces a derivative ∂ . Hence we get the system (4.38) Here ∂ k , denote the R 4k tensor with components sum of ∂ j 1 · · · ∂ j k over all components 0 ≤ j i ≤ 3 for i = 1, . . . , k and k 1 ,...,k m ( P ) are polynomials in P with values in R 4k 1 +...4k m . Here the first dot product is schematic notation to be interpreted as dot products of elements (k) and (k) in some larger dimensional space whose components corresponds to the terms in the sum. We now also want to take d ds derivatives of d X 2 ds = P 2 , d P 2 ds = ( X 1 ) · ( P 1 ), and d X 1 ds = P 1 , Here (k) satisfies the same estimates as ∂ k whereas (k,2) is at least quadratic in (note m ≥ 2 in the summation) and hence satisfies the same estimates as ∂ k−1 multiplied with , which decays s −a faster than another derivative. We have
Using (4.33) we can write this as a linear combination of The difference G I 1 ...I m ( X) − G I 1 ...I m ( X 1 ) can be estimated by the right hand side of (4.44) as for (4.34).
We can now apply the estimates from the previous section with (k) respectively (k,2) in place of . First we obtain the analogue of Proposition 4.15 as follows: Proof. By Proposition 4.15 applied to (k) in place of , and the first part of the proposition follows from Lemma 4.17.
The following corollary of Proposition 4.18 is used at lower orders: The proof of Corollary 4.19 follows from Proposition 4.18 as in the proof of Corollary 4.16.

The Final Conditions
Note that, for any Y (s, t, x, p ), Repeated application (4.45) inductively implies that where Y (k) = d k Y/ds k and Z J i (t) = Z J i (t) are constants times t or x j , for some j.
Applying (4.46) to Y = X and Y = P noting that X (t, t, x, p ) = P(t, t, x, p ) = 0 gives where for the proof of (4.48) we used that X (k+1) = P (k) and for k = 0 we also used (4.47). Hence where everything is evaluated at (s, t, x, p ) where s = t.
We can now apply the estimates from the previous section.

Proposition 4.21.
Suppose that t ≥ t 0 + 1, |x| ≤ ct, (t, x, p ) ∈ supp( f ) and the bounds (4.1) hold. Then, for |I | ≤ N /2 + 2, Proof. We will use induction to prove this. Assuming that we have that |I | < m ≤ N /2 + 1, the assumptions of Corollary 4.19 hold at s = t, and writing X − X 1 = X + X 2 − X 1 and P − P 1 = P + P 2 − P 1 and using the estimates (4.26) and (4.28) for X 2 − X 1 and P 2 − P 1 , respectively, we get, for |L| + k ≤ m − 1, Assuming that the proposition is true for |I | < m it now follows from this and previous lemma that it is also true for |I | = m.

Proposition 4.22.
Suppose t ≥ t 0 + 1, |x| ≤ ct, (t, x, p ) ∈ supp( f ) and the bounds (4.1) hold. Then, for |I | ≤ N , Proof. We will use induction to prove this. Assuming that the proposition is true for |I | < m ≤ N , the assumptions of Proposition 4.18 hold at s = t and so, for since |Z I X 1 | t + |Z I P 1 | ≤ C for any multi index I , using the form of the vector fields Z . Hence, writing X − X 1 = X + X 2 − X 1 and P − P 1 = P + P 2 − P 1 , we get Using induction for the first sum and the estimates (4.25) and (4.27) for X 2 − X 1 and P 2 − P 1 , respectively, the proposition follows.

L ∞ Estimates for Lower Order Derivatives of Geodesics
The estimates in the previous sections easily lead to pointwise bounds for lower order derivatives of X (s, t, x, p) and P(s, t, x, p). Proposition 4.23. Suppose t ≥ t 0 + 1, |x| ≤ ct, (t, x, p ) ∈ supp( f ) and the bounds (4.1) hold. Then, for i = 1, 2, 3, Proof. The proof proceeds by induction. Clearly the result is true when |I | = 0 by Proposition 2.3. Assume the result is true for all |I | ≤ k, for some k ≤ N 2 . Then I clearly satisfies the assumptions of Corollary 4.16 and so, by the equations (2.14), (2.15) and the pointwise bounds (4.1), where we recall X 1 (s, t, x) = s x t and P 1 (t, x) = x t . Writing X − X 1 = X + X 2 − X 1 and P − P 1 = P + P 2 − P 1 and using the estimates (4.26) and (4.28) for X 2 − X 1 and P 2 − P 1 , respectively, gives Integrating backwards from s = t and using Proposition 4.21, and so, after summing over i = 1, 2, 3 and I , the Grönwall inequality and Lemma 4.25, give that The equation (2.14) and Proposition 4.21 then give, after integrating backwards from s = t again, where the fact that, for any function λ(s), for all t 0 ≤ s ≤ t, for |I | = 0, 1, 2, . . . , N 2 + 1. The following form of the Grönwall inequality was used in the proof of Proposition 4.23 above, and will be used in the proof of Proposition 4.26 below:

Higher Order Estimates for Derivatives of Geodesics
The main result of this section is Proposition 4.26. Suppose that t ≥ t 0 + 1, |x| ≤ ct, (t, x, p ) ∈ supp( f ) and the bounds (4.1) hold. Then, for i = 1, 2, 3, and for all t 0 ≤ s ≤ t, |I | ≤ N .
Proof. Let I be a multi index with |I | ≤ N . Using the equation (2.15) and Proposition 4.15, Writing and using Corollary 4.10, Proposition 4.23 and the pointwise bounds (4.1) for gives Integrating backwards from s = t gives Summing over i = 1, 2, 3 and I , the Grönwall inequality, Lemma 4.25, gives Integrating backwards from s = t again, the equation (2.14) implies where the fact that, for any function λ(s), The bound (4.49) follows from Proposition 4.22, along with the fact that Proof. The corollary is an immediate consequence of Propositions 4.9 and 4.26.

Spacetime Derivatives and Small Time
Since the vector fields Z become singular at time t = t 0 , in this section the spacetime ∂ t and ∂ x i derivatives of X (s, t, x, p ) and P (s, t, x, p ) are estimated for t 0 ≤ t ≤ t 0 + 1. Since the results of this section are local in time they are much simpler than those in previous sections. In particular, it is not necessary to subtract the approximations X 2 , P 2 from X and P respectively. Note that ∂ always denotes the spacetime gradient ∂ = (∂ t , ∂ x 1 , ∂ x 2 , ∂ x 3 ). When applied to functions on P the derivatives are, as usual, taken with respect to the (t, x, p ) X (s)) .
Proof. Recall that (X, P ) = (X ) · ( P ) (the s dependence is omitted for brevity). There exist constants c J 1 ...J k I K , c J 1 ...J k I K such that and and so, by the assumed lower order bounds for X and P , from which the result follows.
In order to use the system (1.17) to estimate ∂ I X (s, t, x, p ) and ∂ I P (s, t, x, p ), it is also necessary to estimate the final conditions (note that this is completely straightforward unless ∂ I contains ∂ t derivatives).

Proposition 4.29.
Let I be a multi index with |I | ≥ 1 and suppose |∂ J (t, x)| ≤ C for all |J | ≤ |I | 2 + 1. Then Proof. Recall the notation P (k) and (k) from Section 4.6. By the formula (4. 46) it follows that t, t, x, p ), t, t, x, p ) for some constants C I J 1 ...J k+1 ,L , C I J 1 ...J k+2 ,L , where the proof of the second uses the first and the fact that dX (k) ds = P (k) . Hence, t, t, x, p ) , t, t, x, p ) .
The proof follows by noting that by an appropriate version of Lemma 4.17.
The proof then proceeds exactly as in Proposition 4.30.

Estimates for Components of the Energy Momentum Tensor
In this section a proof of Theorem 1.3 is given. Recall the discussion in Section 2.2. In order to use the results of Section 4 it will again be assumed throughout most of this section that t ≥ t 0 +1, the bounds (4.1) hold, and π (supp( f )) ⊂ {|x| ≤ ct}, where π : P → M is the natural projection. It is shown how Theorem 1.3 then follows in Section 5.4.

Derivatives of Components of the Energy Momentum Tensor in Terms of Derivatives of f
Recall The main result of this section is Proposition 5.4, which uses the bounds on Z I X and Z I P of Corollary 4.27 to give bounds on Z I T μν . In order to prove the bounds for Z I T μν , it is convenient to first rewrite the above integral in terms of the p i variables. h(t, x)).

Proposition 5.1. There exists a non-zero function , smooth provided
Proof. Define p 0 = 1 and note that, since, it follows that 2g αβ p α ∂ p β ∂ p j = 0, and hence, The proof follows by writing g αβ = m αβ + h αβ , noting that and using the fact that det(A −1 ) = (det A) −1 for any matrix A.
In Minkowski space, that is when h = 0, it is straightforward to compute where p j M = p j p 0 M , and p 0 M is defined by the relation m αβ p α p β = −1. It then follows that since | p i | ≤ c < 1, and hence the change of variables (t, x, p) → (t, x, p ) is well defined if ε is sufficiently small. Moreover, recalling that it follows that, for each μ, ν = 0, 1, 2, 3, for some functions μν , smooth when | p | ≤ c < 1, and so For each vector field Z , recall the corresponding functionsZ k (t, x, p ) defined in Section 4.1. Note that, for each Z , theZ k have the form , x) , This notation will be used below.
where Z I is a product of |I | of the vector fields i j , B i , S.

Proof.
Recall that the components of the energy momentum tensor take the form (5.1). Note that h(t, x)), for any function η(t, x, p ) and Therefore, for |I | = 1 and Z = i j , B i , S, h(t, x)) h(t, x)) and the proof follows from the fact that ∂ p l μν and ∂ h αβ μν are smooth functions of p and h. The proof for |I | ≥ 1 follows from a straightforward induction argument.

Proof. Using the Vlasov equation to write
The proof for |I | ≥ 2 follows from a straightforward induction argument.

so Proposition 4.24 and Corollary 4.27 imply
The proof then follows from Proposition 5.2.

Determinants and Changes of Variables
In the proof of Theorem 1.3, the change of variables (t, x, p ) → (t, x, y) will be used, where for i = 1, 2, 3, along with several other changes of variables; see the proof of Propositions 5.8 and 5.9. A first step towards controlling the determinants of these changes is contained in the following: Lemma 5.5. Suppose t ≥ t 0 , |x| ≤ ct, and (t, x, p ) ∈ supp( f ). Then, if the assumption (4.1) holds and ε is sufficiently small, and for all t 0 ≤ s ≤ t, and i, j = 1, 2, 3.

Proof. From the equations (1.17) and the estimates
it follows that d ds Integrating the second inequality backwards from s = t and using the fact that Dividing by s 1+a , integrating again backwards from s = t, and summing over i, j gives Taking ε sufficiently small then gives Inserting back into the above bound gives Integrating this bound backwards from s = t, and using the fact that ∂ X i ∂ p j (t, t, x, p ) = 0, gives where the fact that, for any function λ(s), Inserting back into the above bound then gives and inserting this into (5.4) gives the second bound of (5.2). In a similar manner, it is straightforward to show that d ds

and, using the final conditions
The properties of these changes are collected in the following proposition: Proposition 5.6. For fixed t, x with t ≥ t 0 + 1, |x| ≤ ct, if the assumptions (4.1) hold and ε is sufficiently small then, for p such that (t, x, p ) ∈ supp( f ), the change of variables p → y := X (t 0 , t, x, p ) satisfies s, t, x, y) := σ s t, x, p (t, x, y)) i for i = 1, 2, 3. If t ≥ t 0 + 1 and t 0 ≤ s ≤ t 0 + 1 2 , then the change of variables (x, y) → (x, z 1 (s, t, x, y)) satisfies If t 0 + 1 2 ≤ s ≤ t then the change of variables (x, y) → (z 1 (s, t, x, y), y) satisfies Finally, if t 0 + 1 2 ≤ s ≤ t and 0 ≤ σ ≤ 1, then the change of variables (x, y) → (z 2 (σ, s, t, x, y), y) satisfies Moreover, for t ≥ t 0 , the determinant of the 6 by 6 matrix κ satisfies Proof. Setting s = t 0 in (5.2), it follows that and, if ε is sufficiently small, Since t ≥ t 0 + 1, the bound (5.5) follows. For the remaining bounds of the proposition, it is necessary to consider p i as a function of t, x, y (using the above bound and the Inverse Function Theorem) and estimate ∂ p i ∂ x j and ∂ p i ∂ y j . Clearly the matrix ∂ p i ∂ y j is the inverse of the matrix ∂ y i ∂ p j , and hence it follows from (5.10) that Setting s = t 0 in (5.3) gives and so, since ε, which, if ε is suitably small, implies and the bound (5.7) follows. Finally, and so from which the bound (5.8) follows. We now prove the bound (5.9).
We have and so, The bound (5.9) then follows from the Grönwall inequality.

L 1 and L 2 Estimates of Components of the Energy Momentum Tensor
The main part of the proof of Theorem 1.3 is contained in Propositions 5.8 and 5.9 below. The following Lemma will be used: Lemma 5.7. Suppose π(supp( f 0 )) ⊂ {|x| ≤ K }, t ≥ t 0 + 1 and the assumptions (4.1) hold. Then where Proof. Since f solves the Vlasov equation, and where the change of variables p → y = X (t 0 ) and the bound (5.5) have been used.
Proposition 5.8. Suppose π(supp( f )) ⊂ {|x| ≤ ct} and consider t ≥ t 0 + 1. If the assumptions (4.1) hold and ε is sufficiently small then, for any multi index I with |I | ≤ N − 1 and each μ, ν = 0, 1, 2, 3, where Z I is a product of |I | of the vector fields i j , B i , S.

Proof. Given any function F(t, x) it follows from Proposition 5.4 that
where the change of variables x i → z i 2 := σ s x i t + (1 − σ )X (s) i has been used, along with the fact that det ∂ x i ∂z j 2 ≤ C t s 3 when 0 ≤ σ ≤ 1, t 0 + 1 2 ≤ s ≤ t, by Proposition 5.6.
It follows that, for any F(t, x), The Proposition 5.9. Suppose π(supp( f )) ⊂ {|x| ≤ ct} and consider t ≥ t 0 + 1. If the assumptions (4.1) hold, and ε is sufficiently small, then for any multi index I with |I | ≤ N and each μ, ν = 0, 1, 2, 3, where Z I is a product of |I | of the vector fields i j , B i , S.
Proof. The proof is very similar to that of Proposition 5.8. Recall that, for any F(t, x), where the change of variables x → X (s, t, x, p ) and the bound (5.3) have been used. Clearly then The L 2 bound then follows by setting F = ∂ I T μν , and the L 1 bound follows by setting F(t, x) = χ {|x|≤ct+K } , and using the fact that supp(T μν ) ⊂ {|x| ≤ ct + K }, and χ {|x|≤ct+K } L 2 ≤ Ct Since, for any function F(t, x) and any multi index I , the vector fields Z satisfy it is clear that Propositions 5.8 and 5.9 in fact hold for t ≥ t 0 . Moreover, it is clear from (2.12) that where Z J is a product of |J | of the vector fields i j , B i , S, for · = · L 1 or · L 2 since supp(T μν ) ⊂ {|x| ≤ ct}, and so spacetime derivatives ∂ I can be included in Propositions 5.8 and 5.9. Suppose now that the assumptions of Theorem 1.3 hold. It follows from Proposition 2.1 that the support of f satisfies π(supp( f )) ⊂ {|x| ≤ ct + K } and so, lettingt = t 0 + t where t 0 = K c , it follows from Proposition 5.8, Proposition 5.9 and the above comments that the bounds of Theorem 1.3 hold with t replaced byt and the vector fields Z replaced byZ fort ≥ t 0 , where the vector fieldsZ are as in Section 2.2. The proof of Theorem 1.3 then follows from noting that and so for some constants C I J .

The Einstein Equations
The main results of this section are the following: Proposition 6.1. Suppose that N ≥ 4, ε ≤ 1. Consider a solution of the reduced Einstein equations (1.11) for t < T * such that, with E N and γ as in (1.14), the weak decay estimates and M ≤ ε (6.1) hold for all t ∈ [0, T * ] for some δ such that Then for some constant C N depending only on C N , the weak decay estimates where q ± = max{±q, 0} and q = r − t, hold for all t ∈ [0, T * ].
Note that the inverse of the metric g μν can be expressed as . Therefore H 1 will satisfy the same estimates as h 1 .
We have the following strong decay estimates from the wave coordinate condition for certain tangential components expressed in the null frame:  Then the following strong decay estimates hold: for any −1 ≤ γ < γ − 2δ, and |I | = k ≤ N − 5 there are constants c k such that In addition we have the following estimates for certain tangential components expressed in a null frame: Here all constants depend only on C N in (6.3), on N and on c, K in (6.7).
In the proof of Theorem 1.2 in Section 7, Proposition 6.1 will first be appealed to for the coupled Einstein-Vlasov system (1.2), (1.3), (1.11). As a consequence the assumptions of Proposition 2.1 will be satisfied, which in turn will ensure that the assumptions of Proposition 6.3 and hence of Theorem 6.4 are satisfied.

Weak L ∞ Decay Estimates
Here we assume the weak energy bounds (6.1) and prove that this implies certain decay estimates.  (1.14). Then For a proof of this, see Proposition 14.1 in [37]. Using this we get Proposition 6.6. Suppose that the weak energy bounds (6.1) hold. Then (6.10) Furthermore, (6.11) The same estimates hold for H 1 in place of h 1 , and for h or H in place of h 1 if γ is replaced by δ.
That H 1 satisfy the same estimates as h 1 follows from the discussion after (6.4). That h and respectively, H , only satisfy these estimates with γ replaced with δ follows from the fact that h 0 and H 0 , given by (1.13) and (6.4), respectively, only satisfy these estimates.

The Improved Weak Decay Estimates for the Metric
To get improved decay estimates in the interior we will use Hörmander's L 1 -L ∞ estimates for the fundamental solution of , see Theorem 3.5 in [31].
For the proof below we will also use the following version of Hardy's inequality, see Corollary 13.3 in [37]: Lemma 6.9. For any −1 ≤ a ≤ 1 and any φ ∈ C ∞ 0 (R 3 ), (6.14) On the other hand we can also estimate derivatives in terms of the vector fields as follows: Lemma 6.10. Let ∂ q = (∂ r − ∂ t )/2 and let∂ μ = ∂ μ − L μ ∂ q be the projection of ∂ μ onto the tangent space of the outgoing light cones. Then |∂φ(t, x)| |∂ q φ(t, x)| + |∂φ(t, x)| (6.15) and where q ± = max{±q, 0} and q = r − t. Moreover, The same estimates hold for H 1 in place of h 1 , and h or H in place of h 1 if γ is replaced by 2δ.
Proof. First (6.18) is a consequence of (6.17) using (6.16) so it only remains to prove (6.17) for r < t. Let h 1 (6.19) and We will prove (6.17) for r < t separately for each of v, u, φ. For v and u the proof follows the proof in section 16 of [37]. The estimate for the homogeneous linear part v follows from using Lemma 6.8 with the estimate for v 0 and v 1 obtained from (6.10)-(6.11) when t = 0. Then using the L ∞ bounds (6.10)-(6.11) for a small number of vector fields we have (6.25) We write h = h 0 + h 1 and estimate (6.26) and by Lemma 6.9, where H (r < 3t/4) = 1 when r < 3t/4 and 0 otherwise. Hence It now follows from Lemma 6.7 that which proves (6.17) for r < t and also for u. It remains to prove the estimate for φ as well, but this also follows from Lemma 6.7: (6.32)

The Support and Weak Decay of Matter
The following Sobolev inequality will be used to obtain pointwise bounds for T from the assumptions (6.1): Proof. The proof proceeds by noting that, with c). The inequality for t ≤ 2K /(1 − c) follows from the standard Sobolev inequality, The proof for t ≥ 2K /(1 − c) follows from the identity (2.12).
Lemma 6.12, together with assumption (6.7) on the support of T and the weak energy bounds (6.1), gives

The Sharp Decay Estimates for the First Order Derivatives
Throughout the rest of this section we will assume that the weak decay estimates (6.17)-(6.18) hold for some 0 < 8δ < γ < 1 − 8δ, M ≤ ε ≤ 1, along with the support condition (6.7) for T . However we will not use the weak energy bounds (6.1) any further.

The Sharp Decay Estimates for First Order Derivative of Certain
(6.35) We first express the divergence in a null frame as follows: Lemma 6.13. Let ∂ q = (∂ r − ∂ t )/2 and ∂ s = (∂ r + ∂ t )/2. Then for any tensor k μν , Proof. The proof follows expressing the divergence in a null frame ∂ μ F μ = L μ ∂ q F μ − L μ ∂ s F μ + A μ ∂ A F μ ; for when ∂ q and ∂ s commute with the frame, see Lemma 1 in [34].
Using this, we get Lemma 6.14. We have where χ (s), is a function supported when 1/4 ≤ s ≤ 1/2. Moreover (6.37) also holds for h in place of H .
Proof. It follows from the previous lemmas that Picking U = T and, respectively, U = L, gives (6.38).
The same estimates hold for H in place of H 1 if γ is replaced by 2δ, Proof. The proof follows as in the proof of Proposition 13 in [34]. (6.40) follows from integrating (6.39) in the t − r direction from initial data. When |t − r | > t/8 the estimates follow from Proposition 6.11 so we may assume that |t − r | < t/8. It then follows from Lemma 6.14 and Proposition 6.11 that

The Leading Order Behaviour of the Inhomogeneous Term Towards
Null Infinity The inhomogeneous term in Einstein's equations can be written as where G μν (h)(∂h, ∂h) is cubic: and Q μν (∂h, ∂h) satisfy the standard null condition, and hence |Q(∂h, ∂k)| |∂h| |∂k| + |∂h| |∂k|. (6.43) The main term P(∂ μ h, ∂ ν h) can be further analyzed by first noting that |∂h| |∂k| + |∂h| |∂k|, (6.44) which follows from expressing ∂ μ in a null frame: Expressing P(h, k) = P N (h, k) in a null frame we have Taking into account the wave coordinate condition this reduces in leading order to In fact, by (6.45) we have

The Leading Order of the Geometric Wave Operator Towards Null Infinity
Expanding in a null frame as in the proof of Lemma 6.13 and using (6.16), we get Lemma 6.18. We have As a consequence, we get Lemma 6.19. We have where the asymptotic Schwarzschild wave operator is given by Proof. We apply (6.52) to H αβ 1 ∂ α ∂ β φ using (6.40) and (6.17) to get In spherical coordinates, (6.54) takes the form It would, however, perhaps have been more natural to define it to be a solution of the homogeneous wave h 0 = 0 (or even better 0 h 0 = 0.) with data coinciding with this function at time 0 in which case h 0 μν = χ(r − t) M r δ μν which is equal to (6.56) in the exterior when r ≥ t + 1. We therefore think of (6.56) as an approximate solution to the homogeneous wave equation. By (6.55), we have for some functions χ i (s) that vanish when s ≥ 1/2 or s ≤ 1/4; this, in particular, means that in the exterior and in the wave zone it is a solution of the wave operator 0 . Moreover, for some smooth function χ αβ (s, ω) supported when s ≥ 1/4. Hence

The Sharp Decay Estimates for First Order Derivatives from the Wave Equation
We will now derive sharp estimates for the first order derivatives. Following [31,37], we have Lemma 6.21. Let D t = {(t, x); t − |x| ≤ c 0 t}, for some constant 0 < c 0 < 1 and let w(q) > 0 be an increasing positive weight w (q) ≥ 0. Then Using (6.55), we get (6.64) Integrating this along the flow lines of the vector field ∂ s − M s ∂ q from the boundary of D = ∪ τ ≥0 D τ to any point inside D, using that w is decreasing along the flow lines, gives that, for any (t, x) ∈ D, The lemma now follows from (6.15), since it is trivially true when |r − t| ≥ c 0 t by (6.16).
From Lemma 6.21 and the estimate (6.17), we get Then (6.66) Using Lemma 6.22 and Proposition 6.20, we obtain Proposition 6.23. If the weak energy bounds and initial bounds hold then we have, for any 0 ≤ γ < γ − 4δ, The same estimates hold for h in place of h 1 if γ = 0.
Proof. We want to apply Lemma 6.22 to the decomposition in Proposition 6.20.
To prove (6.67) we note that / c). From the preceeding lemmas it follows that all the terms in the right hand side of (6.61) are bounded independently of t by a constant times ε when (U, V ) = (U, T ) and this proves (6.67). To prove (6.68) we note that the only new term is / P μν , which is controlled by by the first part and multiplying by (1 + t) and integrating gives a logarithm.

The Commutators and Lie Derivatives
We will use Lie derivatives which will simplify the commutators very much by removing the lower order terms. It was first observed in [34] that one can get bounds from the wave coordinate condition for Lie derivatives. Here we take it further and observe that the Lie derivative unlike vector fields preserve the geometric null structure of not only the wave coordinate condition, but also of the nonlinear inhomogeneous terms of Einstein's equations and the commutators with the geometric wave operator.

Modified Lie Derivatives Applied to the Equations
The Lie derivative applied to a (r, s) tensor K is defined by The equality (6.71) follows directly since ∂ x α Z β is constant for each of the vector fields Z . The equality (6.72) follows directly from (6.71).
Since, for an (r, s) tensor K , the quantity ∂ μ 1 · · · ∂ μ k K α 1 ...α r β 1 ...β s appearing in (6.71) (similar quantities also appear below) is not a geometric object, its Lie derivative is defined formally in the {x μ } coordinate system, using the coordinate expression (6.70). Alternatively, one could note that, in the {x μ } coordinate system, where D denotes the connection of the Minkowski metric, since the Christoffel symbols of D with respect to the Cartesian coordinate system {x μ } vanish, D ∂ x α ∂ x β = 0. One could then give a geometric proof of Proposition 6.24 using the fact that the curvature tensor of D vanishes and D 2 Z = 0 for each of the vector fields Z . Let the modified Lie derivative be defined by With this definition L Z m αβ = 0 and L Z m αβ = 0 for the vector fields in our collection, as the modified Lie derivative is defined so it commutes with contractions with the Minkowski metric. Let h αβ and k αβ be (0, 2) tensors and let S μν (∂h, ∂k) be a (0, 2) tensor which is a quadratic form in the (0, 3) tensors ∂h and ∂k with two contractions with the Minkowski metric (in particular P(∂ μ h, ∂ ν h) or Q μν (∂h, ∂k)). Then Moreover, Let L I Z be a product of |I | Lie derivatives with respect to |I | vector fields Z . It follows that and  (H, ∂h), (6.80) it follows that where χ (s), is a function supported when 1/4 ≤ s ≤ 1/2. Moreover, (6.84) also holds for h in place of H .

L ∞ estimates from the Wave Coordinate Condition
For low derivatives, (6.85) leads to the following: The same estimates hold for H in place of H 1 if γ is replaced by 2δ.
The proof is the same as for Proposition 6.15.   and by (6.57), we have

Estimates for the Inhomogeneous Term
where χ (s) is supported in 1/4 ≤ s ≤ 1.
6.3.6. Estimates of the Wave Commutator Term By (6.52), we have Writing H = H 0 + H 1 , this can be divided up in the commutator with 0 = + H αβ 0 ∂ α ∂ β and with Similarly, by (6.87) and (6.17), and we conclude that (6.95)

The Sharp L ∞ Decay Estimates for Higher Order Low Derivatives
As in section 10 of [37], using the methods in Section 6.2 we can also inductively prove sharp decay estimates for higher order low derivatives. As we have already proven the higher order weak decay estimates in Proposition 6.11 and the higher order sharp decay estimates for components we control with the wave coordinate condition in Proposition 6.26, it only remains to generalize Proposition 6.23 to higher order. In order to do that we will, as before, rely on the crucial Lemma 6.21 to control transversal derivatives in terms of tangential derivatives, which we control by Proposition 6.11, and 0 close to the light cone |t − r | < 1 − c. It therefore only remains to get control of 0 L I Z h 1 μν close to the light cone |t − r | < (1 − c)t, where 0 < c < 1. When |t − r | < (1 − c)t, that we have, by (6.57) and (6.58), where L I Z ( g − 0 )h 0 μν is controlled by (6.92) using Proposition 6.11 as follows: (6.97) and by (6.91) as It remains to estimate the difference g − 0 and the commutators By (6.94) and by, (6.100) μν , which can be estimated in the same way, we obtain Summing up, (6.102) From Lemma 6.21 and the estimate (6.17), we get Then Proposition 6.28. If the weak energy bounds and initial bounds hold, then we have, for any 0 ≤ γ < γ − 2δ and |I | = k ≤ N − 5, that there are constants c k such that The same estimates hold for h in place of h 1 if γ = 0.
. We will prove (6.104) by induction, noting that it is true for k = 0 by (6.68). Then, by Lemma 6.27, we have, for k ≥ 1, where the bounds (6.33) have been used. By the induction hypothesis, Using Grönwall's lemma with G denoting the integral we get G (t) ≤ ε(1 + t) −1 c k G(t) + ε(1 + t) 2c k−1 ε and multiplying with the integrating factor we get Assuming that c k ≥ 4c k−1 , we get G(t) ≤ c k ε(1 + t) c k ε , and hence N k (t) ≤ c k ε(1 + t) c k ε .

The Basic Energy Estimate for the Wave Equation
In [37] (see Proposition 6.2 there), the following energy estimate was proven: Lemma 6.29. Let φ be a solution of the wave equation g φ = F, with the metric g such that, for H αβ = g αβ − m αβ , for some μ > 0. Set Then, for any 0 < γ 1 and 0 < ε ≤ γ /C 1 , we have (6.108)

The Lowest Order Energy Estimate for Einstein's Equations Let
(6.109) By Lemma 6.29, we have for j = 0, 1. We have (6.113) and As far as the energy estimate is concerned, one could have picked h 0 μν to satisfy g h 0 μν = 0 and wouldn't have to do anything further. However, since we didn't do this, we will estimate using (6.58) and (6.57) to set (1 + t + r ) 4 +

Higher Order L 2 Energy Estimates
For this section we have to make the following smallness assumption on ε: c k ε ≤ δ, (6.118) where c k are the constants in Proposition 6.28.

Equivalence of Norms
The inhomogeneous terms contain factors of ∂ L I Z h which we estimate by writing h = h 0 +h 1 , and estimate the L 2 norm factors with h 1 in terms of the energy of h 1 whereas the L 2 norms of h 0 can be estimated directly. The commutator terms will in addition contain factors of L I Z H , where we can also write H = H 0 + H 1 and estimate the factors with H 0 directly since it is explicit, and for the factors with L I Z H 1 we first use Hardy's inequality to estimate them in terms of ∂ L I Z H 1 . However H 1 is only approximately equal to −h 1 . We have that H = −h + K (h), where K (h) = O(h 2 ) and hence H 1 = −h 1 + K (h) − h 0 − H 0 . Differentiating, we see that to conclude that the norms of H 1 are approximately bounded by those of h 1 we have to estimate factors of the form L J Z h ∂ L K Z h, with |J | + |K | ≤ |I |, in L 2 with respect to the measure w. Again this can be estimated by writing h = h 0 + h 1 and estimating the factors with h 1 in terms of the energy (after possibly using Hardy's inequality) and estimating the explicit factors with h 0 directly. The conclusion of this process is that and similarly for the space-time integrals of tangential components.
6.6.3. L 2 Estimate of the Wave Operator Applied to h 0 By (6.92), using Hardy's inequality we have It follows that Summing up, we have C ε (1 + τ ) 3/2 + C ε 2 (1+τ ) 2−γ −2δ M E k (τ ) 1/2 dτ. 6.6.5. Higher Order Energy L 2 Estimates Here we give the proof of Theorem 6.4 using the decay estimates proven in the previous section. We will argue by induction so we assume the energy estimate is true for k − 1 and we will prove it for k. Using the energy inequality Lemma 6.29, we get, from adding up the energy contributions from the inhomogeneous term (6.130), the commutator with the wave equation (6.143) and the Vlasov matter E k (τ ) 1/2 dτ for some universal constant M k . We now choose ε so small that 32ε ≤ 1 so that S k (t) in the right can be absorbed into S k (t) on the left, and so that by c k ε ≤ 2δ and C ε ≤ M k , we obtain

The Continuity Argument and the Proof of Theorem 1.2
The proof of Theorem 1.2 is a direct consequence of Theorem 1.3, Propositions 2.1, 6.1, 6.3 and Theorem 6.4, and the following local existence theorem for the reduced Einstein-Vlasov system (the number of derivatives used in the following local existence theorem is far from sharp): for a given time t, define See the work of Choquet-Bruhat [10], and also the textbook of Ringström [42], for related local existence theorems.
The proof of Theorem 1.2 proceeds as follows: let T * be the supremum of all times T 1 such that a solution of the reduced Einstein-Vlasov system (1.2), (1.3), (1.11) attaining the given data exists for all t ∈ [0, T 1 ] and satisfies E N (t) 1 2 ≤ C N ε(1 + t) δ , |I |≤N −1 Z I T (t, ·) L 1 ≤ C N ε (7.2) for all t ∈ [0, T 1 ], where δ > 0 is such that δ < γ < 1 − 8δ and C N is a fixed large constant, to be determined, depending only on N , δ and supp( f 0 ). Recall that 0 < γ < 1 is fixed in the statement of Theorem 1.2. Recall that T μν = T μν − 1 2 trT g μν . Clearly the set of such T 1 is non-empty, by Theorem 7.1, and so T * > 0. Suppose T * < ∞.

Moreover Theorem 1.3 implies that
where the constant C is independent of C N and the constant D N depends on C N . Inserting the above bound for Q N then gives |I |≤N −1 Z I T (t, ·) L 1 Cε + D N ε 2 , provided ε is sufficiently small, for some new C, D N , as above. Now, as above, where the equality (2.12) was used, it follows that and so, where D N depends on C N and C does not. It follows from the bound (7.6) with k = N and the bound (7.7), provided the constant C N is chosen so that C N ≥ max{2(M 0 + · · · + M N ), 4C} and ε is chosen so that ε < min{ δ 2 (d 1 + · · · + d N + (N + 1)C N ) −1 , C N 4D N }, that the bounds Publisher's Note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.