2007 International Conference in Honor of Claude Lobry

Applications de méthodes d'agrégation de variables à l'analyse de modèles spatiaux de dynamique des populations

Tri Nguyen-Huu a,b,* — Pierre Auger a,b

a U. R. GEODES
IRD
Centre de Recherche d'Ile de France
93143 Bondy cedex, France

b IXXI
Ecole Normale Supérieure de Lyon
Lyon Cedex 7, France
pierre.auger@bondy.ird.fr

* Corresponding author
tnguyenh@ens-lyon.fr

**RÉSUMÉ.** Les modèles de dynamique de populations peuvent prendre en compte un nombre important de paramètres et de variables, ce qui les rend difficiles à analyser. Lorsqu'il existe des processus associés à deux échelles de temps différentes, une lente et une rapide, les méthodes d'agrégation de variables permettent de construire un modèle simplifié qui comporte un nombre plus faible de variables. Elles permettent ainsi d'analyser et de décrire un système de manière globale. Nous présentons ces méthodes dans le cas de modèles discrets, puis nous illustrons leur utilisation à l'aide de modèles hôte-parasitoïdes spatialisés.

**ABSTRACT.** Models in population dynamics can deal with an important number of parameters and variables, which can make them difficult to analyse. Aggregation of variables allow reducing complexity of such models by building simplified models governing fewer variables by use of the existence of different time scales associated to the processes governing the whole system. Those reduced models allows analysing and describing the global dynamics of the system. We present those methods for time discrete models and illustrate their use for the study of spatial host-parasitoids models.

**MOTS-CLÉS** : Agrégation de variables, modèles discrets, modèles spatialisés

**KEYWORDS** : Aggregation of variables, discrete models, spatial models
# 1. Introduction

L'écologie cherche à décrire les interactions des êtres vivants entre eux et avec leur environnement. En particulier, la dynamique des populations cherche à décrire l'évolution numérique des populations d'organismes vivants à l'aide de modèles mathématiques ou informatiques. Les facteurs entrant en compte sont nombreux et agissent à des niveaux différents : à l'échelle moléculaire ou de l'individu, de la population (interaction entre les individus d'une même espèce), de l'écosystème (interactions avec d'autres espèces) ou encore environnementale (nature du terrain, climat). Si la modélisation en écologie a tout d'abord proposé des modèles mathématiques simples, la prise en considération d'un nombre important de ces facteurs, favorisée par le développement des ordinateurs et des techniques de calcul constitue l'un des enjeux actuels. Cette tendance à la multiplication des facteurs pris en compte a abouti à la réalisation de modèles mathématiques et informatiques très complexes, difficilement analysables et gourmands en puissance de calcul, rendant parfois leur utilisation très difficile.

Grâce au développement de l'informatique et des ordinateurs et des logiciels de calcul numérique, de nombreux modèles sont apparus, prenant en compte de plus en plus de détails et de données. Des paramètres tels que la fragmentation de l'habitat ou des mécanismes à l'échelle moléculaire peuvent être pris en compte, donnant plus de réalisme aux modèles et permettant la réalisation de simulations. Certains modèles considèrent même de manière séparée les individus d'une population (modèles individu-centrés).

En contrepartie, les modèles utilisent un nombre important de variables et de paramètres. Il devient alors difficile d'analyser mathématiquement ces modèles complexes. De plus, si de tels modèles peuvent être efficaces dans un but prédictif, il apparaît un effet "boîte noire": il devient très difficile de comprendre les mécanismes de l'évolution des populations et de prédire la sensibilité aux paramètres. C'est le cas des modèles pour lesquels l'habitat est divisé en de nombreux patchs, pour lesquels on se retrouve avec un nombre très important de variables : on doit représenter la densité de chaque population ou sous-groupe de population sur chacun des patchs. De plus, si on peut avoir une représentation visuelle de l'évolution des densités de population, il est difficile d'en donner une description globale simple car l'analyse est souvent subordonnée à l'observation de simulations. En outre, il est difficile de décrire la dynamique asymptotique de la population.

Les méthodes d'agrégation des variables ont été développées aussi bien dans le cas de modèles continus que de modèles discrets. Lorsqu'un modèle présente deux échelles de temps, une lente et une rapide, il est possible de s'en servir pour réduire la dimension du système étudié. Elles sont utilisées dans de nombreux domaines et sont en particulier bien adaptées à l'étude de populations qui sont structurées en sous-populations. De telles méthodes ont d'abord été présentées dans [2] et ont été développées initialement pour les Equations aux Dérivées Ordinaires dans [5], [4] et [21]. Ces travaux ont ensuite été étendus aux modèles discrets ([23] et [8]), et enfin aux Equations aux Dérivées Partielles [1, 12] et Equations Différentielles à Retard [24]. Les idées derrière les méthodes d'agrégation des variables sont souvent intuitives et sont parfois utilisées de manière implicite, par exemple en épidémiologie où la dynamique des populations est ignorée car la dynamique épidémiologique est beaucoup plus rapide. Elles sont également couramment utilisées en physique où l'on représente un phénomène global plutôt que de décrire le comportement de chaque particule. Nous allons présenter ces méthodes dans le cadre
de systèmes en temps discret. Nous nous intéresserons ensuite à des applications de ces méthodes à des modèles de dynamiques de populations hôte-parasitoïdes.

## 2. Méthodes d'agrégation des variables en temps discret

On considère un ensemble de plusieurs populations, décomposé en sous-populations appelées états. L'évolution du système s'effectue en temps discret et on connaît la manière de passer de l'état du système de l'instant $t$ à celui à l'instant $t+1$. Le système est divisé en $q$ populations (proies, prédateurs, hôtes, parasitoïdes, etc.), qui peuvent à leur tour être divisées en sous-groupes (répartition spatiale des individus, différents stades de développement d'une espèce, etc). La population est représentée à l'instant $t$ par un vecteur $\mathbf{n}(t) = (\mathbf{n}_1(t), \dots, \mathbf{n}_q(t))^T$ composé de vecteurs représentant les $q$ groupes de la population. Chaque vecteur $\mathbf{n}_i(t) = (n_{i,1}(t), \dots, n_{i,N_i}(t))^T$ pour $i = 1, \dots, q$, représente donc l'état de la $i^{\text{ème}}$ population, divisée en $N_i$ sous-groupes. Le nombre d'états total est $N = N_1 + \dots + N_q$. On a donc $\mathbf{n}(t) \in \mathbb{R}_+^N$. Les $n_{i,j}$ sont appelées variables locales. Un tel système est représenté en figure 1.

<figure><img src="image_1.png" /><figcaption>**Figure 1.** : Vue schématique d'un système hiérarchique. Le système est composé de trois groupes, eux-mêmes composés de trois sous-groupes. Lorsque les interactions à l'intérieur d'un groupe sont plus fortes que celles entre les groupes, il est possible sous certaines conditions de construire un modèle réduit régissant la dynamique à l'échelle des groupes [2].</figcaption></figure>

### 2.1. Modèles discrets linéaires

Dans les modèles que nous allons considérer, les populations sont soumises à deux processus de dynamique relatifs à deux échelles de temps différentes. L'un des processus
est très rapide en comparaison à l'autre : on parlera de dynamique lente et de dynamique rapide. Nous allons présenter les résultats relatifs au cas autonome développé dans [23, 24].

## Modèle complet

Les processus étant associés à des échelles de temps différentes, les intervalles de temps entre deux états du système vont être également différents. Les processus rapides et lents sont représentés par deux matrices différentes $F$ (Fast) et $S$ (Slow). Si à un instant donné, le système est dans un état $\mathbf{n}$, l'état obtenu avec le processus rapide après un pas de temps rapide est donc $F\mathbf{n}$. De même, l'état obtenu pour le processus lent après un pas de temps lent est $S\mathbf{n}$.

Afin d'obtenir un modèle tenant compte des deux processus, nous allons nous placer à l'échelle de temps lente et donner une description approchée de l'effet du processus rapide pour un intervalle de temps plus long. Une manière fréquente de procéder est de considérer que le processus rapide intervient $k$ fois pendant que le processus lent intervient une fois à chaque pas de temps. Plus précisément, lors d'un pas de temps, on répète $k$ fois le processus rapide puis une fois le processus lent. Ainsi, l'état du système au pas de temps $t+1$ sera donné par le système suivant, que nous appellerons modèle complet :

$$ \mathbf{n}_k(t+1) = SF^k\mathbf{n}_k(t) \quad [1] $$

## Modèle agrégé

Lorsque le nombre de sous-populations est important, la dynamique est régie par un système de taille très grande. Sous certaines conditions, l'étude d'un tel système peut être ramenée à celle d'un système réduit. Nous allons pour cela supposer que le rapport entre l'échelle de temps rapide et l'échelle de temps lente est très important, c'est-à-dire que $k$ est très grand.

Nous allons considérer le cas où le fait de répéter la dynamique rapide seule fait tendre le système vers un équilibre que l'on appellera équilibre rapide. Pour $k$ assez grand, on pourra considérer que le système a atteint l'équilibre rapide, et donc remplacer dans le modèle précédent le calcul $F^k\mathbf{n}(t)$ par celui de cet équilibre rapide associé à l'état $\mathbf{n}(t)$. On espère ainsi pouvoir obtenir un système simplifié prenant en considération un nombre plus petit de variables.

Nous supposons ici que pour chaque groupe $i$ la dynamique est interne, conservative d'une certaine variable globale du groupe et tendant asymptotiquement vers une distribution stable. C'est le cas pour une matrice $F$ composée de matrices de projection $F_i$ de taille $N_i \times N_i$, primitives (positives, irréductibles et avec une unique valeur propre de valeur maximale) et de valeur propre dominante 1, telle que $F = \text{diag}(F_1, \dots, F_q)$. Ces conditions sont vérifiées dans le cas de processus de migration gardant constant le nombre total d'individus, le processus de migration étant souvent représenté par une matrice stochastique primitive et irréductible. Le théorème suivant est fondamental pour la construction du modèle agrégé :

**Théorème de Perron-Frobenius** : Soit $A$ une matrice carrée positive, irréductible. Alors il existe une valeur propre de $A$, notée $\lambda$, qui a les propriétés suivantes :

1) $\lambda > 0$;

