Florian Simatos
email: f.simatos@tue.nl.

STATE SPACE COLLAPSE FOR CRITICAL MULTISTAGE EPIDEMICS

We study a multistage epidemic model which generalizes the SIR model and where infected individuals go through K ≥ 1 stages of the epidemic before being removed. An infected individual in stage k ∈ {1, . . . , K } may infect a susceptible individual, who directly goes to stage k of the epidemic; or it may go to the next stage k+1 of the epidemic. For this model, we identify the critical regime in which we establish diffusion approximations. Surprisingly, the limiting diffusion exhibits an unusual form of state space collapse which we analyze in detail. CONTENTS 1. Introduction 1 2. Main results 6 3. Analysis of the SDE (2.4) 10 4. Scaling limits 12 5. Asymptotic behavior of the outbreak size 18 6. An intriguing conjecture 21 Appendix A. Proof of Lemmas 3.1 and 3.

INTRODUCTION

A multistage model. In this paper we study an epidemic model proposed by Antal and Krapivsky [START_REF] Antal | Outbreak size distributions in epidemics with multiple stages[END_REF] which generalizes the classical SIR model. Similarly as for the SIR model, there is a closed population where each individual is either susceptible, infected or removed. In addition, we assume that the disease which is spread progresses through K ≥ 1 stages. For short, we will say that an (infected) individual is in stage k ∈ {1, . . . , K } if it possesses the disease in stage k. It will also be convenient to say that an individual is in stage 0 if it is susceptible and in stage K +1 if it is removed.

Two types of transitions may occur within our model: either (1) an infected individual in stage k ∈ {1, . . . , K } tries to infect a susceptible individual; or [START_REF] Antal | Outbreak size distributions in epidemics with multiple stages[END_REF] an infected individual in stage k ∈ {1, . . . , K } progresses to the next stage k+1. We consider the mean-field regime where in case (1) above, the infected individual samples an individual uniformly at random from the whole population: if the sampled individual is susceptible, then it becomes infected and starts the infection at the same stage than the individual who infected him or her.

In view of its epidemiological interpretation and because it generalizes the SIR model (which corresponds to the single-stage case K =1), it is natural to interpret this model as an epidemic model. From that viewpoint, our model differs from multistage epidemic models previously proposed. For HIV/AIDS for instance, Hyman et al. [START_REF] Hyman | The differential infectivity and staged progression models for the transmission of HIV[END_REF] proposed a multistage (or staged-progression) model where a newly infected individual starts the epidemic in stage one. We believe that the techniques used in the present paper can be adapted to study this case as well.

Besides its natural epidemiological interpretation, this model may also bear insight into other areas of applied probability. We may for instance think of cell population dynamic, e.g., when a cancer tumor progresses, cells accumulate mutations that are passed on to newly infected cells. In this case the stage of an individual corresponds to the number of mutations. Another potential application of this model is in communication networks: in this case we may think of users accumulating information and, upon meeting with a user with no information, sharing all the information the first user has. Such kind of epidemic spreading of information is actually at the heart of modern peer-to-peer systems.

Scaling limits. Epidemic models can broadly be categorized into discrete and continuous ones, see for instance Kendall [START_REF] Kendall | Deterministic and stochastic epidemics in closed populations[END_REF]. Discrete models are individual-based models that focus on the inherent stochasticity of the dynamic, and are often described by finite state space Markov processes. Continuous models correspond to large population approximations: they may be either deterministic (corresponding to some notion of averaging in the large population regime) or stochastic. In the first case, the continuous model is typically described by ordinary differential equations (ODE) and in the second case, by stochastic differential equations (SDE).

In the present paper, we bridge these two viewpoints by starting from a discrete model which, after a suitable renormalization procedure, leads to a continuous model. This is a classical approach which allows to understand the extent to which a continuous model, usually more tractable, suitably approximates a discrete one. In addition to giving insight into how the epidemic develops over time, such schemes are also useful to understand the final outbreak size of the epidemic. In the SIR case for instance, this makes it possible to express the final outbreak size, properly renormalized, as the hitting time of 0 by a Brownian motion with a parabolic drift; this allows in turn for explicit computation, see for instance Martin-Löf [START_REF] Martin-Löf | The final size of a nearly critical epidemic, and the first passage time of a Wiener process to a parabolic barrier[END_REF]. Note that our forthcoming Theorem 2.3 generalizes this result to the multistage case K ≥ 1.

For the SIR model it is well-known that different scaling limits may be obtained depending on the discrete model's parameters, see for instance von Bahr and Martin-Löf [START_REF] Von | Threshold limit theorems for some epidemic processes[END_REF]. In particular, the limit may be deterministic or retain some stochasticity of the original discrete model, and the intermediate regime at which transition occurs is usually referred to as critical regime: this is the regime in which we are interested in in the present paper. In this regime, scaling limits of the SIR model have for instance been obtained by Martin-Löf et al. [START_REF] Martin-Löf | Symmetric sampling procedures, general epidemic processes and their threshold limit theorems[END_REF][START_REF] Martin-Löf | The final size of a nearly critical epidemic, and the first passage time of a Wiener process to a parabolic barrier[END_REF][START_REF] Von | Threshold limit theorems for some epidemic processes[END_REF] and by Aldous [START_REF] Aldous | Brownian excursions, critical random graphs and the multiplicative coalescent[END_REF] in the closely related context of the Erdös-Rényi random graph.

A popular method to establish such scaling limits relies on semi-group techniques, see for instance Ethier and Kurtz [START_REF] Stewart | Markov processes[END_REF]. This approach has for instance recently been followed by Dolgoarshinnykh and Lalley [START_REF] Dolgoarshinnykh | Critical scaling for the SIS stochastic epidemic[END_REF] to study the related SIS model, who also sketch a proof for the SIR model. However, because of the different boundary conditions between the SIS and the SIR models, the authors only mention that different analytical results from the theory of parabolic partial differential equations are needed in the SIR case. In the present paper we use a different and more probabilistic approach, closer to Aldous' approach [START_REF] Aldous | Brownian excursions, critical random graphs and the multiplicative coalescent[END_REF], that relies on semimartingale arguments. These arguments rely on elementary calculations made on the infinitesimal generator of the Markov process, and require almost no analytical results. This is the reason why we believe that this method is robust and can be adapted to study other models such as the HIV/AIDS model of Hyman et al. [START_REF] Hyman | The differential infectivity and staged progression models for the transmission of HIV[END_REF].

Presentation and discussion of the main results. Let us present our results in the particular case where all individual transition rates are equal to one: in general these will be allowed to be close to one. For n ≥ 1 we consider the K +2-dimensional Markov process a n (t ) = (a n,k (t ), k = 0, . . . , K + 1) that represents a multistage epidemic process in a population of size n. It is defined by an initial state a n (0) = (a n,k (0), k = 0, . . . , K + 1) satisfying the relation a n,0 (0) + • • • + a n,K +1 (0) = n and, for a = (a 0 , a 1 , . . . , a K +1 ) ∈ N K +2 , by the following transition rates:

(1.1) a -→ a -e k + e k+1 at rate a k , 1 ≤ k ≤ K , a -e 0 + e k at rate a k a 0 /n, 1 ≤ k ≤ K ,
where e 0 = (1, 0, . . . , 0), e 1 = (0, 1, 0, . . . , 0), etc. For t ≥ 0, a n,k (t ) is the number of individuals in stage k at time t , i.e., a n,0 (t ) is the number of susceptible individuals; a n,k (t ) for k = 1, . . . , K is the number of infected individuals in stage k; and a n,K +1 (t ) is the number of removed individuals. The transition rates (1.1) preserve the total population size, so that the relation a n,0 (t )

+ • • • + a n,K +1 (t ) = n is satisfied for any t ≥ 0.
These rates formalize the description of the dynamic given above: each infected individual progresses to the next stage of the epidemic at rate one; each infected individual makes an infection attempt at rate one, which is successful with probability proportional to the number of susceptible individuals (which gives the factor a 0 /n).

Consider the initial condition where a n,1 (0) = n β for some β ∈ (0, 1) and a n,k (0) = 0 for k = 2, . . . , K + 1. The asymptotic behavior of the K +1-dimensional process (a n,k , k = 1, . . . , K + 1) as n → +∞ is then governed by the precise value of β, and there are three different regimes:

Small initial condition β < 1/(K + 2): on the time scale n β , a n,k lives on the space scale n kβ and the scaling limit is governed by the following SDE:

(1.2) dA 1 (t ) = (2A 1 (t )) 1/2 dB (t ), dA k (t ) = A k-1 (t )dt for k = 2, . . . , K + 1;
Intermediate initial condition β = 1/(K + 2): on the time scale n 1/(K +2) , a n,k lives on the space scale n k/(K +2) , the outbreak size is of the order of n (K +1)/(K +2) and the scaling limit is governed by the following SDE:

(1.3)      dA 1 (t ) = -A K +1 (t )A 1 (t )dt + (2A 1 (t )) 1/2 dB (t ), dA k (t ) = (A k-1 (t ) -A K +1 (t )A k (t ))dt for k = 2, . . . , K , dA K +1 (t ) = A K (t )dt ;
Large initial condition β > 1/(K + 2): on the time scale n γ with γ = (1-β)/(K +1), a n,k lives on the space scale n β+(k-1)γ and the scaling limit is governed by the following ODE:

(1.4)

     dA 1 (t ) = -A K +1 (t )A 1 (t )dt , dA k (t ) = (A k-1 (t ) -A K +1 (t )A k (t ))dt for k = 2, . . . , K , dA K +1 (t ) = A K (t )dt .
In the above, B is a standard Brownian motion. Moreover, since the total population size is fixed, the asymptotic behavior of a n,0 is recovered from the behavior of the a n,k 's: since na n,0 = a n,1 + • • • + a n,K +1 and since in all three regimes a n,K +1 a n,k for k = 1, . . . , K , na n,0 obeys the same scaling and has the same limit than a n,K +1 .

It is also interesting to note that the intermediate regime interpolates between the small and large ones in two ways: (1) the scaling at play there can be obtained by letting β ↑ 1/(K +2) in the small regime or β ↓ 1/(K +2) in the large regime, and (2) the evolution of A 1 in (1.3) is a mixture of (1.2) and (1.4). However, β does not appear in the asymptotic dynamics (1.2)-(1.4) and so (1.3) appears somehow discontinuously.

When specializing the results of the intermediate regime to K = 1 we recover the classical result on the SIR model [START_REF] Von | Threshold limit theorems for some epidemic processes[END_REF]: starting with of the order of n 1/3 infected individuals, we end up with an outbreak size of the order of n 2/3 . When K = 2, these results show that starting with n 1/4 individuals in the first stage, the number of individuals in the second stage during the outbreak is of the order n 1/2 , while the total outbreak size is of the order of n 3/4 .

