Petit cours d'épidémiologie mathématique
Le modèle SIS

Julien Arino

Department of Mathematics & Data Science Nexus
University of Manitoba*

Centre canadien de modélisation des maladies (CCDM/CCMM)
NSERC-PHAC EID Modelling Consortium (CANMOD, MfPH, OMNI/RÉUNIS)

Plan de ce cours

• Modèles en compartiments
• Modèle SIS endémique
• Analyse mathématique (version 1)
• Le nombre de reproduction élémentaire
• Analyse mathématique (version 2)
• Un peu de computationnel
• Un peu plus au sujet de la stabilité

Modèles en compartiments

Modèles en compartiments

• Devenus synonymes des modèles épidémiologiques
• Bien des modèles épidémiologiques sont en compartiments, mais le développement des modèles en compartiments dans les années 1970s-1980s n'était pas spécifique à l'épidémiologie
• Voir en particulier les travaux de John Jacquez, Carl Simon, GG Walter
• Un peu oubliés, de façon injuste: il y a de très jolis résultats dans le domaine

Compartiment (Jacquez 1979)

A compartment is an amount of some material which acts kinetically like a distinct, homogeneous, well-mixed amount of material. A compartmental system consists of one or more compartments which interact by exchanging the material. There may be inputs into one or more compartments from outside the system and there may be excretions from the compartments of the system.

Un compartiment

• taille du compartiment, i.e., quantité de matériau cynetiquement homogène présent dans ;
• et coefficients/fonctions de transfert
• coefficient/fonction d'excrétion
• entrées d'en dehors du système

Modèle SIS endémique