2) $\lambda$ est associée à un vecteur propre à gauche $g$ et à un vecteur propre à droite $d$ strictement positifs;
3) pour toute valeur propre $\lambda'$ de A, $|\lambda'| \le \lambda$;

4) l'espace propre associé à $\lambda$ est de dimension 1;

D'après le théorème de Perron-Frobenius, chaque matrice $F_i$ a deux vecteurs propres à gauche et à droite $u_i$ (vecteur ligne) et $v_i$ (vecteur colonne) strictement positifs associés à la valeur propre 1 tels que $F_i v_i = v_i$ et $u_i F_i = u_i$. On a de plus $\lim_{k \to \infty} F^k = v_i u_i$. Posons $\bar{F}_i = \lim_{k \to +\infty} F^k$ et $\bar{F} = \text{diag}(\bar{F}_1, \dots, \bar{F}_q)$, et notons $V = \text{diag}(v_1, \dots, v_q)$ et $U = \text{diag}(u_1, \dots, u_q)$. On a donc

$$ \bar{F} = \lim_{k \to +\infty} F^k = VU \quad [2] $$

Dans le cas où l'échelle de temps rapide devient très rapide par rapport à l'échelle de temps lente, c'est-à-dire $k \to +\infty$, on peut approximer le système complet (1) par le modèle suivant, appelé modèle auxiliaire :

$$ \mathbf{n}(t+1) = S\bar{F}\mathbf{n}(t) \quad [3] $$

Ce modèle auxiliaire comprend autant de variables que le modèle initial. Il est cependant possible d'aboutir à un système réduit en le réécrivant de la manière suivante :

$$ \mathbf{n}(t+1) = SVU\mathbf{n}(t) \quad [4] $$

Il est alors possible de réduire le nombre de variables du système en posant $X(t) = Un(t)$. $X(t)$ vérifie $X(t+1) = USVX(t)$. En posant $\bar{S} = SVU$, et comme $X(t) \in \mathbb{R}^q$, on se ramène au système à $q$ variables suivant, que nous appelerons modèle agrégé :

$$ X(t+1) = \bar{S}X(t) \quad [5] $$