Note that for the above comparison we have referred to the work by von Bahr and Martin-Löf [START_REF] Von | Threshold limit theorems for some epidemic processes[END_REF] which actually considers a discrete-time version of the SIR model, the so-called Reed-Frost model. In this model, individuals stay infected for a deterministic amount of time and, at the end of this period, try to infect every other individual with a fixed probability. In the case K = 1 the two models are indeed equivalent, at least with respect to the final outbreak size: heuristically, it does not matter when and by whom a given individual gets infected. However, in the truly multistage case K ≥ 2 this is no longer the case: the stage of an individual that spreads the disease matters. In that respect, individuals from different stages "compete" for the pool of susceptible individuals and it is natural to consider a continuous-time model. We suspect that, compared to the case K = 1, this additional difficulty is causing the discrepancy discussed in Section 6.

In order to get insight into the scalings at play, let us consider the continuous-time Markov branching process a B = (a B k , k = 1, . . . , K + 1) with K + 1 types, whose transition rates for a ∈ N K +1 are given by

(1.5) a -→ a -e k + e k+1 at rate a k , 1 ≤ k ≤ K , a + e k at rate a k , 1 ≤ k ≤ K .
These rates are obtained from (1.1) by setting a 0 = n: this corresponds to an infinite population limit, where an infection attempt is always successful. Then a B 1 is a continuous-time Markov branching process (also known as Bellman-Harris branching process), and it can be seen that if a B 1 (0) = n β and time is sped up by n β , then it converges to an interesting object (Feller diffusion). Thus if we are to have an interesting limit for a B 1 , the correct time scale (for all processes) needs to be n β . But then, the total number individuals it ever begets is of the order of ∞ 0 a B 1 (t )dt which is of the order of n 2β . These individuals form the initial condition for a B 2 , which can thus be seen as a Bellman-Harris branching process, started at n 2β but for which time is sped up by n β . However, we know from the above discussion that in order for a B 2 to evolve on the space scale n 2β , it needs to be sped up in time by n 2β : thus on the time scale n β , a B 2 remains essentially constant. However, on this time scale it begets of the order of n β 0 a B 2 (t )dt ∝ n 3β individuals, which form the basis for a B 3 . Iterating this argument, we see the geometric progression observed in all three regimes appearing, where a n,k+1 is of the order of the time scale times a n,k . But actually, more can be learned from this simple approximation.

First of all, we see that the approximation of a n by a B is asymptotically exact in the small regime, in the sense that (1.2) is also the scaling limit of a B subject to the same scaling. This approximation begins to break down in the intermediate regime, which is therefore the regime where finite-size population effects begin to kick in. This phenomenon is also sometimes called the depletion of points effect: in the intermediate regime the epidemic begins to feel the decreasing number of susceptible individuals.

Moreover, we have just argued that we are considering a time scale which is suited for a n,1 but not for a n,2 . More precisely, on the time scale where a n,1 evolves, a n,2 has not enough time to evolve on its own and all its randomness comes from a n,1 . This problem of time scales leads to a classical phenomenon in queueing theory, called state space collapse, but which is more surprising to find in the context of epidemic.

An unusual form of state space collapse. Let K = K -1 and T n be the first time at which a n,1 hits zero, i.e., T n = inf{t ≥ 0 : a n,1 (t ) = 0}. It is then important to observe that the K +2-dimensional process ān = ( ān,k , k = 0, . . . , K + 1) obtained by keeping track of the individuals in stages 0, 2, . . . , K + 1 after time T n , i.e., with ān,0 (t ) = a n,0 (T n + t ) and ān,k (t ) = a n,k+1 (T n + t ) for k = 1, . . . , K + 1, is a multistage epidemic process with K stages. However, we know that in the intermediate regime, a n,2 , and thus ān,1 (0), is of the order of n 2/(K +2) = n 2/( K +3) . As such, this corresponds to a multistage epidemic process with a large initial condition (with β = 2/( K +3) > 1/( K +2)), and we know from the above discussion that its scaling limit is deterministic. In summary, once the first stage has got extinct, the subsequent evolution of the process is deterministic, which supports the previous claim that all the randomness lies in the first stage. This property can actually be directly read off the SDE (1.3).

Indeed, a striking feature of this SDE is that the diffusion coefficients governing A k for k = 1 are equal to zero: in particular, only the coordinate A 1 is truly random, the other ones being obtained deterministically from it via an ODE. More precisely, let F : R × R K → R K be the function such that, if Ā = (A 2 , . . . , A K +1 ), then (1.3) rewrites d Ā(t ) = F (A 1 (t ), Ā(t ))dt : then Ā satisfies an ODE of the kind (1.6) dx(t ) dt = F y(t ), x(t ) , t ≥ 0.

It will be seen that if y : [0, ∞) → [0, ∞) is continuous, then this ODE has a unique solution. This naturally defines a one-dimensional manifold M (y) ⊂ R K , namely, if x y is the unique solution to (1.6), the space M (y) = x y (t ) : t ≥ 0 .

If A = (A 1 , . . . , A K +1 ) is the solution to (1.3), then by definition Ā is the unique solution to (1.6) with y = A 1 , which justifies seeing Ā as deterministically obtained from A 1 . Thus we can see A K +1 as a function of A 1 and the initial condition ā = Ā(0), and write A K +1 (t ) = ϕ(A 1 , ā; t ) for some deterministic map ϕ. With this notation, we see that A 1 on its own satisfies the following SDE

(1.7) dA 1 (t ) = γ 1 -ϕ(A 1 , ā; t ) A 1 (t )dt + (2A 1 (t )) 1/2 dB (t ).
Such an SDE is sometimes referred to as an Itô process and is, in the terminology of Rogers and Williams [START_REF] Rogers | Diffusions, Markov processes, and martingales[END_REF], not of diffusion type since ϕ(A 1 , ā; t ) depends on A 1 through the entire path (A 1 (s), 0 ≤ s ≤ t ) and not only through the value of A 1 at time t .

In summary, we start with a sequence of K +2-dimensional stochastic processes (or K +1-dimensional, given that the total population size is fixed), but the diffusion A = (A 1 , . . . , A K +1 ) appearing in the limit is truly one-dimensional: one coordinate, A 1 , is random and given by (1.7), while the other coordinates A 2 , . . . , A K +1 are obtained deterministically from it via (1.6).

Such a phenomenon of reduction of the dimension of the state space in the critical (or near-critical) regime is well-known in queueing theory, where it is usually referred to as state space collapse. It has first been systematically investigated in [START_REF] Bramson | State space collapse with application to heavy traffic limits for multiclass queueing networks[END_REF][START_REF] Williams | Diffusion approximations for open multiclass queueing networks: sufficient conditions involving state space collapse[END_REF] in the context of multiclass queueing networks, and has since then been observed in various settings [START_REF] Borst | A stochastic network with mobile users in heavy traffic[END_REF][START_REF] Reiman | Some diffusion approximations with state space collapse[END_REF][START_REF] Stolyar | Maxweight scheduling in a generalized switch: state space collapse and workload minimization in heavy traffic[END_REF][START_REF] Verloop | Heavy-traffic analysis of a multiple-phase network with discriminatory processor sharing[END_REF]. In all these examples and similarly as what we observe in our model, the randomness of the limiting diffusion is contained in one coordinate (the workload process for queueing systems), while the other coordinates are obtained deterministically from it. In queueing theory typically, the workload process converges to a reflected Brownian motion W and the queue length processes Q are obtained via a deterministic map, e.g., Q = G • W .

However, there are two main differences between these results from queueing theory and our present form of state space collapse. First of all, in all these queueing examples, the map G is linear, i.e., we can write Q k (t ) = α k W (t ) for some deterministic coefficients α k > 0. This makes the queue length processes Q = (Q k ) live in a deterministic one-dimensional manifold, namely, the space M Q = {(α k w) : w ≥ 0}. In contrast, the coordinates A 2 , . . . , A K +1 live in our case in the random one-dimensional manifold M (A 1 ).

Moreover, the reasons leading to this state space collapse phenomenon are also distinct. As explained above, in our case this comes from a matter of time scales while in queueing theory, this essentially comes from the law of large numbers since the workload is the sum of the residual service times over the number of customers.

Organization of the paper. In the next section we introduce the general model considered in the present paper, which generalizes (1.1) to the near-critical case, and state the main results, which formalize the results described previously. In Section 3 we analyze the limiting SDE that generalizes (1.3), where we prove existence and uniqueness of solutions as well as a useful sample path property. These results rely on the analysis of the ODE that generalizes (1.6), and the deterministic results that are needed are proved in the Appendix A. Sections 4 and 5 are devoted to the proofs of the main results: Section 4 gives the proof for the scaling limits, i.e., the asymptotic behavior at the process level, and Section 5 is concerned with the proof of the asymptotic behavior of the outbreak size. Finally, we conclude the paper by discussing in Section 6 the conjecture formulated in [START_REF] Antal | Outbreak size distributions in epidemics with multiple stages[END_REF], and how it relates to our results.

Acknowledgements. I am extremely grateful to Remco van der Hofstad and Johan van

Leeuwaarden for introducing me to the topic of this paper and for very stimulating discussions on this particular model. I would also like to thank the two anonymous reviewers whose comments greatly helped to improve the exposition of the paper.

MAIN RESULTS

2.1. Notation. Let in the sequel D be the set of real-valued, càdlàg functions. For f ∈ D and ε > 0, we define

T 0 ( f ) = inf{t ≥ 0 : f (t ) = 0}, T ↑ ε ( f ) = inf{t ≥ 0 : f (t ) ≥ ε} and T ↓ ε ( f ) = inf{t ≥ 0 : f (t )
≤ ε}, with the convention inf = +∞. If t ≥ 0 and f is càdlàg (real-or vector-valued) we consider θ t the shift operator, defined by

θ t ( f ) = f (t + • ).
The space of càdlàg functions with values in R d is endowed with the topology of uniform convergence on compact sets.

In the sequel we fix some integer K ≥ 1 and we consider the map π : R K +2 → R K +1 defined by π(a 0 , . . . , a K +1 ) = (a 0 , a 2 , . . . , a K +1 ) (this unusual indexing of vectors will be convenient for our purposes).

Model and main results.

We now present the full model investigated in the rest of the paper, which generalizes (1.1) by allowing individual transition rates to be close to one. For each n ≥ 1, let a n (t ) = (a n,k (t ), k = 0, . . . , K +1) be the K +2-dimensional Markov process corresponding to a finite population of size n, i.e., a n,0 (0)

