C. R. Biol. 341 (2018) p. 145–151
https://doi.org/10.1016/j.crvi.2018.01.009
Nicolas Bacaër
Institut de recherche pour le développement
Unité de modélisation mathématique et informatique des systèmes complexes
Les Cordeliers, Paris, France
nicolas.bacaer@ird.fr
On étudie le taux d'extinction d'un processus de branchement avec plusieurs types en temps continu dans un environnement aléatoire. Des calculs numériques dans un exemple particulier inspiré d'un modèle épidémique suggèrent une formule explicite pour ce taux d'extinction mais seulement pour certaines valeurs des paramètres.
On suppose que l'environnement, noté k, oscille de manière aléatoire entre un nombre fini d'états 1,…,K selon une chaîne de Markov en temps continu. Pour k≠h, la probabilité pour que l'environnement bascule de h vers k est \(Q_{k,h}\, dt\) pendant chaque intervalle de temps infinitésimal dt, avec \(Q_{k,h}\geq 0\). On définit la matrice \(\mathbf{Q}=(Q_{k,h})\) avec \[Q_{h,h}=-\sum_{k\neq h} Q_{k,h}\quad \forall h.\] C'est la matrice transposée du générateur infinitésimal de la chaîne (Sericola, 2013).
On considère une population d'individus qui peuvent être de plusieurs types : 1, 2, …, J. Cette population évolue dans l'environnement aléatoire que l'on vient de décrire. On suppose qu'il y a au moins un individu dans la population au temps \(t=0\). Un individu de type i dans l'environnement k a une probabilité \(c_i^{(k)}\, dt\) de subir quelque événement pendant chaque intervalle de temps infinitésimal dt, avec \(c_i^{(k)} > 0\). Si l'événement se produit, on trouve à la place de cet individu \(n_j\) individus de type j pour 1≤j≤J avec une probabilité \(\pi_i^{(k)}(n_1,\ldots,n_J)\). Autrement dit, entre deux sauts de l'environnement, chaque individu a un « temps de vie » aléatoire qui suit une loi exponentielle de paramètre \(c_i^{(k)}\). Au bout de ce temps, la loi de reproduction est donnée par \((\pi_i^{(k)}(n_1,\ldots,n_J))\), indépendamment des autres individus. C'est donc un processus de branchement avec plusieurs types en temps continu (Méléard, 2016).
On définit \(p^{(k)}(t,n_1,\ldots,n_J)\), la probabilité que la population soit composée de \(n_i\) individus de type i pour 1 ≤ i ≤ J et que l'environnement soit k au temps t. On range les états \((k,n_1,\ldots,n_J)\) du système par groupes selon le nombre total d'individus \(n=n_1+\cdots+n_J\), de manière à avoir un vecteur colonne infini \(\mathbf{p}(t)\). On observe dans la section 2 que \(\mathbf{p}(t)\) est solution d'un système linéaire d'équations différentielles \(\frac{d\mathbf{p}}{dt}=\mathbf{Z} \, \mathbf{p}(t)\). \(\mathbf{Z}\) est une matrice infinie de la forme \begin{equation} \left (\begin{array}{c|ccccc} \mathbf{Z}_{0,0} & \mathbf{Z}_{0,1} & \mathbf{0} & \cdots&\mathbf{0}&\cdots\\ \hline \mathbf{0} & \mathbf{Z}_{1,1} & \mathbf{Z}_{1,2}& \ddots &\vdots &\\ \mathbf{0} & \mathbf{Z}_{2,1} & \mathbf{Z}_{2,2} & \ddots&\mathbf{0}&\\ \vdots & \vdots & \vdots & \ddots & \mathbf{Z}_{n-1,n}&\ddots\\ \mathbf{0} & \mathbf{Z}_{n,1} & \mathbf{Z}_{n,2} & \ddots & \mathbf{Z}_{n,n}&\ddots\\ \vdots & \vdots& \vdots& \ddots& \vdots&\ddots \end{array} \right ) \tag{1} \end{equation} et les blocs \(\mathbf{Z}_{m,n}\) sont eux-mêmes des matrices de tailles différentes. On voit sur la structure de cette matrice que lorsqu'un individu subit un événement, le nombre total d'individus ne peut décroître que d'une unité mais il peut croître de plusieurs unités. De plus, la classe des états avec zéro individu est absorbante : elle correspond à l'extinction de la population. On se restreint dans la suite au cas sous-critique où la population s'éteint presque sûrement. Un résultat de (Athreya et Karlin, 1971) relatif aux modèles en temps discret permet de déterminer à quelle condition il y a extinction dans notre modèle en temps continu. L'objectif est alors d'essayer de déterminer le taux d'extinction de la population.
Les articles (Dyakonova, 2008 ; Dyakonova, 2013 ; Vatutin et Wachtel, 2017) ont calculé ce taux d'extinction dans un modèle analogue mais en temps discret, où les environnements successifs sont aléatoires, indépendants et identiquements distribués. Les deux premières références donnent une formule simple pour le taux d'extinction mais avec des conditions assez restrictives sur les différents environnements (les matrices moyennes doivent avoir un vecteur propre commun). (Vatutin et Wachtel, 2017) donne une formule moins explicite mais avec des hypothèses plus générales, en restant néanmoins dans le cas « fortement sous-critique ».
Le taux d'extinction dépend des propriétés spectrales de la sous-matrice infinie extraite de (1) avec \(\mathbf{Z}_{1,1}\) dans son coin supérieur gauche. On définit \begin{equation} M_{i,j}^{(k)}=c_j^{(k)} \left ( \sum_{n_1,\ldots,n_J\geq 0} n_i\, \pi_j^{(k)}(n_1,\ldots,n_J) -\delta_{i,j} \right ), \tag{2} \end{equation} avec \(\delta_{i,j}=1\) si \(i=j\), \(\delta_{i,j}=0\) si \(i\neq j\). On définit
Dans la section 3, on considère tout d'abord le cas particulier des processus de naissance et de mort à plusieurs types, puis on se restreint aux populations avec seulement deux types d'individus. On présente un exemple où l'on compare la valeur numérique de \(\omega_1\) avec la limite \(\alpha_1\) si \(n\to \infty\) de la borne spectrale de la sous-matrice finie de (1) avec \(\mathbf{Z}_{1,1}\) dans son coin supérieur gauche et \(\mathbf{Z}_{n,n}\) dans son coin inférieur droit. Les résultats numériques suggèrent l'égalité de ces deux nombres pour certaines valeurs des paramètres mais pas pour toutes, comme dans le cas où il n'y a qu'un seul type d'invidu (Bacaër, 2017c). On conjecture par ailleurs que la limite \(\alpha_1\) est bien le taux d'extinction de la population, défini par exemple comme étant la limite \[ \lim_{t \to +\infty} \frac{1}{t} \log p^{(k)}(t,n_1,\ldots,n_J)\, .\] D'après (Collet et coll., 2013, section 4.5), cette limite (appelée paramètre de Kingman) ne dépend ni de k, ni de \((n_1,\ldots,n_J)\) pourvu que \(n_1+\cdots+n_J\geq 1\), ni des conditions initiales (environnement et nombre d'individus des différents types).
Notations :
On définit \(\mathbf{1}=(1,\ldots,1)\) et \(\mathbf{x}=(x_1,\ldots,x_J)\). On définit les fonctions génératrices \[f^{(k)}(t,\mathbf{x})=\sum_{\mathbf{n}\geq \mathbf{0}} p^{(k)}(t,\mathbf{n})\, x_1^{n_1}\ldots x_J^{n_J},\] où les \(x_i\) sont des nombres complexes et les indices \(n_i\) des entiers. Parce que \[\sum_k \sum_{n_i\geq 0} p^{(k)}(t,n_1,\ldots,n_J)=1,\] le domaine de convergence de ces séries inclut l'ensemble \(\{(x_1,\ldots,x_J);\ |x_i| < 1 \ \forall i\}\) (Cartan, 1961, chapitre IV). On définit \[g_j^{(k)}(\mathbf{x})=\sum_{\mathbf{n}\geq \mathbf{0}} \pi^{(k)}_j(\mathbf{n})\ x_1^{n_1}\ldots x_J^{n_J}\, . \] On a \(g_j^{(k)}(\mathbf{1})=1\). On a alors \[\frac{\partial f^{(k)}}{\partial x_j}(t,\mathbf{x})=\sum_{\mathbf{r}\geq \mathbf{0}} (r_j+1)\ p^{(k)}(t,\mathbf{r}+\mathbf{u}_j)\ x_1^{r_1}\ldots x_J^{r_J}.\] et \[\frac{\partial f^{(k)}}{\partial t}(t,\mathbf{x}) =\sum_{\mathbf{n}\geq \mathbf{0}} \frac{dp^{(k)}}{dt}(t,\mathbf{n})\, x_1^{n_1}\ldots x_J^{n_J}\, .\] Avec le système (3), on obtient \begin{equation} \frac{\partial f^{(k)}}{\partial t}=\sum_{h=1}^K Q_{k,h} \, f^{(h)}+\sum_{j=1}^J c_j^{(k)} (g_j^{(k)}(x_1,\ldots,x_J)-x_j) \frac{\partial f^{(k)}}{\partial x_j}\, . \tag{4} \end{equation} C'est une généralisation du système (3) dans (Bacaër, 2017c), qui correspond à \(J=1\).
On définit les valeurs moyennes \[E_i^{(k)}(t)=\frac{\partial f^{(k)}}{\partial x_i}(t,\mathbf{1})=\sum_{\mathbf{n}\geq \mathbf{0}} n_i\, p^{(k)}(t,\mathbf{n})\, .\] Avec \(g_j^{(k)}(\mathbf{1})=1\), on déduit de l'équation (4), en prenant sa dérivée partielle par rapport à \(x_i\) puis en prenant \(\mathbf{x}=\mathbf{1}\), \[\frac{dE_i^{(k)}}{dt} = \sum_{h=1}^K Q_{k,h} E_i^{(h)} + \sum_{j=1}^J M_{i,j}^{(k)} \, E_j^{(k)}\, ,\] avec \begin{equation} M_{i,j}^{(k)}=c_j^{(k)} \left ( \frac{\partial g_j^{(k)}}{\partial x_i}(\mathbf{1})-\delta_{i,j} \right ) \tag{5} \end{equation} comme dans l'introduction, et où l'on suppose que \(M_{i,j}^{(k)} < +\infty\) pour toutes les valeurs des indices. C'est la généralisation de l'équation (4) de (Bacaër, 2017c), qui correspondait à J=1. On définit E, le vecteur colonne \[(E_1^{(1)},\ldots,E_J^{(1)},\ldots,E_1^{(K)},\ldots,E_J^{(K)}).\] On a alors \(\frac{d\mathbf{E}}{dt} = \mathbf{\Omega}^{[1]} \mathbf{E}\) avec \[\mathbf{\Omega}^{[1]}= \left ( \begin{array}{cccc} Q_{1,1} \mathbf{I}+\mathbf{M}^{(1)} & Q_{1,2} \mathbf{I} & \cdots & Q_{1,K} \mathbf{I}\\ Q_{2,1} \mathbf{I} & Q_{2,2} \mathbf{I}+\mathbf{M}^{(2)} &\ddots & \vdots\\ \vdots & \ddots & \ddots & Q_{K-1,K}\mathbf{I}\\ Q_{K,1} \mathbf{I} & \cdots & Q_{K,K-1} \mathbf{I}& Q_{K,K} \mathbf{I} + \mathbf{M}^{(K)} \end{array} \right )\] comme dans l'introduction. On rappelle que \(\omega_1\) est la borne spectrale de cette matrice, c'est-à-dire la valeur propre de plus grande partie réelle. Parce que \(Q_{k,h}\geq 0\ \forall k\neq h\) et \(M^{(k)}_{i,j}\geq 0\ \forall i\neq j\), on remarque en effet que \(\mathbf{\Omega}^{[1]}\) est une matrice dont les coefficients en dehors de la diagonale sont tous ≥ 0. D'après un corollaire du théorème de Perron et Frobenius, \(\mathbf{\Omega}^{[1]}\) a bien une valeur propre réelle dominante, c'est-à-dire supérieure à la partie réelle de toutes les autres valeurs propres.
On suppose pour simplifier que la matrice Q est irréductible: \(\forall k\neq h\), \(\exists (k_0,\ldots,k_N)\), \(k_0=k\), \(k_N=h\), \(Q_{k_n,k_{n+1}} > 0\ \forall n\). On suppose de plus que les matrices \(\mathbf{M}^{(k)}\) sont toutes irréductibles. Alors la matrice \(\mathbf{\Omega}^{[1]}\) est elle aussi irréductible. Avec \(\mathbf{E}(0)\neq 0\), on a \[\frac{1}{t} \log \| \mathbf{E}(t) \| \mathop{\longrightarrow}_{t \to +\infty} \omega_1\, \] et \(\omega_1\) est le taux de croissance ou de décroissance du vecteur des espérances \(\mathbf{E}(t)\).
Considérons une suite fixée d'environnements engendrés par la chaîne de Markov en temps continu de matrice Q : l'environnement est d'abord \(k_0\) pour \(t_0 < t < t_1\) avec \(t_0=0\), puis \(k_1\) pour \(t_1 < t < t_2\), etc. Entre deux sauts de l'environnement, la population évolue selon un processus de branchement en temps continu avec plusieurs types dans un environnement constant. À notre processus en temps continu, on peut donc associer un processus en temps discret qui ne considère que l'état de la population aux instants \(t_n\) où l'environnement bascule. Les deux processus sont simultanément surcritiques, critiques ou sous-critiques.
Le vecteur des espérances des populations de chaque type, \(\mathbf{e}(t)=(e_1(t),\ldots,e_J(t))\), sachant que l'environnement est \(k_n\) pour \(t_n < t < t_{n+1}\), est solution de \(\frac{d\mathbf{e}}{dt}=\mathbf{M}^{(k_n)} \mathbf{e}(t)\ \) pendant cet intervalle de temps. Avec \(\tau_n=t_{n+1}-t_n\), on a \(\mathbf{e}(t_{n+1})=\exp(\tau_{n}\, \mathbf{M}^{(k_n)})\, \mathbf{e}(t_n)\). D'après (Athreya et Karlin, 1971, section 4), la suite \begin{equation} \frac{1}{n} \log \| \exp(\tau_n \mathbf{M}^{(k_n)}) \cdots \exp(\tau_0 \mathbf{M}^{(k_0)}) \| \tag{6} \end{equation} converge presque sûrement vers une limite \(\lambda_1\) indépendante de la suite particulière d'environnements. De plus, la population dans le processus en temps discret s'éteint presque sûrement dans le cas sous-critique où \(\lambda_1 < 0\) et il n'y a pas extinction avec une probabilité strictement positive si \(\lambda_1 > 0\). Ainsi le processus en temps continu dans notre modèle de départ est aussi sous-critique si \(\lambda_1 < 0\). Si T est la limite de \(t_n/n\) si \(n\to +\infty\), on remarque que \(\lambda_1/T\) est l'exposant de Lyapounoff du système différentiel pour \(\mathbf{e}(t)\).
Notons enfin que la borne spectrale \(\omega_1\) de la section précédente peut être positive alors que \(\lambda_1\) est négatif. C'est déjà possible lorsqu'il n'y a qu'un seul type d'individus (Bacaër et Ed-Darraz, 2014).
Pour le système (4), on cherche des solutions \[f^{(k)}(t,\mathbf{x})=e^{\omega t} F^{(k)}(\mathbf{x})\] avec des fonctions \(F^{(k)}(x)\) qui ne sont pas toutes identiquement nulles. On a ainsi \begin{equation} \omega\, F^{(k)}(\mathbf{x})=\sum_{h=1}^K Q_{k,h}\, F^{(h)}(\mathbf{x})+\sum_{j=1}^J c_j^{(k)} (g_j^{(k)}(\mathbf{x})-x_j) \frac{\partial F^{(k)}}{\partial x_j}(\mathbf{x})\, . \tag{7} \end{equation} En prenant \(\mathbf{x}=\mathbf{1}\), on voit que \[\omega\, F^{(k)}(\mathbf{1})=\sum_{h=1}^K Q_{k,h}\, F^{(h)}(\mathbf{1}).\] Il y a deux cas possibles :
Si les fonctions \(F^{(k)}(\mathbf{x})\) sont analytiques dans un voisinage de \(\mathbf{x}=\mathbf{1}\) (c'est le « cas régulier»), alors on voit comme dans la section 2.3, en dérivant l'équation (7) par rapport à \(x_i\) et en prenant \(\mathbf{x}=\mathbf{1}\), \[\omega\, \phi_i^{(k)} = \sum_{h=1}^K Q_{k,h}\, \phi_i^{(h)} + \sum_{j=1}^J M_{i,j}^{(k)} \phi_j^{(k)}\, ,\] avec \(\phi_i^{(k)}=\frac{\partial F^{(k)}}{\partial x_i}(\mathbf{1})\). Donc il y a aussi deux cas possibles :
On définit \[\forall n\geq 1,\quad \forall (i_1,\ldots,i_n)\in \{1,\ldots,J\}^n,\quad \phi_{i_1,i_2,\ldots,i_n}^{(k)}=\frac{\partial^{n} F^{(k)}}{\partial x_{i_1} \partial x_{i_2} \cdots \partial x_{i_n}}(\mathbf{1}).\] On choisit \(n\geq 2\). On suppose que \(\ \phi_{i_1,i_2,\ldots,i_m}^{(k)}=0\ \forall \, 1\leq k\leq K\), \(1\leq m < n\) et \((i_1,\ldots,i_m)\in \{1,\ldots,J\}^m\). On dérive l'équation (7) par rapport à \(x_{i_1},\ldots,x_{i_n}\) et on prend \(\mathbf{x}=\mathbf{1}\). À cause de l'hypothèse sur les dérivées partielles d'ordre < n et parce que \(\ g_j^{(k)}(\mathbf{1})=1\), il ne reste que les termes suivants: \begin{equation} \omega\, \phi_{i_1,i_2,\ldots,i_n}^{(k)} = \sum_{j=1}^J \left [M^{(k)}_{i_1,j}\, \phi_{j,i_2,\ldots,i_n}^{(k)} + \cdots + M^{(k)}_{i_n,j}\, \phi_{i_1,i_2,\ldots,j}^{(k)} \right ]+\sum_{h=1}^K Q_{k,h}\, \phi_{i_1,i_2,\ldots,i_n}^{(h)} \quad \forall (i_1,\ldots,i_n) \in \{1,\ldots,J\}^n,\quad 1\leq k\leq K. \tag{8} \end{equation} Donc il y a deux cas possibles :
En résumé, on conclut que si les fonctions propres \(F^{(k)}(\mathbf{x})\) associées à la valeur propre ω sont analytiques dans un voisinage de \(\mathbf{x}=\mathbf{1}\), alors ω est une valeur propre de la matrice Q, ou bien ω est une valeur propre d'une matrice \(\mathbf{\Omega}^{[n]}\) pour un certain n ≥ 1. En effet, si tel n'était pas le cas, on aurait \[F^{(k)}(\mathbf{1})=0,\quad \phi_i^{(k)}=0\quad \forall i,\ \forall k.\] D'après ce qui précède, on en déduirait par récurrence que \(\phi^{(k)}_{i_1,\ldots,i_n}=0\) pour tous les indices avec \(n\geq 2\). D'après le principe du prolongement analytique (Cartan, 1961, chapitre IV), \(F(x)\) serait identiquement nulle, ce qui n'est pas possible.
Si J=1, on écrit plus simplement \[\psi^{(k)}_n=\frac{\partial^{n} F^{(k)}}{\partial x_1^n}(1),\quad M^{(k)}=c_1^{(k)} \left ( \frac{\partial g_1^{(k)}}{\partial x_1}(1)-1 \right ).\] L'équation (8) s'écrit alors \[\omega\, \psi^{(k)}_n = \sum_{h=1}^K Q_{k,h}\, \psi^{(h)}_n + n\, M^{(k)} \psi^{(k)}_n\, .\] On en déduit que \(\omega\) est une valeur propre de la matrice \(\mathbf{Q}\) ou d'une matrice \(\mathbf{Q}+n \, \mathrm{diag} ( M^{(1)},\ldots,M^{(K)})\) avec \(n\geq 1\). La section 4.2 de (Bacaër, 2017a) avait déjà remarqué ce cas particulier pour les processus linéaires de naissances et de mort avec un seul type.
On définit \(\mathbf{Y}_n\), la sous-matrice finie de la matrice (1) avec \(\mathbf{Z}_{1,1}\) dans son coin supérieur gauche et \(\mathbf{Z}_{n,n}\) dans son coin inférieur droit. La borne spectrale de cette matrice est telle que \(\mu_n \leq \mu_{n+1} \leq 0\). Comme dans la proposition 2 de (Bacaër, 2017a), on définit \(\mathbf{1}\) le vecteur ligne \((1,\ldots,1)\) de taille convenable. On a alors \[\mathbf{1} \mathbf{Y}_n \leq \mathbf{0} = 0\cdot \mathbf{1}.\] On a donc \(\mu_n \leq 0\ \) (voir par exemple (Nkague Nkamba, 2012, théorème 30.1)). Il existe un vecteur colonne \(\mathbf{v}_n\neq 0\) avec \[\mathbf{Y}_n \mathbf{v}_n = \mu_n \, \mathbf{v}_n\] et \(\mathbf{v}_n\geq 0\). On définit \(\mathbf{w}_n\) le vecteur colonne \((\mathbf{v}_n, \mathbf{0})\). On a alors \[\mathbf{Y}_{n+1} \mathbf{w}_n \geq \mu_n \mathbf{w}_n\] parce que les coefficients des matrices \(\mathbf{Z}_{n+1,1}\), \(\ldots\), \(\mathbf{Z}_{n+1,n}\) sont \(\geq 0\). On en déduit que \(\ \mu_{n+1}\geq \mu_n\ \) d'après (Nkague Nkamba, 2012, théorème 30.3). De ceci, il résulte aussi que la limite \(\alpha_1\) de \(\mu_n\) si \(n\to \infty\) existe.
Pour chaque environnement k, on se donne comme dans (Bacaër, 2017b) trois matrices de même taille: une matrice de naissance \(\mathbf{A}^{(k)}=(A_{i,j}^{(k)})\) à coefficients \(\geq 0\), une matrice de transfert \(\mathbf{T}^{(k)}=(T^{(k)}_{i,j})\) avec \[\sum_i T^{(k)}_{i,j}=0\quad \forall j,\quad T^{(k)}_{i,j}\leq 0 \quad \forall i\neq j,\] et une matrice diagonale de sortie \(\mathbf{S}^{(k)}=(S_{i,j}^{(k)})\) avec \(S^{(k)}_{j,j}\geq 0\ \forall j\). Autrement dit, pendant chaque intervalle de temps infinitésimal dt, chaque individu de type j qui se trouve dans l'environnement k a
Prenons le cas d'un processus de naissance et de mort où il n'y a que J=2 types, ce qui permet d'ordonner facilement les différents états de la population. Introduisons les matrices diagonales \[\mathbf{C}_j=\mathrm{diag}(c_j^{(1)},\ldots,c_j^{(K)}),\quad \mathbf{A}_{i,j}=\mathrm{diag}(A_{i,j}^{(1)},\ldots,A_{i,j}^{(K)}),\] \[\mathbf{T}_{i,j}=\mathrm{diag}(T_{i,j}^{(1)},\ldots,T_{i,j}^{(K)}),\quad \mathbf{S}_{j,j}=\mathrm{diag}(S_{j,j}^{(1)},\ldots,S_{j,j}^{(K)}).\] On range les fonctions \(p^{(k)}(t,n_1,n_2)\) selon le nombre total \(n_1+n_2\) d'individus, et ce nombre fixé, par le nombre d'individus de type 1 puis par l'environnement. Avec cet ordre, on considère le vecteur colonne infini \begin{eqnarray*} \mathbf{p}(t)&=&(p^{(1)}(t,0,0),\ldots,p^{(K)}(t,0,0),\\ &&p^{(1)}(t,1,0),\ldots,p^{(K)}(t,1,0),\\ &&p^{(1)}(t,0,1),\ldots,p^{(K)}(t,0,1),\\ &&p^{(1)}(t,2,0),\ldots,p^{(K)}(t,2,0),\\ &&p^{(1)}(t,1,1),\ldots,p^{(K)}(t,1,1),\\ &&p^{(1)}(t,0,2),\ldots,p^{(K)}(t,0,2),\ldots). \end{eqnarray*} Alors le système (3) s'écrit aussi \(d\mathbf{p}/dt=\mathbf{Z} \, \mathbf{p}(t)\).
Comme exemple de processus avec deux types d'individus, considérons le cas du modèle épidémique linéaire de (Bacaër et Ait Dads, 2014), où les individus de type 1 sont les personnes infectées mais pas encore infectieuses (c'est-à-dire dans la phase latente) et les individus de type 2 ceux qui sont infectieux. Le cas linéaire sous-critique correspond par exemple à la situation où la maladie est importée dans un environnement aléatoire défavorable à sa propagation. On a alors \[\mathbf{T}^{(k)}=\left ( \begin{array}{cc} \tau & 0 \\ -\tau & 0 \end{array} \right ),\quad \mathbf{A}^{(k)}=\left ( \begin{array}{cc} 0 & \beta^{(k)} \\ 0 & 0 \end{array} \right ),\quad \mathbf{S}^{(k)}=\left ( \begin{array}{cc} 0 & 0 \\ 0 & \gamma \end{array} \right ).\] Le paramètre τ est le taux auquel les personnes dans la phase latente deviennent infectieuses, indépendamment de l'environnement. Le paramètre \(\beta^{(k)}\) est le taux auquel les personnes infectieuses infectent de nouvelles personnes au début d'une épidémie; il dépend de l'environnement à cause de l'influence du climat sur la probabilité de transmission. Le paramètre γ est le taux de guérison des personnes infectieuses. On a ainsi \[\mathbf{M}^{(k)}=\left ( \begin{array}{cc} -\tau & \beta^{(k)} \\ \tau & -\gamma \end{array} \right ).\] On suppose de plus qu'il n'y a que K=2 environnements différents et que \[\mathbf{Q}=\left ( \begin{array}{cc} -q_1 & q_2 \\ q_1 & -q_2 \end{array} \right ).\] Le système (9) s'écrit \[ \frac{\partial f^{(1)}}{\partial t} = q_2 f^{(2)}-q_1 f^{(1)} + \tau (x_2-x_1) \frac{\partial f^{(1)}}{\partial x_1}+ \left [\beta^{(1)}(x_1-1)x_2-\gamma(x_2-1)\right ] \frac{\partial f^{(1)}}{\partial x_2}\] et \[ \frac{\partial f^{(2)}}{\partial t} = q_1 f^{(1)} - q_2 f^{(2)} + \tau (x_2-x_1) \frac{\partial f^{(2)}}{\partial x_1}+ \left [\beta^{(2)}(x_1-1)x_2-\gamma(x_2-1)\right ] \frac{\partial f^{(2)}}{\partial x_2}\, .\]
On a choisi les valeurs numériques suivantes, utilisées dans (Bacaër, Ait Dads, 2014) pour la rougeole: \(1/\tau=8\) jours, \(1/\gamma=5\) jours. Pour \(\beta^{(k)}\), on suppose que \(\beta^{(1)}=4\, \varepsilon\) par mois (avec un mois de 30 jours) et que \(\beta^{(2)}=8\, \varepsilon\) par mois. Dans (Bacaër et Ait Dads, 2014), le coefficient β variait de manière périodique entre 4 et 8 par mois pour avoir un bon ajustement avec la courbe épidémique. Le paramètre ε est destiné à varier. Supposons enfin que \(q_1=q_2=1\), de sorte que l'environnement passe en moyenne la moitié du temps dans chacun des deux états.
Avec une méthode itérative qui tire avantage de la structure tridiagonale par blocs (Ciarlet, 2006), on estime la borne spectrale \(\mu_n\) de la sous-matrice finie \(\mathbf{Y}_n\) de \({\bf Z}\) si n vaut successivement 25, 50, 100 et 200. Cette sous-matrice carrée est de taille \[2K+3K+\cdots+(n+1)K=n(n+3)K/2.\] On estime par ailleurs le paramètre critique \(\lambda_1\) lié à l'exposant de Lyapounov en utilisant par exemple 5000 sauts de l'environnement. Les résultats sont exposés dans la figure 1. Notons que si l'environnement était constant, le processus serait sous-critique pour \(\beta < \gamma\).
La figure suggère qu'on a \(\alpha_1=\omega_1\) si \(\varepsilon\) est petit, en particulier tant que \(\beta^{(1)} < \gamma\) et \(\beta^{(2)} < \gamma\), autrement dit \(\frac{8\, \varepsilon}{30} < \frac{1}{5}\), c'est-à-dire \(\varepsilon < \mbox{0,75}\). Mais \(\alpha_1 < \omega_1\) dans une zone où \(\lambda_1\) reste strictement négatif. Malheureusement, on n'est pas parvenu à déterminer \(\alpha_1\) plus explicitement. On s'attend à ce que \(\alpha_1=0\) si \(\lambda_1=0\). Rappelons que lorsqu'il n'y a qu'un seul type d'individus (\(J=1\)), les matrices \(\mathbf{M}^{(k)}\) sont en fait des nombres scalaires \(M^{(k)}\) et (Bacaër, 2017c) tendait à montrer dans ce cas que \[\alpha_1 = \min_{0\leq \theta\leq 1} s(\mathbf{Q}+\theta\, \mathrm{diag}(M^{(1)},\ldots,M^{(K)})).\] \(s(\cdot)\) désigne la borne spectrale d'une matrice. L'analogue de cette formule lorsqu'il y a plusieurs types d'individus reste à déterminer. (Dyakonova, 2008 ; Dyakonova, 2013 ; Vatutin et Wachtel, 2017) n'y étaient pas non plus parvenus dans le cadre des modèles en temps discret. Sans doute une meilleure compréhension du comportement des fonctions \(F^{(k)}(\mathbf{x})\) près du point singulier \(\mathbf{x}=\mathbf{1}\) dans le système (7) permettrait de progresser.