Etant donné une condition initiale $\mathbf{n}_0$, on peut définir une condition initiale $X_0 = U\mathbf{n}_0$ pour le système (5) et trouver une solution $\{X(t)\}_{t \in \mathbb{N}}$ relative à cette condition initiale. D'après (4), il est possible de retrouver $\mathbf{n}(t+1)$ à partir de $X(t)$ grâce à la relation $\mathbf{n}(t+1) = SVX(t)$. Il est donc possible de trouver une solution à (3) pour la condition initiale $\mathbf{n}_0$ à partir du modèle (5). Le passage du modèle auxiliaire au modèle agrégé est un cas d'agrégation parfaite dans le sens de Iwasa [15].

Une fois le modèle agrégé (5) réalisé, il faut voir si la dynamique de celui-ci permet l'étude du système complet (1). Dans [23], il est prouvé que les éléments du système complet définissant son comportement asymptotique (la valeur propre dominante et les vecteurs propres associés) peuvent être approximés par ceux du modèle agrégé lorsque la matrice associée au modèle agrégé est primitive.

En effet, supposons que $\bar{S}$ est une matrice primitive. Considérons $\lambda > 0$ son unique valeur propre dominante, et $\bar{w}_l$ et $\bar{w}_r$ ses vecteurs propres associés respectivement à gauche (vecteur ligne) et à droite (vecteur colonne). Dans [23] Sánchez et al. démontrent que pour le système auxiliaire (3); $\lambda$ est l'unique valeur dominante de $\bar{S}\bar{F}$ et ses vecteurs propres associés respectivement à gauche et à droite sont $\bar{w}_l U$ et $\bar{S}V\bar{w}_r$.

Le système initial (1) peut être vu comme une perturbation du système auxiliaire (3) en considérant que $SF^k$ est une perturbation de $\bar{S}\bar{F}$. Ordonnons les valeurs propres de $F$
par module décroissant : $\lambda_1 = \dots = \lambda_q = 1 > |\lambda_{q+1}| \ge \dots \ge |\lambda_N|$. Si $\|\cdot\|$ est une norme consistante sur l'espace des matrices de taille $N \times N$, alors pour tout $\alpha > |\lambda_{q+1}|$ nous avons $\|SF^k - SF\|_{k \to +\infty} = o(\alpha^k)$. On en déduit que $SF^k$ a une unique valeur propre dominante de la forme $\mu_k = \lambda + O(\alpha^k)$ et des vecteurs propres associés respectivement à gauche et à droite $\bar{w}_l U + O(\alpha^k)$ et $SV \bar{w}_r + O(\alpha^k)$ pouvant être choisis arbitrairement, on peut prendre $\alpha < 1$. Les éléments définissant le comportement asymptotique de chaque modèle peuvent ainsi être liés précisément, leur valeur devenant proche quand $k$ tend vers l'infini.

Une application du théorème de Perron-Frobenius nous donne pour toute condition initiale positive $X_0$ la description du comportement asymptotique pour le modèle agrégé.

$$ \lim_{t \to +\infty} \frac{X(t)}{\lambda^t} = \frac{\bar{w}_l X_0}{\bar{w}_l \bar{w}_r} \bar{w}_r $$

Le comportement asymptotique du modèle initial vérifie alors

$$ \lim_{t \to +\infty} \frac{\mathbf{n}_k(t)}{\mu_k^t} = \frac{\bar{w}_l U \mathbf{n}_0}{\bar{w}_l \bar{w}_r} \frac{1}{\lambda} \bar{w}_r + O(\alpha^k) $$

On peut déduire plusieurs choses de ces résultats. La population croit asymptotiquement à une vitesse $\lambda$, valeur propre dominante de $SF$. De plus, la distribution de la population converge vers une distribution stable indépendante des conditions initiales, correspondant au vecteur propre à droite de $SF$ associé à la valeur propre $\lambda$. On dit que le système est fortement ergodique. Cette convergence est résumée :

$$ \lim_{t \to +\infty} \frac{X(t)}{\|X(t)\|} = \bar{w}_r $$

Plusieurs extensions de ces résultats ont été réalisées. Dans [26] Sanz et Bravo de la Parra se sont intéressés à un cas linéaire plus général, lorsque les matrices de projections $F_i$ ne sont plus nécessairement primitives.

D'autres cas ont été étudiés. Dans [7], les processus rapides sont toujours linéaires mais varient en fonction du temps à l'échelle rapide. La dynamique rapide est alors décrite par la composition de $k$ matrices différentes. Si la suite donnée par la composition de ces matrices converge, alors sous certaines conditions, la limite de la suite de matrices joue le même rôle que $\bar{F}$ pour la fabrication d'un modèle agrégé. On peut alors prouver des résultats similaires aux précédents.

Il est également possible de construire le système initial non plus sur la base de l'échelle de temps correspondant au processus lent mais sur l'échelle correspondant au processus rapide. Dans ce cas le processus rapide n'est réalisé qu'une seule fois par intervalle de temps et l'effet de la matrice $S$ durant cet intervalle de temps doit être approximé en utilisant une matrice $S_\varepsilon = \varepsilon S + (1 - \varepsilon)I$ où $I$ est la matrice identité. Une étude similaire au cas précédent peut alors être réalisée [23, 8, 9, 10].

Les méthodes d'agrégation de modèles linéaires en temps discret ont été étendues aux cas non autonomes et stochastiques. Dans [25], le cas d'environnement changeant avec le temps pour un modèle non autonome est étudié. La propriété de faible ergodicité, qui est en rapport avec la capacité d'un système à devenir asymptotiquement indépendant des
conditions initiales, est comparée pour les modèles complets et agrégés. Dans [28] Sanz et Bravo de la Parra ont obtenu des majorations de l'erreur obtenue lors de l'approximation de la dynamique du système complet par celle du système agrégé. Dans [7] Blasco et al. étendent également les résultats obtenus au cas d'un environnement changeant rapidement.

Enfin Sanz et Bravo de la Parra [27, 29] ont étendu les résultats précédents à l'étude de modèles discrets dans un environnement qui change de manière aléatoire avec le temps. Les variables d'état du modèle complet et les variables globales du modèle agrégé sont liées dans le cas où le processus de variation temporelle est Markovien.

## 2.2. Modèles discrets non linéaires