+ • • • + a n,K +1 (0) = n,
and with non-zero transition rates given for a = (a 0 , a 1 , . . . ,

a K +1 ) ∈ N K +2 by (2.1) a -→ a -e k + e k+1 at rate (1 + δ n,k )a k , 1 ≤ k ≤ K , a -e 0 + e k at rate (1 + ε n,k )a k a 0 /n, 1 ≤ k ≤ K , where δ n,k , ε n,k > -1 for k = 1, .
. . , K , and we will always assume that

δ n,k , ε n,k → 0 as n → +∞ for each k = 1, . . . , K : the critical case (1.1) is recovered by setting δ n,k = ε n,k = 0.
Moreover, it will be convenient to define

δ n,0 = δ n,K +1 = ε n,K +1 = -1.
We first state the results concerning the intermediate regime, which constitute the main results of the paper.

Theorem 2.1 (Intermediate initial condition). Let β = 1/(K +2) and assume that for each k = 1, . . . , K , there exists

γ k ∈ R such that (2.2) n β ε n,k -δ n,k -→ n→+∞ γ k . Let A n = (A n,k , k = 0, . . . , K + 1)
be the K +2-dimensional process defined as follows:

(2.3) A n,0 (t ) = n -a n,0 (n β t ) n (K +1)β and A n,k (t ) = a n,k (n β t ) n kβ , k = 1, . . . , K + 1, t ≥ 0. If A n (0) → a ∈ [0, ∞) K +2 ,

then the sequence of processes (A n , n ≥ 1) converges weakly as n → +∞ to the unique solution

A = (A k , 0 ≤ k ≤ K + 1) of the SDE (2.4) A(t ) = a + t 0 b(A(u))du + t 0 σ(A(u))dB (u),
where B is a standard K +2-dimensional Brownian motion and for a = (a 0 , . . . , a K +1 ) ∈ [0, ∞) K +2 we have σ(a) = diag(0, (2a 1 ) 1/2 , 0, . . . , 0) and

   b 0 (a) = b K +1 (a) = a K , b 1 (a) = (γ 1 -a 0 )a 1 , b k (a) = a k-1 + (γ k -a 0 )a k , k = 2, . . . , K .
Since all the diffusion coefficients but one are equal to 0, when convenient we will identify the (K +2) × (K +2) matrix σ(a) and its only non-zero entry (2a 1 ) 1/2 . For the same reason, if we write B = (B 0 , . . . , B K +1 ) then only the coordinate B 1 matters and we will identify B and B 1 , i.e., whenever convenient we will consider that B is a onedimensional standard Brownian motion; note that this is in line with the notation used in (1.2) and (1.3). We now state some properties of solutions to (2.4).

Proposition 2.2. Consider the assumptions and notation of Theorem 2.1. Then A almost surely satisfies the following properties:

(i) the two processes A 0 and A K +1 are equal; (ii) if K ≥ 2, then A 0 is strictly increasing, and its terminal value A 0 (∞) is finite and

satisfies A 0 (∞) > γ k for every k = 2, . . . , K ; (iii) T 0 (A 1 ) < +∞; (iv) for every k = 2, . . . , K , A k (t ) → 0 as t → +∞ but T 0 (A k ) = +∞.
Note that the first property of A above comes from the relation K +1 k=0 a n,k (t ) = n, which translates after scaling to K +1 k=1 n kβ A n,k (t ) = n (K +1)β A n,0 (t ). Next, we turn to the asymptotic behavior of the outbreak size A n,K +1 (∞): the following result essentially states that we can interchange the limits n → +∞ and t → +∞. Theorem 2.3. Under the assumptions and notation of Theorem 2.1, the renormalized outbreak size A n,K +1 (∞) converges weakly as n → +∞ toward A K +1 (∞).

Let us conclude the presentation of results in the intermediate regime by comparing with the SIR case K = 1. In the case K = 1 the natural notion of criticality is that every infected individual tries to infect in average (close to) one other individual. The notion of criticality that we consider here is at first sight different, since every infected individual tries to infect in average one other individual in each stage. If an individual starts in stage k, it will therefore make in average K -k+1 infection attempts. But the two notions actually coincide. Indeed, Theorem 2.1 shows that individuals in the last stage K of the epidemic dominate: there are n K /(K +2) such individuals and n k/(K +2) n K /(K +2) individuals at stage k = 1, . . . , K -1. Thus with overwhelming probability, an individual eventually infected has actually started the epidemic in stage K and has made in average one infection attempt. Beware however: this does not mean that only stage K matters, since a n,K becomes of the order of n K /(K +2) thanks to the help of the individuals in the previous stages! Moreover, according to Proposition 2.2, the case K = 1 is the only case where A K +1 (∞) is reached in finite time. Indeed, for K = 1 the dynamic is frozen after the time T 0 (A 1 ) which is finite, while if K ≥ 2 then A K +1 remains strictly increasing at all times. As will be discussed in Section 5 when proving Theorem 2.3, this difference creates an additional difficulty in order to control the asymptotic behavior of the outbreak size.

We now state the result for the large regime. Proposition 2.4 (Large initial condition). Fix some sequence n 1/(K +2)

α n,1 n and for n ≥ 1 define

(2.5) τ n = n α n,1 1/(K +1) , α n,0 = n τ n and α n,k = τ k-1 n α n,1 for k = 2, . . . , K + 1.
Assume that for each k = 1, . . . , K , there exists γ k ∈ R such that

(2.6) τ n ε n,k -δ n,k -→ n→+∞ γ k .
Let A n = (A n,k , k = 0, . . . , K + 1) be the K +2-dimensional renormalized process defined as follows:

(2.7) A n,0 (t ) = n -a n,0 (τ n t ) α n,0 and A n,k (t ) = a n,k (τ n t ) α n,k , k = 1, . . . , K + 1, t ≥ 0. If A n (0) → a ∈ [0, ∞) K +2
, then the sequence of processes (A n , n ≥ 1) converges weakly as n → +∞ to the unique solution of the ODE (2.8)

A(t ) = a + t 0 b(A(u))du, with b defined in Theorem 2.1.
Note that this result is coherent with Theorem 2.1, in the sense that the limit of

θ T 0 (A n,1 ) • π • A n = π • θ T 0 (A n,1 )
• A n is the same whether we consider this process as the process π • A n shifted at time T 0 (A n,1 ) and then use Theorem 2.1, or whether we consider this process as a multistage epidemic process with K = K -1 stages started from a large initial condition and then use Proposition 2.4.

We complete the results of Theorem 2.1 and Proposition 2.4 by studying the case of a small initial condition a n,1 (0) n 1/(K +2) . When a n,1 (0) converges to some finite number, the sequence of processes (a n,k , k = 1, . . . , K + 1) converges to the multitype branching process given by (1.5). As the next result shows, such a branching approximation continues to be valid in the small regime. Proposition 2.5 (Small initial condition). Fix some sequence 1 α n,1 n 1/(K +2) and assume that for each k = 1, . . . , K , there exists γ k ∈ R such that

(2.9) α n,1 ε n,k -δ n,k -→ n→+∞ γ k . Let A n = (A n,k , k = 0, . . . , K +1)
for n ≥ 1 be the K +2-dimensional renormalized process defined as follows:

(2.10) A n,0 (t ) = n -a n,0 (α n,1 t ) α K +1 n,1 and A n,k (t ) = a n,k (α n,1 t ) α k n,1 , k = 1, . . . , K + 1, t ≥ 0. If A n (0) → a ∈ [0, ∞) K +2
, then the sequence of processes (A n , n ≥ 1) converges weakly as n → +∞ to the unique solution of the SDE (2.11)

A(t ) = a + t 0 b S (A(u))du + t 0 σ(A(u))dB (u),
where σ is as in Theorem 2.1 and b S is given by

     b S 0 (a) = b S K +1 (a) = a K , b S 1 (a) = γ 1 a 1 , b S k (a) = a k-1 + γ k a k , k = 2, . . . , K .
As discussed in the introduction for the strictly critical case δ n,k = ε n,k = 0, the limiting diffusion (2.11) obtained in the small regime is the same as the limit of the branching process corresponding to the infinite population setting, i.e., where a 0 /n = 1 in (2.1) (this could be proved using the techniques of the present paper). Thus a n,1 (0) ∝ n 1/(K +2) is the threshold at which finite-size population effects (or depletion of points effects) begin to kick in: this is the threshold at which the branching approximation ceases to be valid.

Note also that for A = (A k , k = 0, . . . , K + 1) satisfying (2.11), A 1 is Feller diffusion with drift γ 1 , i.e., is the unique solution to the SDE dA 1 = γ 1 A 1 dt + (2A 1 ) 1/2 dB (see (3.2) below), and then A k for k = 2, . . . , K + 1 is given recursively by A k = γ k A k + A k-1 (here and in the sequel, prime denotes differentiation with respect to the time variable). In particular, existence and uniqueness of solutions to (2.11) follow immediately.

The next three sections are devoted to proving these results. In the next section we prove uniqueness and existence to solutions of (2.4), and we also prove Proposition 2.2. The scaling limits results of Theorem 2.1 and Propositions 2.4 and 2.5 are proved in Section 4, and the proof of Theorem 2.3 on the asymptotic behavior of the outbreak size is given in Section 5.

ANALYSIS OF THE SDE (2.4)

As mentioned in the introduction, the fact that σ(a) = diag(0, (2a 1 ) 1/2 , 0, . . . , 0) is a manifestation of the state space collapse property, and it makes the process π • A deterministically obtained from A 1 by an ODE. More precisely, we consider in the sequel F : R × R K +1 → R K +1 the function defined by F (a 1 , π(a)) = π(b(a)) for a ∈ R K +2 and with b as in Theorem 2.1 (when all parameters γ k = 0, this coincides with the function F in (1.6)). With this notation, we see that if A satisfies (2.4), then π • A is a solution of the ODE

(3.1) dx(t ) dt = F (y(t ), x(t )), x(0) = ā, with y = A 1 .
The following properties of this ODE will be needed: their proof is postponed to the Appendix A. Then the function x 0 is non-decreasing and its limit x 0 (∞) as t → +∞ exists in [0, ∞]. Moreover:

(i) if y(0) > 0 or a 2 > 0, then x k (t ) > 0 for every k = 0, 2, . . . , K and t > 0;