• On considère une population close (pas d'entrées ni de sorties de la population)
• Supposons que les individus dans la population sont dans l'un de deux états:
  • susceptibles (à la maladie) s'ils ne sont pas en ce moment porteurs du pathogène
  • infectieux (et infectés) s'ils ont été infectés par la maladie et participent à sa propagation

il y a deux compartiments et le but de la modélisation est de décrire l'évolution du nombre d'individus dans chacun des compartiments

• nombre de susceptibles au temps
• nombre d'infectieux au temps
• la population totale

Les hypothèses qui suivent décrivent une maladie pour laquelle la durée de la période d'incubation est très courte voire inexistante

On suppose aussi que l'infection ne persiste pas chez un individu et qu'aussitôt après s'être remis, les individus sont de nouveau susceptibles

Type de compartiments

Individus susceptibles

• Naissent au taux per capita proportionnel à la population totale
• Meurent au taux per capita , proportionnel à la population susceptible
• Les nouveaux nés sont susceptibles (on ignore la transmission verticale)

Individus infectieux

• Meurent au taux per capita , proportionnel à la population infectieuse
• Guérissent au taux per capita
• On ne considère pas de mortalité induite par la maladie

Fonction d'incidence

• On revient sur les fonctions d'incidence dans le Cours 05
• Pour le moment, on suppose juste que l'incidence est une incidence proportionnelle (ou incidence standard) de la forme

Diagramme de flot du modèle

Le modèle

En balançant les flux entrant et sortant dans les deux compartiments, on obtient

Comme souvent en EDO, on omet la dépendance en la variable indépendante (le temps ), et on note , ce qui donne ..

Le modèle

En balançant les flux entrant et sortant dans les deux compartiments, on obtient

On considère le problème de Cauchy consistant en - auquel on adjoint les conditions initiales et

Remarques

• - est un modèle SIS (Susceptible-Infectieux-Susceptible)
• Si (pas de guérison), le modèle est un modèle SI
  • Dans ce cas, un individu infectieux reste infectieux toute sa vie
• Des maladies présentant ce type de caractéristiques sont des maladies bactériennes telles que celles causées par staphylococcus aureus, streptococcus pyogenes, chlamydia pneumoniae ou neisseria gonorrhoeae

La naissance et la mort sont relatives

Les notions de naissance et mort se réfèrent à la population que l'on considère

Par exemple, supposons un modèle pour l'immunodéficience humaine (VIH) dans une population à risque d'utilisateurs de drogue administrées de façon intraveineuse. Dans ce cas

• la naissance correspond au moment du début du comportement à risque
• la mort survient au moment où le comportement à risque cesse, que ce soit parce que l'individu meurt ou parce qu'il/elle cesse d'utiliser de la drogue

Analyse mathématique (version 1)

Analyse du système

Le système - est planaire non-linéaire

En principe, on devrait appliquer les méthodes usuelles dans le plan. Et on le fera plus tard

Toutefois, il est ici possible de trouver une solution explicite, sous certaines conditions

NB: Ceci est une illustration utile, mais est une exception!! Pratiquement aucun autre modèle que nous rencontrerons ne sera intégrable de cette façon

Le système est-il bien posé?

• On note que l'on a une fonction de la forme , donc cela peut poser problème si
• Normalement, on commence par vérifier que le système - est bien posé
• Cela fait partie de l'arsenal classique

Dynamique de

On a

Supposons que , alors et par conséquent, pour tout ,

On suppose dorénavant que et que (le cas n'est pas très captivant)

Les variables en proportions

Puisque , on peut considérer

Remarquons que . La dérivée de est

On a

mais, puisque ,

En substituant le coté droit de dans cette équation, on obtient

Le système en proportions

Puisque , on peut utiliser dans l'équation , ce qui donne

Par conséquent, le système en proportions est

Puisque est constant, les solutions de - sont déduites directement de celles de - et on considère maintenant -

Réécrivons comme

Ceci est une équation de Bernoulli et le changement de variables donne l'équation linéaire

Donc enfin

Un facteur intégrant est

et par conséquent

avec , donc enfin

La condition initiale prend la forme . Par conséquent,

ce qui implique que

Par conséquent, la solution de est

et celle de est

La solution du système en proportions est donc donnée par

et l'on voit que le terme joue un rôle important, et l'on identifie deux cas selon le signe de cette quantité

• Si , alors

donc et

• Si alors

donc et

Le nombre de reproduction élémentaire

Le nombre de reproduction élémentaire

Reformulons le résultat en termes épidémiologiques en utilisant le nombre de reproduction élémentaire

On a alors les équivalences suivantes

Aussi,

On a prouvé le résultat suivant

On dit que le modèle SIS est un modèle endémique du fait de la possibilité d'observer ce second équilibre

Quelques remarques additionnelles au sujet de

• determine la propension d'une maladie à s'établir dans une population
• Les politiques de contrôle ont par conséquent pour objectif de réduire à des valeurs plus petites que 1
• La définition "verbale" de est le nombre moyen de cas secondaires d'infection produits par l'introduction d'un unique individu infectieux dans une population complètement naïve
• Dans notre modèle, est le temps moyen de séjour dans le compartiment avant de guérir ou de mourir et est la probabilité d'infection

Quelques valeurs de (estimées)

Analyse mathématique (version 2)

Méthode classique d'analyse du système SIS

On considère le système

Comportement de la population totale

On a déjà vu que

La solution de cette EDO scalaire est facile:

Donc il y a 3 possibilités:

• si , , la population totale explose
• si , , la population totale reste constante
• si , , la population s'effondre

Donc on suppose que

• On veut un cas raisonnable, on suppose donc que

• Le système est donc

• Puisque est constant, on pourrait simplifier et étudier le système avec seulement l'une des deux équations. On choisit toutefois de garder le système en l'état

Le système est-il bien posé ?

Dans le cadre d'un modèle épidémiologique:

• les solutions de - existent-elles et sont-elles uniques ?
• le cône positif (quadrant ici) est-il invariant sous le flot de - ?
• les solutions de - sont-elles bornées ?
  • Il peut se trouver des modèles sans bornitude, mais ils sont plutôt rares et devront être considérés à part

• On évite le cas , qui ne sert à rien
• et on évite également le cas
• Donc le champ de vecteur est toujours , ce qui entraîne que les solutions existent et sont uniques
• Vérifions maintenant que le cône positif est invariant sous le flot de -

Invariance de sous le flot (1)

Supposons pour commencer que . Alors le système se réduit à l'équation scalaire

puisque ici

Il suit que l'axe est invariant et toute solution initiée dans cet ensemble est constante

Cela implique qu'une solution avec et ne peut atteindre l'axe

Invariance de sous le flot (2)

Cela implique qu'une solution avec et ne peut atteindre l'axe

Supposons que et , et qu'il existe tel que et

Invariance de sous le flot (3)

Mais en et , il passe une autre solution, celle telle que et , puisqu'avec ces conditions initiales, on a et pour tout

Cela contredit l'unicité des solutions si

Invariance de sous le flot (4)

On a vu que si

Supposons alors que . L'équation est alors

Donc si , alors pour tout . Si, d'autre part , alors pour petit; par ce qui précède, cela est aussi vrai pour tout

On dit que le champ est entrant

ne peut pas devenir nul

En résumé, pour l'invariance

• Si , alors ,

• Si , alors ,

Le modèle est donc satisfaisant, en ce qu'il n'autorise pas des solutions à changer de signe

Remarque

Ce genre de raisonnement a toute sa place dans une thèse de MSc ou de PhD: il faut montrer que vous savez faire

Dans un papier de recherche, cela n'est pas vraiment nécessaire, c'est souvent même superflu

Bornitude

Il suit de ce qui précède que le quadrant positif est (positivement) invariant sous le flot de -

De l'invariance et du fait que la population totale est bornée (constante, en fait), on déduit que les solutions de - sont bornées

Recherche des équilibres

On cherche les équilibres du système en supposant . - devient

De , on déduit que

ou alors . Substituant dans , il vient que , i.e, . C'est l'équilibre sans maladie (ESM)

De la relation

on déduit l'équilibre endémique : substituant cette valeur de dans ,

d'où il suit que

Par conséquent, l'équilibre endémique est

Notant que , il vient

Méthode classique de calcul de

est la courbe dans l'espace des paramètres où ESM perd sa stabilité locale. Pour trouver , on étudie donc la stabilité asymptotique locale de ESM

En un point arbitraire , la matrice Jacobienne de - prend la forme

La stabilité asymptotique locale de ESM dépend du signe de la partie réelle des valeurs propres de en cet équilibre, donc on évalue

Matrice triangulaire valeurs propres sont et stabilité asymptotique locale ESM déterminée par le signe de , donnant

Méthode plus efficace de calcul de

Diekmann & Heesterbeek, caractérisé dans le cas ODE par PvdD & Watmough (2002)

On considère seulement les compartiments avec des individus infectés et écrivons

• entrée au sein des compartiments infectés du fait de nouvelles infections
• contient tous les autres flux (attention au signe )

On calcule les dérivées (de Fréchet) et par rapport aux variables infectées (les Jacobiennes) et on évalue en l'ESM

Alors

où est le rayon spectral

Résultat de PvdD & Watmough (2002)

(Les conditions (A1)-(A5) sont explicitées dans le Cours 07 et nous discutons de la raison pour il est important de vérifier ces conditions dans le Cours 11)

Calcul de par matrice de prochaine génération

Ce calcul se substitue à celui de la matrice Jacobienne. Lorsque l'on a trouvé l'ESM, on ne considère que les variable infectées de -, c'est à dire

On écrit cette equation sous la forme

avec , les nouvelles infections, i.e.,

et tous les autres flux, en prenant garde au signe :

On calcule les matrices Jacobiennes de et (ici, des scalaires puisque et sont scalaires)

et

On obtient et en évaluant ces dérivées en l'ESM,

Enfin, on inverse , soit, ici,

Notons que si les conditions (A1)-(A5) du théorème sont satisfaites, on récupère aussi le fait que ESM est LAS

On trouve donc comme

Intérêt pas forcément évident ici, mais on verra plus loin à quel point cette méthode peut simplifier les calculs

Un peu de computationel

library(deSolve)
rhs_SIS <- function(t, x, p) {
  with(as.list(c(x, p)), {
    N = S + I
    dS = d*(N-S) + gamma * I - d * S - beta * S * I / N
    dI = beta * S * I / N - (d + gamma) * I
    return(list(c(dS, dI)))
  })
}
# Paramètres "connus"
d = 1/(80 * 365.25)
gamma = 1/14
# On règle beta pour que R_0 = 1.5
R_0 = 1.5
beta = R_0 * (d + gamma)
params = list(d = d, gamma = gamma, beta = beta)
IC = c(S = 1000, I = 1)
times = seq(0, 365, 1)
# On appelle l'intégrateur numérique
sol <- ode(y = IC, times = times, func = rhs_SIS, 
           parms = params, method = "ode45")

Je viens de faire ...

ce que je préconise de ne pas faire: illustrer un résultat mathématique sans rien vraiment ajouter au résultat lui même

Améliorons les choses un peu. Voir le code

On pourrait continuer, mais pour un modèle aussi simple, il y a peu de choses à faire: les 3 paramètres du système se combinent dans et ce dernier résume bien la dynamique

On va toutefois représenter une chose importante: le schéma de bifurcation

On a vu que lorsque , tandis que lorsque , . Représentons ceci (code)

Un peu plus au sujet de la stabilité

Ici, j'ai triché !

Dans cette 2ème analyse, tout ce qu'on a montré, c'est la stabilité/instabilité locale de l'ÉSM. (On a tout montré avec la première analyse...)

Ce qui reste à faire (2ème analyse)

• Stabilité asymptotique locale de l'ÉE lorsqu'il est relevant biologiquement
• Stabilité asymptotique globale de l'ÉSM lorsque
• Stability asymptotique globale de l'ÉE lorsque

Équilibre relevant biologiquement

On se souvient que l'on avait trouvé

Ici, la valeur de n'a de sens comme équilibre d'un système en compartiments que si , i.e., si . On dit alors que l'ÉE est relevant biologiquement

SAL de l'ÉE lorsqu'il est relevant biologiquement

Outre l'ÉE

on avait aussi calculé la Jacobienne en un point arbitraire

Évaluons en :