Nous allons nous intéresser au cas non-linéaire, en étendant le cadre général de l'étude précédente. Les matrices représentant les processus rapides et lents sont cependant remplaçées par des fonctions :

$$S, F : \Omega_N \to \Omega_N; S, F \in C^1(\Omega_N)$$

où $\Omega_N$ est un ensemble non vide de $\mathbb{R}^N$.

Comme pour le modèle linéaire, nous choisissons un pas de temps qui correspond à la dynamique lente, le processus rapide pouvant être répété $k$ fois durant cet intervalle de temps. En notant $\mathbf{n}_k(t) \in \mathbb{R}^N$ le vecteur des variables d'état à l'instant $t$, le système complet est défini de manière semblable au modèle (1) de la manière suivante :

$$\mathbf{n}_k(t+1) = S(F^k(\mathbf{n}_k(t))) \quad [6]$$

où $F^k$ représente la composée $k$ fois de la fonction $F$ par elle-même.

On suppose que pour toute condition initiale $X \in \Omega_N$, la dynamique rapide tend vers un équilibre, c'est-à-dire qu'il existe une fonction $\bar{F} : \Omega_N \to \Omega_N; F \in C^1(\Omega_N)$ telle que pour tout $X \in \Omega_N$, $\lim_{k \to +\infty} F^k(X) = \bar{F}(X)$. On suppose également que cet équilibre peut dépendre d'un nombre moins important de variables : il existe un ensemble non vide $\Omega_q \in \mathbb{R}^N$ avec $q < N$, et deux fonctions $G : \Omega_N \to \Omega_q; G \in C^1(\Omega_N)$ et $E : \Omega_q \to \Omega_N; E \in C^1(\Omega_q)$ telles que $\bar{F} = E \circ G$.

Comme pour le cas linéaire, on peut définir un système auxiliaire qui est une approximation de (6) lorsque $k$ tend vers l'infini, c'est-à-dire quand le processus rapide arrive à l'état d'équilibre. Ce système devient :

$$\mathbf{n}(t+1) = S(\bar{F}\mathbf{n}(t)) \quad [7]$$

Ce système peut se réécrire également $\mathbf{n}(t+1) = S \circ E \circ G(\mathbf{n}(t))$.

On peut alors définir le vecteur $Y(t) \in \mathbb{R}^q$ des variables globales de la manière suivante :

$$Y(t) = G(\mathbf{n}(t)) \quad [8]$$

En appliquant $G$ aux 2 côtés de (7), on obtient
$$
\begin{align*}
Y(t+1) &= G(\mathbf{n}(t+1)) \\
&= G \circ S \circ E \circ G(\mathbf{n}(t)) \\
&= G \circ S \circ E(Y(t))
\end{align*}
$$

En posant $\bar{S} = G \circ S \circ E$, on obtient un système autonome qui est un système approché de (6). C'est le système agrégé suivant :

$$
Y(t + 1) = \bar{S}(Y(t)) \quad [9]
$$

Comme dans le cas linéaire, nous allons nous intéresser à la comparaison des compos-
tements du modèle complet et du modèle agrégé pour de grandes valeurs du paramètre
k. Les résultats présentés dans cette partie sont exposés de manière plus générale dans
[6, 16]. Nous allons dans un premier temps comparer les solutions des deux systèmes
à un instant t fixé. Nous allons énoncer un théorème qui affirme que la dynamique du
système auxiliaire est entièrement déterminée par celle du système réduit, et que la solu-
tion du système complet peut être approchée par celle du modèle agrégé pour t fixé sous
quelques hypothèses supplémentaires :

**Théorème 2.1** : Soit **n**₀ ∈ Ω_N et soit Y₀ = G(**n**₀) ∈ Ω_q. On a alors :

1) La solution {$\mathbf{n}(t)}_{t \in \mathbb{N}}$ de (7) correspondant à la condition initiale $\mathbf{n}_0 \in \Omega_N$ et la solution {$Y(t)}_{t \in \mathbb{N}}$ de (9) correspondant à la condition initiale $Y_0$ sont liées par les relations suivantes :

$$
\begin{cases}
Y(t) = G(\mathbf{n}(t)) \\
\mathbf{n}(t) = S \circ E(Y(t-1))
\end{cases}
$$

2) Soit $t \in \mathbb{N}$ fixé, et supposons qu'il existe un ensemble borné $\Omega$ tel que $\bar{\Omega} \subset \Omega_N$,
$\mathbf{n}(i) \in \Omega_N$ pour $i \in \mathbb{N}$, et tel que $F^k$ converge uniformément vers $F$ sur $\Omega$. Alors la
solution de (6) correspondant à la condition initiale $\mathbf{n}_0$ et la solution $Y(t)$ de (9) sont liées
par les relations suivantes :

$$
Y(t) = \lim_{k \to +\infty} G(\mathbf{n}_k(t)) \quad \text{et} \quad \lim_{k \to +\infty} \mathbf{n}_k(t) = S \circ E(Y(t-1)) \quad [10]
$$

Nous allons ensuite étudier les relations entre les points d'équilibre du modèle complet et du modèle agrégé. La relation (10) du théorème 2.1 permet d'obtenir une relation entre les points d'équilibre du système auxiliaire et ceux du modèle agrégé :

- Si $\mathbf{n}^* \in \Omega_N$ est un équilibre de (7), alors $Y^* = G(\mathbf{n}^*) \in \Omega_q$ est un équilibre de (9).

- Réciproquement, si $Y^*$ est un équilibre de (9), alors $\mathbf{n}^* = S \circ E(Y^*)$ est un équilibre de (7).

Les points d'équilibre correspondants de (7) et (9) sont simultanément stables ou instables.

Le résultat suivant montre que sous certaines conditions, si il existe un équilibre $Y^*$
pour le système agrégé, alors il existe pour des valeurs de $k$ suffisamment grandes un point
d'équilibre $\mathbf{n}_k^*$ pour le modèle complet. On peut donner une approximation de cet équilibre à partir de $Y^*$. De plus, dans le cas hyperbolique, la stabilité de $Y^*$ est équivalente à celle de $\mathbf{n}_k^*$. Dans le cas d'un équilibre asymptotiquement stable, le bassin d'attraction de $\mathbf{n}_k^*$ peut être approché par celui de $Y^*$.