(ii) if x k (t ) > 0 for every k = 2, . . . , K and t > 0 and x 0 is bounded, then ∞ 0 y < +∞ and moreover ∞ 0 x k < +∞ and x 0 (∞) > γ k for every k = 2, . . . , K ; (iii) if y(t ) = 0 for all t ≥ 0 and a 2 > 0, then x 0 is bounded and x k (t ) → 0 as t → +∞ for every k = 2, . . . , K .

We now turn to the study of the SDE (2.4). In the sequel we will call Feller diffusion with drift γ ∈ R the unique solution to the SDE (3.2)

Z (t ) = Z (0) + γ t 0 Z (u)du + t 0 (2Z (u)) 1/2 dB (u).
It is well-known that there is a unique strong solution to (3.2). If Z is this solution, then it does not explode in finite time, P(T 0 (Z ) < +∞) > 0 and P(T 0 (Z ) [START_REF] Rogers | Diffusions, Markov processes, and martingales[END_REF] implies that A 1 (t ) ≤ Z (t ) for all t ≥ 0 almost surely, where Z is Feller diffusion with drift γ 1 started at Z (0) ≥ A 1 (0). This comparison argument will be used several times. The proof of the following result uses standard arguments, and we only sketch the proof.

< +∞) = 1 if γ ≤ 0. Moreover, since A 1 satisfies A 1 (t ) = A 1 (0) + t 0 b 1 (A(u))du + t 0 (2A 1 (u)) 1/2 dB (u) with b 1 (a) = (γ 1 -a 0 )a 1 ≤ γ 1 a 1 , Theorem V.43.1

in Rogers and Williams

Lemma 3.3. Uniqueness in law holds for the SDE (2.4).

Proof. The problem to invoke classical results is that some of the coefficients b k grow quadratically. However, standard localization and change of drift arguments (to go back to the case of Feller diffusion, since σ(a) = diag(0, (2a 1 ) 1/2 , 0, . . . , 0)) show that for each N ≥ 1, the law of a solution A to (2.4) stopped at inf{t ≥ 0 : A(t ) ≥ N } is uniquely determined (say, with a = a 0 + • • • + a K +1 ). By successively patching up these solutions, we obtain uniqueness to (2.4) until the time of explosion, and so it only remains to show that solutions to (2.4) Before proving this result, we first explain briefly how it yields Proposition 2.2 in combination with Lemma 3.2.

Proof of Proposition 2.2. The fact that T 0 (A 1 ) is almost surely finite is precisely the content of Lemma 3.4. Then by shifting A at time T 0 (A 1 ), we see that θ T 0 (A 1 ) • A satisfies the ODE (3.1) with y = 0, and so the results of Lemma 3.2 are precisely those that we need to prove for Proposition 2.2.

Proof of Lemma 3.4. To prove P(T 0 (A 1 ) < +∞) = 1, it is enough to prove that

(3.3) P T 0 (A 1 ) < +∞, T ↑ γ 1 (A 0 ) < +∞ = P T ↑ γ 1 (A 0 ) < +∞ and (3.4) P T 0 (A 1 ) = +∞, T ↑ γ 1 (A 0 ) = +∞ = 0
which we do now.

In the event {T ↑ γ 1 (A 0 ) < +∞}, we have dA

↑ (t ) = b(A ↑ (t ))dt + σ(A ↑ (t ))dB ↑ (t )
, where A ↑ and B ↑ are the processes A and B shifted at time T ↑ γ 1 (A 0 ), so that, since T ↑ γ 1 (A 0 ) is a stopping time, B ↑ is a Brownian motion according to the strong Markov property. Moreover, we have by definition b 1 (A ↑ (t )) = (γ 1 -A ↑ 0 (t ))A ↑ 1 (t ) ≤ 0 and so A ↑ 1 is dominated by Feller diffusion without drift (see, e.g., Theorem V.43.1 in Rogers and Williams [START_REF] Rogers | Diffusions, Markov processes, and martingales[END_REF]), which almost surely hits 0 in finite time. This proves (3.3).

We now prove (3.4). Assume that T ↑ γ 1 (A 0 ) = +∞: then A 0 is bounded (by γ 1 , and in particular γ 1 > 0) and according to (i)+(ii) of Lemma 3.2, we obtain that ∞ 0 A 1 is finite. In particular L = 0, where we define L = lim inf t →∞ A 1 (t ), and so

P T 0 (A 1 ) = +∞, T ↑ γ 1 (A 0 ) = +∞ ≤ P (T 0 (A 1 ) = +∞, L = 0) .
Thus to prove (3.4), we only have to show that P(T 0 (A 1 ) = +∞, L = 0) = 0. Let ϑ i be defined recursively by ϑ 0 = 0 and ϑ i +1 = inf{t ≥ 1 + ϑ i : A 1 (t ) ≤ 1}. Then in {L = 0}, ϑ i is almost surely finite for every i ≥ 0 and so we can define Z i as the solution to the SDE

dZ i = γ 1 Z i dt + (2Z i ) 1/2 dB i with initial condition Z i (0) = 1
, where B i is the process B shifted at time ϑ i . Note that conditionally on {ϑ i < +∞}, Z i is a Feller diffusion with drift γ 1 > 0 started at 1, and that the strong Markov property and the comparison theorem V.43.1 in Rogers and Williams [START_REF] Rogers | Diffusions, Markov processes, and martingales[END_REF] show that A 1 (t + ϑ i ) ≤ Z i (t ) for t ≥ 0 and i ≥ 1. In particular,

P (T 0 (A 1 ) = +∞, L = 0) ≤ P (∀i ≥ 1 : ϑ i < +∞ and T 0 (Z i ) = +∞)
and so we only have to show that this last probability is equal to 0. For I ≥ 1 and in the event {ϑ I < +∞}, the strong Markov property at time ϑ I gives P ϑ I < +∞ and T 0 (Z i ) = +∞ for i = 1, . . . , I = P T 0 (Z 1 ) = +∞ P ϑ I < +∞ and T 0 (Z i ) = +∞ for i = 1, . . . , I -1 ≤ P T 0 (Z 1 ) = +∞ P ϑ I -1 < +∞ and T 0 (Z i ) = +∞ for i = 1, . . . , I -1 and so we obtain by induction P ϑ I < +∞ and T 0 (Z i ) = +∞ for i = 1, . . . , I ≤ P T 0 (Z 1 ) = +∞ I .

Since P T 0 (Z 1 ) = +∞ < 1, letting I → +∞ achieves the proof of (3.4).

SCALING LIMITS

We now prove the convergence results of Theorem 2.1 and Propositions 2.4 and 2.5. The proofs of these three results can be cast into the same framework by defining A n as in (2.7) with:

Intermediate initial condition (Theorem 2.1): α

n,1 = n 1/(K +2) , τ n = α n,1 , α n,0 = α n,K +1
and α n,k = α k n,1 for k = 2, . . . , K + 1; Large initial condition (Proposition 2.4): n 1/(K +2) α n,1 n and τ n and α n,k for k = 2, . . . , K + 1 are as in (2.5); Small initial condition (Proposition 2.5):

1 α n,1 n 1/(K +2) , τ n = α n,1 , α n,0 = α n,K +1
and α n,k = α k n,1 for k = 2, . . . , K + 1; and by assuming, with this notation, that τ n (ε n,k -δ n,k ) → γ k , which is consistent with the assumptions (2.2), (2.6) and (2.9).

In order to prove Theorem 2.1 and Propositions 2.4 and 2.5 with this notation, we have to prove that the sequence (A n , n ≥ 1) converges weakly toward: the solution of (2.4) in the intermediate regime; the solution of (2.8) in the large regime; and the solution of (2.11) in the small regime. In the sequel we will use the fact that in all three regimes, the algebraic relations α n,0 = α n,K +1 and α n,k = τ n α n,k-1 for k = 2, . . . , K + 1 hold.

The proof relies on the standard machinery: we first prove tightness and then identify accumulation points. Both steps rely on semimartingale arguments based on the explicit form of the generator of A n . Indeed, A n is by definition a Markov process with generator Ω n given, for any function f : R K +2 → R and any a ∈ R K +2 , by

(4.1) Ω n ( f )(a) = τ n K k=1 f a - e k α n,k + e k+1 α n,k+1 -f (a) (1 + δ n,k )α n,k a k + τ n 1 -α n,0 a 0 /n K k=1 f a + e 0 α n,0 + e k α n,k -f (a) (1 + ε n,k )α n,k a k
(note that a n,0 /n = 1 -α n,0 A n,0 /n, which gives the factor 1 -α n,0 a 0 /n in the above expression). Since A n lives on a finite state space, for any function f the process

M f n (t ) = f (A n (t )) -f (A n (0)) -V f n (t ) with V f n (t ) = t 0 Ω n ( f )(A n (u))du
is a martingale whose quadratic variation process 〈M

f n 〉 is equal to 〈M f n 〉(t ) = t 0 Θ n ( f )(A n (u))du with Θ n ( f ) = Ω n ( f 2 ) -2 f Ω n ( f ),
see for instance Lemma VII.3.68 in Jacod and Shiryaev [START_REF] Jacod | Grundlehren der Mathematischen Wissenschaften[END_REF]. Recall that δ n,0 = δ n,K +1 = ε n,K +1 = -1, and let in the rest of this section π k for k = 0, . . . , K + 1 be the projection on the kth coordinate, i.e., π k (a) = a k for any a = (a 0 , . . . , a K +1 ) ∈ R K +2 .

Lemma 4.1. For each n ≥ 1 and t ≥ 0, (4.2)

V π 0 n (t ) = τ n α n,0 K k=1 (1 + ε n,k )α n,k t 0 (1 -α n,0 A n,0 (u)/n)A n,k (u)du, (4.3) V π k n (t ) = t 0 (1 + δ n,k-1 )A n,k-1 (u) + τ n (ε n,k -δ n,k )A n,k (u) - τ n α n,0 n (1 + ε n,k )A n,0 (u)A n,k (u) du for k = 1, . . . , K + 1 and (4.4) 〈M π k n 〉(t ) = 1 α n,k V π k n (t ) + 2 τ n α n,k (1 + δ n,k ) t 0 A n,k (u)du for k = 0, . . . , K + 1. Proof. Recall that V π k n (t ) = t 0 Ω n (π k )(A n (u))du and that 〈M π k n 〉(t ) = t 0 Θ n (π k )(A n (u)
)du, so that we only have to compute Ω n (π k ) and Θ n (π k ) for k = 0, . . . , K + 1. By writing

Ω n ( f 2 )(a) -2 f (a)Ω n ( f )(a) in the form j (x 2 j -y 2 )β j -2y j (x j -y)β j = j (x 2 j -y 2 ) -2y(x j -y) β j = j (x j -y) 2 β j
for some x j , β j and y, we see that Θ n ( f ) can alternatively be written as follows:

