Home | About | Journals | Submit | Contact Us | Français |

**|**Int J Environ Res Public Health**|**v.7(3); 2010 March**|**PMC2872325

Formats

Article sections

- Abstract
- 1. Introduction
- 2. Single-type Bienaymé-Galton-Watson Branching Processes
- 3. Single-type Branching Processes with Population-Dependent Offsprings
- 4. Set of Single-type Branching Processes with Population-Dependent Offsprings
- 5. Multitype Branching Processes with Age and Population-Dependent Offsprings
- 6. Discussion
- References

Authors

Related links

Int J Environ Res Public Health. 2010 March; 7(3): 1186–1204.

Published online 2010 March 19. doi: 10.3390/ijerph7031204

PMCID: PMC2872325

National Agricultural Research Institute, UR341, Department of Applied Mathematics and Informatics, F-78352 Jouy-en-Josas, France; E-Mail: rf.arni.yuoj@bocaj.enitsirhc

Received 2010 January 11; Revised 2010 March 1; Accepted 2010 March 16.

Copyright © 2010 by the authors; licensee Molecular Diversity Preservation International, Basel, Switzerland.

This article is an open-access article distributed under the terms and conditions of the Creative Commons Attribution license (http://creativecommons.org/licenses/by/3.0/).

Branching processes are stochastic individual-based processes leading consequently to a bottom-up approach. In addition, since the state variables are random integer variables (representing population sizes), the extinction occurs at random finite time on the extinction set, thus leading to fine and realistic predictions. Starting from the simplest and well-known single-type Bienaymé-Galton-Watson branching process that was used by several authors for approximating the beginning of an epidemic, we then present a general branching model with age and population dependent individual transitions. However contrary to the classical Bienaymé-Galton-Watson or asymptotically Bienaymé-Galton-Watson setting, where the asymptotic behavior of the process, as time tends to infinity, is well understood, the asymptotic behavior of this general process is a new question. Here we give some solutions for dealing with this problem depending on whether the initial population size is large or small, and whether the disease is rare or non-rare when the initial population size is large.

Mathematical models of propagation of a disease in given populations play a central role for understanding this propagation, for predicting the future extension of the outbreak, its extinction time, and for evaluating the efficiency of control measures. Of course the validity and the richness of results of a model strongly depend on the reliability and the accuracy of the model. A fine predictive model should be built as far as possible in a rigorous mechanistic way starting from the mechanism of exposure/infection of each individual and taking into account their variability. The population dynamic described by births, deaths, migrations should also be taken into account, especially when the incubation time is relatively long in comparison with this dynamic. Of course, the time unit of the model should be chosen in keeping with the respective durations of each health state and the data.

Nevertheless, since the asymptotic behavior, as time tends to infinity, of deterministic models in continuous time are more easily studied than that of discrete time models or stochastic models, a large literature in applied mathematical journals is devoted to such theoretical studies [19], while for statistical purposes, either descriptive non dynamical models or very simple stochastic models are used, such as a chain binomial model for the propagation of *SIS* or *SIR* diseases in closed populations [5, 6], or BGW (Bienaymé-Galton-Watson) branching processes on the clinical (or diagnosed) cases as an approximate model of the beginning of an outbreak [2, 3, 4, 6, 10, 11, 20]. These two approaches are both individual-based (each individual transition is given) and stochastic (individual variability is taken into account) but the first one takes into account the individual health state evolution *S* → *I* → *S* or *S* → *I* → *R*, where *S* means susceptible, *I* means infective or clinical cases (according to authors and purposes) and *R* means removed from the susceptible population, either by immunization or death; while the second approach directly models the process of *I* individuals, without analyzing this mechanism of exposure/infection. The behavior of these two types of processes are well-known. The chain binomial models are just Markov chains in a finite state space, and since there is no immigration and no incubation period, the state “0 infected individual” is an absorbing state. Concerning BGW processes, a huge literature exists for describing their properties. In fact these two different types of models may be written as particular cases of a large class of branching processes with individual transitions that may be population and age dependent, and which can take into account both the population dynamics and the disease dynamic.

Branching processes were initiated in the nineteen century by Bienaymé and then by Galton and Watson, for studying the extinction of some family names. Since this time, the complexity of these processes continues to increase allowing to describe more and more realistic population dynamics. These processes are based on the simple property that the population size of each considered type (such as clinical cases here) at each time is calculated as the sum of all the new individuals (“offspring”) of this type who are generated by the individuals of the population at the previous times. Since the modelled variables are integers, then the extinction time of the population of each type is finite on the set of trajectories which extinct. This is a finer and more realistic property than the asymptotic extinction time given by a deterministic model. Since the population dynamic may influence the disease dynamic, and conversely, these two dynamics should be explicitly taken into account in the model, thus leading to rigorous, but not simple, multitype models where each type should represent an health state crossed with influence factors levels. Typical examples of such influence factors are age, geographical locations. However, some authors model directly the time evolution of the incidence of clinical cases, without explicating all the intermediary steps.

In the following subsections, we present such direct models starting from the simplest one, the single-type Bienaymé-Galton-Watson process, and finishing by a general and rigorous approach that takes into account the intermediary steps and the population dynamic, and is based on age-dependent and population-dependent individual transitions. We focus on models in discrete time since they have the double advantage to be easily written as recursive models and to easily correspond to the time unit of observation, which offers a pedagogic framework and moreover facilitates the model validation and the estimation of unknown parameters. We study here the behavior of these models. From now on, all results are given conditionally to the initial value *F*_{0} of the process and the notation “|*F*_{0}” will be therefore omitted for the sake of simplicity of formulas. The proofs are given in detail in [23] and will be omitted here.

Let *I _{n}* be the incidence of clinical (or diagnosed) cases at time

$${I}_{n}=\sum _{i=1}^{{I}_{n-1}}{Y}_{n,i},$$

(1)

where the variables {*Y _{n,i}*}

The behavior of this process has been deeply analysed for a long time (see for example [1, 13, 16, 24]). We recall here the main properties. The trajectories of this process either die out or explode in an exponential way: *P*({lim_{n}*I _{n}* =

In addition, the probability of extinction, *q*, is the solution of *q* = *f*(*q*) := *E*(*q*^{Y1,1}) and in the subcritical/critical cases *m* ≤ 1, the “epidemic size”
$N\hspace{0.17em}:={\sum}_{k=0}^{\mathit{\text{Text}}.}{I}_{k}$ (total number of cases generated until the extinction time *T _{ext.}*), follows a power series distribution when

$$P(N=k)=\text{exp}(-mk)\frac{{I}_{0}}{k}\frac{{(mk)}^{k-{I}_{0}}}{(k-{I}_{0})!},\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}k\ge {I}_{0},$$

(2)

and *E*(*N*) = *I*_{0}(1 − *m*)^{−1}, *Var*(*N*) = *I*_{0}*m*(1 − *m*)^{−3} are increasing functions of *m*, *I*_{0}.

When the population is small or when the individual migrations are slow compared to the infection process, the depletion of susceptible individuals due to the infection should be taken into account, which is not the case in the BGW process. The typical bell curve form of outbreaks is due to this depletion. The simplest such model is just an extension of the BGW process, that is,
${I}_{n}={\sum}_{i=1}^{{I}_{n-1}}{Y}_{n,i}$, where the {*Y _{n,i}*}

A simple example is when
$\mathcal{L}({F}_{n-1})=\mathcal{L}({\sum}_{l=1}^{{d}_{n}}{\mu}^{l}{I}_{n-l})$, *μ* [0, 1], where
$m({\sum}_{l=1}^{{d}_{n}}{\mu}^{l}{I}_{n-l})\hspace{0.17em}(\text{resp}.\hspace{0.17em}P({Y}_{n,1}=0|{F}_{n-1})=p({\sum}_{l=1}^{{d}_{n}}{\mu}^{l}{I}_{n-l}))$ is a decreasing (resp. increasing) function of
${\sum}_{l=1}^{{d}_{n}}{\mu}^{l}{I}_{n-l}$. For example, *E*(*Y*_{n,1}|*Y*_{n,1} ≠ 0) = *α* and
$P({Y}_{n,1}=0|{F}_{n-1})=1-K{[K+{\sum}_{l=1}^{{d}_{n}}{\mu}^{l}{I}_{n-l}]}^{-1}$, 0 < μ ≤ 1, *K* > 0. We will see in this section that this kind of models may exhibit a random cycle and therefore may be used for modeling recurrent epidemics. The extremal case *μ* = 0 means no depletion in *S*, which is got by an immediate healing with no immunization, after one time unit in the state *I* (*SIS* disease) and this model is reduced to the previous BGW process.

The other extremal case *μ* = 1 with *d _{n}* =

The intermediate case 0 < *μ* ≤ 1 with *d _{n}* =

Let *A*1: there exists a positive function *f*(.) such that *P*(*Y _{n,i}* = 0

This assumption is checked as soon as *P*(*Y _{n,i}* = 0

**Proposition 1** *Let us assume A*1*. Then P*(lim_{n}*I _{n}* = 0 lim

Let *A*2: there exists *m _{*}* < ∞ such that

Let us define the bifurcation parameter *R _{∞}* by

$${R}_{\infty}=\text{sup}\{R:0<R<1\Rightarrow P(\underset{n}{\text{lim}}\hspace{0.17em}{I}_{n}=0)=1\},$$

where *R* is any real quantity depending on the parameters of the distribution of {*I _{n}*}.

**Proposition 2** *The bifurcation parameter R*_{∞} *is equal to*
lim_{|F|→∞}*m*(*F*)*. Let us assume A*1 *and A*2*. Then R*_{∞} < 1*, i.e., the process dies out a.s.. Moreover for the deterministic trajectory* {*X _{n}*}

**Remark 1** *When d _{n}* =

$$\tilde{\mathbf{\text{M}}}({F}_{n-1})=\left(\begin{array}{lllll}m({F}_{n-1})\hfill & 1\hfill & 0\hfill & \dots \hfill & 0\hfill \\ 0\hfill & 0\hfill & 1\hfill & \dots \hfill & 0\hfill \\ \dots \hfill & \dots \hfill & \hfill & \hfill & \hfill \\ 0\hfill & 0\hfill & 0\hfill & \dots \hfill & 1\hfill \\ 0\hfill & 0\hfill & 0\hfill & \dots \hfill & 0\hfill \end{array}\right)$$

*which implies that, for all n* ≥ *d* − 1, **M͂*** ^{n}*[1

**Remark 2** *If m*((1, 0, ..., 0)) > 1 *(supercritical assumption at the beginning of the disease spread), since m*(*F*_{n−1}) < 1 *under A*2*, for |F*_{n−1}*| and n sufficiently large, then we may observe oscillations until the a.s. extinction, the process increasing when |***Ĩ**_{n−1}|*is small enough, and decreasing when it is too large. It is the case of the logistic model*
$m({F}_{n-1})=\alpha K{[K+{\sum}_{l=1}^{{d}_{n}}{\mu}^{l}{I}_{n-l}]}^{-1}$, *when d _{n}* =

We generalize here the model of the previous section to a set of *J* similar diseases in competition:

$${I}_{n}^{(j)}=\sum _{i=1}^{{I}_{n-1}^{(j)}}{Y}_{n,i}^{(j)},\hspace{0.17em}\hspace{0.17em}j=1,\hspace{0.17em}\dots ,\hspace{0.17em}J,$$

where, for each *j*, the
${\left\{{Y}_{n,i}^{(j)}\right\}}_{i}$ are i.i.d. according to
${\mathcal{L}}^{(j)}({F}_{n-1})=:\mathcal{L}({m}^{(j)}({F}_{n-1}),{\sigma}^{2(j)}({F}_{n-1}))$ given
${F}_{n-1}\hspace{0.17em}:\hspace{0.17em}{\left\{{\left\{{I}_{k}^{(j)}\right\}}_{j=1,\dots J}\right\}}_{k\le n-1}$. We assume that *m*(*F*_{n−1}) is non increasing in
${I}_{k}^{(j)}$, for each *j*, each *k* ≤ *n* − 1.

A simple example is when as previously
${\mathcal{L}}^{(j)}({F}_{n-1})={\mathcal{L}}^{(j)}({\sum}_{l=1}^{{d}_{n}}{\mu}^{l}{\sum}_{{j}^{\prime}}{I}_{n-l}^{({j}^{\prime})})$, *μ* [0, 1], and
${m}^{(j)}({\sum}_{l=1}^{{d}_{n}}{\mu}^{l}{\sum}_{{j}^{\prime}}{I}_{n-l}^{({j}^{\prime})})\hspace{0.17em}(\text{resp}\text{.}\hspace{0.17em}P({Y}_{n,1}^{(j)}=0|{F}_{n-1})={p}^{(j)}({\sum}_{l=1}^{{d}_{n}}{\mu}^{l}{\sum}_{{j}^{\prime}}{I}_{n-l}^{({j}^{\prime})}))$ is a decreasing (resp. increasing) function of
${\sum}_{l=1}^{{d}_{n}}{\mu}^{l}{\sum}_{{j}^{\prime}}{I}_{n-l}^{({j}^{\prime})}$. A typical example is given by different influenza viruses.

As in the single-type case, assuming *m*^{(j)}(*F*) < *m _{*}* < ∞, for all

Two populations of infectives from similar diseases in competition following the same logistic Poisson model
${I}_{n}^{(j)}={\sum}_{i=1}^{{I}_{n-1}^{(j)}}{Y}_{n,i}^{(j)}$,
${Y}_{n,i}^{(j)}|{F}_{n-1}\sim \mathit{\text{Poisson}}(\alpha K{[K+{\sum}_{l\le d}{\mu}^{l}({I}_{n-l}^{(1)}+}^{}$ **...**

Let the deterministic associate trajectory
$\left\{{X}_{n}^{(j)}\right\}$ be defined by
${X}_{n}^{(j)}\hspace{0.17em}:={m}^{(j)}({F}_{n-1}^{X}){X}_{n-1}^{(j)}$, *n*≥ 1,
${X}_{0}^{(j)}\hspace{0.17em}:={I}_{0}^{(j)}$.

**Proposition 3** *The reproductive number*
${R}_{\infty}^{(j)}\hspace{0.17em}:={\overline{\text{lim}}}_{|F|\to \infty}{m}^{(j)}(F)$ *is the bifurcation parameter for*
$\left\{{I}_{n}^{(j)}\right\}$*, that is, if*
${R}_{\infty}^{(j)}<1$*, then*
${\text{lim}}_{n}{I}_{n}^{(j)}\stackrel{a.s.}{=}0$*, and the reproductive number*
${R}_{0}^{(j)}\hspace{0.17em}:={\underset{\_}{\text{lim}}}_{|{F}^{X}|\to 0}{m}^{(j)}({F}^{X})$ *is the bifurcation parameter of*
$\left\{{X}_{n}^{(j)}\right\}$*, that is if*
${R}_{0}^{(j)}<1$*, then*
${\text{lim}}_{n}{X}_{n}^{(j)}=0$*. Moreover if*
${R}_{0}^{(j)}\hspace{0.17em}>\hspace{0.17em}1$*, then*
${\overline{\text{lim}}}_{n}{X}_{n}^{(j)}\hspace{0.17em}>\hspace{0.17em}0$.

A possible extension of the previous models consists in adding an immigration. Then:

$${I}_{n}=\sum _{i=1}^{{I}_{n-1}}{Y}_{n,i}+{J}_{n}{\delta}_{n},$$

(3)

where *J _{n}* is an immigration at time

This type of processes has been deeply studied in the subcritical case *m* < 1, when ({*Y _{n,i}*}

We may also write (3) in the following way:

$${I}_{n}={1}_{\{{I}_{n-1}\ne 0\}}\sum _{i=1}^{{I}_{n-1}}{\tilde{Y}}_{n,i}+{1}_{\{{I}_{n-1}=0\}}{J}_{n}{\delta}_{n},\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{\tilde{Y}}_{n,i}={Y}_{n,i}+\frac{{J}_{n}}{{I}_{n-1}}{\delta}_{n}.$$

Consequently, when ({*Y _{n,i}*}

Let ${\mathbf{\text{N}}}_{n}=({N}_{n}^{1},\hspace{0.17em}\dots {N}_{n}^{D})$ satisfying

$${N}_{n}^{k}=\sum _{h=1}^{D}\sum _{l=1}^{d}\sum _{i=1}^{{N}_{n-l}^{h}}{Y}_{n-l,n,i}^{(h),k},\hspace{0.17em}\hspace{0.17em}k=1,\hspace{0.17em}\dots ,\hspace{0.17em}D,$$

(4)

where *h*, *k* represent types and *l* is the maturation time for getting the “offsprings”
${Y}_{n-l,n,i}^{(h),k}$ of the individual *i* (number of *k* individuals at time *n* generated by *i* belonging to the type *h* at time *n − l*). We assume that the
${\left\{{Y}_{n-l,n,i}^{(h),k}\right\}}_{i,l}$ are mutually independent given *F _{n−1}* = {

The models of the previous sections (Sections 2. and 4.) belong to this class when *d _{n}* =

When rigorously modeling the propagation of a disease in a given population, taking explicitly into account the population dynamic and the disease dynamic, then a type should be a health state (*S*, *E* (incubation) or *I*, for example) crossed with an influence factors level, and the “offsprings”
$\left\{{Y}_{n-l,n,i}^{(h),k}\right\}$ should be expressed according to the number of newborn individuals of *i*, who are of the *k* type, and to the proper transition of *i* from the state *h* at *n − l* to the state *k* at *n*. For example assuming that there is no influence factors and that the set of individual health states is {*S, E, I*} and assuming that the life duration in the state *I* is at most one time unit, then, for *k* = *I*,

$${Y}_{n-l,n,i}^{(h),I}={1}_{\{h=S\}}{\delta}_{2,n-l+1,n,i}^{E,I|h}+\sum _{j=1}^{{Y}_{n-l+1,i}^{(h)}}{\delta}_{1,n-l+1,n,i,j}^{E,I|(h)},$$

where the
${\left\{{\delta}_{2,n-l+1,n,i}^{E,I|S}\right\}}_{i}$ are i.i.d. Bernoulli variables given *F _{n−1}*, with
${\delta}_{2,n-l+1,n,i}^{E,I|S}$ equal to 1 if the individual

Let us write * _{n}* := (

$$\mathbb{\text{M}}({F}_{n-1})=\left(\begin{array}{lllll}{\mathbf{\text{M}}}_{n-1,n}({F}_{n-1})\hfill & \mathbf{\text{I}}\hfill & \mathbf{0}\hfill & \dots \hfill & \mathbf{0}\hfill \\ {\mathbf{\text{M}}}_{n-2,n}({F}_{n-1})\hfill & \mathbf{0}\hfill & \mathbf{\text{I}}\hfill & \dots \hfill & \mathbf{0}\hfill \\ \dots \hfill & \dots \hfill & \hfill & \hfill & \hfill \\ {\mathbf{\text{M}}}_{n-(d-1),n}({F}_{n-1})\hfill & \mathbf{0}\hfill & \mathbf{0}\hfill & \dots \hfill & \mathbf{\text{I}}\hfill \\ {\mathbf{\text{M}}}_{n-d,n}({F}_{n-1})\hfill & \mathbf{0}\hfill & \mathbf{0}\hfill & \dots \hfill & \mathbf{0}\hfill \end{array}\right)$$

**I** is the *D × D* identity matrix, and
${\mathbf{\text{M}}}_{n-l,n}({F}_{n-1})[h,\hspace{0.17em}k]\hspace{0.17em}:=E({Y}_{n-l,n,i}^{(h),k}|{F}_{n-1})$, for *h* = 1, ..., *D*, *k* = 1, ..., *D*. Then, assuming that
$\mathbb{\text{M}}$(*F*) is invertible, for any value *F* of *F _{n−1}*,

In Section 5.1., we give some approaches for studying the behavior of model (4) when the total population size remains bounded. In Sections 5.2. and 5.3., we study the asymptotic behavior of limit models derived from (4) as the initial population size increases to infinity.

Let
${N}_{n}\hspace{0.17em}:={\sum}_{k}{N}_{n}^{k}$ be the total population size at time *n*. Let us first assume that this population size remains bounded by some control. For example the population is a herd of farm animals. We assume that the number of newborns *Y _{n,i}* at time

$$\sum _{{N}_{1}}Q({N}_{1},{N}_{2})\pi ({N}_{1})=\pi ({N}_{2}),\hspace{0.17em}\forall {N}_{2}\in {\mathbb{\text{N}}}^{dD}.$$

But due to combinatorial aspects, the transition probabilities {*Q*(*N*_{1}*, N*_{2})} may be difficult to compute. A solution is then to work in continuous time: in [26], we determine the distribution of the population process from the individual transitions that may be age and population-dependent and that may be all different, and in [40], we use this type of process for studying the propagation of the BVD (Bovine Diarrhoea Virus) within a dairy herd. For that, a renewal approach is used, which generalizes to a population the semi-Markov process theory relative to one individual. This approach also allows to build a rigorous and general simulation algorithm for this individual based model, but when used for non bounded branching populations, it cannot lead to fine stochastic behaviors such as
${\text{lim}}_{n}{N}_{n}{[{N}_{0}{\rho}^{n}]}^{-1}\stackrel{a.s.}{=}\hspace{0.17em}W$ determined by martingale theory in the frame of BGW processes.

We assume here that the initial population size *N*_{0} is large which allows to study the limit, as *N*_{0} → ∞, of the following quantities: **N**_{n}/N_{0} =: * _{n}* (empirical densities) or

A simple theoretical example that we (artificially) apply here on epidemics is the following model on densities with *D* = 1, *d* = 1 [34, 35, 36, 37, 38]:
${I}_{n}={\sum}_{i=1}^{{I}_{n-1}}{Y}_{n,i}$, where the {*Y _{n,i}*}

Let us study now the empirical probabilities in the general case *D* ≥ 1, *d* ≥ 1 under a density-dependence assumption. We prove in Proposition 5 that
${\text{lim}}_{n}{\text{lim}}_{{N}_{0}}{\widehat{\text{P}}}_{n}|{N}_{n}\ne 0\stackrel{P}{=}{\text{lim}}_{{N}_{0}}{\text{lim}}_{n}{\widehat{\text{P}}}_{n}|{N}_{n}\ne 0$ (and similarly for {**W*** _{n}*}), which allows to approximate lim

Let us first study the behavior of the total population size, under the assumption that the number of newborns is independent of the mother health state.

**Proposition 4** *Let us assume that the distribution of*
${\sum}_{k=1}^{D}{Y}_{n-l,n,1}^{(h),k}$ *only depends on l and is denoted*
${\sum}_{k=1}^{D}{Y}_{n-l,n,i}^{(h),k}={Y}_{n-l,n,i}$*. Let us denote*
${\Phi}_{l}\hspace{0.17em}:=E({Y}_{n-l,n,1}|{F}_{n-1}),\hspace{0.17em}{\sigma}_{l}^{2}\hspace{0.17em}:=Var({Y}_{n-l,n,1}|{F}_{n-1})$*. Then the total size of the alive population at time n,*
${N}_{n}={\sum}_{k=1}^{D}{N}_{n}^{k}$*, may be expressed as a multitype branching process* **Ñ*** _{n}* := (

$${N}_{n,j}=\sum _{l=1}^{d}\sum _{i=1}^{{N}_{n-1,l}}{Y}_{n,i}^{(l,j)},$$

*with, for j* = 1,
${Y}_{n,i}^{(l,1)}={Y}_{n-l,n,i}$, *for j* > 1,
${Y}_{n,i}^{(j-1,j)}=1$ *and*
${Y}_{n,i}^{(l,j)}=0$, *for l* ≠ *j* − 1, *and the*
${\left\{{Y}_{n,i}^{(l,j)}\right\}}_{i}$ *are independent given F _{n−1}*.

*Moreover, assuming* 0 < Φ_{l} < ∞, *l* = 1, …, *d*,
$0<{\sum}_{l=1}^{d}{\sigma}_{l}^{2}<\infty $*, and P*(*Y _{n−l,n,i}* ≥ 2

- {
*N*}_{n}*is positively regular, nonsingular, and checks the xlogx property;* *P*(lim_{n}*N*= 0 lim_{n}_{n}*N*= ∞) = 1_{n}*Let ρ be the**first eigenvalue of***M͂***defined by E*(**Ñ**|_{n}**Ñ**_{n}_{−1}) =:**Ñ**_{n}_{−1}**M͂**,*that is,*$$\tilde{\mathbf{\text{M}}}=\left(\begin{array}{lllll}{\Phi}_{1}\hfill & 1\hfill & 0\hfill & \dots \hfill & 0\hfill \\ {\Phi}_{2}\hfill & 0\hfill & 1\hfill & \dots \hfill & 0\hfill \\ \dots \hfill & \dots \hfill & \hfill & \hfill & \hfill \\ {\Phi}_{d-1}\hfill & 0\hfill & 0\hfill & \dots \hfill & 1\hfill \\ {\Phi}_{d}\hfill & 0\hfill & 0\hfill & \dots \hfill & 0\hfill \end{array}\right)$$*Then E*(**Ñ**_{n}*ξ*) =^{t}**Ñ**_{0}*ρ*^{n}*ξ*,^{t}*where***M͂***ξ*=^{t}*ρ**ξ*≤ 1^{t}. Moreover ρ*(subcritical and critical cases) implies that P*(lim_{n}*N*= 0) = 1_{n}*(a.s. extinction), and ρ*> 1*(supercritical case) implies the existence of an integrable random variable W such that*${\text{lim}}_{n\to \infty}{W}_{n}\stackrel{a.s.}{=}W$*, where W*:=_{n}*N*[_{n}*N*_{0}*ρ*]^{n}^{−1}*, and P*(lim_{n}*N*= ∞) =_{n}*P*(*W*> 0)*. Finally let us assume that N*_{−}=_{l}*ρ*^{−l}*N*_{0},*l*= 1, ...,*d*− 1.*Then if ρ*$\mathbb{\text{Q}}$, ${\text{lim}}_{{N}_{0}\to \infty}{W}_{n}\stackrel{a.s.}{=}{\text{lim}}_{{N}_{0}\to \infty}{\text{lim}}_{n\to \infty}{W}_{n}\stackrel{a.s.}{=}1$.*ρ is solution of*${\sum}_{l=1}^{d}{\Phi}_{l{\rho}^{-l}}=1$*, and ρ*≤ 1*R*_{∞}≤ 1*, where*${R}_{\infty}={\sum}_{l=1}^{d}{\Phi}_{l}$*(total mean number of offspring generated by an individual).*

Let us notice that the process {*N _{n}*} is the counterpart in discrete time of a single-type age-dependent Crump-Mode-Jagers branching process in continuous-time [29].

Let us write *α _{l}* := Φ

*A*3:
${\text{lim}}_{n}{\rho}^{-n/2}{\sum}_{k={n}_{d}+1}^{n}{\sum}_{{l}_{1}\hspace{0.17em},\dots ,{l}_{k}}[{\alpha}_{{l}_{1}}{\rho}^{{l}_{1}/2}\dots {\alpha}_{{l}_{k}}{\rho}^{{l}_{k}/2}]{1}_{\left\{{\sum}_{j=1}^{k}{l}_{j}\le n\right\}}<\infty $.

For *d* = 1, or *d* > 1 with * _{l}* = 0, for all

$$\sum _{k={n}_{d}+1}^{n}\sum _{{l}_{1}\hspace{0.17em},\dots ,{l}_{k}}[{\alpha}_{{l}_{1}}\dots {\alpha}_{{l}_{k}}]{1}_{\left\{{\sum}_{j=1}^{k}{l}_{j}\le n\right\}}<\infty ,$$

where, we notice that, for all *k* ≤ *n*,
${\alpha}_{1}^{k}\le {\sum}_{{l}_{1}\hspace{0.17em},\dots ,{l}_{k}}[{\alpha}_{{l}_{1}}\dots {\alpha}_{{l}_{k}}]{1}_{\left\{{\sum}_{j=1}^{k}{l}_{j}\le n\right\}}\le 1$ and
${\sum}_{{l}_{1}\hspace{0.17em},\dots ,{l}_{k}}[{\alpha}_{{l}_{1}}\dots {\alpha}_{{l}_{k}}]{1}_{\left\{{\sum}_{j=1}^{k}{l}_{j}\le n\right\}}$ is decreasing in *k*.

**Proposition 5** *Let us assume, as in Proposition 4, that the distribution of*
${\sum}_{k=1}^{D}{Y}_{n-l,n,1}^{(h),k}$ *depends only on l. Let us moreover assume ρ* > 1, *A*3, **M*** _{n−l,n}*(

$$\begin{array}{l}\underset{{N}_{0}}{\text{lim}}\underset{n}{\text{lim}}P(|{\widehat{\mathbf{\text{P}}}}_{n}-{\mathbf{\text{P}}}_{n}|>\eta |{N}_{n}\ne 0)=\underset{n}{\text{lim}}\underset{{N}_{0}}{\text{lim}}P(|{\widehat{\mathbf{\text{P}}}}_{n}-{\mathbf{\text{P}}}_{n}|>\eta |{N}_{n}\ne 0)=0,\\ \underset{{N}_{0}}{\text{lim}}\underset{n}{\text{lim}}P(|{\mathbf{\text{W}}}_{n}-{\mathbf{\text{P}}}_{n}|>\eta |{N}_{n}\ne 0)=\underset{n}{\text{lim}}\underset{{N}_{0}}{\text{lim}}P(|{\mathbf{\text{W}}}_{n}-{\mathbf{\text{P}}}_{n}|>\eta |{N}_{n}\ne 0)=0,\end{array}$$

*where* **P**_{n}*is the vector of probabilities defined by:*

$${\mathbf{\text{P}}}_{n}\hspace{0.17em}:=({\sum}_{l}{\sum}_{h}{P}_{n-l}^{h}{\rho}^{-l}{\mathbf{\text{M}}}_{n-l,n}({F}_{n-1})[h,\hspace{0.17em}1],\hspace{0.17em}\dots ,{\sum}_{l}{\sum}_{h}{P}_{n-l}^{h}{\rho}^{-l}{\mathbf{\text{M}}}_{n-l,n}({F}_{n-1})[h,D]),\hspace{0.17em}n\ge 1.$$

Proposition 5 allows to use the attractors of the dynamical model {**P*** _{n}*} as an approximation of the asymptotic behavior of the empirical frequencies if the initial population size is very large, and moreover generalizes the result of the BGW process:
${\text{lim}}_{{N}_{0}}{\text{lim}}_{n}{N}_{n}{[{N}_{0}{\rho}^{n}]}^{-1}\stackrel{a.s.}{=}1$ deduced from
${\text{lim}}_{n}{N}_{n}{[{N}_{0}{\rho}^{n}]}^{-1}\stackrel{a.s.}{=}W$ (Section 2.). But the main drawbacks of this approach consist in the loss of the population variability, since the quantities are normalized by

Another approach for studying the asymptotic behavior of the process {**N*** _{n}*} is to study the asymptotic behavior of the limit model assuming now that the disease is rare at the initial time, which forbids the use of densities or probabilities as in Section 5.2. We present here such an approach generalizing the case of the BSE propagation at the scale of a country [28].

Let us recall model (4):
${N}_{n}^{k}={\sum}_{l=1}^{d}{\sum}_{h=1}^{D}{\sum}_{i=1}^{{N}_{n-l}^{h}}{Y}_{n-l,n,i}^{(h),k}$, where the
${\left\{{Y}_{n-l,n,i}^{(h),k}\right\}}_{i,l}$ are mutually independent given *F _{n−}*

$${Y}_{n-l,n,i}^{(h),I\times a}:={1}_{\{a>l\}}{1}_{\{h=S\times a-l\}}{\delta}_{2,n-l+1,n,i}^{E,I|h}+{1}_{\{a=l\}}\sum _{j=1}^{{Y}_{n-l+1,i}^{(h)}}{\delta}_{1,n-l+1,n,i,j}^{E,I|(h)},$$

where we assume that the
${\left\{{\delta}_{2,n-l+1,n,i}^{E,I|h}\right\}}_{i}$ are i.i.d. Bernoulli variables given *F _{n−1}*, the
${\left\{{\delta}_{1,n-l+1,n,i,j}^{E,I|(h)}\right\}}_{i,j}$ are i.i.d. Bernoulli variables given

$${N}_{n}^{I\times a}=\sum _{l}\left[{1}_{\{a>l\}}\sum _{h}\sum _{i=1}^{{N}_{n-l}^{h}}[{1}_{\{h=S\times a-l\}}{\delta}_{2,n-l+1,n,i}^{E,I|h}+{1}_{\{a=l\}}\sum _{h}\sum _{i=1}^{{N}_{Y,n-l+1}^{(h)}}{\delta}_{1,n-l+1,n,i,j}^{E,I|(h)}\right],$$

(5)

where ${N}_{Y,k}^{(h)}\hspace{0.17em}:={\sum}_{i=1}^{{N}_{k-1}^{h}}{Y}_{k,i}^{(h)}$.

Moreover, if a *I* individual may survive in this state a longer time than a time unit, the infection process will depend on the total number of infectives including all the *I* individuals. In this case, let *k* = *I ^{tot}* ×

$${N}_{n}^{{I}^{tot}\times a}=\sum _{l}\sum _{i=1}^{{N}_{n-l}^{I\times a-l}}{\delta}_{n-l,n,i}^{I,{I}^{\mathit{\text{tot}}}|I\times a-l}+{N}_{n}^{I\times a}$$

(6)

where
${\delta}_{n-l,n,i}^{I,{I}^{tot}|I\times a-l}$ is a Bernoulli variable, equal to 1 if *i* survives in the state *I* from age *a − l* at *n − l* to age *a* at *n* at least. Since the distribution of
${\left\{{N}_{n}^{{I}^{tot}\times a}\right\}}_{a,n}$ is easily calculated from the distribution of
${\left\{{N}_{n}^{I\times a}\right\}}_{a,n}$ and since clinical cases are generally rapidly isolated and then removed from the infection process, and since moreover only the new cases are counted by surveillance systems, we will study here the limit of
${N}_{n}^{I\times a}$ given
${F}_{n-1}^{*}\hspace{0.17em}:=\{{\left\{{N}_{k}^{I}\right\}}_{k\le n-1},{\left\{{N}_{k}\right\}}_{k\le n-1},{\left\{{N}_{Y,k}\right\}}_{k\le n}\}$, where
${N}_{Y,k}={\sum}_{h}{N}_{Y,k}^{(h)}$. Of course, expressions similar to (5) and (6) can be written for
${N}_{n}^{E\times a}$ and
${N}_{n}^{{E}^{\mathit{\text{tot}}}\times a}$.

Let *δ _{n−l,h}*(

**Proposition 6** *Let us assume that*
${\text{lim}}_{{N}_{0}}\hspace{0.17em}{N}_{n}^{H\times a}$ *exist, for H * {*I, E, I ^{tot}, E^{tot}*}

$$\begin{array}{l}E({\delta}_{n-l,S\times a-l}(i){\delta}_{2,n-l+1,n,i}^{E,I|S\times a-l}|{F}_{n-1}^{*})=E({\delta}_{n-l,S\times a-l}(i){\delta}_{2,n-l+1,n,i}^{E,I|S\times a-l}|{F}_{n-l}^{*})=:\hspace{0.17em}{p}_{2,a,l|n-l}({F}_{n-l}^{*})\\ \sum _{h}E({\delta}_{n-l,h}(i){\delta}_{1,n-l+1,n,i,j}^{E,I|(h)}|{F}_{n-1}^{*})=\sum _{h}E({\delta}_{n-l,h}(i){\delta}_{1,n-l+1,n,i,j}^{E,I|(h)}|{F}_{n-l}^{*})=:{p}_{1,a,l|n-l}({F}_{n-l}^{*})\\ \underset{{N}_{0}}{\text{lim}}{N}_{n-l}{p}_{2,a,l|n-l}({F}_{n-l}^{*})={\Psi}_{2,a,l|n-l}{N}_{n-l}^{I},\hspace{0.17em}\hspace{0.17em}\underset{{N}_{0}}{\text{lim}}{N}_{Y,n-l+1}{p}_{1,a,l|n-l}({F}_{n-l}^{*})={\Psi}_{1,a,l|n-l}{N}_{n-l}^{I},\end{array}$$

*where* Ψ_{1,a,l|n−l} *and* Ψ_{2,a,l|n−l} *are independent of*
${F}_{n-l}^{*}$*. Then the limit process*
${I}_{a,n}\hspace{0.17em}:\stackrel{\mathcal{D}}{=}{\text{lim}}_{{N}_{0}}{N}_{n}^{I\times a}$ *of incidence of cases aged a at time n is a single-type Markovian process of order d with a Poissonian transition law:*

$${I}_{a,n}|{I}_{n-1},\hspace{0.17em}{I}_{n-2},\hspace{0.17em}\dots \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\sim \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\mathit{\text{Poisson}}(\sum _{l=1}^{d}{\Psi}_{a,l|n-l}{I}_{n-l}),$$

(7)

$${\Psi}_{a,l|n-l}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}:=\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{1}_{\{a>l\}}{\Psi}_{2,a,l|n-l}+{1}_{\{a=l\}}{\Psi}_{1,a,l|n-l},$$

(8)

*where I _{n}* := ∑

As an example of such an approach, we studied the propagation of the BSE in Great-Britain by a general model of type (4). Since the time unit is large (one year), we used Proposition 6 with, for *l* ≥ 2,
${\tilde{F}}_{n-l}^{*}\hspace{0.17em}:=\left\{{\left\{{N}_{k}^{I}\right\}}_{k\le n-l+1,}{\left\{{N}_{k}\right\}}_{k\le n-l,}\hspace{0.17em}{\left\{{N}_{Y,k}\right\}}_{k\le n-l+1}\right\}$ instead of
${F}_{n-l}^{*}$. Then the assumptions of this proposition 6 were satisfied under the following assumptions [28]:

- The
*S*and*E*individuals have the same time-homogeneous survival law {*S*}_{a};_{a} - There is no over-contamination during the incubation period or the clinical state;
- The number of newborn animals ${Y}_{l,i}^{(h)}$ at each birth per individual at time
*l*, is independent of*l*,*i*, and of the health state*h*of*i*(but the health state of each newborn and his survival during the first time unit may depend on*h*); - The population is roughly stable: ${\text{lim}}_{{N}_{0}}{N}_{1}{N}_{0}^{-1}\stackrel{a.s.}{=}1$;
- The disease is rare at the initial time: ${\text{lim}}_{{N}_{0}}{({N}_{0}^{E}+{N}_{0}^{I})}^{2}{N}_{0}^{-1}\stackrel{a.s.}{=}0$;
- The probability for a given
*S*to be infected at time*k*+ 1 via the horizontal route of excretion, follows a Reed-Frost’s type model, that iswhere ${N}_{k}^{{I}_{*}^{\mathit{\text{tot}}}}$ is the number of infectious animals at time$$\begin{array}{l}P(i\hspace{0.17em}\mathit{\text{aged}}\hspace{0.17em}a\hspace{0.17em}\mathit{\text{at}}\hspace{0.17em}k,\hspace{0.17em}\mathit{\text{is}}\hspace{0.17em}\mathit{\text{in}}\hspace{0.17em}\mathit{\text{fected}}\hspace{0.17em}\mathit{\text{at}}\hspace{0.17em}k+1\hspace{0.17em}\mathit{\text{by}}\hspace{0.17em}\mathit{\text{excretion}}|{F}_{k}^{*},\hspace{0.17em}i\hspace{0.17em}\mathit{\text{survives}}\hspace{0.17em}\mathit{\text{at}}\hspace{0.17em}k+1)=\\ E(1-{(1-\frac{{\theta}^{a,\hspace{0.17em}{R}^{c}}}{{N}_{k}})}^{{N}_{k}^{{I}_{*}^{\mathit{\text{tot}}}}}|{F}_{k}^{*})\stackrel{{N}_{k}\hspace{0.17em}\mathit{\text{large}}}{\simeq}\hspace{0.17em}{\theta}^{a,\hspace{0.17em}{R}^{c}}\frac{E({N}_{k}^{{I}_{*}^{\mathit{\text{tot}}}}|{F}_{k}^{*})}{{N}_{k}},\end{array}$$*k*(including those in the latest stage of their incubation period). We assume a similar expression for the infection via the horizontal route concerning contaminated meat and bone meal produced from dead infectious animals.

Under these assumptions, then
${\Psi}_{l|n-l}={\sum}_{a=l}^{d}\left[({\theta}^{a-l,{R}^{c}}+{\theta}^{a-l,R}{\phi}_{n-l})+{1}_{\{a=l\}}{p}_{mat}\right]P(a){P}_{inc.}(l)$, where
${\theta}_{l}^{a-l,{R}^{c}}$ and
${\theta}_{l}^{a-l,R}{(\text{exp}\hspace{0.17em}\lambda -1)}^{-1}{\phi}_{n-l}$ represent the mean numbers of new infected animals aged *a − l* + 1 per infective at time *n − l* + 1 respectively via a horizontal route (excretion of alive animals and slaughtered animals respectively), * _{n−l}* [0, 1] represents the efficiency at time

Let us assume that Φ* _{l|n−l}* depends only on

Let **M͂** be defined by *E*(**Ĩ**_{n}|F_{n−}_{1}) =: **Ĩ**_{n−}_{1}**M͂**, where

$$\tilde{\mathbf{\text{M}}}=\left(\begin{array}{lllll}{\Psi}_{1}\hfill & 1\hfill & 0\hfill & \dots \hfill & 0\hfill \\ {\Psi}_{2}\hfill & 0\hfill & 1\hfill & \dots \hfill & 0\hfill \\ \dots \hfill & \dots \hfill & \hfill & \hfill & \hfill \\ {\Psi}_{d-1}\hfill & 0\hfill & 0\hfill & \dots \hfill & 1\hfill \\ {\Psi}_{d}\hfill & 0\hfill & 0\hfill & \dots \hfill & 0\hfill \end{array}\right)$$

**Proposition 7** {**Ĩ**_{n}_{}} *is a multitype BGW branching process:*

$${I}_{n,j}=\sum _{l=1}^{d}\sum _{i=1}^{{I}_{n-1,l}}{Y}_{n,i}^{(l,j)},$$

*where, for j* = 1,
${Y}_{n,i}^{(l,1)}\hspace{0.17em}:={Y}_{n-l,n,i}$, *for j* > 1,
${Y}_{n,i}^{(j-1,j)}:=1$ *and*
${Y}_{n,i}^{(l,j)}\hspace{0.17em}:=0$, *for l* ≠ *j* − 1*, and the*
${\{{Y}_{n,i}^{(l,j)}\}}_{i,l,j}$ *are independent given F _{n}*

*Moreover E*(**Ĩ**_{n}*ξ** ^{t}*) =

Let s := (*s*_{1}, .., *s _{d}*),

**Proposition 8** *The generating function of* {**Ĩ**_{n}},
${F}_{n}(s)\hspace{0.17em}:=E({s}_{1}^{{I}_{n,1}}\dots {s}_{d}^{{I}_{n,d}})$ *is equal to*

$${F}_{n}(\mathbf{\text{s}})={F}_{n-1}(\mathbf{\text{f}}(\mathbf{\text{s}}))={F}_{n-2}(\mathbf{\text{f}}(\mathbf{\text{f}}(\mathbf{\text{s}})))=\dots ={F}_{0}({\mathbf{\text{f}}}_{n}(\mathbf{\text{s}}))\hspace{0.17em}:={[{f}_{n,1}(\mathbf{\text{s}})]}^{{I}_{0,1}}\dots {[{f}_{n,d}(\mathbf{\text{s}})]}^{{I}_{0,d}},$$

*where* **f*** _{n}*(

As in Proposition 4, the bifurcation parameter of the process is given by the first eigenvalue *ρ* of **M͂**.

**Proposition 9** *ρ is solution of*
${\sum}_{l=1}^{d}{\Psi}_{l}{\rho}^{-1}=1$*, and ρ* ≤ 1 *is equivalent to R _{∞}* ≤ 1

Let us assume from now on that Ψ_{1} > 0,..., Ψ* _{d}* > 0.

**Proposition 10**

*P*(lim_{n}*I*= 0 lim_{n}_{n}*I*= ∞) = 1;_{n}*ρ*> 1*is equivalent to R*_{∞}> 1*(supercritical case) which itself implies the existence of a positive integrable random variable W such that*${\text{lim}}_{n\to \infty}{I}_{n}{[{I}_{0}{\rho}^{n}]}^{-1}\stackrel{a.s.}{=}W$,*where P*(lim_{n}*I*= ∞) =_{n}*P*(*W*> 0)*. Moreover let us assume that I*_{−}=_{l}*c*,_{l}I_{0}*l*= 1, ...,*d*− 1, where*c*is independent of_{l}*I*_{0}. Then ${\text{lim}}_{{I}_{0}\to \infty}{I}_{n}{[{I}_{0}{\rho}^{n}]}^{-1}\stackrel{a.s.}{=}{\text{lim}}_{{I}_{0}\to \infty}{\text{lim}}_{n\to \infty}{I}_{n}{[{I}_{0}{\rho}^{n}]}^{-1}\stackrel{a.s.}{=}1$.*ρ*≤ 1*is equivalent to R*_{∞}≤ 1*(subcritical and critical cases) which implies P*(lim_{n}*I*= 0) = 1_{n}*(a.s. extinction).*

Let us point out that the deterministic model derived from *E*(**Ĩ**_{n}|**Ĩ**_{n−}_{1)} = **Ĩ**_{n−}_{1} **M͂** is defined by: **X͂*** _{n}* :=

Let *T _{ext.}* be the (random) extinction time of {

**Proposition 11** *We have*
$q={q}_{1}^{{I}_{0}}\dots {q}_{d}^{{I}_{-(d-1)}}$*, where*
${q}_{h}=\text{exp}(-{\sum}_{l=h}^{d}{\Psi}_{l}(1-{q}_{1}))$*, is the extinction probability starting from* **Ĩ**_{0} = (0*, ..,* 0*,* 1*,* 0*, ...,* 0) *(one individual of the h type), h* = 1*, ..., d.*

As a consequence, in the particular case **Ĩ**_{0} = (*I*_{0}*,* 0*, ...,* 0), then
$q={q}_{1}^{{I}_{0}}$, where *q*_{1} is solution of the equation: 1 *− q*_{1} = 1 − exp(*−R _{∞}*(1

Let *G _{n}* :=

**Proposition 12** *For any value of R _{∞}, the extinction time distribution satisfies, for n* ≥

$${G}_{n}=E(E({1}_{\{{\tilde{\mathbf{\text{I}}}}_{n}=0\}}|{F}_{n-d}))=E(\text{exp}(-\sum _{l=0}^{d-1}{I}_{n-d-l}\sum _{h=l+1}^{d}{\Psi}_{h})).$$

Let us notice that *G _{n}* decreases (and

In addition, using Proposition 8, we get the following result, allowing an exact iterative computation of {*G _{n}*}:

**Proposition 13** *For any value of R _{∞}, the distribution* {

$${G}_{n}=\prod _{l=2}^{d}{[{f}_{n-1,l}({\mathbf{\text{s}}}_{0})]}^{{I}_{0,l-1}}\text{exp}[-\sum _{l=1}^{d}{\Psi}_{l}{I}_{0,l}(1-{(\frac{{G}_{n-1}}{{\prod}_{l=2}^{d}{[{f}_{n-1,l}({\mathbf{\text{s}}}_{0})]}^{{I}_{0,l}}})}^{1/{I}_{0,1}})],$$

*where* **s**_{0} = (0*, ...,* 0) *and* **f*** _{n}*(

Let
$N:={\sum}_{l=-(d-1)}^{Text.-1}{I}_{l}$ be the size of the tree generated by **Ĩ**_{0} (epidemic size). Let us recall that the Borel-Tanner distribution with parameter (*λ, l*) is a (*f _{λ}, s^{l}*) GLPD (Generalized Lagrange Probability Distribution), where

Let *g*(*s*) := *E*(*s ^{N}*) be the generating function of

**Proposition 14** *Let us assume that*
${\sum}_{l=1}^{d}{\Psi}_{l}<1$ *(subcritical case).*

*Let***Ĩ**_{0}:= (*I*_{0}, 0, ..., 0)*. Then*$g(s)={s}^{{I}_{0}}{f}_{{I}_{0}{\sum}_{l=1}^{d}{\Psi}_{l}}(g(s))$*, that is*$N\sim \mathit{\text{Borel}}-\mathit{\text{Tanner}}({\sum}_{l=1}^{d}{\Psi}_{l},{I}_{0})$:$$P(N=k)=\text{exp}(-k\sum _{l=1}^{d}{\Psi}_{l})\frac{{I}_{0}}{k}\frac{{(k{\sum}_{l=1}^{d}{\Psi}_{l})}^{k-{I}_{0}}}{(k-{I}_{0})!},\hspace{0.17em}k\ge {I}_{0}$$*and*$E(N)={I}_{0}{(1-{\sum}_{l=1}^{d}{\Psi}_{l})}^{-1}$, $\mathit{\text{Var}}(N)={I}_{0}({\sum}_{l=1}^{d}{\Psi}_{l}){(1-{\sum}_{l=1}^{d}{\Psi}_{l})}^{-3}$.*Let***Ĩ**_{0}:= (*I*_{0},*I*_{−1}, , ...,*I*_{−(}_{d}_{−1)})*. Then*$N\stackrel{\mathcal{D}}{=}{\oplus}_{h=1}^{d}{\oplus}_{i=1}^{{I}_{-(h-1)}}({\oplus}_{j=1}^{{Y}_{0,h,i}}{N}_{i,j}+1)$,*where*${Y}_{0,h,i}\sim \mathit{\text{Poisson}}({\sum}_{l=h}^{d}{\Psi}_{l})$,*the*{*N*}_{i,j}*are i.i.d. with*${N}_{i,j}\sim \mathit{\text{Borel}}-\mathit{\text{Tanner}}({\sum}_{l=1}^{d}{\Psi}_{l},\hspace{0.17em}1)$*and**means the mutual independence, that is*$$E({s}^{N})=\prod _{h=1}^{d}{s}^{{I}_{-(h-1)}}{\left[{f}_{{\sum}_{l=h}^{d}{\Psi}_{l}}(g(s))\right]}^{{I}_{-(h-1)}},$$*and*$$\begin{array}{lll}E(N)& =& \sum _{h=1}^{d}{I}_{-(h-1)}[\sum _{l=h}^{d}{\Psi}_{l}{(1-\sum _{l=1}^{d}{\Psi}_{l})}^{-1}+1]\\ Var(N)& =& \sum _{h=1}^{d}{I}_{-(h-1)}\sum _{l=h}^{d}{\Psi}_{l}{(1-\sum _{l=1}^{d}{\Psi}_{l})}^{-3}\end{array}$$

We presented a general class of branching processes in discrete time for modeling in a stochastic way some diseases propagation when the infected period is long respectively to the time frequency of births. However when the transitions are population dependent, the long-term prediction of these processes is an open problem in the general case. We indirectly solved this problem by studying the behavior of the limit models, as the total initial population size increases to infinity, assuming at the initial time, either a non-rare disease with a density-dependence assumption, or a rare disease.

Under the first assumption, since
${\text{lim}}_{n}{\text{lim}}_{{N}_{0}}{\widehat{\mathbf{\text{P}}}}_{n}|{N}_{n}\ne 0\stackrel{P}{=}{\text{lim}}_{{N}_{0}}{\text{lim}}_{n}{\widehat{\mathbf{\text{P}}}}_{n}|{N}_{n}\ne 0$, we proved, for a large *N*_{0}, that the proportion of infected individuals in the whole population,
${\widehat{P}}_{n}^{E}+{\widehat{P}}_{n}^{I}$, behaves, as *n* → ∞, as the probability
${P}_{n}^{E}+{P}_{n}^{I}$ for an individual in an infinite population to be infected. Under the second assumption, we got a branching process on the incidence of cases which has the advantage to correspond to the usual observations, to be rigorously built, starting from a detailed multitype branching process {**N*** _{n}*} taking into account the different disease steps together with the population dynamic, to keep the population variability, and to belong to the simple class of multitype BGW processes for which many analytical results exist and that we used or generalized to this model. We calculated the distribution of the extinction time and the distribution of the epidemic size. Moreover we proved that the bifurcation parameter of this process was the total mean number of secondary cases who began to be generated by one case during his first time unit.

This result would validate and generalize the current use in epidemiology of the reproductive number as a bifurcation parameter [8, 17, 18], if we could establish that
${\text{lim}}_{n}{\text{lim}}_{{N}_{0}}{\mathbf{\text{N}}}_{n}^{I}={\text{lim}}_{{N}_{0}}{\text{lim}}_{n}{\mathbf{\text{N}}}_{n}^{I}$. We studied the left quantity, while the valid realistic quantity is in reality the right one. In fact this equality is in general not true because the bell form of an epidemic can be modelled only by a population dependent process. But for large *N*_{0}, the stochastic limit model is a good approximation of the epidemic growth (in the supercritical case), or of the epidemic decay (in the subcritical case). In any case, if *N*_{0} is too small, the limit in *N*_{0} cannot be used and therefore the limit models developed here cannot be used.

Let us finally notice that we proved that in a size-dependent model on the clinical cases, then the quantity determining the extinction of the process was the total mean number of secondary cases that will be produced in the future by a case as the whole current population is infected (*R _{∞}*), while the corresponding quantity in the associate deterministic model derived from the conditional expectation of the process is the total mean number of secondary cases that will be produced by a case as the whole current population is susceptible (

1. Athreya KB, Ney PE. Branching Processes. Springer-Verlag; Berlin, Germany: 1972.

2. Becker NG. On parametric estimation for mortal branching processes. Biometrika. 1974;61:393–399.

3. Becker NG. Estimation for an epidemic model. Biometrics. 1976;32:769–777. [PubMed]

4. Becker NG. Estimation for discrete time branching processes with applications to epidemics. Biometrics. 1977;33:515–522. [PubMed]

5. Becker NG. A general chain binomial model for infectious diseases. Biometrics. 1981;37:251–258. [PubMed]

6. Becker NG. Analysis of Infectious Disease Data. Chapman and Hall/CRC Press; Boca Raton, FL, USA: 1989.

7. Devroye L. The branching process method in Lagrange random variate generation. Commun. Statistics-Simulat. Comput. 1992;21:1–14.

8. Diekmann O, Heesterbeek JAP. Mathematical Epidemiology of Infectious Diseases Model Building, Analysis and Interpretation Wiley Series in Mathematical and Computational Biology. John Wiley Ltd; Chichester, UK: 2000

9. Elaydi S. An Introduction to Difference Equations 3rd ed Undergraduate Texts in Mathematics. Springer; New York, NY, USA: 2005

10. Farrington CP, Grant AD. The distribution of time to extinction in subcritical branching processes: Applications to outbreaks of infectious disease. J. Appl. Probab. 1999;36:771–779.

11. Farrington CP, Kanaan MN, Gay NJ. Branching process models for surveillance of infectious diseases controlled by mass vaccination. Biostatistics. 2003;4:279–295. [PubMed]

12. Good IJ. The Lagrange distributions and branching processes. SIAM J. Appl. Math. 1975;28:270–275.

13. Haccou P, Jagers P, Vatutin VA. Branching Processes Variation, Growth, and Extinction of Populations. Cambridge Studies in Adaptative Dynamics; Cambridge, UK: 2005.

14. Haight FA, Breuer AB. The Borel-Tanner distribution. Biometrika. 1960;47:143–150.

15. Hall P, Heyde CC. Martingale Limit Theory and its Application. Probability and Mathematical Statistics; New York, NY, USA: 1980.

16. Harris TE. The Theory of Branching Processes Die Grundlehren der Mathematischen Wissenschaften, Bd 119 Springer-Verlag; Berlin, Germany: Prentice-Hall, Inc; Englewood Cliffs, N.J., USA: 1963

17. Heesterbeek JAP, Dietz K. The concept of Ro in epidemic theory. Stat. Neerl. 1996;50:89–110.

18. Hefferman JM, Smith RJ, Wahl LM. Perspectives on the basic reproductive ratio. J. R. Soc. Interface. 2005;2:281–293. [PMC free article] [PubMed]

19. Hethcote HW. The Mathematics of Infectious Diseases. SIAM Rev. 2000;42:599–653.

20. Heyde CC. On assessing the potential severity of an outbreak of a rare infectious disease: a Bayesian approach. Austral. J. Statist. 1979;21:282–292.

21. Jacob C, Viet AF. Epidemiological modeling in a branching population. Particular case of a general SIS model with two age classes. Math. Biosci. 2003;182:93–111. [PubMed]

22. Jacob C, Magal P. Influence of the routine slaughtering on the evolution of BSE. Example of the British and French slaughterings. Risk Anal. 2007;27:1151–1167. [PubMed]

23. Jacob C. Branching Processes: Their rôle in Epidemiology Technical report, UR341, MIA department, INRA; Jouy-en-Josas, France: 2009

24. Jagers P. Branching Processes with Biological Applications. Wiley Series in Probability and Mathematical Statistics, John Wiley; London, UK: 1975.

25. Jacob C, Lalam N, Yanev N. Statistical inference for processes depending on environments and application in regenerative processes. Pliska Stud. Math. Bulgar. 2005;17:109–136.

26. Jacob C, Viet AF. A new class of processes for formalizing and generalizing individual-based models: the semi-semi-Markov processes. Pliska Stud Math Bulgar. 2007:121–144.

27. Jacob C. Saturation Effects in Population Dynamics: Use Branching Processes or Dynamical Systems? In: Deutsch A, Bravo de la Parra R, de Boer R, Dieckmann O, Jagers P, Kisdi E, Kretzschmar M, Lansky P, Metz H, editors. Mathematical Modeling of Biological Systems. Birkhäuser; Boston, MA, USA: 2008. pp. 339–352.

28. Jacob C, Maillard-Teyssier L, Denis JB, Bidot C. A branching process approach for the propagation of the Bovine Spongiform Encephalopathy in Great-Britain Accepted in *Branching processes and their Applications*, Lecture Notes in Statistics, Springer-Verlag; Berlin, Germany: 2010. 227–242.242

29. Jagers P. A general stochastic model for population development. Skan. Aktuarietidskr. 1969;52:84–103.

30. Klebaner FC. Population-size-dependent branching process with linear rate of growth. J. Appl. Probab. 1983;20:242–250.

31. Klebaner FC. On population-size-dependent branching processes. Adv. Appl. Probab. 1984;16:30–55.

32. Klebaner FC. Linear growth in near-critical population-size-dependent multitype Galton-Watson processes. J. Appl. Probab. 1989;26:431–445.

33. Klebaner FC. Correction: “Linear growth in near-critical population-size-dependent multitype Galton-Watson processes” [J. Appl. Probab. 26 (1989), no. 3, 431–445] and “Asymptotic behavior of near-critical multitype branching processes” [ibid. 28 (1991), no. 3, 512–519] J. Appl. Probab. 1992;29:246.

34. Klebaner FC. Population-dependent branching processes with a threshold. Stochastic Process. Appl. 1993;46:115–127.

35. Klebaner FC, Nerman O. Autoregressive approximation in branching processes with a threshold. Stochastic Process. Appl. 1994;51:1–7.

36. Klebaner FC, Zeitouni O. The exit problem for a class of density-dependent branching systems. Ann. Appl. Probab. 1994;4:1188–1205.

37. Klebaner FC. Population and density dependent branching processes Classical and Modern Branching Processes IMA Vol Math Appl, 84, Springer; New York, NY, USA: 1997. 165–169.169

38. Klebaner FC, Lazar J, Zeitouni O. On the quasi-stationary distribution for some randomly perturbed transformations of an interval. Ann. Appl. Probab. 1998;8:300–315.

39. Yanev NM, Mitov KV. Regenerative Branching Processes Records and Branching processes Ahsanullah M, Yanev GP, editors. Nova Science Publishers, Inc; New York, NY, USA: 2008. Chapter 3; 37–62.62

40. Viet AF, Fourichon C, Seegers H, Jacob C, Guihenneuc-Jouyaux C. A model of the spread of the bovine viral-diarrhoea virus within a dairy herd. Prev. Vet. Med. 2004;63:211–236. [PubMed]

Articles from International Journal of Environmental Research and Public Health are provided here courtesy of **Multidisciplinary Digital Publishing Institute (MDPI)**