**Théorème 2.2** : Supposons que $\bar{F} \in C^1(\Omega_N)$ et que $F^k$ et $DF^k$ convergent uniformément respectivement vers $\bar{F}$ et $DF$ sur tout compact $K \subset \Omega_N$.

Si $Y^* \in \mathbb{R}^q$ est un équilibre hyperbolique asymptotiquement stable (resp. asymptotiquement instable) de (9), alors il existe $k_0 \in \mathbb{N}$ tel que pour $\forall k > k_0$, $k \in \mathbb{N}$, il existe un point fixe hyperbolique $\mathbf{n}_k^*$ de (6) qui est asymptotiquement stable (resp. asymptotiquement instable) et qui vérifie $\lim_{k \to +\infty} \mathbf{n}_k(t) = S \circ E(Y^*)$.

De plus, si pour $\mathbf{n}_0 \in \Omega_N$ la solution $\{Y(t)\}_{t \in \mathbb{N}}$ de (6) correspondant à la condition initiale $Y_0 = G(\mathbf{n}_0)$ est telle que $\lim_{k \to +\infty} Y(t) = Y^*$ alors pour $\forall k > k_0$, $k \in \mathbb{N}$, la solution $\{\mathbf{n}_k(t)\}_{t \in \mathbb{N}}$ de (6) correspondant à $\mathbf{n}_0$ vérifie $\lim_{k \to +\infty} \mathbf{n}_k(t) = \mathbf{n}_k^*$.

Ce résultat est utilisé dans des modèles présentés dans [11] et [3]. Comme dans le cas linéaire, il est également possible de construire un modèle complet en prenant comme pas de temps celui de la dynamique rapide. Dans [9] et [10], un système avec une dynamique rapide linéaire et une dynamique lente non linéaire est réduit à l'aide d'un théorème de variété centrale.

### 3. Applications des méthodes d'agrégation des variables

Les modèles hôtes-parasitoïdes spatialisés sur des grilles bidimensionnelles carrées sont des systèmes faisant appel à un nombre très important de variables. Les différences d'échelles de temps entre les processus démographiques et les processus de dispersion permettent d'utiliser les méthodes d'agrégation des variables vues précédemment pour décrire la dynamique de tels systèmes et l'analyser plus facilement.

Pour ce type de modèle et dans un environnement homogène, la dynamique d'interaction hôte-parasitoïde avec des générations discrètes peut être décrite par un système d'équations différentielles du premier ordre (11)

$$ \begin{cases} n_{t+1} = \lambda(n_t)n_t e^{-ap_t} \\ p_{t+1} = cn_t(1 - e^{-ap_t}) \end{cases} \quad [11] $$

qui correspond au modèle de Nicholson-Bailey avec un taux de croissance des hôtes $\lambda(n_t)$. Ici $n_t$ et $p_t$ sont les densités respectives des hôtes et des parasitoïdes à la génération $t$. Le paramètre $a$ est l'efficacité de recherche du parasitoïde et $c$ est le nombre moyen de parasitoïdes femelles émergeant de chaque hôte infecté.

Pour le taux de croissance des hôtes, nous adoptons soit la forme $\lambda(n_t) = rn_t$, soit la forme logistique $\lambda(n_t) = e^{(r(1-\frac{n_t}{K}))}$ où $r$ est un paramètre positif et $K$ est la capacité limite. A ces deux formes correspondent deux dynamiques différentes. En effet, dans le premier cas la dynamique est instable et nous assistons à des variations d'amplitudes croissantes du nombre d'individus, alors que dans le deuxième cas elle est persistante. Nous considérons un environnement composé d'une grille carrée avec $A \times A$ sites. La dynamique implique deux phases : La première phase correspond aux interactions hôte-
parasitoïde sur chaque site. Les équations (11) nous permettent de calculer les populations locales sur chaque site d'une génération à la suivante. La deuxième phase est la dispersion entre les sites voisins. Nous supposons qu'une proportion constante $\mu_n$ (resp. $\mu_p$) de migrants pour la population d'hôtes (resp. parasitoïdes) quitte chaque site pour aller vers les huit plus proches voisins. Pour chaque patch, la dynamique de dispersion correspond au système (12) :

$$ egin{cases} n_{t+1} = (1 - \mu_n)n_t + \frac{\mu_n}{8} \sum_{\text{voisins}} n_t \\ p_{t+1} = (1 - \mu_p)p_t + \frac{\mu_p}{8} \sum_{\text{voisins}} p_t \end{cases} \quad [12] $$

Le processus complet de dispersion peut être décrit par une matrice $M$ non négative. Le modèle complet combine la dispersion et les interactions locales selon l'équation (13).

$$ N_{t+1} = F(M^k N_t) \quad [13] $$

Les composantes du vecteur $N_t$ de dimension $2A^2$ sont les densités des hôtes et des parasitoïdes sur chaque site. Le vecteur des fonctions $F$ est de la même dimension et incorpore les deux équations précédentes (11) d'interactions locales pour tous les sites de la grille. Le paramètre $k$ est un nombre entier que nous appelons fréquence de dispersion. Quand le paramètre $k = 1$, il y a un seul événement de dispersion entre deux générations. Quand $k > 1$, le processus de dispersion est répété plusieurs fois. Dans ce cas, au lieu d'aller seulement vers les huit sites voisins les plus proches, les migrants peuvent atteindre des sites plus lointains. Cette situation est susceptible de se produire pour la plupart des associations hôte-parasitoïde concernant des insectes adultes ailés (comme des guêpes, des papillons, des mouches ou des moucherons). Quand $k$ devient assez grand, le processus de dispersion est rapide par rapport aux interactions locales. Il est alors possible de construire un modèle agrégé qui régit les densités totales de populations d'hôtes et de parasitoïdes de la grille.

Le modèle agrégé est obtenu en deux étapes. À la première étape, nous négligeons le processus lent de l'interaction locale et nous étudions le processus rapide de dispersion qui est caractérisé par un équilibre rapide stable, qui correspond à la même proportion $\frac{1}{A^2}$ d'hôtes et de parasitoïdes sur chaque site. Cet équilibre rapide correspond à une distribution homogène spatiale des hôtes et des parasitoïdes sur la grille, et donc à un chronisme spatial. Dans la deuxième étape, nous faisons l'approximation suivante : Pour toute génération $t$, nous supposons que le système est à l'équilibre rapide. Puis en substituant l'équilibre rapide dans les équations (11) du système complet, et par sommation des densités de population d'insectes de tous les sites de la grille, nous obtenons un modèle agrégé régissant les densités totales d'hôtes et de parasitoïdes. Nous obtenons alors le modèle agrégé (14) pour le modèle de Nicholson-Bailey classique :