(4.5) Θ n ( f )(a) = τ n K k=1 f a - e k α n,k + e k+1 α n,k+1 -f (a) 2 (1 + δ n,k )α n,k a k + τ n 1 -α n,0 a 0 /n K k=1 f a + e 0 α n,0 + e k α n,k -f (a) 2 (1 + ε n,k )α n,k a k .
For k = 0,

Ω n (π 0 )(a) = τ n α n,0 (1 -α n,0 a 0 /n) K k=1 (1 + ε n,k )α n,k a k and Θ n (π 0 )(a) = 1 α n,0 Ω n (π 0 )(a)
which proves the result for k = 0. For k = K + 1,

Ω n (π K +1 )(a) = τ n α n,K +1 (1 + δ n,K )α n,K a K = (1 + δ n,K )a K using α n,K +1 = τ n α n,K
to get the last equality, and Θ n (π K +1 ) = Ω n (π K +1 )/α n,K +1 , which proves the result for k = K + 1. Consider now k = 1, . . . , K :

Ω n (π k )(a) = τ n - 1 α n,k (1 + δ n,k )α n,k a k + 1 α n,k (1 + δ n,k-1 )α n,k-1 a k-1 + τ n α n,k (1 -α n,0 a 0 /n)(1 + ε n,k )α n,k a k = (1 + δ n,k-1 )a k-1 + τ n (ε n,k -δ n,k )a k - τ n α n,0 n (1 + ε n,k )a 0 a k , using δ n,0 = -1 for k = 1, and α n,k-1 τ n = α n,k for k ≥ 2. This proves the result for V π k n , while writing Θ n (π k )(a) = τ n 1 α 2 n,k (1 + δ n,k )α n,k a k + 1 α 2 n,k (1 + δ n,k-1 )α n,k-1 a k-1 + τ n α 2 n,k (1 -α n,0 a 0 /n)(1 + ε n,k )α n,k a k = 1 α n,k Ω n (π k )(a) + 2 τ n α n,k (1 + δ n,k )a k
proves the result for M π k n , which concludes the proof.

4.1. Tightness. We now prove that the sequence (A n , n ≥ 1) is tight. Since A n makes jumps of vanishing size (at most 1/α n,1 ), in order to show that (A n , n ≥ 1) is tight it is sufficient to show that for every T ≥ 0, (4.6) lim

η→0 lim sup n→+∞ ∆ n (η) = 0 where ∆ n (η) = sup Υ sup 0≤t ≤η E K +1 k=0 A n,k (Υ + t ) -A n,k (Υ)
and where the first supremum in the definition of ∆ n (η) is taken over all the random variables Υ ≤ T that are stopping times relatively to the filtration generated by A n (see Corollary on page 179 in Billingsley [START_REF] Billingsley | Convergence of probability measures[END_REF]). Fix in the rest of the proof some T ≥ 0. Because of the strong Markov property and the fact that we consider stopping times Υ ≤ T , we have

(4.7) ∆ n (η) ≤ sup 0≤t ≤η E sup 0≤y≤T Φ t n (A n (y))
where

(4.8) Φ t n (a) = K +1 k=0 E a A n,k (t ) -A n,k (0)
and the subscript refers to the initial state of the process A n (when there is no subscript, this refers to an initial condition as in the statement of the theorem or the propositions). By definition we have

E a (|A n,k (t )-A n,k (0)|) = E a (|M π k n (t )+V π k n (t )|)
, and so combining the triangular inequality with Cauchy-Schwarz inequality and summing over k = 0, . . . , K + 1 leads to (4.9)

Φ t n (a) ≤ K +1 k=0 E a 〈M π k n 〉(t ) + K +1 k=0 E a V π k n (t ) .
In order to do control the first moment of 〈M π k n 〉(t ) and of V π k n (t ), we introduce the functions ψ = π 0 + • • • + π K , i.e., ψ(a) = a 0 + • • • + a K for a = (a 0 , . . . , a K +1 ) ∈ R K +2 , and Ψ = ψ + ψ 2 . Since α n,1 → +∞ as n → +∞, α n,k+1 = τ n α n,k for k = 1, . . . , K and τ n → +∞, we will assume for convenience that 1 ≤ α n,1 ≤ • • • ≤ α n,K +1 . Moreover, it can be checked that the constant

C 0 = sup n≥1,1≤k≤K τ n |ε n,k -δ n,k |, 1 + |δ n,k |, 1 + |ε n,k |, τ n α n,0 n , 2τ n α n,k , τ n α n,k α n,0
is finite. From these various definitions and (4.2)-(4.4), and also recalling that

δ n,0 = δ n,K +1 = ε n,K +1 = -1, it follows that (4.10) max V π k n (t ) , 〈M π k n 〉(t ) ≤ 2C 2 0 t 0 Ψ(A n (u))du, k = 0, . . . , K + 1, t ≥ 0.
Combined with (4.9) this leads to

Φ t n (a) ≤ C 1 t 0 E a [Ψ(A n (u))] du +C 1 t 0 E a [Ψ(A n (u))] du with C 1 = 2C 2 0 (K + 2). For t ≤ η ≤ 1 we obtain (4.11) Φ t n (a) ≤ C 1 η 1/2 sup 0≤u≤1 E a [Ψ(A n (u))] + sup 0≤u≤1 E a [Ψ(A n (u))] .

Lemma 4.2.

There exists a finite constant C 2 such that for every initial state a,

(4.12) sup 0≤u≤1 E a ψ i (A n (u)) ≤ C 2 i j =0 ψ j (a), i = 1, 2.
In particular,

(4.13) sup n≥1 E sup 0≤y≤T Ψ(A n (y)) < +∞.
Let us quickly finish the proof of the tightness of (A n ) based on this lemma: first, plugging in (4.12) into (4.11), we obtain the existence of a finite constant C 3 such that

Φ t n (a) ≤ C 3 η 1/2 Ψ(a) + Ψ(a) , 0 ≤ t ≤ η ≤ 1,
and so in view of (4.7), we see that for every η ≤ 1,

∆ n (η) ≤ C 3 η 1/2 E sup 0≤y≤T Ψ(A n (y)) + Ψ(A n (y)) .
Thus, (4.13) implies the existence of a finite constant C 4 such that ∆ n (η) ≤ C 4 η 1/2 which achieves to prove that (A n ) is tight. It remains to prove Lemma 4.2.

Proof of Lemma 4.2. Let us first prove (4.12). Defining η n,k = 1/α n,k -1 {k =K } /α n,k+1 ≥ 0 and writing

ψ i (a) -ψ i a - e k α n,k + e k+1 α n,k+1 = ψ(a) i -ψ(a) -η n,k i = i ψ(a)
ψ(a)-η n,k

x i -1 dx, we get

ψ i (a) -ψ i a - e k α n,k + e k+1 α n,k+1 ≥ i η n,k ψ(a) -η n,k i -1 .
Defining µ n,k = 1/α n,0 + 1/α n,k similarly leads to

ψ i a + e 0 α n,0 + e k α n,k -ψ i (a) ≤ i µ n,k ψ(a) + µ n,k i -1 .
Plugging these inequalities in the definition (4.1) of Ω n , we obtain

Ω n (ψ i )(a) ≤ i τ n (1 -α n,0 a 0 /n) K k=1 µ n,k ψ(a) + µ n,k i -1 (1 + ε n,k )α n,k a k -i τ n K k=1 η n,k ψ(a) -η n,k i -1 (1 + δ n,k )α n,k a k .
Using 1 -α n,0 a 0 /n ≤ 1, expanding the terms raised to the power i -1 and changing the order of summation, we end up with

Ω n (ψ i )(a) ≤ i i -1 j =0 i -1 j ψ(a) i -1-j × τ n K k=1 µ j +1 n,k (1 + ε n,k )α n,k a k -(-1) j τ n K k=1 η j +1 n,k (1 + δ n,k )α n,k a k .
Plugging in the definitions of η n,k and µ n,k , we see that for j = 0 the term between brackets in the previous display is given by

τ n K k=1 (1 + ε n,k )(1 + α n,k /α n,0 )a k -τ n K k=1 (1 + δ n,k )(1 -1 {k =K } α n,k /α n,k+1 )α n,k a k = τ n K k=1 (ε n,k -δ n,k )a k + K k=1 τ n α n,k α n,0 (1 + ε n,k )a k + K -1 k=1 (1 + δ n,k )a k ≤ 3C 2 0 ψ(a).
For j ≥ 1, we obtain similarly, using also µ n,k ≤ 2/α n,k ≤ 2 and η n,k ≤ 1/α n,k ≤ 1, that the term between brackets is upper bounded by

τ n K k=1 4 α 2 n,k (1 + ε n,k )α n,k a k + τ n K k=1 1 α 2 n,k (1 + δ n,k )α n,k a k ≤ 3C 2 0 ψ(a).
We thus obtain, for some finite constant C (i ) and every a,

(4.14) Ω n (ψ i )(a) ≤ C (i ) i j =1 ψ(a) j .
In particular,

E a ψ i (A n (t )) = ψ i (A n (0)) + t 0 E a Ω n (ψ i )(A n (u)) du ≤ C (i ) t 0 E a i -1 j =0 ψ j (A n (u)) du +C (i ) t 0 E a ψ i (A n (u)) du
assuming without loss of generality that C (i ) ≥ 1 for the last inequality. Thus Gronwall's lemma implies

E a ψ i (A n (t )) ≤ C (i ) t 0 i -1 j =0 E a ψ j (A n (u)) du × e C (i )t .
Then (4.12) follows from this inequality by induction on i . We now derive (4.13): since Ψ = ψ + ψ 2 it is enough to prove (4.13) with ψ 2 in place of Ψ. First of all, note that the previous reasoning shows the existence of a finite constant C 2 such that (4.15) sup 0≤u≤T

E ψ i (A n (u)) ≤ C 2 i j =0 ψ j (A n (0)), i = 1, 2, 3, 4.
By definition,

ψ 2 • A n = ψ 2 (A n (0)) + V n + M n , defining V n = V ψ 2 n and M n = M ψ 2
n , and so

E sup 0≤y≤T ψ 2 (A n (y)) ≤ ψ 2 (A n (0)) + E sup 0≤y≤T V n (y) + E sup 0≤y≤T M n (y) .
The first term of the above upper bound is bounded (in n) since the sequence A n (0) converges by assumption. For the second term, we write

V n (y) = y 0 Ω n (ψ 2 )(A n (u))du ≤ C (2) y 0 Ψ(A n (u))du
where the inequality comes from (4.14). Thus

E sup 0≤y≤T V n (y) ≤ C (2) T 0 E Ψ(A n (u)) du ≤ C (2)T sup 0≤u≤T E Ψ(A n (u))
which is finite by (4.15). We now control the last martingale term. For any real-valued random variable X we have E(X ) ≤ 1 + E(X 2 ) and so Doob's inequality gives

E sup 0≤t ≤T M n (t ) ≤ 1 + 4E [〈M n (T )〉] = 1 + 4 T 0 E Θ n (ψ 2 )(A n (u)) du.
Defining Θ n as in (4.5) but with the term -α n,0 a 0 /n taken equal to 0, we have

Θ n (ψ 2 )(a) ≤ Θ n (ψ 2 )(a) = Ω n (ψ 4 )(a) -2ψ 2 (a)Ω n (ψ 2 )(a)
defining Ω n similarly as Ω n in (4.1) but with the term -α n,0 a 0 /n taken equal to 0. Since the first step in the derivation of (4.14) was to use 1-α n,0 a 0 /n ≤ 1, the reasoning leading to (4.14) also leads to an upper bound on |Ω n (ψ 2 )| of the same kind, i.e., it leads to the existence of a finite constant

C (i ) such that |Ω n (ψ i )(a)| ≤ C (i ) i j =1 ψ j (a)
. This finally proves that Θ n (ψ 2 )(a) ≤ C 4 j =1 ψ j (a) for some finite constant C , and in particular

E sup 0≤t ≤T M n (t ) ≤ 1 + 4C T 4 j =1 sup 0≤u≤T E ψ j (A n (u)) .
Since the supremum over n ≥ 1 of the right-hand side is finite by (4.15), the proof is complete.

Characterization of accumulation points.

Let A = (A k , k = 0, . . . , K + 1) be any accumulation point of (A n ) and assume without loss of generality that A n ⇒ A, i.e., A n converges weakly to A. Remember that we have to prove that A satisfies (2.4) in the intermediate regime, (2.8) in the large regime and (2.11) in the small regime. We treat the three regimes separately.

In each regime, if f ∈ D we denote by f the function ( t 0 f (u)du, t ≥ 0), and we will use the following result in conjunction with the continuous mapping theorem. If ( f n ), (g n ) are two sequences of functions with f n → f and g n → g for some continuous functions f and g (where → denotes uniform convergence on compact sets), then f n g n → f g and f n → f .

In particular, since in all three regimes it holds that τ n α n,k /α n,0 = α n,k+1 /α n,K +1 → 1 {k=K } and α n,0 /n → 0, we have from (4.2) that V π 0 n ⇒ A K . Since α n,0 → +∞ and δ n,0 = -1 we have 〈M π 0 n 〉 ⇒ 0 by (4.4), which implies by Doob's inequality that M

π 0 n ⇒ 0. Similarly, recalling that δ n,K +1 = ε n,K +1 = -1 we have V π K +1 n ⇒ A K and M π K +1 n ⇒ 0. 4.2.1. Large initial condition. In this regime, we have τ n α n,0 = n. Since moreover δ n,0 = δ n,K +1 = ε n,K +1 = -1, we obtain from (4.3) that V π 1 n ⇒ (γ 1 -A 0 )A 1 and V π k n ⇒ (A k-1 + (γ k -A 0 )A k ) for k = 2, . . . , K . In other words, V π k n ⇒ b k • A for k = 0, . . . , K + 1. Moreover, since α n,k → 0 and τ n /α n,k → 0 for k = 1, . . . , K , we obtain from (4.4) that 〈M π k n 〉 ⇒ 0 for any k = 0, . . . , K + 1. By Doob's inequality, this implies M π k n ⇒ 0 and since A n,k = A n,k (0) + M π k n + V π k n we finally get A n,k ⇒ a k + b k • A.
Since all the above convergences hold jointly as a consequence of the continuous mapping theorem, we get on the one hand that A n ⇒ a+ b•A, while on the other hand, since A n ⇒ A by assumption, we get A = a + b • A. Thus A solves (2.8) as desired, and since uniqueness of solutions to this ODE is guaranteed by Lemma 3.1, this uniquely characterizes A.

Intermediate initial condition.

We still have τ n α n,0 = n, so that as in the large regime we have V

π k n ⇒ b k • A.
Moreover, in this regime we have τ n /α n,k → 0 for k = 2, . . . , K which implies as in the previous regime M π k n ⇒ 0 for k = 0, 2, . . . , K + 1. The difference with the large regime is that since τ n = α n,1 , we have 〈M

π 1 n 〉 ⇒ 2 A 1 . Note that 2 A 1 = 〈M 〉 where M (t ) = t 0 (2A 1 (u)) 1/2 dB (u). Moreover, since M π k n ⇒ 0 for k = 1 while 〈M π 1
n 〉 ⇒ 〈M 〉, we get that the quadratic co-variation processes 〈M π k n , M π n 〉 vanish for any k = by polarization. Since all these convergences hold jointly, Theorem IX.2.4 in Jacod and Shiryaev [START_REF] Jacod | Grundlehren der Mathematischen Wissenschaften[END_REF] shows that A is the semimartingale with characteristics (b • A, M ) in the sense that A = a + b • A + M , i.e., A solves (2.4) (and thus is uniquely determined by Lemma 3.3). 4.2.3. Small initial condition. It follows similarly as in the two previous regimes, noting that in this regime we have τ n α n,0 /n → 0 (which leads to the drift term b S instead of b) and τ n /α n,k → 1 {k=1} (which leads to the non-vanishing diffusion term as in the intermediate regime).

ASYMPTOTIC BEHAVIOR OF THE OUTBREAK SIZE

This section is devoted to proving Theorem 2.3. For f ∈ D we define the operator T0 ( f ) = sup{t ≥ 0 : f (t ) > 0}; recall also the various operators defined in Section 2.1.

Then we have A n,0 (∞) = A n,0 ( T0 (A n,K )) and, in order to compute T0 (A n,K ), we will use the relation

T0 (A n,2 ) = T 0 (A n,1 ) + T 0 • θ T 0 (A n,1 ) (A n,2 ) which, iterated, leads to T0 (A n,K ) = K -1 k=1 T 0 • θ T 0 (A n,k-1 ) • • • • • θ T 0 (A n,1 ) (A n,k ).
There are two difficulties to solve in order to prove Theorem 2.3: the first one is that hitting times are in general not continuous functional, i.e., we may have

f n → f but T 0 ( f n ) → T 0 ( f ).
The second difficulty is that, as Proposition 2.1 shows, T0 (A K ) = +∞ for K ≥ 2 while the convergence A n ⇒ A holds uniformly on compact sets (see the discussion following Theorem 2.3). We address these two difficulties in two steps. 5.1. First step. The goal of this first step is to prove that (A n , T 0 (A n,1 )) ⇒ (A, T 0 (A 1 )). Assume first that for every δ > 0, (5.1) lim sup

n→+∞ P T 0 (A n,1 ) -T ↓ ε (A n,1 ) ≥ δ -→ ε→0 0.
We now argue that this implies (A n , T 0 (A n,1 )) ⇒ (A, T 0 (A 1 )). First of all, note that if T 0 (A n,1 ) ⇒ T 0 (A 1 ) then the joint convergence automatically holds, see for instance Corollary 2.2 in Lambert et al. [START_REF] Lambert | Scaling limits via excursion theory: Interplay between Crump-Mode-Jagers branching processes and Processor-Sharing queues[END_REF].

To see that (5.1) implies T 0 (A n,1 ) ⇒ T 0 (A 1 ), let us say that A 1 goes across ε if for every η > 0 we have inf 0≤t ≤η A 1 (T ↓ ε (A 1 ) + t ) < ε, and define the random set G through G = {ε > 0 : A 1 goes across ε}. Then it is known (and actually easy to show) that if

T ↓ ε (A 1 ) is almost surely finite and P(ε ∈ G ) = 1, then T ↓ ε (A n,1 ) ⇒ T ↓ ε (A 1 )
, see for instance Proposition VI.2.11 in Jacod and Shiryaev [START_REF] Jacod | Grundlehren der Mathematischen Wissenschaften[END_REF] or Lemma 3.1 in Lambert and Simatos [START_REF] Lambert | Asymptotic behavior of local times of compound Poisson processes with drift in the infinite variance case[END_REF]. Note that in our case, T ↓ ε (A 1 ) is finite by Proposition 2.2. On the other hand, the complement G c of G is precisely the set of discontinuities of the process (T ↓ ε (A 1 ), ε > 0). Since (T ↓ ε (A), ε > 0) is càglàd, as the left-continuous inverse of the process (inf [0,t ] A 1 , t ≥ 0), the set {ε > 0 : P(ε ∈ G c ) > 0} (sometimes called set of fixed times of discontinuities) is at most countable, see for instance Billingsley [4, Section 13]. Gathering these two observations, we see that

T ↓ ε (A n,1 ) ⇒ T ↓ ε (A 1 )
for all ε > 0 outside a countable set. Then, writing

P T 0 (A n,1 ) ≥ x = P T 0 (A n,1 ) ≥ x, T 0 (A n,1 ) -T ↓ ε (A n,1 ) ≥ δ + P T 0 (A n,1 ) ≥ x, T 0 (A n,1 ) -T ↓ ε (A n,1
) < δ , using (5.1) and playing with quantifiers gives the convergence of T 0 (A n,1 ) toward T 0 (A 1 ). Indeed, we can for instance write

P T 0 (A n,1 ) ≥ x ≤ P T 0 (A n,1 ) -T ↓ ε (A n,1 ) ≥ δ + P T ↓ ε (A n,1 ) ≥ x -δ , then choose ε such that T ↓ ε (A n,1 ) → T ↓ ε (A 1 )
to get by the portmanteau theorem, for any δ > 0, lim sup

n→+∞ P T 0 (A n,1 ) ≥ x ≤ lim sup n→+∞ P T 0 (A n,1 ) -T ↓ ε (A n,1 ) ≥ δ + P T ↓ ε (A 1 ) ≥ x -δ .
Since T ↓ ε (A 1 ) → T 0 (A 1 ) as ε → 0, we get by letting first ε → 0 and then δ → 0, and using (5.1), lim sup n→+∞ P T 0 (A n,1 ) ≥ x ≤ P (T 0 (A 1 ) ≥ x) .

Since x was arbitrary, this shows that T 0 (A n,1 ) ⇒ T 0 (A 1 ) by the portmanteau theorem. In conclusion, (5.1) indeed implies T 0 (A n,1 ) ⇒ T 0 (A 1 ).