$$ egin{cases} N_{t+1} = \lambda N_t e^{-a \frac{P_t}{A^2}} \\ P_{t+1} = c N_t (1 - e^{-a \frac{P_t}{A^2}}) \end{cases} \quad [14] $$

Dans le cas du modèle de Nicholson-Bailey avec croissance logistique des hôtes, on obtient le modèle agrégé (15) :
$$
\begin{cases}
n_{t+1} = n_t e^{\left(r\left(1 - \frac{n_t}{kA^2}\right) - \frac{a}{A^2}p_t\right)} \\
p_{t+1} = c n_t \left(1 - e^{-\frac{a}{A^2}p_t}\right)
\end{cases}
\quad [15]
$$

Nous étudions maintenant l'adéquation entre les modèles complets et les modèles
agrégés correspondants en mettant plus particulièrement en avant le rôle du paramètre *k*.

**3.1. Modèle de Nicholson-Bailey classique**

Dans le cas du modèle de Nicholson-Bailey classique, le système non spatialisé (11) peut être réécrit sous la forme du système (16) :

$$
\begin{cases}
n_{t+1} = \lambda n_t e^{-ap_t} \\
p_{t+1} = c n_t (1 - e^{-ap_t})
\end{cases}
\quad [16]
$$

La dynamique de ce système est bien connue. Il existe un unique équilibre dans positif donné par :

$$
N^* = \frac{\lambda \ln(\lambda)}{ac(\lambda - 1)}, \quad P^* = \frac{\ln(\lambda)}{a}
$$

Cet équilibre est instable. Les densités d'individus présentent des oscillations d'ampli-
tudes croissantes, aboutissant à l'explosion ou l'extinction des espèces.

La spatialisation de ce modèle dans le cas d’un déplacement par génération ($k=1$) peut conduire à une dynamique globale du système persistante [14] au-delà d’une taille minimale de la grille. Les densités totales d’hôtes et de parasitoïdes restent alors bornées. Différentes structures spatiales peuvent apparaître, comme des spirales, des réseaux cristallins ou du chaos spatial en fonction des valeurs des paramètres $\mu_n$ et $\mu_p$ [14, 13, 22].

Nous avons étudié l'influence du paramètre *k* sur la dynamique du modèle spatialisé.
A paramètres constants, lorsque *k* augmente, on observe une augmentation de la taille des
structures spatiales observées. La figure 2 présente des spirales obtenues pour différentes
valeurs de *k*, ainsi que l'évolution temporelle des densités d'hôtes et de parasitoïdes cor-
respondantes. Lorsque *k* devient suffisamment important, la taille des structures dépasse
celle de la grille. Il n'est alors plus possible d'obtenir la persistance du système, et on
assiste à la réapparition d'une dynamique d'oscillations croissantes des densités d'hôtes
et de parasitoïdes.

Le modèle agrégé (14) qui correspond à ce modèle spatialisé est lui-même un mo-
dèle de Nicholson-Bailey classique, dont la dynamique est instable. Les dynamiques des
modèles complet et agrégé sont représentées en figure (3.1).

Pour des valeurs faibles de *k*, on observe des structures spatiales et le système complet est persistant, alors que le système agrégé n’est pas persistant. Le modèle agrégé ne permet pas de décrire qualitativement la dynamique du système complet. En revanche, pour des valeurs élevées de *k*, on observe pour les deux systèmes des oscillations d’amplitudes croissantes non bornées. La dynamique du système agrégé permet alors de décrire qualitativement celle du modèle complet [17].

Les conditions initiales ont une influence sur la dynamique du système, certaines étant plus favorables que d'autres pour la mise en place d'une dynamique persistante. Lorsque les densités d'hôtes et de parasitoïdes sont uniformes dans l'espace, il est impossible d'obtenir la persistance du système, la dynamique du modèle agrégé correspond alors
<figure><img src="image_2.png" /><figcaption>**Figure 2.** *Structures spatiales de la densité d'hôtes obtenues avec le modèle complet de Nicholson-Bailey pour les paramètres a = 0.2, λ = 2, c = 0.4, μn = μp = 1 sur une grille 50 × 50 pour différentes valeurs de k.*</figcaption></figure>

exactement à la dynamique globale du système. A l'opposé, on observe que si les individus sont répartis selon un motif en spiral, la structure spatiale est plus difficile à casser et on obtient plus facilement une dynamique persistante. Ainsi, pour toute condition initiale, il existe une valeur du paramètre *k* au-delà de laquelle le système n'est plus persistant et le modèle agrégé permet de décrire qualitativement la dynamique du modèle complet. On peut définir un seuil global *k̃* au-delà duquel la dynamique du système complet est non persistante, indépendamment de la condition initiale. Au-delà de cette valeur, le modèle complet se comporte alors comme le modèle agrégé, quelle que soit la condition initiale. Les résultats obtenus pour différentes tailles de grilles sont représentés dans le tableau (1). Par exemple, un *k̃* de 3 correspond à trois déplacements élémentaires par génération, les individus ne pouvant atteindre des patchs localisés uniquement à moins de 3 patchs de distance du patch de départ. Ce seuil est déterminé numériquement à partir d'un grand nombre de simulations avec des conditions initiales différentes.