The proof of (5.1) relies on a simple coupling between A n,1 and a continuous-time branching process (more precisely, a Bellman-Harris branching process). Looking at the transition rates of the process a n,1 in (2.1), we see that A n,1 decreases by 1/α n,1 at rate (1 + δ n,1 )A n,1 α 2 n,1 and increases by 1/α n,1 at rate

(1 + ε n,1 )A n,1 (1 -A n,0 /α n,1 )α 2 n,1 ≤ (1 + ε n,1 )A n,1 α 2 n,1 .
In particular, shifting the origin of time at T ↓ ε (A n,1 ) and using the strong Markov property, one sees that A n,1 can be coupled with a Markov process Z n,1 in such a way that Z n,1 (0) = α n,1 ε/α n,1 , A n,1 (T ↓ ε (A 1 ) + t ) ≤ Z n,1 (t ) for t ≥ 0, and Z n,1 decreases by 1/α n,1 at rate (1 + δ n,1 )Z n,1 α 2 n,1 and increases by 1/α n,1 at rate (1 + ε n,1 )Z n,1 α 2 n,1 (note that the law of Z n,1 depends on ε, but we omit this dependency in order to ease the notation). More concretely, this coupling can for instance be realized by adding a "ghost" individual in the population each time an individual in stage one makes an unsuccessful infection attempt (alternatively, we could also invoke the comparison result of Rogers and Williams [START_REF] Rogers | Diffusions, Markov processes, and martingales[END_REF]Theorem V.43.1]). Using the strong Markov property at time T ↓ ε (A n,1 ) and this coupling, we get

P T 0 (A n,1 ) -T ↓ ε (A n,1 ) ≥ δ ≤ P T 0 (Z n,1 ) ≥ δ .
It is well-known that Z n,1 ⇒ Z 1 , where Z 1 is Feller diffusion with drift γ 1 started at ε (see for instance Ethier and Kurtz [START_REF] Stewart | Markov processes[END_REF]Chapter 9]). Moreover, standard arguments can be used to show that T 0 (Z n,1 ) ⇒ T 0 (Z 1 ), for instance by using the fact that Z n,1 and Z 1 are time-change of Lévy processes killed at 0, say Y n and Y , so that T 0 (Z n,1 ) = ∞ 0 Y n ⇒ ∞ 0 Y = T 0 (Z 1 ) (this time-change transformation is usually called Lamperti transformation, see for instance Lamperti [START_REF] Lamperti | Continuous state branching processes[END_REF]). Thus we have (making clear the role of the initial condition) lim sup

n→+∞ P T 0 (A n,1 ) -T ↓ ε (A n,1 ) ≥ δ ≤ P (T 0 (Z 1 ) ≥ δ | Z 1 (0) = ε) .
Since T 0 (Z 1 ) ⇒ 0 as Z 1 (0) → 0, we have finally proved (5.1), which concludes the first step. 5.2. Second step. The first step shows, by using the strong Markov property at time T 0 (A n,1 ), that we only need to prove Theorem 2.3 when A n,1 (0) = 0. In this case, Theorem 2.1 shows that A n ⇒ A, where A is a solution to the ODE (3.1) with y = 0.

With this initial condition, we have A n,1 (t ) = 0 for all t ≥ 0, and π • A n is a Markov process: actually, it is a multistage epidemic process with K -1 stages. The problem to iterate the arguments of the first step is that Proposition 2.2 shows that, although T ↓ ε (A 2 ) < +∞ for every ε > 0, we have T 0 (A 2 ) = +∞. To get round this problem, we will use a time-change argument. Such an idea is classical in the SIR case K = 1, see for instance von Bahr and Martin-Löf [START_REF] Von | Threshold limit theorems for some epidemic processes[END_REF].

Let us set up a similar coupling as in the first step. The process A n,2 decreases by 1/α 2 n,1 at rate (1 + δ n,2 )α 3 n,1 A n,2 and increases by 1/α 2 n,1 at rate

(1 + ε n,2 )α 3 n,1 A n,2 (1 -A n,0 /α n,1 ) ≤ (1 + ε n,2 )α 3 n,1 A n,2 (1 -A n,0 (0)/α n,
+ ε n,2 )α 3 n,1 Z n,2 (1 -A n,0 (0)/α n,1 ). In particular, Z n,2 ⇒ z 2 with z 2 (t ) = z 2 (0) exp(-(A 0 (0) -γ 2 )t ).
Since by Proposition 2.2, A 0 is strictly increasing with A 0 (∞) > γ 2 , we can assume without loss of generality by shifting the processes at time T ↑ γ 2 (A n,0 ) + 1 that A 0 (0) > γ 2 , so that each z 2 vanishes exponentially fast. The problem, as mentioned earlier, is that T 0 (z 2 ) = +∞: we now introduce the time-change argument.

Let C n,2 be the right-continuous inverse of t → t 0 Z n,2 and c 2 be the right-continuous inverse of t → t 0 z 2 , in the sense that

A n • C n,2 ⇒ A • c 2 and Z n,2 • C n,2 ⇒ z 2 • c 2 ,
see for instance Helland [START_REF] Helland | Continuity of a class of random time transformations[END_REF]. Time-changing Z n,2 with C n,2 actually corresponds to the Lamperti transformation mentioned above: Z n,2 •C n,2 is a continuous-time random walk (killed at 0), z 2 •c 2 starts at A 2 (0) and decays linearly at rate A 0 (0) -γ 2 , i.e., z 2 (c 2 (t )) = A 2 (0) -(A 0 (0) -γ 2 )t for t ≤ ∞ 0 z 2 , and T 0 (Z n,2 • C n,2 ) ⇒ T 0 (z 2 • c 2 ). In particular, since

P T 0 (A n,2 • C n,2 ) -T ↓ ε (A n,2 • C n,2 ) ≥ δ ≤ P T 0 (Z n,2 • C n,2 ) ≥ δ | Z n,2 (0) = ε we obtain lim sup n→+∞ P T 0 (A n,2 • C n,2 ) -T ↓ ε (A n,2 • C n,2 ) ≥ δ -→ ε→0 0
and the arguments of the first step imply that

A n • C n,2 , T 0 (A n,2 • C n,2 ) ⇒ (A • c 2 , T 0 (A 2 • c 2 )) .
In particular, the strong Markov property at time T

0 (A n,2 • C n,2 ) shows that (A n • C n,2 ) T 0 (A n,2 • C n,2 ) = A n (T 0 (A n,2 )) ⇒ (x • c 2 ) (T 0 (A 2 • c 2 )) = x(∞) = 0.
By using the strong Markov property at time T 0 (A n,2 ) and iterating this argument, we finally end up with the desired result that A n,0 (∞) ⇒ A 0 (∞), which concludes the proof of Theorem 2.3.

AN INTRIGUING CONJECTURE

We conclude this paper by discussing a conjecture formulated in [START_REF] Antal | Outbreak size distributions in epidemics with multiple stages[END_REF] which initially motivated the present work.

Conjecture (Antal and Krapivsky [START_REF] Antal | Outbreak size distributions in epidemics with multiple stages[END_REF]). Assume that ε n,k = δ n,k = 0 and let N n,k be the number of individuals ever being of type k = 1, . . . , K over the course of the epidemic, starting from the initial condition a n,1 (0) = 1 and a n,k (0) = 0 for k = 2, . . . , K + 1. Then E(N n,k ) grows as n → +∞ like n kλ K , where

λ K = 2 K -1 (K + 1)2 K -1 .
Note that, with our notation, N n,K = a n,K +1 (∞), since any individual ever removed must have been in stage K of the epidemic at some point (and vice-versa). On the other hand, Theorem 2.3 shows that starting with of the order of n 1/(K +2) individuals in stage one (instead of just one as in the above conjecture), a n,K +1 (∞) is, in distribution, of the order of n (K +1)/(K +2) . When K = 1, there is a well-known argument that links these two objects: the connection goes through a random partitioning of {1, . . . , n}.

Consider indeed the following model, where n individuals are assigned a unique label from the set {1, . . . , n} and which results in a random partition Π 1 , . . . , Π S of the set {1, . . . , n}. Imagine that Π 1 , . . . , Π s have been generated and that the set F s = {1, . . . , n} \ (Π 1 ∪• • •∪Π s ) is not empty: then the iteration proceeds as follows. Choose an individual v uniformly at random from F s , and run the epidemic with the following initial condition: at time 0 the individuals in F s \{v} are susceptible and the individuals in Π 1 ∪• • •∪Π s are removed, so that only v is infected (and, in the case K ≥ 2, is in the first stage of the epidemic). Eventually, this epidemic will die out and we define Π s+1 as the set of individuals infected over the course of this epidemic.

This connection between random partition and epidemic processes is well-known, see for instance Barbour and Mollison [3]. In the case K =1 this links clusters of the Erdös-Rényi random graph to the SIR process. This construction leads to several interesting by-products, one of them being that it makes it possible to compute the mean size of a typical cluster in terms of the mean size of the largest ones. More precisely, choose an individual v uniformly at random in {1, . . . , n} and let Π * be the cluster to which it belongs. Let moreover (Π (i ) , i ≥ 1) be the clusters ordered in decreasing size, i.e., {Π (i ) } = {Π i } and

|Π (1) | ≥ |Π (2) | ≥ • • • with |E | the size of a set E ⊂ N. Then the proba- bility that v belongs to Π (i ) is exactly |Π (i ) |/n, in which case Π * = Π (i ) and so E(|Π * |) = i ≥1 E |Π (i ) | n |Π (i ) | .
Assuming that the largest term in this sum dominates, we get the approximation (6.1)

E |Π * | ≈ 1 n E |Π 2 (1) | . When K = 1, |Π * | is, by exchangeability, equal in distribution to N n,1 and so E(|Π * |) = E(N n,1 ). Moreover, |Π (1)
| is of the same order than an outbreak size started from an intermediate initial condition. Roughly speaking, this comes from the fact that the intermediate regime is precisely the one where the finite-size population effects begin to kick in. In particular, |Π (1) | is of the order of n (K +1)/(K +2) = n 2/3 . Gathering these two observations, we end up thanks to (6.1) with the relation E(N n,1 ) ≈ n 4/3-1 = n 1/3 , when K = 1. This answer coincides with Antal and Krapivsky's conjecture and the above reasoning through random partitioning can be made rigorous, see, e.g., Janson and Luczak [START_REF] Janson | Susceptibility in subcritical random graphs[END_REF].

It is tempting to also use this reasoning for K ≥ 2, and this initially motivated us to identify the intermediate regime and the scalings at play there. However, this reasoning would lead to the estimate

E(N n,K ) ≈ 1 n n 2(K +1)/(K +2) = n K /(K +2)
which is different (for K ≥ 2) from the n K λ K predicted by Antal and Krapivsky [START_REF] Antal | Outbreak size distributions in epidemics with multiple stages[END_REF]. We find this discrepancy, and the related fact that the above reasoning through random partitioning seems to fail, very intriguing. We believe that the temporal aspects intrinsic to this multistage epidemic (see the discussion of the main results in the introduction) play a major role in this discrepancy, although it is challenging to obtain rigorous results in that direction.

APPENDIX A. PROOF OF LEMMAS 3.1 AND 3.2

Let x be any solution to (3.1) defined on the interval J = [0, t * ) for some t * ∈ (0, ∞]. It will be convenient to define x 1 = y and to index the R K +1 -valued function x by the set {0, 2, . . . , K + 1}, i.e., to write x = (x 0 , x 2 , . . . , x K +1 ). Let also in the sequel x k (t ) = x k (0) + t 0

I k (t ) = t 0 x 0 (u) -γ k du, t ∈ J , k = 1, . . . , K . Then (e I k x k ) e -I k = (x 0 -γ k )x k + x k is equal to x k-1 for k = 2, .
x k-1 (u)e I k (u) du e -I k (t ) , k = 2, . . . , K , t ∈ J .

A.1. Proof of Lemma 3.1. Note first that the representation (A.1) implies that x(t ) ∈ [0, ∞) K +1 for every t ∈ J . Indeed, since x 1 (t ) = y(t ) ≥ 0, this implies that x 2 (t ) ≥ 0 and by induction on k, this implies that x k (t ) ≥ 0 for every k = 2, . . . , K and t ∈ J . Since finally x 0 (0), x K +1 (0) ∈ [0, ∞) and x 0 = x K +1 = x K which has just been showed to stay non-negative, we obtain that x 0 and x K +1 also stay non-negative.

Let us now prove Lemma 3.1, i.e., existence and uniqueness of solutions to (3.1). Since F is locally Lipschitz, the Picard-Lindelhöf theorem implies local existence and uniqueness to (3.1). To show this globally, we only have to show that local solutions do not explode. Since x 2 = y + (γ 2x 0 )x 2 and x 2 , x 0 ≥ 0, we obtain x 2 ≤ y + γ + 2 x 2 where γ + = max(0, γ) for any γ ∈ R. Gronwall's lemma thus shows that x 2 does not explode. Similarly, for k = 3, . . . , K we have x k ≤ x k-1 + γ + k x k and so x k also does not explode. Finally, since x 0 = x K +1 = x K , solutions stay locally bounded which proves the global existence and uniqueness on [0, ∞).

A.2. Proof of Lemma 3.2. We prove each property separately. As mentioned earlier, in the rest of the proof we define x 1 = y.

Proof of (i). Assume that y(0) > 0 or a 2 > 0: then it is clear from (A.1) that x 2 (t ) > 0 for all t > 0. By induction, we see that x k (t ) > 0 for all k = 2, . . . , K and t > 0.

Proof of (ii). Assume that x k (t ) > 0 for k = 2, . . . , K and t > 0, and that x 0 is bounded: we prove by backwards induction on k that ∞ 0 x k is finite for k = 1, . . . , K and that x 0 (∞) > γ k for k = 2, . . . , K . For k = K , the fact that ∞ bounded and non-decreasing, and so its derivative x K is integrable on [0, ∞). Consider now any 2 ≤ k ≤ K and assume that ∞ 0 x k is finite: we prove that ∞ 0 x k-1 is finite and that x 0 (∞) > γ k . Since ∞ 0 x k is finite, there must exist a sequence t n → +∞ such that x k (t n ) → 0. Moreover, we have by definition x k = x k-1 +(γ k -x 0 )x k and so integrating between times 0 and t n we obtain

t n 0 x k-1 = x k (t n ) -x k (0) + t n 0 (x 0 (u) -γ k )x k (u)du ≤ x k (t n ) + (x 0 (∞) -γ k ) t n 0 x k .
Letting n → +∞, we obtain the inequality ∞ 0 x k-1 ≤ (x 0 (∞)-γ k ) ∞ 0 x k which shows, by induction, that ∞ 0 x k-1 is finite and also that x 0 (∞) > γ k (since

∞ 0 x k-1 > 0).
Proof of (iii). Assume that y(t ) = 0 for all t ≥ 0 and that a 2 > 0: in particular, according to (i), we see that x k (t ) > 0 for every t > 0 and k = 0, 2, . . . , K , a fact that will repeatedly be used in the sequel. We begin by proving the following formula:

(A.2) x k (t ) = k-2 i =0
x k-i (0)φ k,i (t ) × exp - We prove (A.2) by induction on k. For k = 2, we have by (A.1), and since x 1 = 0, x 2 (t ) = x 2 (0)e -I 2 (t ) which is precisely (A.2). Assume now that (A.2) holds for k ≥ 2 and let us prove it for k + 1: plugging in (A.2) into (A.1), we obtain x k+1 (t ) = x k+1 (0)

+ t 0 k-2 i =0
x k-i (0)φ k,i (u) e -I k (u) × e I k+1 (u) du e -I k+1 (t ) .

Using I k+1 (u) -I k (u) = (γ k -γ k+1 )u = η k+1 u, exchanging the integral and the sum and changing variables in the sum, we obtain

x k+1 (t ) = x k+1 (0) + k-1 i =1
x k+1-i (0) t 0 φ k,i -1 (u)e η k+1 u du e -I k+1 (t ) from which we get (A.2) by (A.3). We now prove that x 0 (∞) is finite by contradiction, so assume that x 0 (∞) = +∞. Starting from the definition (A.3) of the φ k,i 's, we get by induction that φ k,i (t ) ≤ e i η * t with η * = 1+max 3≤ j ≤K |η j |. In particular, we get from (A.2) for k = K that

x K (t ) ≤ k-2 i =0
x k-i (0) e K η * t × exp -t 0 (x 0 (u) -γ K )du . Since x 0 (t ) -γ K > K η * for t large enough (since we are assuming x 0 (∞) = +∞), the last display implies that x K converges to 0 exponentially fast, and in particular ∞ 0 x K < +∞. Since x 0 = x K , x 0 (∞) is finite, which yields the contradiction. Thus x 0 (∞) must be finite and so the conclusions of (ii) apply, in particular x 0 (∞) > γ k for every 2 ≤ k ≤ K .

We now complete the proof and show that x k (t ) → 0 for every k = 2, . . . , K : in view of (A.2) we only have to show that (A.4) lim t →+∞ φ k,i (t )e -t 0 (x 0 -γ k ) = 0 for every i = 0, . . . , K -2 and every k = i + 2, . . . , K . We prove this by induction on i : for i = 0 this comes immediately from the facts that φ k,0 (t ) = 1 and x 0 (∞) > γ k . So assume that (A.4) holds for some i = 0, . . . , K -2 and every k = i + 2, . . . , K : we show that it also holds for i + 1 and k = i + 3, . . . , K . By definition (A.3) we have φ k,i +1 (t )e -t 0 (x 0 -γ k ) = t 0 φ k-1,i (u)e η k u du e -t 0 (x 0 -γ k ) .

Let ε > 0 and, by induction hypothesis, t * such that φ k-1,i (t ) ≤ εe t 0 (x 0 -γ k-1 ) for every t ≥ t * . Then for such t , Since the first term of the above upper bound vanishes as t → +∞ and we can rewrite the second term as Thus to achieve the proof we only have to show that this last supremum is finite. Let κ > 0 and s * < +∞ be such that x 0 (t ) -γ k ≥ κ for t ≥ s * and k = 2, . . . , K . Then for t ≥ s *

Lemma 3 . 1 .Lemma 3 . 2 .

 3132 For any ā ∈ [0, ∞) K +1 and any continuous function y : [0, ∞) → [0, ∞), the ODE (3.1) has a unique solution defined on [0, ∞). Fix ā = (a 0 , a 2 , . . . , a K +1 ) ∈ [0, ∞) K +1 and y : [0, ∞) → [0, ∞) a continuous function, and let x = (x 0 , x 2 , . . . , x K +1 ) be the unique solution to (3.1) given by Lemma 3.1.

1 . 3 . 4 .

 134 do not explode. But A 1 cannot explode since it is dominated by Feller diffusion, and since π • A satisfies (3.1), no other coordinate can explode by Lemma 3.Lemma If A satisfies (2.4), then T 0 (A 1 ) is almost surely finite.

2 <

 2 +∞, c 2 blows up at time ∞ 0 z 2 . Moreover, such random timechange transformations induce continuous mappings, so that

  . . , K by(3.1). Thus for these k we have (e I k x k ) e -I k = x k-1 , which can be rewritten as (A.1)

  u) -γ k )du , k = 2, . . . , K , t ≥ 0,where the functions φ k,i for k = 2, . . . , K and i = 0, . . . , k -2 are defined recursively by φ k,0 (t ) = 1 and for i = 1, . . . , k -2,(A.3) φ k,i (t ) = t 0 φ k-1,i -1 (u)e η k u du, with η k = γ k-1 -γ k (k = 3, . . . , K ).

  φ k,i +1 (t )e -t 0 (x 0 -γ k ) ≤ t * 0 φ k-1,i (u)e η k u du e -t 0 (x 0 -γ k ) 0 -γ k-1 ) e η k u du e -t 0 (x 0 -γ k ) .

  0 -γ k-1 ) e η k u du e -t 0 (x 0 -γ k ) = t 0 e u 0 (x 0 -γ k ) du e -t 0 (x 0 -γ k ) = t 0 e -t u (x 0 -γ k ) du, we obtain lim sup t →+∞ φ k,i +1 (t )e -t 0 (x 0 -γ k ) ≤ ε sup t ≥0 t 0 e -t u (x 0 -γ k ) du .

  t u (x 0 -γ k ) du ≤ s * 0 e -s * u (x 0 -γ k ) du so that writing t 0 e -t u (x 0 -γ k ) du = s * 0 e -t u (x 0 -γ k ) du + t s * e -t u (x 0 -γ k ) duachieves the proof.

  [START_REF] Aldous | Brownian excursions, critical random graphs and the multiplicative coalescent[END_REF] ), where the inequality follows from the monotonicity of A n,0 . Thus similarly as in the first step, we can couple A n,2 with a continuous-time Markovian branching process Z n,2 such that Z n,2 (0) = A n,2 (0), A n,2 (t ) ≤ Z n,2 (t ) for t ≥ 0, and Z n,2 decreases by 1/α 2

	n,1 at rate
	(1 + δ n,2 )α 3 n,1 Z n,2 and increases by 1/α 2 n,1 at rate (1