Il est intéressant de remarquer que les valeurs obtenues pour ce seuil sont relativement faibles, et correspondent à un déplacement assez localisé à chaque génération (environ 10% de la taille d'un côté de la grille). Cela signifie qu'il n'est pas nécessaire que les individus se déplacent sur de grandes distances pour que le modèle agrégé puisse être utilisé pour décrire la dynamique du modèle complet.

<figure><img src="image_3.png" /><figcaption>**Figure 3.** *Evolution de la densité d'hôtes et de parasitoïdes pour k = 1 (a) et k = 4 (b) sur une grille 50 × 50 pour les paramètres a = 0.2, λ = 2, c = 0.4, μn = μp = 1.*</figcaption></figure>
<table><caption>**Tableau 1.** : valeur du seuil $\tilde{k}$ pour différentes tailles de la grille, pour les paramètres $a = 0.2$, $\lambda = 2$, $c = 0.4$, $\mu_n = \mu_p = 1$.</caption><tr><td>A</td><td>30</td><td>50</td><td>100</td></tr><tr><td>k&#770;</td><td>3</td><td>5</td><td>11</td></tr></table>

## 3.2. Modèle de Nicholson-Bailey avec croissance logistique des hôtes

Dans le cas du modèle de Nicholson-Bailey classique, le système non spatialisé (11) peut être réécrit sous la forme du système (17) :

$$ \begin{cases} n_{t+1} = \exp\left(r\left(1 - \frac{n_t}{K}\right)\right) n_t e^{-ap_t} \\ p_{t+1} = cn_t(1 - e^{-ap_t}) \end{cases} \quad [17] $$

où $r$ est un paramètre strictement positif et $K$ est une capacité limite. En l'absence de parasitoïdes, la population d'hôtes tend vers la capacité limite $K$. Contrairement au cas précédent, le modèle non spatialisé ne présente pas d'oscillations d'amplitude croissante des populations d'hôtes et de parasitoïdes, mais donne une dynamique persistante. Le modèle agrégé (15) obtenu a la même forme mathématique que le modèle complet (17), avec des paramètres différents. Les deux modèles présentant une dynamique persistante, le modèle agrégé permet donc de décrire qualitativement la dynamique du modèle complet. La précision de cette approximation de la dynamique du modèle complet dépend du paramètre $k$.

On observe également dans ce modèle l'apparition de structures spatiales (spirales). Lorsque $k$ augmente, la taille de ces spirales augmente. Au-delà d'un certain seuil, cette taille dépasse celle de la grille, et on aboutit finalement à une répartition homogène des individus dans l'espace. Dans le cas d'une répartition homogène des individus, la dynamique du modèle complet est exactement régie par le modèle agrégé (15). La figure (4) représente les structures spatiales obtenues avec le modèle complet pour différentes valeurs de $k$.

<figure><img src="image_4.png" /><figcaption>**Figure 4.** Structures spatiales de la densité d'hôtes obtenue avec le modèle complet de Nicholson-Bailey avec croissance logistique des hôtes pour les paramètres $a = 0.2$, $r = 2.6$, $K = 50$, $c = 0.4$, $\mu_n = \mu_p = 1$ sur une grille $100 \times 100$ pour différentes valeurs de $k$.</figcaption></figure>

La figure (5) représente la dynamique du modèle complet et du modèle agrégé pour différentes valeurs de $k$ dans l'espace des phases (hôte, parasitoïde) pour $A= 100$. Les
trajectoires des modèles complet et agrégé tendent chacune vers un attracteur. Lorsque $k$ est faible, les deux attracteurs sont distincts. Lorsque $k$ augmente, les deux attracteurs se rapprochent jusqu'à se confondre. On observe en même temps une diminution du coefficient $V$, somme des variances des distributions des hôtes et des parasitoïdes. Lorsque ce coefficient atteint la valeur 0, on se retrouve dans le cas d'homogénéité spatiale.

<figure><img src="image_5.png" /><figcaption>**Figure 5.** Simulations numériques du modèle agrégé (noir) et du modèle complet (gris) avec les paramètres $a = 0.2$, $r = 2.6$, $c = 0.4$, $K = 50$, $\mu_n = \mu_p = 1$ pour $k = 1$ (a) et $k = 3$ (c). (b) et (c) représentent les variations de $V$, somme des variances de la distribution des hôtes et des parasitoïdes.</figcaption></figure>

Il est également possible de définir un seuil $\tilde{k}$, indépendant des conditions initiales, en dessous duquel les structures spatiales peuvent exister, et au-dessus duquel elles disparaissent systématiquement, c'est-à-dire lorsque $V$ tend vers 0 [18]. Ce seuil est déterminé numériquement à partir d'un grand nombre de simulations avec différentes conditions initiales. Les valeurs du seuil obtenues pour différentes tailles de grilles sont représentées dans le tableau (2) :

<table><caption>**Tableau 2.** Valeur du seuil $\tilde{k}$ pour différentes tailles de la grille, pour les paramètres $a = 0.2$, $r = 2.6$, $c = 0.4$, $K = 50$, $\mu_n = \mu_p = 1$.</caption><tr><td>A</td><td>30</td><td>50</td><td>100</td><td>200</td></tr><tr><td>k&#771;</td><td>2</td><td>3</td><td>8</td><td>10</td></tr></table>

Le paramètre $\tilde{k}$ est ici aussi relativement faible comparé à la taille de la grille. Cela correspond à des déplacements d'un nombre de patchs d'au plus 10% de la taille d'un côté de la grille à chaque génération. Cela signifie que pour des déplacements assez locaux, le modèle agrégé permet de décrire fidèlement le modèle complet.
## 4. Conclusion

Nous avons montré que dans ces deux cas il existe un seuil pour la fréquence de dispersion $k$ au-delà duquel les structures spatiales disparaissent et les dynamiques des patchs se synchronisent. Dans ce cas, le modèle agrégé gouvernant les densités totales d'hôtes et de parasitoïdes permet de prédire fidèlement la dynamique du modèle complet. L'existence d'un seuil $\tilde{k}$ n'est pas surprenante. La dispersion a en effet tendance à promouvoir l'homogénéité spatiale. Lorsque celle-ci devient plus fréquente, il devient plus facile d'obtenir un état d'homogénéité spatiale. Ainsi, pour une valeur élevée de $k$, la proportion d'hôtes et de parasitoïdes devient similaire entre les patchs, ce qui fait disparaître les structures spatiales. Cependant, il est assez inattendu que cette valeur du seuil $\tilde{k}$ soit relativement faible. Pour une dispersion dans un disque de 5 à 10% de la taille d'un côté de la grille, c'est-à-dire relativement peu en comparaison à la taille de l'habitat, on aboutit à l'homogénéité spatiale. Ainsi, pour des rayons de dispersion réalistes, il est possible d'utiliser le modèle agrégé pour décrire la dynamique globale du modèle complet.

Les méthodes d'agrégation des variables commencent à être utilisés dans des cas concrets, comme l'étude de la dynamique d'une population de carabes [20]. De plus, ces méthodes ne se limitent pas à des cas où la dispersion est uniforme et aboutit à l'homogénéité spatiale. C'est également le cas pour des modèles avec dispersion selon un gradient ou dispersion densité-dépendante [19].

## 5. Bibliographie

[1] O. ARINO, E. SÁNCHEZ, R. BRAVO DE LA PARRA, P. AUGER, « A singular perturbation in an age-structured population model », *SIAM Journal of Applied Mathematics*, vol. 60, 1999.

[2] P. AUGER, « Dynamics and thermodynamics in hierarchically organized systems », *Pergamon Press, Oxford*, 1989.

[3] P. AUGER, R. BRAVO DE LA PARRA, « Methods of aggregation of variables in population dynamics », *Comptes Rendus de l'Académie des Sciences Paris, Sciences de la Vie*, vol. 323, 2000.

[4] P. AUGER, J.C. POGGIALE, « Aggregation and emergence in systems of ordinary differential Equations », *Mathematical and Computer Modelling*, vol. 27, n° 4, 1998.

[5] P. AUGER, R. ROUSSARIE, « Complex ecological models with simple dynamics : from individuals to populations », *Acta Biotheoretica*, vol. 42, 1994.

[6] P. AUGER, R. BRAVO DE LA PARRA, J.C. POGGIALE, E. SÁNCHEZ, T. NGUYEN-HUU, « Aggregation of variables and applications to population dynamics », in *Structured population models in biology and epidemiology*, Eds Magal, P. and Ruan, S., *Lecture notes in Mathematics*, Springer, 2008.

[7] A. BLASCO, L. SANZ, P. AUGER, R. BRAVO DE LA PARRA « Linear discrete population models with two time scales in fast changing environments I : autonomous case », *Acta Biotheoretica*, vol. 49, 2001.

[8] R. BRAVO DE LA PARRA, P. AUGER, E. SÁNCHEZ, « Aggregation methods in discrete models », *Journal of Biological Systems*, vol. 3, 1995.

[9] R. BRAVO DE LA PARRA, E., SÁNCHEZ, P. AUGER, « Time scales in density dependent discrete models », *Journal of Biological Systems*, vol. 5, 1997.
[10] R. BRAVO DE LA PARRA, E. SÁNCHEZ, « Aggregation methods in population dynamics discrete models », *Mathematical and Computer Modelling*, vol. 27, n° 4, 1998.

[11] R. BRAVO DE LA PARRA, E. SÁNCHEZ, O., ARINO P. AUGER, « A discrete model with density dependent fast migration », *Mathematical Biosciences*, vol. 157, 1999.

[12] R. BRAVO DE LA PARRA, O. ARINO, « A model of an age structured population with two time scales », *Mathematical and Computer Modelling*, vol. 31, 2000.

[13] H.N. COMINS, M.P. HASSELL, R.M. MAY, « The spatial dynamics of host-parasitoid systems », *Journal of Animal Ecology*, vol. 61, 1992.

[14] M.P. HASSELL, H.N. COMINS, R.M. MAY, « Spatial structure and chaos in insect population dynamics », *Nature*, vol. 353, 1991.

[15] Y. IWASA, V. ANDREASEN, S. LEVIN, « Aggregation in model ecosystems. I. Perfect Aggregation », *Ecological Modelling*, vol. 37, 1987.

[16] M. MARVÁ, E. SÁNCHEZ, R. BRAVO DE LA PARRA, L. SANZ, « Reduction of slow-fast discrete models coupling migration and demography », *soumis*.

[17] T. NGUYEN-HUU, C. LETT, J.C. POGGIALE, P. AUGER, « Effect of movement frequency on global host-parasitoid spatial dynamics with unstable local dynamics », *Ecological Modelling*, vol. 197, n° 3-4, 2006.

[18] T. NGUYEN-HUU, C. LETT, P. AUGER, J.C. POGGIALE, « Spatial synchrony in host-parasitoid models using aggregation of variables », *Mathematical Biosciences*, vol. 203, n° 2, 2006.

[19] T. NGUYEN-HUU, C. LETT, P. AUGER, « Effects of density-dependent dispersal on the stability of host-parasitoid associations in a chain of patches », *Ecological Complexity*, vol. 5, n° 1, 2008.

[20] J.B. PICHANCOURT, F. BUREL, P. AUGER, « Assessing the effect of habitat fragmentation on population dynamics : an implicit modeling approach », *Ecological Modelling*, vol. 192, 2006.

[21] J.C. POGGIALE, « Applications des variétés invariantes à la modélisation de l'hétérogénéité en dynamique de populations », *thèse de doctorat, université de Bourgogne*, 1994.

[22] P. ROHANI, O. MIRAMONTES, « Host-parasitoid metapopulations : the consequences of parasitoid aggregation on spatial dynamics and searching efficiency », *Proceedings of the Royal Society of London Series B*, vol. 260, 1995.

[23] E. SÁNCHEZ, R. BRAVO DE LA PARRA, P. AUGER, « Discrete models with different timescales », *Acta Biotheoretica*, vol. 43, 1995.

[24] E. SÁNCHEZ, R. BRAVO DE LA PARRA, P. AUGER, P. GÓMEZ-MOURELO, « Time scales in linear delayed differential equations », *Journal of Mathematical Analysis and Applications*, vol. 323, 2006.

[25] L. SANZ, R. BRAVO DE LA PARRA, « variables aggregation in time varying discrete systems », *Acta Biotheoretica*, vol. 46, 1998.

[26] L. SANZ, R. BRAVO DE LA PARRA, « Variables aggregation in a time discrete linear model », *Mathematical Biosciences*, vol. 157, 1999.

[27] L. SANZ, R. BRAVO DE LA PARRA, « Time scales in stochastic multiregional models », *Nonlinear Analysis - Real World Applications*, vol. 1, 2000.

[28] L. SANZ, R. BRAVO DE LA PARRA, « Approximate reduction techniques in population models with two time scales : study of the approximation », *Acta Biotheoretica*, vol. 50, n° 4, 2002.

[29] L. SANZ, R. BRAVO DE LA PARRA, « Approximate reduction of multiregional models with environmental stochasticity », *Mathematical Biosciences*, vol. 206, n° 1, 2007.