# Modélisation probabiliste de la dynamique de transmission de la fièvre typhoïde à Mayotte avec étude de risques épidémiques

Ibrahim Bouzalmat

► To cite this version:

Ibrahim Bouzalmat. Modélisation probabiliste de la dynamique de transmission de la fièvre typhoïde à Mayotte avec étude de risques épidémiques. Statistiques [math.ST]. Université de Montpellier, 2023. Français. NNT: . tel-04439679

HAL Id: tel-04439679

https://hal.science/tel-04439679

Submitted on 5 Feb 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# THÈSE POUR OBTENIR LE GRADE DE DOCTEUR DE L'UNIVERSITÉ DE MONTPELLIER

En Biostatistiques

École doctorale : Information, Structures, Systèmes

Unité de recherche : Institut Montpelliérain Alexander Grothendieck

## Modélisation probabiliste de la dynamique de transmission de la fièvre typhoïde à Mayotte avec étude de risques épidémiques

Présentée par Ibrahim BOUZALMAT
le 15/11/2023

Sous la direction de Benoîte DE SAPORTA
et Solym MANOU-ABI

Devant le jury composé de

Benoîte De Saporta
Marie-Pierre Etienne
Viet Chi Tran
Nathalie Peyrard
Benoit Liquet
Manon Costa
Bertrand Cloez
Solym Mawaki Manou-Abi

Professeur, Université de Montpellier
Maître de conférences, Institut Agro Rennes-Angers
Professeur, Université Gustave Eiffel
Directeur de recherche, INRAE
Professor, UPPA
Maître de conférences, Université Toulouse 3 - Paul Sabatier
Chargé de recherche, INRAE Montpellier
Maître de conférences, CUFR de Mayotte

Directrice de thèse
Rapporteur
Rapporteur
Présidente
Examinateur
Examinatrice
Examinateur
Co-encadrant de thèse
<figure><img src="image_2.png" /></figure>
À l'âme de mon père.

"L'âme de mon père réside en moi, tissant un fil invisible qui me guide, m'inspire et m'accompagne tout au long de ma vie."
<figure><img src="image_3.png" /></figure>
Remerciements

Je tiens à remercier toutes les personnes qui ont contribué peu ou prou au succès de cette thèse, et je tiens à saisir cette occasion unique pour les remercier.

En premier lieu, je tiens à remercier profondément ma directrice de thèse, Benoîte De Saporta. Son remarquable sens de l'organisation a été une source d'inspiration. Elle a guidé de près mon travail de recherche sur trois années, apportant constamment clarté et rigueur à nos échanges. Sa capacité à maintenir un haut niveau d'organisation dans notre collaboration a été essentielle pour le succès de ma thèse. De plus, je tiens à souligner que ma directrice est non seulement bien organisée, mais aussi d'une gentillesse exceptionnelle. Elle a toujours pris le temps de m'écouter, de répondre à mes questions et de m'apporter son soutien, peu importe la situation. Sa gentillesse naturelle a créé un environnement de travail agréable et motivant, favorisant ma progression dans mes recherches. Je voudrais également remercier mon co-encadrant de thèse, Solym Manou-Abi, qu'a fait preuve d'un dévouement remarquable tout au long de cette collaboration, en trouvant toujours du temps pour nos réunions malgré la distance et les contraintes qui se sont présentées. Malgré les barrières géographiques qui nous séparaient, mon co-encadrant a su maintenir un suivi régulier de notre travail en organisant des réunions à distance. Sa disponibilité et son engagement ont été essentiels pour le succès de notre projet de recherche.

Je tiens ensuite à remercier Marie-Pierre Etienne et Viet Chi Tran qui me font l'honneur de rapporter cette thèse ainsi que Nathalie Peyrard, Benoit Liquet, Manon Costa et Bertrand Cloez qui ont gracieusement accepté de faire partie du jury.

Je tiens également à exprimer ma gratitude envers toutes les personnes qui ont contribué à ce travail, avec une mention spéciale à Julien Balicchi et à l'Agence régionale de santé (ARS) de Mayotte pour leur précieuse collaboration dans la collecte des données. De plus, je souhaite souligner que ce projet a bénéficié d'un soutien financier partiel du LabEx NUMEV. Je suis reconnaissant d'avoir eu l'opportunité de participer aux journées scientifiques organisées par les membres de ce laboratoire d'excellence. J'ai également beaucoup apprécié les échanges constructifs avec le directeur, Andrea Parmeggiani de L2C, ainsi que Nabil Zemiti du LIRMM. Enfin, je tiens à exprimer ma gratitude envers Geneviève Carrière, la gestionnaire qui a veillé à la coordination de mes missions scientifiques. Toutes ces contributions ont été essentielles à la réussite de ce travail.

Intégrer l'équipe des chercheur/es et enseignant/es de l'IMAG à l'Université de Montpellier m'a ouvert la porte à de précieuses rencontres et à l'acquisition de connaissances. J'ai eu le privilège de faire la connaissance de personnes exceptionnelles, des amis et des collègues, chacun/une ayant sa place spéciale dans mon esprit. Je tiens à exprimer ma profonde gratitude à Thiziri Moulla, Florian Miralles, Marwa Salah, Pascal Azerad, Ali Gannoun, ainsi qu'à ceux qui ont depuis quitté le laboratoire, à savoir Ali Haidar,
Chayma Daayeb, Franschesko, Hanine Awada, Zaineb Smida et Antonia Jabbour. Un remerciement particulier revient à Sophie Cazanave, Nathalie Collain et Brigitte Fabric, dont le dévouement exemplaire envers les membres de l'équipe a été remarquable.

Je ne saurais oublier de mentionner mes co-bureaux, Marien Hanno, notre maestro du laboratoire, un collègue hors pair et un ami toujours présent pour des discussions mathématiques, et Paul-Aimé-Latsouck Faye, un ami doté de nombreuses vertus et de valeurs admirables. Comme l'a si bien dit le chef du restaurant administratif, il est un véritable "prêtre sur terre".

Et bien sûr, je tiens également à exprimer ma gratitude envers ma famille : ma mère et ma sœur, qui m'ont apporté un soutien inestimable tout au long de cette aventure. Leur présence et leurs encouragements ont été une source de réconfort et de motivation, et je leur suis profondément reconnaissant pour leur soutien indéfectible.

Ibrahim BOUZALMAT
# Résumé

Ce manuscrit de thèse vise à étudier la transmission de la fièvre typhoïde à Mayotte en utilisant des approches de modélisation mathématique. Nous introduisons tout d'abord le contexte de notre étude, les problématiques associées et les objectifs de la thèse. Un état de l'art sur la modélisation mathématique de la transmission de la fièvre typhoïde est présenté, mettant en évidence la spécificité de notre approche.

Nous proposons un modèle initial en deux versions, déterministe et stochastique, pour décrire la dynamique de transmission de la maladie à Mayotte. Nous explorons le comportement du modèle à travers des simulations numériques dans différents scénarios, en mettant en évidence les facteurs clés de la transmission.

Cependant, en raison des limitations du jeu de données disponibles, nous proposons un modèle stochastique simplifié et une méthode d'estimation paramétrique. Cette approche nous permet d'ajuster le modèle aux données disponibles et d'estimer les caractéristiques clés de la transmission de la fièvre typhoïde à Mayotte.

En enrichissant notre modèle, nous introduisons de nouvelles extensions. Nous incluons un compartiment pour les individus exposés, prenant en compte la durée d'incubation de la maladie. Les propriétés théoriques de ce modèle sont étudiées et illustrées par des simulations numériques. De plus, nous proposons une méthodologie d'estimation des paramètres adaptée à ce nouveau modèle, et des simulations numériques ont été réalisées pour évaluer la performance de notre approche d'estimation.

Nous examinons ensuite l'impact de la pluviométrie sur la transmission de la fièvre typhoïde à Mayotte, en utilisant des données publiques de précipitations. Nous identifions une saisonnalité de la pluie et estimons les paramètres du modèle sous différents régimes. Les résultats soulignent l'importance de cette variable météorologique dans la propagation de l'épidémie.

Ce manuscrit ouvre des perspectives de recherche, telles que l'extension du modèle à d'autres maladies infectieuses présentes à Mayotte et sa généralisation à d'autres territoires. Ces travaux contribueront à une meilleure compréhension et gestion des maladies infectieuses à Mayotte et dans d'autres régions similaires.
Abstract

The aim of this thesis manuscript is to study the transmission of typhoid fever in Mayotte using mathematical modelling approaches. We first introduce the context of our study, the associated issues, and the objectives of the thesis. A state-of-the-art review on mathematical modeling of typhoid fever transmission is presented, highlighting the specificity of our approach.

We propose an initial model in two versions, deterministic and stochastic, to describe the transmission dynamics of the disease in Mayotte. We explore the behavior of the model through numerical simulations in different scenarios, highlighting key factors of transmission.

However, due to the limitations of the available dataset, we propose a simplified sto-
chastic model and a parametric estimation method. This approach enables us to fit the
model to the available data and to estimate the key characteristics of typhoid fever trans-
mission in Mayotte.

In enriching our model, we are introducing new extensions. We include a compart-
ment for individuals exposed, taking into account the incubation period of the disease.
The theoretical properties of this model are studied and illustrated by numerical simula-
tions. In addition, we propose a parameter estimation methodology adapted to this new
model, and numerical simulations have been carried out to evaluate the performance of
our estimation approach.

We then examine the impact of rainfall on the transmission of typhoid fever in Mayotte, using publicly available precipitation data. We identify rainfall seasonality and estimate model parameters under different regimes. The results highlight the importance of this meteorological variable in the spread of the epidemic.

This manuscript opens up research perspectives, such as the extension of the model to other infectious diseases present in Mayotte and its generalisation to other territories. This work will contribute to a better understanding and management of infectious diseases in Mayotte and other similar regions.
<table>
  <tbody>
    <tr>
      <td>3.3.1</td>
      <td>Expression analytique des paramètres</td>
      <td>61</td>
    </tr>
    <tr>
      <td>3.3.2</td>
      <td>Modèle de Markov caché</td>
      <td>69</td>
    </tr>
    <tr>
      <td>3.4</td>
      <td>Étude numérique</td>
      <td>74</td>
    </tr>
    <tr>
      <td>3.4.1</td>
      <td>Données synthétiques</td>
      <td>74</td>
    </tr>
    <tr>
      <td>3.4.2</td>
      <td>Données de la typhoïde à Mayotte</td>
      <td>79</td>
    </tr>
    <tr>
      <td>3.4.3</td>
      <td>Application à l'estimation de la probabilité d'atteindre un seuil épi-démique dans une période de temps</td>
      <td>80</td>
    </tr>
    <tr>
      <td>3.5</td>
      <td>Conclusion</td>
      <td>84</td>
    </tr>
    <tr>
      <td>4</td>
      <td>Modélisation et estimation des paramètres d'un modèle multi-chaîne cachée</td>
      <td>85</td>
    </tr>
    <tr>
      <td>4.1</td>
      <td>Introduction</td>
      <td>86</td>
    </tr>
    <tr>
      <td>4.2</td>
      <td>Modèle et propriété</td>
      <td>87</td>
    </tr>
    <tr>
      <td>4.2.1</td>
      <td>Définition du modèle</td>
      <td>87</td>
    </tr>
    <tr>
      <td>4.2.2</td>
      <td>Propriétés mathématiques du modèle</td>
      <td>88</td>
    </tr>
    <tr>
      <td>4.3</td>
      <td>Estimation des paramètres</td>
      <td>92</td>
    </tr>
    <tr>
      <td>4.3.1</td>
      <td>Observations</td>
      <td>93</td>
    </tr>
    <tr>
      <td>4.3.2</td>
      <td>Stratégie d'estimation</td>
      <td>94</td>
    </tr>
    <tr>
      <td>4.3.3</td>
      <td>Expression analytique des paramètres</td>
      <td>94</td>
    </tr>
    <tr>
      <td>4.3.4</td>
      <td>Modèle de Markov Multi-Chaîne Cachée</td>
      <td>97</td>
    </tr>
    <tr>
      <td>4.4</td>
      <td>Simulation numérique</td>
      <td>102</td>
    </tr>
    <tr>
      <td>4.4.1</td>
      <td>Observations complètes</td>
      <td>102</td>
    </tr>
    <tr>
      <td>4.4.2</td>
      <td>Observations du cumul des nouveaux cas déclarés</td>
      <td>103</td>
    </tr>
    <tr>
      <td>4.5</td>
      <td>Conclusion</td>
      <td>107</td>
    </tr>
    <tr>
      <td>5</td>
      <td>Effet de la pluviométrie sur l'évolution de la fièvre typhoïde à Mayotte</td>
      <td>108</td>
    </tr>
    <tr>
      <td>5.1</td>
      <td>Introduction</td>
      <td>109</td>
    </tr>
    <tr>
      <td>5.2</td>
      <td>Description des données de précipitations</td>
      <td>109</td>
    </tr>
    <tr>
      <td>5.3</td>
      <td>Détection de rupture</td>
      <td>110</td>
    </tr>
    <tr>
      <td>5.3.1</td>
      <td>Méthode CUSUM</td>
      <td>111</td>
    </tr>
    <tr>
      <td>5.3.2</td>
      <td>Application aux données de précipitations à Mayotte</td>
      <td>111</td>
    </tr>
    <tr>
      <td>5.4</td>
      <td>Estimation des paramètres d'un modèle de Markov caché à deux régimes</td>
      <td>114</td>
    </tr>
    <tr>
      <td>5.5</td>
      <td>Conclusion</td>
      <td>118</td>
    </tr>
    <tr>
      <td>6</td>
      <td>Perspectives</td>
      <td>119</td>
    </tr>
    <tr>
      <td>Annexe</td>
      <td></td>
      <td>123</td>
    </tr>
    <tr>
      <td>A</td>
      <td>Formulaire de déclaration de la fièvre typhoïde</td>
      <td>123</td>
    </tr>
    <tr>
      <td>B</td>
      <td>Matrice Jacobienne Dg</td>
      <td>124</td>
    </tr>
    <tr>
      <td>C</td>
      <td>Log-vraisemblance l(O|M<sup>n</sup>)</td>
      <td>127</td>
    </tr>
    <tr>
      <td>D</td>
      <td>Fonction générateur des moments du processus (E<sub>t</sub>, I<sub>t</sub>)<sub>t≥0</sub></td>
      <td>128</td>
    </tr>
    <tr>
      <td>E</td>
      <td>Séries chronologiques</td>
      <td>130</td>
    </tr>
    <tr>
      <td>F</td>
      <td>Détection de rupture</td>
      <td>135</td>
    </tr>
    <tr>
      <td>Bibliographie</td>
      <td></td>
      <td>144</td>
    </tr>
  </tbody>
</table>
# Chapitre 1

## Introduction générale

### Sommaire

<table><tr><td>1.1</td><td>Introduction</td><td>9</td></tr><tr><td>1.2</td><td>Contexte</td><td>9</td></tr><tr><td>1.2.1</td><td>La fièvre typhoïde</td><td>9</td></tr><tr><td>1.2.2</td><td>Description de Mayotte</td><td>10</td></tr><tr><td>1.2.3</td><td>Accès à l'eau à Mayotte</td><td>12</td></tr><tr><td>1.2.4</td><td>La Typhoïde à Mayotte</td><td>13</td></tr><tr><td>1.3</td><td>Objectifs de la thèse</td><td>15</td></tr><tr><td>1.4</td><td>Outils mathématiques</td><td>15</td></tr><tr><td>1.4.1</td><td>Modèles à compartiments</td><td>15</td></tr><tr><td>1.4.2</td><td>Processus Markovien Déterministe par Morceaux</td><td>17</td></tr><tr><td>1.4.3</td><td>Processus de naissance et mort linéaire avec immigration</td><td>18</td></tr><tr><td>1.4.4</td><td>Modèle de Markov caché</td><td>24</td></tr><tr><td>1.5</td><td>État de l'art sur la modélisation mathématique de la fièvre typhoïde</td><td>29</td></tr><tr><td>1.6</td><td>Organisation de la thèse</td><td>30</td></tr></table>
# 1.1 Introduction

L'épidémiologie est une science qui cherche à comprendre la dynamique régissant la propagation de maladies infectieuses, et étudier les facteurs de risques qui contribuent à cette propagation, leur prévalence et leurs conséquences afin d'établir des stratégies de prévention et d'intervention permettant de diminuer leur impact sur la santé publique. Le premier modèle utilisé en épidémiologie est celui qui décrit les infections par la vérole (ou variole). Il a été publié par Daniel Bernoulli [1] en 1760. C'est la première fois que l'ensemble des observations épidémiologiques sur la variole ont été collectées, combinées avec les données démographiques de l'époque et insérées dans un modèle mathématique afin de tester les hypothèses et d'évaluer le rapport bénéfice-risque de la variolisation. Depuis lors et au début du XXe siècle, un grand nombre de modèles épidémiologiques ont été développés [2]-[3], avec en général l'un des objectifs principaux suivants :

* Traduire schématiquement la propagation de certaines maladies infectieuses, ce qui devrait permettre de mieux comprendre les phénomènes épidémiques et donc de mieux les contrôler.

* Décrire théoriquement et numériquement la dynamique épidémique des agents pathogènes d'intérêt.

* Quantifier l'impact et l'efficacité des différentes méthodes de prévention, évaluer également les risques associés à cette épidémie et aux interventions mise en place.

* Planifier, tester et optimiser des plans d'expérience, en complément de l'expérimentation par exemple.

Le but des modèles en épidémiologie est de décrire, comprendre, prévoir et/ou contrôler la dynamique des pathogènes. Dans ce contexte, la modélisation permet d'isoler les effets de certains facteurs et de sélectionner les interactions dont on tiendra compte. Par ailleurs, l'étude des effets de différents facteurs, ainsi que de leurs interactions est grandement enrichie par la possibilité de réplication des mesures bien supérieures à ce que l'expérimentation pourrait apporter. De plus, par le jeu des changements des valeurs des paramètres, les modèles permettent d'étudier rapidement un grand nombre de situations. La littérature et le développement de l'épidémiologie mathématique sont bien documentés et peuvent être trouvés par exemple dans [4], [5] et [6].

Dans ce travail, nous nous intéressons à modéliser et étudier la dynamique de transmission de la maladie de la fièvre typhoïde à Mayotte dans le but de quantifier son influence et son impact sur la santé globale de la population.

# 1.2 Contexte

## 1.2.1 La fièvre typhoïde

La fièvre typhoïde (du grec tuphos, torpeur) ou typhus abdominal est une maladie infectieuse décrite en 1818 par Pierre Bretonneau, causée par une bactérie de la famille Entérobactérie, du genre des salmonelles, une bactérie qui se multiplie dans l'estomac et les intestins de personnes infectées, et dont les espèces responsables sont Salmonella enterica -Typhi ou Paratyphi A, B, C.

La principale voie de transmission de la fièvre typhoïde est la voie oro-fécale. Cela signifie qu'elle peut se propager directement à partir de porteurs infectés par le biais de matières
fécales contaminées, d'urines et de vomissements, ou indirectement par la consommation d'eau ou d'aliments contaminés (Fig.1.1).

Les symptômes de la fièvre typhoïde comprennent de la fièvre élevée, des maux de tête, une perte d'appétit, des douleurs abdominales et une faiblesse générale. Dans les cas graves, elle peut causer des complications sérieuses, y compris des saignements internes, une perforation intestinale et même la mort. Une fois diagnostiquée, la fièvre typhoïde peut être traitée avec des antibiotiques. Il est important de prévenir la fièvre typhoïde en adoptant des mesures d'hygiène alimentaire, en se vaccinant contre cette maladie et en traitant une infection présente dès que possible.

<figure><img src="image_4.png" /><figcaption>FIGURE 1.1 – Modes de transmission de la Typhoide.</figcaption></figure>

1.2.2 Description de Mayotte

Mayotte est une île de l'océan Indien, située dans l'archipel des Comores. C'est à
la fois une région insulaire française et un département de France d'outre-mer qui sont
administrés dans le cadre d'une collectivité territoriale unique dirigée par le conseil dé-
partemental de Mayotte [7]. Mayotte se compose de deux îles principales, Grande-Terre
et Petite-Terre, entourées d'un lagon turquoise. Les reliefs sont formés de collines et de
montagnes, couvertes de forêts tropicales luxuriantes. Elle est célèbre pour ses plages de
sable blanc, ses eaux cristallines et sa faune marine diversifiée. Les principales activités
économiques de l'île sont l'agriculture, la pêche et le tourisme.

Mayotte est le plus petit département d'outre-mer français, et le moins peuplé, cependant c'est aussi le territoire avec la plus forte croissance démographique, et la population de l'île devrait bientôt dépasser celle de la Guyane. Mayotte est cependant déjà largement en tête en terme de densité de population, avec 682 hab./km² (devant la Réunion, avec seulement 340 hab./km² [8]).

Elle est en particulier marquée par une expansion démographique spectaculaire : la population est passée de 11 000 habitants en 1911 à 256 518 en 2017 avec un taux de fécondité qui s'élève à 4,68 enfants par femme en 2019 [8], très au-dessus des autres départements français, et similaire à celui de l'Union des Comores voisine, l'hôpital de Mamoudzou est ainsi la première maternité d'Europe (Fig.1.2).
<figure><img src="image_5.png" /></figure>

Source : Insee, recensements de la population.

<caption>FIGURE 1.2 – Croissance démographique de la population de Mayotte entre 1958 et 2017.</caption>

Sources : Insee, recensements de la population.

Entre 2012 et 2017, la répartition de la population de Mayotte par âge a très peu changé. L'âge moyen des résidents en 2017 était de 23 ans. Du fait de l'immigration et du taux de natalité élevé, le vieillissement de la population est encore limité. La population est encore beaucoup plus jeune que dans les autres régions de France : l'âge moyen est de 41 ans en métropole, de 35 ans à La Réunion et de 28 ans en Guyane [9].

En 2012, la moitié de la population avait moins de 18 ans et les trois dixièmes des Mahorais avaient moins de 10 ans. 11 000 personnes de plus de 60 ans ne représentent que 4% de la population totale, soit six fois moins qu'en France métropolitaine [9].

La pyramide des âges présente un creux marqué autour de 20 ans, du fait du déficit migratoire à ces âges (Fig.1.3).
<figure><img src="image_6.png" /><figcaption>FIGURE 1.3 – Pyramide des âges à Mayotte (2012 et 2017) et en France métropolitaine (2017).</figcaption></figure>

### 1.2.3 Accès à l'eau à Mayotte

Mayotte, en tant que département qui représente un contexte particulier de précarisation, présente une forte prévalence des maladies hydriques. Les difficultés pour accéder à l'eau courante et les conditions d'hygiène insuffisantes posent évidemment des problèmes sanitaires à Mayotte. Un habitant sur trois n'a pas à accès à l'eau à l'intérieur de son logement et six logements sur dix ne proposent toujours pas le confort sanitaire de base (eau courante, toilettes, douche). A cette situation de très forte précarité s'ajoute l'absence d'un système d'assainissement de qualité.

Dans ce contexte, l'agence régionale de santé (ARS) est chargée d'organiser et de mettre en oeuvre le contrôle sanitaire des eaux destinées à la consommation. Outre cette mission, ses actions portent également sur l'amélioration de l'accès à l'eau potable pour tous puisque les habitants iraient s'alimenter chez un voisin, à une borne fontaine monétique (BFM), au puits ou encore à la rivière et peuvent recourir à un stockage de l'eau dans des conditions douteuses. Ainsi, l'installation de BFM (Fig.1.4) sur les quartiers d'habitat insalubre où résident les populations les plus démunies a été identifiée comme action de prévention prioritaire, depuis 2000, avec l'apparition de cas de choléra à Mayotte.

Les BFM permettent de fournir de l'eau potable aux habitants des quartiers qui ne sont pas raccordés au réseau de distribution d'eau potable. Financées par l'ARS, installées par la Société Mahoraise des Eaux (SMAE), elles permettent à la population de récupérer le volume d'eau souhaité, à l'aide d'une carte (en vente à la SMAE). Depuis 2010, de
nouvelles implantations et des déplacements de BFM réalisés en faveur des quartiers de
grande précarité portent le nombre de BFM en fonctionnement à 52, dont 33 se situent
sur les communes de Mamoudzou et de Koungou. Au cours de la période 2015-2020, une
demande de crédit de 1 million d'euros sur six ans a été déposée dans le cadre du contrat
de plan État-région, ce qui équivaut à l'installation de 10 BFM/an pendant six ans en
sachant qu'au regard des données du recensement de 2012, le nombre de logements sans
accès à l'eau nécessiterait actuellement l'installation de 43 BFM (considérant qu'un BFM
dessert environ 100 maisons) [10].

<figure><img src="image_7.png" /><figcaption>FIGURE 1.4 – Des fontaines d'eau potable à Mayotte [10].</figcaption></figure>

Par conséquent, la procédure d'installation de BFM à Mayotte devrait permettre
notamment de lutter contre la propagation des maladies liées à la consommation d'eau
non conforme comme la fièvre typhoïde.

1.2.4 La Typhoïde à Mayotte

La fièvre typhoïde est une maladie endémique à Mayotte, et il est essentiel de souligner
qu'elle est soumise à une déclaration obligatoire (voir la formulaire de déclaration en An-
nexe A). Au cours de l'année 2022, 100 cas de fièvre typhoïde confirmés par hémoculture
et/ou coproculture ont été signalés par le laboratoire du Centre Hospitalier de Mayotte
(CHM) à la Cellule de Veille, d'Alerte et de Gestion Sanitaire (CVAGS) de l'ARS. Le
taux d'incidence de la fièvre typhoïde est en hausse au cours des trois dernières années,
passant de 14 cas pour 100 000 habitants en 2020 à 35.7 cas pour 100 000 habitants en
2022. De 2016 à 2021, 39 cas de fièvre typhoïde par an, en moyenne, ont été signalés par
l'ARS de Mayotte. [11] (Fig.1.5).

Parmi les 100 cas de fièvre typhoïde enregistrés en 2022, 71 cas ont été hospitalisés entre 1 et 7 jours (71.2 %) dont 9 dans le service de réanimation. 74 cas (73,7 %) sont survenus entre les semaines 17 et 25. Parmi les cas examinés, 53 cas avaient accès à l'eau courante ou à une borne-fontaine (BFM) à proximité et 9 cas ont déclaré ne pas avoir accès à l'eau courante et utiliser de l'eau de rivière ou de puits. Les cas de fièvre typhoïde sont survenus
<figure><img src="image_8.png" /><figcaption>FIGURE 1.5 – Nombre de cas confirmés de fièvre typhoïde par année de prélèvement, janvier 2016-août 2022, Mayotte.</figcaption></figure>

Source : enquêtes réalisées par la CVAGS de Mayotte, ARS.

tout au long de l'année à Mayotte, néanmoins, une recrudescence de cas a été observée entre les mois de février et juin de chaque année, avec des pics épidémiques variés d'une année à l'autre (Fig.1.6). Les raisons de ces variations peuvent être multiples, notamment des facteurs environnementaux, socio-économiques et comportementaux. Par exemple, des conditions météorologiques particulières, une densité de population élevée dans certaines zones, une insalubrité des installations sanitaires ou une faible hygiène personnelle peuvent favoriser la propagation de la bactérie Salmonella Typhi, qui est responsable de la fièvre typhoïde.

<figure><img src="image_9.png" /><figcaption>FIGURE 1.6 – Évolution hebdomadaire des nouveaux cas confirmés de fièvre typhoïde par semaine de prélèvement, janvier 2018-décembre 2022.</figcaption></figure>
## 1.3 Objectifs de la thèse

Cette thèse vise à proposer et étudier de nouveaux modèles de propagation de la fièvre typhoïde adaptés à Mayotte. L'objectif principal est de prendre en compte l'accès à l'eau potable des populations, ainsi que la durée d'incubation de la maladie et l'impact de la pluviométrie sur la dynamique de transmission. Plus précisément, nous visons à :

* Définir un modèle global de propagation de la fièvre typhoïde pour Mayotte, en tenant compte des spécificités de l'île comme l'accès ou non à l'eau potable.

* Étudier théoriquement et numériquement le comportement des modèles proposés.

* Estimer les paramètres des modèles à partir des données fournies par l'Agence Régionale de Santé de Mayotte correspondant au nombre cumulé de nouveaux cas déclarés par jour.

* Prendre en compte la durée d'incubation de la maladie dans notre modèle afin de mieux comprendre la dynamique de cette maladie infectieuse.

* Étudier l'impact de la pluviométrie sur la transmission de la fièvre typhoïde à Mayotte en utilisant des données de précipitations et en identifiant des ruptures saisonnières.

* Explorer les perspectives de recherche futures dans ce domaine.

La thèse contribuera également à développer l'axe III du volet recherche du projet d'établissement du Centre Universitaire de Formation et de Recherche de Mayotte, à savoir : *Modélisations et analyse épistémologique des systèmes complexes pour comprendre le fonctionnement et simuler le devenir des systèmes insulaires : développer des outils de modélisation statistique, descriptive ou dynamique dans une approche interdisciplinaire.*
Cette démarche correspond aux priorités de recherche identifiées par la commission scientifique du Centre Universitaire de Mayotte.

## 1.4 Outils mathématiques

Afin de garantir la rigueur et l'analyse nécessaires pour atteindre les objectifs de cette thèse, il est essentiel de s'appuyer sur des concepts et des outils mathématiques appropriés. Cette section fournit un aperçu détaillé des outils mathématiques clés qui sont utilisés tout au long de cette thèse.

### 1.4.1 Modèles à compartiments

Les modèles compartimentaux, également connus sous le nom de modèles à compartiments, sont des outils mathématiques puissants pour décrire la dynamique d'une population où les individus traversent successivement plusieurs états, représentés par des compartiments distincts. Dans ces modèles, la population est répartie dans les différents compartiments en fonction de leur statut infectieux par exemple. Chaque compartiment est alors associé à une dynamique propre et des échanges avec les autres compartiments.
Les modèles compartimentaux sont largement utilisés dans le domaine de l'épidémiologie pour étudier la propagation des maladies infectieuses au sein d'une population. Ils permettent de modéliser comment les individus passent d'un compartiment à un autre en fonction de leur statut infectieux et des interactions qu'ils ont les uns avec les autres. Par exemple, les modèles SI, SIS, SIR, SIRS, etc., sont couramment utilisés pour représenter différentes classes épidémiologiques, tels que les Susceptibles, les Infectés, les Retirés (morts ou immunisés), etc. La dynamique de ces modèles, dans le cas déterministe, est généralement décrite par un ensemble d'équations différentielles [2, 12]. Elles peuvent être basées sur diverses hypothèses, telles que la transmission de la maladie à travers des contacts directs ou indirects entre individus, la distribution de la population dans différents groupes d'âge ou de susceptibilité, et la variation des taux d'infection, de guérison et de mortalité. Ces modèles mathématiques permettent de prédire l'évolution temporelle de l'épidémie, d'estimer les paramètres épidémiologiques clés, et d'évaluer l'impact des interventions de santé publique sur la propagation de la maladie.

La **Fig.1.7** représente un schéma d'un modèle à compartiment de type **SIR** dont les retirés ne redeviennent jamais susceptibles, où les bacs représentent les compartiments dans lesquels sont regroupés les individus de même statut épidémiologique (susceptibles, infectés, retirés) et les robinets représentent les échanges entre ces compartiments.

Dans le chapitre 2 (2.2), nous utiliserons un modèle à compartiments pour étudier la dynamique de propagation de la fièvre typhoïde et analyser l'interaction entre les différentes populations impliquées dans le processus épidémique.

<figure><img src="image_10.png" /><figcaption>FIGURE 1.7 – Représentation schématique d'un modèle à compartiment.</figcaption></figure>
### 1.4.2 Processus Markovien Déterministe par Morceaux

Selon une définition de Davis [13], la classe des Processus Markoviens Déterministes par Morceaux (en anglais : Piecewise Deterministic Markov Process (PDMP)) est une famille générale de modèles stochastiques hybrides couvrant pratiquement toutes les applications de non diffusion.

Un PDMP ($X_t$) = ($M_t, Y_t$) est un processus de Markov en temps continu à valeurs dans un espace d'état $E = M \times \bigcup_{m \in M} F_m$, où M est un ensemble fini ou dénombrable et $F_m$ est un ouvert de $\mathbb{R}^d$ pour chaque $m \in M$. Il est constitué de deux variables. La première $M_t$, appelée *régime* ou *mode*, est discrète et permet de caractériser l'environnement et/ou le régime dans lequel le processus évolue. La seconde variable, notée $Y_t$, est appelée *variable euclidienne* et décrit les grandeurs physiques du modèle.

Un exemple concret pourrait être un modèle météorologique où la variable discrète *mode* pourrait représenter les différents types de temps (ensoleillé, nuageux, pluvieux, etc.), tandis que la variable continue *variable euclidienne* pourrait représenter des grandeurs physiques telles que la température, la pression atmosphérique, l'humidité, etc. Le modèle pourrait alors utiliser la variable discrète pour déterminer le régime météorologique dans lequel il se trouve, et la variable continue pour représenter les grandeurs physiques associées à ce régime.

Les PDMP sont une grande famille de différents processus stochastiques qui comprend les processus de Markov en temps discret, les chaînes de Markov en temps continu, les processus déterministes avec sauts, les systèmes dynamiques avec commutation aléatoire et certains processus ponctuels. Ce type de modèle peut être utilisé dans de nombreux domaines, tels que la physique, la biologie, la dynamique des populations, la fiabilité des systèmes, etc. [14, 15, 16, 17].

L'évolution d'un PDMP est déterminée par trois caractéristiques locales (dépendant du mode $m \in M$) :

* Le flot $\phi_m : \mathbb{R}^d \times \mathbb{R}^+ \to \mathbb{R}^d$ qui décrit la trajectoire déterministe du processus entre deux sauts. Il peut s'interpréter comme une solution au temps $t$ d'un système d'équations différentielles partant de $y$ au temps 0. Il vérifié également la propriété de semi-groupe,

$$ \forall s, t \ge 0, \quad \phi_m(\cdot, t+s) = \phi_m(\phi_m(\cdot, s), t). $$

On appelle $t^*(m, y)$ le temps déterministe d'atteinte de la frontière de l'espace $F_m$ en partant du point $x = (m, y)$ :

$$ t^*(m, y) = \inf \{t > 0 : \phi_m(y, t) \in \partial F\}. $$

Si $t^*(m, y) = \infty$, cela signifie que le processus n'atteint pas la frontière de $F$ en temps fini.

* L'intensité de saut dans le mode $m$, $\lambda_m : \bar{F} \to \mathbb{R}^+$ une fonction mesurable. Elle vérifie la propriété d'intégrabilité locale suivante :

$$ \forall x = (m, y) \in \mathbb{E}, \quad \exists \varepsilon > 0, \quad \int_0^\varepsilon \lambda_m(\phi_m(y, t)) dt < +\infty. $$

On suppose que le temps de saut aléatoire $T$ est distribué comme suit :

$$ \mathbb{P}_{(m,y)}(T > t) = e^{-\int_0^t \lambda_m(\phi_m(y,s))ds}. $$
Ainsi, le premier temps du saut $T_1$ est défini par :

$$T_1 = \min\{t^*(m, y), T\}.$$

C'est le minimum entre le temps déterministe d'atteinte de la frontière et le temps de saut aléatoire, sa loi est donnée par :

$$\mathbb{P}_{(m,y)}(T_1 > t) = e^{-\int_0^t \lambda_m(\phi_m(y,s))ds} \mathbb{1}_{\{t<t^*(m,y)\}}. \quad (1.1)$$

* Le noyau de transition $Q_m : \bar{F} \times \mathcal{B}(\bar{F}) \to [0, 1]$, qui permet de déterminer la position après le saut. Le nouveau mode et/ou position $(M_{T_1}, Y_{T_1})$ sélectionnés au temps $T_1$ par le noyau markovien $Q_m$ :

$$\mathbb{P}_{(m,y)}((M_{T_1}, Y_{T_1}) \in A | T_1) = \int_A Q_m(\phi_m(y, T_1), dx).$$

La dynamique d'un PDMP peut être définie de manière itérative suivante : le processus démarre d'un point $x_0 = (m, y)$ dans l'espace d'état $\mathbb{E}$. Le premier temps de saut $T_1 = S_1$ est déterminé par sa loi **1.1**, il peut se produire soit aléatoirement avec un taux $\lambda_m$ ou bien lorsqu'on touche la frontière $\partial F$. Sur l'intervalle $[0, T_1[$, la variable euclidienne $Y_t$ suit le flot déterministe $t \to \phi_m(y, t)$. La position du processus à l'instant de saut $T_1$ est sélectionnée selon la distribution $Q_m(\phi_m(y, T_1), \cdot)$, et le processus repart à partir d'une nouvelle position $X_{T_1} = (m_1, Y_1)$. Ensuite, et de même manière, un temps inter-saut $S_2 = T_2 - T_1$ est déterminée à partir de la fonction de survie **1.1**, et sur l'intervalle $[T_1, T_2[$, le processus $X_t$ suit le flot déterministe $\phi_{m_1}(Y_1, t - T_1)$. La nouvelle position après le saut $T_2$ est $X_{T_2}$ et ainsi de suite. L'exemple d'un PDMP représenté jusqu'à son deuxième temps de saut est illustré dans la **Fig.1.8**.

<figure><img src="image_11.png" /><figcaption>FIGURE 1.8 – Exemple d'une trajectoire d'un PDMP.</figcaption></figure>

Dans le chapitre **2.3**, nous utiliserons ce type de modèle afin de proposer une modélisation stochastique de la dynamique de transmission de la fièvre typhoïde à Mayotte.

### 1.4.3 Processus de naissance et mort linéaire avec immigration

Le chapitre **3** de ce manuscrit se concentre sur l'estimation des paramètres d'un processus de Markov de saut pur appartenant à une famille particulière. Cette famille représente un cas particulier de PDMP et une généralisation des processus de comptage, dans lesquels la population peut augmenter ou diminuer d'une unité à chaque saut. Dans cette
partie, nous présentons les notions de base pour caractériser ce processus et décrire ses
propriétés mathématiques.

Un processus de naissance et de mort linéaire avec immigration (en anglais : linear Birth-Death Process with Immigration (LBDPI)), noté $(X_t)_{t \ge 0}$, est un processus de Markov à temps continu à valeurs dans $\mathbb{N}$, comptant le nombre d'individus dans une population à l'instant $t$, son générateur infinitésimal $Q$ est donné par :

$$b_i = Q_{i,i+1} = \lambda i + \nu \qquad \text{pour } i \ge 0,$$

$$d_i = Q_{i,i-1} = \mu i \qquad \text{pour } i \ge 1,$$

$$Q_{i,i} = -(\lambda + \mu)i - \nu \qquad \text{pour } i \ge 0,$$

où $\lambda, \mu$ et $\nu$ sont positifs.

$(X_t)_{t \ge 0}$ est un processus qui décrit une population qui reçoit des immigrants au taux $\nu$ et dans laquelle chaque individu donne naissance à un autre individu au taux $\lambda$ et meurt au taux $\mu$, indépendamment des autres.

Les fonctions $b_i$ et $d_i$ sont respectivement appelées taux de naissance (birth) et taux de mort (death) du processus $(X_t)_{t \ge 0}$ en présence de $i$ individus. Le paramètre $\lambda$ s'interprète comme le taux de naissance individuel, et $\mu$ est comme le taux de mort (retrait) individuel.

La dynamique du processus $X_t$ est défini comme suit : une fois à l'état $i \ge 0$, on y reste un temps de loi $\text{Exp}(b_i + d_i)$ avant d'aller soit à l'état $i + 1$ avec probabilité $b_i/(b_i + d_i)$, soit à l'état $i - 1$ avec probabilité $d_i/(b_i + d_i)$ (voir Fig.1.9 et Fig.1.10).

<figure><img src="image_12.png" /><figcaption>FIGURE 1.9 – Processus de naissance et de mort.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>FIGURE 1.10 – Diagramme de Markov du processus $(X_t)_{t \ge 0}$ d'intensités de naissance $b_i$ et de mort $d_i$, si $b_i > 0 \ \forall i \ge 0$ et $d_i > 0 \ \forall i \ge 1$.</figcaption></figure>
Le semi-groupe de transition défini par $(P(t) = (p_{i,j}(t), i, j \in \mathbb{N}))_{t \ge 0}$ où

$$p_{i,j}(t) = \mathbb{P}(X_t = j | X_0 = i), i, j \in \mathbb{N}, t \ge 0,$$

satisfait l'équation de Kolmogorov progressive suivante :

$$\frac{dP(t)}{dt} = P(t)Q, \qquad (1.2)$$

pour $t > 0$, qui peut être réécrite pour $i, j \in \mathbb{N}$ avec $j \ge 1$ comme suit

$$\frac{dp_{i,0}(t)}{dt} = \mu p_{i,1}(t) - \nu p_{i,0}(t),$$

$$\frac{dp_{i,j}(t)}{dt} = (\lambda(j-1) + \nu) p_{i,j-1}(t) + \mu(j+1)p_{i,j+1}(t) - ((\lambda + \mu)j + \nu) p_{i,j}(t),$$

et avec la condition initiale $p_{i,i}(0) = 1$ et $p_{i,j}(0) = 0$ pour $j \ne i$. L'équation de Kolmogorov correspond à la relation infinitésimale : pour $h > 0$,

$$\mathbb{P}(X_{t+h} = j | X_t = i) = Q_{i,j}h + o(h) \quad \text{si } j \neq i,$$

$$\mathbb{P}(X_{t+h} = i | X_t = i) = 1 + Q_{i,i}h + o(h),$$

uniformément en $t$. Les trajectoires de $(X_t)_{t \ge 0}$ sont constantes par morceaux et ponctuées de sauts d'amplitude +1 ou -1.

Par la suite, nous présentons certaines propriétés théoriques bien connues du processus $X_t$ (voir par exemple [18]), qui permettent de décrire son comportement moyen et asymptotique. Nous complétons ces résultats par des simulations numériques qui viennent illustrer ces propriétés.

**Théorème 1.1.** (Comportement moyen de $X_t$). Pour tout $i \ge 0$ et $t \ge 0$,

$$\mathbb{E}_i[X_t] = \begin{cases} \nu t + i & \text{si } \lambda = \mu, \\ ie^{(\lambda-\mu)t} + \frac{\nu}{\lambda-\mu} (e^{(\lambda-\mu)t} - 1) & \text{si } \lambda \neq \mu. \end{cases}$$

*Preuve.* En utilisant l'équation de Kolmogorov progressive **1.2**, on obtient

$$\begin{aligned} \frac{d}{dt} \mathbb{E}_i[X_t] &= \sum_{j \ge 0} j p'_{i,j}(t) \\ &= \sum_{j \ge 0} j \sum_{k \ge 0} p_{i,k}(t) Q_{k,j} \\ &= \sum_{k \ge 0} (\lambda - \mu) k p_{i,k}(t) + \nu. \end{aligned}$$

On a alors

$$\frac{d}{dt} \mathbb{E}_i[X_t] = (\lambda - \mu) \mathbb{E}_i[X_t] + \nu,$$

avec la condition initiale $\mathbb{E}_i[X_0] = i$.

On résout cette équation différentielle par la méthode de variation de la constante, on trouve

$$\mathbb{E}_i[X_t] = \begin{cases} \nu t + i & \text{si } \lambda = \mu, \\ ie^{(\lambda-\mu)t} + \frac{\nu}{\lambda-\mu} (e^{(\lambda-\mu)t} - 1) & \text{si } \lambda \neq \mu. \end{cases}$$

□
**Corollaire 1.1.** Asymptotiquement et conditionnellement à $X_0 = i$, si $\lambda \ge \mu$, le processus moyen tend vers l'infini ($\lim_{t \to \infty} \mathbb{E}_i[X_t] = \infty$) de manière linéaire lorsque $\lambda = \mu$ (Fig.1.12), et de manière exponentielle lorsque $\lambda > \mu$ (Fig.1.11). En revanche, si $\lambda < \mu$, l'espérance de $X_t$ se stabilise et dans ce cas $\lim_{t \to \infty} \mathbb{E}_i[X_t] = \nu/(\mu-\lambda)$ (Fig.1.13).

**Théorème 1.2.** (Comportement asymptotique de $X_t$). Le processus $(X_t)_{t \ge 0}$ admet une unique probabilité invariante si et seulement $\mu > \lambda$.

De plus, si on note $\pi = (\pi_k, k \in \mathbb{N})$ cette probabilité invariante, alors pour tout $k \ge 0$ :

$$ \pi_k = \left( \frac{\nu}{\lambda} + k - 1 \middle| k \right) \left( \frac{\lambda}{\mu} \right)^k \left( 1 - \frac{\lambda}{\mu} \right)^{\nu/\lambda}, \quad (1.3) $$

où

$$ \binom{\frac{\nu}{\lambda} + k - 1}{k} = \frac{(\frac{\nu}{\lambda} + k - 1)(\frac{\nu}{\lambda} + k - 2) \cdots \frac{\nu}{\lambda}}{k!}. $$

Preuve. $(X_t)_{t \ge 0}$ est un processus de naissance et mort linéaire, tel que :

$$ \forall i \neq 0, \quad b_i \neq 0, d_i \neq 0 \quad \text{et} \quad b_0 > 0. $$

Donc, il est irréductible.

• Existence de la loi invariante :

D'après [19, Théorème 1.4], le processus $(X_t)_{t \ge 0}$ admet une loi invariante si et seulement si :

$$ C = \sum_{k \ge 1} \frac{b_0 \dots b_{k-1}}{d_1 \dots d_k} = \sum_{k \ge 1} \frac{\nu(\lambda + \nu) \dots ((k-1)\lambda + \nu)}{\mu 2\mu \dots k\mu} < +\infty. $$

Ainsi, $(X_t)_{t \ge 0}$ admet une loi invariante $(\pi_k)_{k \in \mathbb{N}}$ si et seulement si $\mu > \lambda$.

• Expression de $(\pi_k)_{k \in \mathbb{N}}$ :

D'après [19, Théorème 1.4], l'unique loi invariante est donnée par :

$$ \pi_0 = \frac{1}{1+C}, \quad \pi_k = \frac{b_0 \dots b_{k-1}}{d_1 \dots d_k} \pi_0 \quad \forall k \ge 1. $$

Nous commençons par calculer $\pi_0$. On a :

$$ C = \sum_{k \ge 1} \frac{(\frac{\nu}{\lambda} + k - 1)(\frac{\nu}{\lambda} + k - 2) \cdots \frac{\nu}{\lambda}}{k!} \frac{\lambda^k}{\mu^k}, $$

donc :

$$ 1 + C = \sum_{k \ge 0} \binom{\frac{\nu}{\lambda} + k - 1}{k} \left(\frac{\lambda}{\mu}\right)^k. $$

Or, on sait que : pour tout $x \in [0, 1[$ et pour tout $n \in \mathbb{R}_+$, on a :

$$ \frac{1}{(1-x)^n} = \sum_{k=0}^{\infty} \binom{n+k-1}{k} x^k. $$

Ainsi :

$$ \pi_0 = \left(1 - \frac{\lambda}{\mu}\right)^{\nu/\lambda}. $$
Nous déduisons maintenant l'expression de $\pi_k$ pour tout $k \ge 1$,

$$ \begin{aligned} \pi_k &= \frac{\nu(\lambda + \nu) \dots ((k-1)\lambda + \nu)}{\mu 2\mu \dots k\mu} \left(1 - \frac{\lambda}{\mu}\right)^{\nu/\lambda} \\ &= \frac{(\frac{\nu}{\lambda} + k - 1)(\frac{\nu}{\lambda} + k - 2) \dots \frac{\nu}{\lambda}}{k!} \frac{\lambda^k}{\mu^k} \left(1 - \frac{\lambda}{\mu}\right)^{\nu/\lambda}. \end{aligned} $$

Par conséquence, pour tout $k \ge 1$ :

$$ \pi_k = \left( \frac{\nu}{\lambda} + k - 1 \middle| k \right) \left( \frac{\lambda}{\mu} \right)^k \left( 1 - \frac{\lambda}{\mu} \right)^{\nu/\lambda}. $$

□

**Théorème 1.3.** Le LBDPI $(X_t)_{t \ge 0}$ (avec $\nu > 0$) a trois régimes possibles de comportement à long terme : croissance exponentielle (transitoire) si $\lambda > \mu$ (Fig.1.11). Tous les états sont visités infiniment souvent avec un temps de retour moyen infini (récurrence nulle) si $\lambda = \mu$ (Fig.1.12). Tous les états sont visités infiniment souvent avec un temps de retour moyen fini (récurrence positive) si $\lambda < \mu$ (Fig.1.13). Dans ce dernier cas, il existe une unique loi invariante donnée par 1.3.

<figure><img src="image_14.png" /><figcaption>FIGURE 1.11 – Exemple des trajectoires de $X_t$ et son espérance $\mathbb{E}[X_t]$ pour $\lambda = 0.1, \mu = 0.08$ et $\nu = 0.02$ (régime transitoire).</figcaption></figure>
<figure><img src="image_15.png" /><figcaption>FIGURE 1.12 – Exemple des trajectoires de $X_t$ et son espérance $E[X_t]$ pour $\lambda = \mu = 0.1$ et $\nu = 0.02$ (régime récurrent nul).</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>FIGURE 1.13 – Exemple des trajectoire de $X_t$ et son espérance $E[X_t]$ pour $\lambda = 0.1, \mu = 0.15$ et $\nu = 0.02$ (régime récurrent positif).</figcaption></figure>
### 1.4.4 Modèle de Markov caché

Cette section introductive de la thèse est dédiée aux problèmes d'information cachée, qui constituent une grande partie des problèmes abordés dans ce travail. Nous y présentons les principes fondamentaux, les notations clés et les algorithmes principaux qui seront développés et adaptés pour notre contexte.

Dans la pratique, de nombreux processus suivent un modèle de Markov, mais ne sont pas observés directement. Dans ce cas, on parle de chaînes de Markov cachées, où seule une partie de l'état du processus est observée. L'objectif est alors d'obtenir le maximum d'informations sur la partie cachée de la chaîne à partir des observations partielles disponibles.

Les modèles de Markov cachés (en anglais : Hidden Markov Model (HMM)) sont large- ment utilisés pour modéliser des données séquentielles dans de nombreux domaines, tels que la reconnaissance de la parole [20], la reconnaissance d'écriture manuscrite [21], la bioinformatique [22], la finance et la surveillance [23]. Les HMM sont des modèles sta- tistiques qui permettent de représenter des données observées qui sont la résultante de processus stochastiques non observables, également appelés processus cachés. L'objectif de cette modélisation est de découvrir les états cachés à partir des données observées, et d'inférer les paramètres du modèle à partir des données observées. Les HMM sont par- ticulièrement utiles pour les problèmes de classification et de prédiction de séquences, et constituent un outil puissant pour l'analyse de données séquentielles. Nous allons présen- ter les principes et les notions de base des HMM.

#### Caractéristiques des HMM à temps discret

On considère une chaîne de Markov homogène $X := X_{0:T} = (X_0, X_1, \dots, X_T)$ à valeurs dans un espace d'états fini ou dénombrable $\mathcal{S}$. La chaîne $X$ n'est pas directement observable, elle est dite cachée. Considérons une autre suite des observations $O := O_{0:T} = (O_0, O_1, \dots, O_T)$ prenant des valeurs dans un ensemble $\mathcal{O}$ appelé l'espace d'observation. Cet espace peut être au plus dénombrable. Nous supposons que les variables aléatoires $O_{0:T}$ sont indépendantes conditionnellement à la suite d'états cachée $X_{0:T}$ et que chaque observation $O_k$ ne dépend que de l'état $X_k$ (l'hypothèse du canal d'observation sans mémoire : HMM standard). C'est-à-dire :

$$ \mathbb{P}(O_{0:T} = o_{0:T} | X_{0:T} = x_{0:T}) = \prod_{t=0}^{T} \mathbb{P}(O_t = o_t | X_t = x_t), \quad x_{0:T} \in \mathcal{S}^{T+1}, o_{0:T} \in \mathcal{O}^{T+1}. $$

Un HMM $(X, O)$ à espace d'observation dénombrable est caractérisé par :

* Une matrice de transition $P$ du processus d'état caché $X$ :

$$ p_{x,x'} = \mathbb{P}(X_{t+1} = x'|X_t = x), \quad x, x' \in \mathcal{S}. $$

* Une matrice d'émission $\psi$ définie par :

$$ \psi_x(o) = \mathbb{P}(O_t = o | X_t = x), \quad x \in \mathcal{S}, o \in \mathcal{O}. $$

* Une loi initiale $\rho$ du processus $X$ :

$$ \rho_x = \mathbb{P}(X_0 = x), \quad x \in \mathcal{S}. $$
On note par le triplet $M = (P, \psi, \rho)$ les paramètres du HMM.

La probabilité conjointe des états cachés et des observations se factorise comme suit en utilisant les propriétés de Markov et d'indépendance conditionnelle du HMM :

$$ \mathbb{P}(X_{0:T} = x_{0:T}, O_{0:T} = o_{0:T}) = \rho_{x_0} \psi_{x_0}(o_0) \prod_{t=1}^{T} \psi_{x_t}(o_t) p_{x_{t-1},x_t}. $$

La figure 1.14 illustre une représentation graphique du HMM standard.

<figure><img src="image_17.png" /><figcaption>FIGURE 1.14 – Schéma d'émission du HMM standard.</figcaption></figure>

## Algorithme Forward-Backward de Baum [20, 24]

Dans les problèmes du HMM, il est essentiel de pouvoir évaluer efficacement la probabilité $\mathbb{P}_M(O_{0:T})$ d'observer la suite d'observations $O_{0:T}$, étant donné le modèle $M$ (la vraisemblance des observations).

Il est important de remarquer que cette probabilité peut être exprimée sous la forme suivante :

$$ \begin{aligned} \mathbb{P}_M (O_{0:T} = o_{0:T}) &= \sum_{x_{0:T} \in \mathcal{S}^{T+1}} \mathbb{P}_M (X_{0:T} = x_{0:T}, O_{0:T} = o_{0:T}) \\ &= \sum_{x_{0:T} \in \mathcal{S}^{T+1}} \rho_{x_0} \psi_{x_0}(o_0) \prod_{t=1}^{T} \psi_{x_t}(o_t) p_{x_{t-1},x_t}, \end{aligned} \quad (1.4) $$

pour tout $o_{0:T} \in \mathcal{O}^{T+1}$.

Il convient de noter que le calcul direct de la probabilité $\mathbb{P}_M(O_{0:T})$ est possible, mais requiert des calculs d'ordre exponentiel ($T \times N^T$, où $N$ représente le nombre d'états cachés dans le modèle), ce qui est très coûteux en temps de calcul. C'est pourquoi, une approche récursive (Algorithme Forward-Backward) plus astucieuse et plus rapide a été développée pour réduire le temps de calcul nécessaire à l'évaluation de cette probabilité.

### • Algorithme Forward :

On s'intéresse à évaluer la probabilité conjointe des $t$ premières observations $O_{0:t}$ et $X_t = x \in \mathcal{S}$, étant donné le modèle $M$ :

$$ \alpha_x(t) = \mathbb{P}_M(O_{0:t} = o_{0:t}, X_t = x), \quad \text{pour } t \ge 0. $$

Cette probabilité satisfait une formule récursive donnée par le théorème suivant.
**Théorème 1.4.** Pour tout $1 \le t \le T$, $x \in S$,

$$ \left\{ \begin{array}{l} \alpha_x(0) = \rho_x \psi_x(o_0), \\ \alpha_x(t) = \psi_x(o_t) \sum_{x' \in S} p_{x',x} \alpha_{x'}(t-1). \end{array} \right. \qquad (1.5) $$

*Preuve.* Soit $x \in S$, d'abord, pour $t=0$, on a :

$$ \alpha_x(0) = \mathbb{P}_M(O_0 = o_0, X_0 = x) = \rho_x \psi_x(o_0). $$

Ensuite, on peut remarquer que : pour tout $t \ge 2$,

$$ \begin{aligned} \alpha_x(t-1) &= \sum_{x_0:t-2 \in S^{t-1}} \mathbb{P}_M(O_{0:t} = o_{0:t}, X_{0:t-2} = x_{0:t-2}, X_{t-1} = x) \\ &= \sum_{x_0:t-2 \in S^{t-1}} \rho_{x_0} \psi_{x_0}(o_0) \prod_{l=1}^{t-1} \psi_{x_l}(o_l) p_{x_{l-1},x_l}. \end{aligned} $$

Ainsi,

$$ \begin{aligned} \alpha_x(t) &= \sum_{x_0:t-1 \in S^t} \rho_{x_0} \psi_{x_0}(o_0) \prod_{l=1}^{t} \psi_{x_l}(o_l) p_{x_{l-1},x_l} \\ &= \psi_{x_t}(o_t) \sum_{x_{t-1} \in S} p_{x_{t-1},x_t} \sum_{x_0:t-2 \in S^{t-1}} \rho_{x_0} \psi_{x_0}(o_0) \prod_{l=1}^{t-1} \psi_{x_l}(o_l) p_{x_{l-1},x_l} \\ &= \psi_{x_t}(o_t) \sum_{x_{t-1} \in S} p_{x_{t-1},x_t} \alpha_x(t-1). \end{aligned} $$

Le résultat souhaité est donc obtenu. □

**Remarque 1.1.** La fonction de vraisemblance $\mathbb{P}_M(O_{0:T})$ peut être calculée à partir de l'algorithme Forward en utilisant la formule suivante :

$$ \mathbb{P}_M(O_{0:T}) = \sum_{x \in S} \alpha_x(T). $$

• **Algorithme Backward :**

Maintenait, nous définissons la probabilité Backward présentant la probabilité d'observer la suite partielle ultérieure $O_{t+1:T}$, étant donné le modèle $M$ et $X_t = x \in S$ :

$$ \beta_x(t) = \mathbb{P}_M(O_{t+1:T} = o_{t+1:T} | X_t = x), \quad \text{pour } 0 \le t < T. $$

Contrairement à la méthode Forward, la valeur de $\beta_x(T)$ est arbitrairement choisie :

$$ \beta_x(T) = 1. $$

De même, nous pouvons calculer la probabilité Backward à partir d'une formule de récurrence, présentée dans le théorème suivant.

**Théorème 1.5.** Pour tout $0 \le t < T$, $x \in S$,

$$ \left\{ \begin{array}{l} \beta_x(T) = 1, \\ \beta_x(t) = \sum_{x' \in S} \psi_{x'}(o_{t+1}) p_{x,x'} \beta_{x'}(t+1). \end{array} \right. \qquad (1.6) $$
Preuve. Soit $x \in S$, nous commençons par $t = T - 1$, on a

$$
\begin{align*}
\beta_x(T-1) &= \mathbb{P}_M(O_T = o_T \mid X_{T-1} = x) \\
&= \sum_{x_T \in S} \mathbb{P}_M(O_T = o_T, X_T = x_T \mid X_{T-1} = x) \\
&= \sum_{x_T \in S} \mathbb{P}_M(O_T = o_T \mid X_T = x_T) \mathbb{P}_M(X_T = x_T \mid X_{T-1} = x) \\
&= \sum_{x_T \in S} \psi_{x_T}(o_T) p_{x,T}.
\end{align*}
$$

Ensuite, pour tout $0 \le t < T - 1$,

$$
\begin{align*}
\beta_x(t) &= \mathbb{P}_M(O_{t+1:T} = o_{t+1:T} \mid X_t = x) \\
&= \sum_{x_{t+1:T} \in S^{T-t}} \mathbb{P}_M(O_{t+1:T} = o_{t+1:T}, X_{t+1:T} = x_{t+1:T} \mid X_t = x) \\
&= \sum_{x_{t+1:T} \in S^{T-t}} \mathbb{P}_M(O_{t+1:T} = o_{t+1:T} \mid X_{t+1:T} = x_{t+1:T}) \mathbb{P}_M(X_{t+1:T} = x_{t+1:T} \mid X_t = x) \\
&= \sum_{x_{t+1:T} \in S^{T-t}} \prod_{l=t+1}^{T} \psi_{x_l}(o_l) p_{x,T} \prod_{l=t+2}^{T} p_{x_{l-1},T} \\
&= \sum_{x_{t+1:T} \in S} \psi_{x_{t+1}}(o_{t+1}) p_{x,T} \sum_{x_{t+2:T} \in S^{T-t-1}} \prod_{l=t+2}^{T} \psi_{x_l}(o_l) p_{x_{l-1},T} \\
&= \sum_{x_{t+1:T} \in S} \psi_{x_{t+1}}(o_{t+1}) p_{x,T} \sum_{x_{t+2:T} \in S^{T-t-1}} \mathbb{P}_M(O_{t+2:T} = o_{t+2:T} \mid X_{t+2:T} = x_{t+2:T}) \\
&\qquad \mathbb{P}_M(X_{t+2:T} = x_{t+2:T} \mid X_{t+1} = x_{t+1}) \\
&= \sum_{x_{t+1:T} \in S} \psi_{x_{t+1}}(o_{t+1}) p_{x,T} \mathbb{P}_M(O_{t+2:T} = o_{t+2:T} \mid X_{t+1} = x_{t+1}) \\
&= \sum_{x_{t+1:T} \in S} \psi_{x_{t+1}}(o_{t+1}) p_{x,T} \beta_{x_{t+1}}(t+1)
\end{align*}
$$

Ainsi, la preuve est complète.

<figure><img src="image_18.png" /></figure>

**Proposition 1.1.** La fonction de vraisemblance des observations $\mathbb{P}_M(O_{0:T})$ est donnée par :

$$
\mathbb{P}_M(O_{0:T}) = \sum_{x \in S} \alpha_x(t) \beta_x(t),
$$

où $\alpha_x(t)$ et $\beta_x(t)$ sont les probabilités Forward et Backward calculées à partir de **1.5** et **1.6**.

Preuve. On a :

$$
\mathbb{P}_M(O_{0:T}) = \sum_{x \in S} \mathbb{P}_M(O_{0:T}, X_t = x),
$$

où, pour tout $x \in S$,

$$
\begin{align*}
\mathbb{P}_M(O_{0:T}, X_t = x) &= \sum_{x_{0:t-1} \in S^t} \sum_{x_{t+1:T} \in S^{T-t}} \mathbb{P}_M(O_{0:T}, X_{0:T} = x_{0:T}) \\
&= \sum_{x_{0:t-1} \in S^t} \sum_{x_{t+1:T} \in S^{T-t}} \rho_{x_0} \psi_{x_0}(o_0) \prod_{t=1}^{T} \psi_{x_t}(o_t) p_{x_{t-1},x_t} \\
&= \sum_{x_{0:t-1} \in S^t} \rho_{x_0} \psi_{x_0}(o_0) \prod_{l=1}^{t} \psi_{x_l}(o_l) p_{x_{l-1},x_l} \sum_{x_{t+1:T} \in S^{T-t}} \prod_{l=t+1}^{T} \psi_{x_l}(o_l) p_{x_{l-1},x_l} \\
&= \alpha_x(t) \beta_x(t). \tag{1.7}
\end{align*}
$$
Par conséquence,

$$ \mathbb{P}_M(O_{0:T}) = \sum_{x \in S} \alpha_x(t)\beta_x(t). $$

□

**Remarque 1.2.** Le schéma Forward-Backward permet de calculer l'évaluation de la probabilité $\mathbb{P}_M(O_{0:T})$ de manière plus efficace que la méthode naïve donnée par **1.4**, en réduisant l'ordre de calcul de $T \times N^T$ à $T \times N^2$. Cette amélioration est obtenue grâce à la combinaison des probabilités Forward et Backward (Pour $N = 4$ et $T = 100$, l'évaluation de $\mathbb{P}(O_{0:100})$ grâce à la procédure forward-backward passe de l'ordre de calcul $T \times N^T = 100 \times 4^{100}$ à $T \times N^2 = 100 \times 4^2 = 1600$. Cela représente une réduction massive du temps de calcul nécessaire pour évaluer la vraisemblance de la séquence d'observations pour ce modèle de HMM).

## Ré-estimation des paramètres du HMM

Une fois que nous avons déterminé la vraisemblance de notre modèle HMM pour une séquence d'observations donnée, nous pouvons ré-estimer les paramètres du modèle pour améliorer sa performance. Cette étape est appelée la ré-estimation des paramètres du HMM et se base sur l'algorithme de Baum-Welch [25].

L'idée générale de l'algorithme de Baum-Welch est d'ajuster les paramètres du modèle de telle sorte que la probabilité d'observer la séquence d'observations soit maximisée. Il s'agit d'un processus itératif qui commence par une initialisation des paramètres, puis ré-estime à chaque itération les paramètres en utilisant l'algorithme Expectation-Maximization (EM).

La première étape de l'algorithme EM consiste à calculer une fonction objectif appelée fonction de vraisemblance attendue définie par :

$$ Q(M, M^{(n)}) = \sum_{x_{0:T}} P(X_{0:T} = x_{0:T} | O_{0:T}, M^{(n)}) \log P(O_{0:T}, X_{0:T} = x_{0:T} | M), \quad (1.8) $$

où $M$ représente les paramètres du HMM, $M^{(n)}$ est une estimation courante de ces paramètres.

La deuxième étape consiste à maximiser cette fonction de vraisemblance attendue par rapport aux paramètres du HMM $M$, ce qui donne la nouvelle estimation $M^{(n+1)}$ :

$$ M^{(n+1)} = \arg \max_{M} Q(M, M^{(n)}). $$

Ces deux étapes sont répétées jusqu'à convergence, c'est-à-dire lorsque la variation entre les estimations successives est suffisamment petite. Ainsi, l'algorithme EM permet de trouver les paramètres optimaux du HMM qui maximisent la vraisemblance des observations.

Pour calculer explicitement les ré-estimations de Baum-Welsh, nous définissons les deux probabilités suivantes :

$$ \xi_{x,x'}(t) = \mathbb{P}_M(X_{t+1} = x', X_t = x | O_{0:T} = o_{0:T}), $$

et

$$ \gamma_x(t) = \sum_{x' \in S} \xi_{x,x'}(t) = \mathbb{P}_M(X_t = x | O_{0:T} = o_{0:T}). $$
En utilisant les probabilités calculées par l'algorithme Forward-Backward $\alpha_x(t)$ et $\beta_x(t)$, Ces deux probabilités peuvent être exprimées de la manière suivante :

$$ 
\begin{aligned}
\xi_{x,x'}(t) &= \frac{\mathbb{P}_M(X_{t+1} = x', X_t = x, O_{0:T} = o_{0:T})}{\mathbb{P}_M(O_{0:T} = o_{0:T})} \\
&= \frac{\alpha_x(t)p_{x,x'}\psi_{x'}(o_{t+1})\beta_{x'}(t+1)}{\sum_{x,x' \in S} \alpha_x(t)p_{x,x'}\psi_{x'}(o_{t+1})\beta_{x'}(t+1)}.
\end{aligned}
$$

Pour la deuxième probabilité, d'après **1.7**, on obtient :

$$ 
\begin{aligned}
\gamma_x(t) &= \frac{\mathbb{P}_M(X_t = x, O_{0:T} = o_{0:T})}{\mathbb{P}_M(O_{0:T} = o_{0:T})} \\
&= \frac{\alpha_x(t)\beta_x(t)}{\sum_{x \in S} \alpha_x(t)\beta_x(t)}.
\end{aligned}
$$

En sommant les coefficients $\xi_{x,x'}(t)$ et $\gamma_x(t)$ sur les instants $t = 0, \dots, T-1$, on obtient des estimateurs du nombre de transitions de l'état $x$ à l'état $x'$ et du nombre de transitions de l'état $x$ dans la suite $X_{0:T}$, respectivement. Ces estimateurs sont basés sur les observations de la suite $O_{0:T} = o_{0:T}$.

Par conséquence, les ré-estimateurs de Baum-Welch qui maximise la fonction **1.8**, étant donné le modèle $M^{(n)}$ sont donnés par

$$ 
\begin{aligned}
p_{x,x'}^{(n+1)} &= \frac{\sum_{t=0}^{T-1} \xi_{x,x'}^{(n)}(t)}{\sum_{t=0}^{T-1} \gamma_x^{(n)}(t)}, \\
\psi_x^{(n+1)}(o) &= \frac{\sum_{t=0}^{T} 1_{o_t=o} \gamma_x^{(n)}(t)}{\sum_{t=0}^{T} \gamma_x^{(n)}(t)}, \\
\rho_x^{(n+1)} &= \gamma_x^{(n)}(0).
\end{aligned}
$$

L'algorithme de Baum-Welch se termine lorsque la vraisemblance des observations converge vers un maximum local, ce qui indique que les paramètres ont été ajustés de manière optimale pour la suite d'observations donnée. La qualité de la solution obtenue dépend de la qualité de l'initialisation des paramètres et de la qualité de la suite d'observations utilisée pour l'entraînement.

Cet algorithme sera étendu et appliqué à différents cas d'étude dans les chapitres 3 et 4.

## 1.5 État de l'art sur la modélisation mathématique de la fièvre typhoïde

La modélisation de la transmission de la fièvre typhoïde est un sujet de recherche important, abordé par de nombreux chercheurs et mathématiciens tels que Lauria et al.[26], Kalajdzievska et al. [27], Mushayabasa et al.[28], Abboubakar et al. [29] et Pitzer et al.[30]. Cependant, malgré les nombreux travaux existants, l'influence de la qualité de l'eau sur la transmission de la maladie n'a pas été suffisamment considérée [28, 29, 30] et aucun modèle ne prend en compte deux populations en interaction avec accès ou non à l'eau potable. De plus, les modèles classiques utilisés sont souvent des modèles compartimentaux déterministes (1.4.1) du type S, I, R et C (le compartiment des porteurs
chroniques), qui supposent que toutes les variables d'entrée sont des fonctions déterministes du temps. Ces modèles ne permettent pas de prendre en compte les variations aléatoires de certains facteurs qui peuvent influencer la dynamique de transmission de la maladie.

Plusieurs travaux ont été réalisés dans ce cadre. Par exemple, dans [27], les auteurs ont développé un modèle mathématique pour évaluer les effets des porteurs chronique (porteurs de germes, demeurant asymptomatique, pendant une très longue période, parfois jusqu'à la mort) sur la dynamique de transmission des maladies infectieuses telles que la typhoïde. Ils ont conclu que les porteurs jouent un rôle important dans la transmission des maladies infectieuses. **Mushayabasa et al.** [28] ont étudié un modèle épidémiologique pour la transmission directe (personne à personne) et indirecte (par l'eau) de la fièvre typhoïde. Ils ont montré à partir d'une analyse de sensibilité que la transmission indirecte de la typhoïde a plus d'impact sur la détermination de la prévalence de la maladie que la transmission directe. Ils ont également souligné l'importance de prendre en compte la qualité de l'eau dans les modèles de transmission de la maladie. **Abboubakar et al.** [29] présentent un modèle mathématique déterministe qui prend en compte l'effet de l'environnement sur la transmission de la fièvre typhoïde, et une stratégie de contrôle optimale pour réduire la propagation de la maladie. Les auteurs ont utilisé des données épidémiologiques pour paramétrer le modèle et une analyse de sensibilité pour identifier les paramètres les plus importants pour la transmission de la maladie. Ils ont également établi un problème de contrôle optimal pour étudier l'impact des différentes stratégies de vaccination sur la transmission de la maladie et ont conclu que la vaccination généralisée était plus efficace pour réduire la prévalence de la maladie que la vaccination ciblée sur des groupes à risque élevé.

**Pitzer et al.** [30] ont pris en compte l'impact de la qualité de l'eau sur la propagation de la typhoïde en élaborant un modèle parcimonieux structuré par âge qui décrit l'histoire naturelle et l'immunité à l'infection typhoïde. Ils ont montré que l'amélioration de la qualité de l'eau peut avoir un impact significatif sur la prévention de la transmission de la typhoïde et ont évalué l'efficacité de différentes stratégies de vaccination en milieu scolaire utilisant des vaccins vivants oraux, des vaccins polysaccharidiques et des vaccins conjugués.

Aucun de ces modèles présentés ci-dessus ne peut décrire de manière adéquate la propagation de la typhoïde dans le contexte particulier de Mayotte, vu qu'ils ne considèrent pas les spécificités de Mayotte, comme l'accès aux BFM limité à certaines zones géographiques ainsi que de la pression migratoire que subit l'île. Cependant, ces modèles restent pertinents pour servir de références permettant d'identifier les ordres de grandeur de chaque paramètre de notre modèle décrit dans la section 2.2.

## 1.6 Organisation de la thèse

Dans cette section, nous allons présenter l'organisation des chapitres de cette thèse.

Le chapitre 2 propose un modèle initial à compartiments pour caractériser et décrire la dynamique de transmission de la fièvre typhoïde à Mayotte, en prenant en compte le contexte spécifique de cette maladie hydrique ainsi que les particularités de la population mahoraise. Le modèle est présenté dans ses versions déterministe et stochastique (PDMP), et est exploré tant sur le plan théorique que numérique, avec des simulations
de trajectoires dans différents scénarios, ainsi que l’analyse de leur comportement qualitatif. En outre, nous évaluons la virulence de la maladie en utilisant le paramètre clé R₀ représentant le nombre de reproduction de base.

Le jeu de données fourni par l'ARS de Mayotte est restreint et ne permet pas une calibration précise de notre modèle initial, ni une quantification précise de la propagation épidémique à Mayotte. Dans le chapitre 3, nous proposons donc un modèle stochastique simplifié à un compartiment de type naissance et mort avec immigration 1.4.3, qui permet de décrire les effectifs de la population infectée. Nous utilisons cette approche afin d'obtenir un modèle ajusté aux données disponibles et d'estimer les paramètres clés de la dynamique de la fièvre typhoïde à Mayotte en adaptant un modèle de Markov caché 1.4.4 à notre contexte particulier. La procédure d'estimation et les principales propriétés des estimateurs sont établies théoriquement et illustrées numériquement sur des données synthétiques ainsi que sur des données réelles provenant de Mayotte.

Dans les chapitres 4 et 5, nous enrichissons notre modèle décrit dans le chapitre 3 en proposant deux extensions. Tout d'abord, dans le chapitre 4, nous introduisons un nouveau compartiment dans le modèle, qui correspond aux individus exposés à la bactérie Salmonella Typhi. Cette nouvelle extension nous permet de prendre en compte la durée d'incubation de la maladie. Les propriétés théoriques de ce nouveau modèle, appelé modèle infecté-exposé, sont étudiées et illustrées par des simulations numériques. Nous présentons également une nouvelle méthodologie d'estimation des paramètres du modèle bivarié en utilisant un modèle de Markov multi-chaîne cachée. Les résultats obtenus sont illustrés par des études numériques.

La deuxième extension de notre modèle présentée dans le chapitre 5 vise à étudier l'impact de la pluviométrie sur la transmission de la maladie à Mayotte, en exploitant des données publiques de précipitations. Pour cela, nous avons utilisé une méthode de détection de rupture afin de sélectionner et caractériser la saisonnalité de la pluie. Les résultats d'estimation des paramètres de notre modèle sous les deux régimes mettent en évidence l'importance de cette variable météorologique sur la propagation de l'épidémie, soulignant la nécessité de prendre en compte cet effet dans la planification des interventions de santé publique pour lutter contre la fièvre typhoïde à Mayotte.

Le chapitre 6 conclut ce manuscrit en explorant les pistes de recherche futures.
# Chapitre 2

## Modèle à compartiments pour la transmission de la fièvre typhoïde à Mayotte

### Sommaire

<table><tr><td>2.1</td><td>Introduction</td><td>33</td></tr><tr><td>2.2</td><td>Modélisation déterministe</td><td>33</td></tr><tr><td>2.2.1</td><td>Formulation du modèle</td><td>34</td></tr><tr><td>2.2.2</td><td>Calibration des paramètres</td><td>36</td></tr><tr><td>2.2.3</td><td>Discussion et simulations numériques</td><td>38</td></tr><tr><td>2.2.4</td><td>Le nombre de reproduction de base</td><td>45</td></tr><tr><td>2.3</td><td>Modélisation par un processus markovien déterministe par morceaux</td><td>48</td></tr><tr><td>2.3.1</td><td>Formulation du modèle</td><td>49</td></tr><tr><td>2.3.2</td><td>Discussion et simulations numériques</td><td>50</td></tr><tr><td>2.4</td><td>Conclusion</td><td>54</td></tr></table>
## 2.1 Introduction

La fièvre typhoïde est une maladie infectieuse d'origine bactérienne, causée par la Salmonella typhi, qui se propage principalement par la consommation d'eau ou d'aliments contaminés. À Mayotte, territoire français d'outre-mer, la fièvre typhoïde constitue un problème de santé publique majeur, classant l'île parmi les foyers endémiques de la maladie en France. La population mahoraise, en particulier les individus les plus vulnérables, sont confrontés à des défis importants liés à l'accès limité à une eau potable de qualité et à des systèmes d'assainissement efficaces 1.2.3.

Ce chapitre constitue une première étape importante dans notre travail de modélisation de la fièvre typhoïde à Mayotte, visant à mieux comprendre les mécanismes sousjacents de propagation de la maladie dans ce contexte spécifique. Dans ce chapitre, nous proposons un nouveau modèle à compartiments qui prend en compte l'impact de la qualité de l'eau sur la dynamique de la fièvre typhoïde à Mayotte.

Le modèle que nous présentons comprend deux populations distinctes, chacune ayant un accès à une eau de qualité différente. Nous explorons les variations de l'évolution de la maladie dans ces deux populations, en tenant compte des caractéristiques sanitaires spécifiques à Mayotte. Nous examinons le comportement du modèle dans deux versions : déterministe et stochastique, afin d'analyser la dynamique de la maladie dans différents scénarios.

L'organisation du chapitre est comme suit : dans la section 2.2, nous présentons la version déterministe de notre modèle et décrivons la formulation du modèle à compartiments. Nous effectuons ensuite des simulations pour caractériser le comportement qualitatif du modèle dans différents scénarios. Nous évaluons également la virulence de la maladie en utilisant le paramètre du nombre de reproduction de base, qui nous permet de quantifier la capacité de la maladie à se propager dans la population.

Dans la section 2.3, nous présentons la version stochastique du modèle, basée sur un processus markovien déterministe par morceaux (PDMP), qui permet de prendre en compte la variabilité intrinsèque des systèmes biologiques et d'étudier la dynamique de la maladie dans des conditions plus réalistes. Nous effectuons ensuite des simulations des trajectoires du PDMP pour différents valeurs des paramètres du modèle.

Enfin, dans la section 2.4, nous concluons ce chapitre en résumant les résultats obtenus.

## 2.2 Modélisation déterministe

Comme évoqué dans l'état de l'art sur la modélisation mathématique de la fièvre typhoïde 1.5, de nombreux modèles déterministes ont été proposés pour décrire la dynamique de propagation de cette maladie hydrique (voir par exemple [28]). Cependant, très peu de ces modèles prennent explicitement en compte la qualité de l'eau [30], et aucun ne considère spécifiquement deux populations interagissant avec un accès à des eaux de qualité différente.

Notre modèle déterministe, développé dans cette section, tient compte des particularités épidémiologiques de la fièvre typhoïde à Mayotte, notamment de deux populations présentant des caractéristiques sanitaires distinctes et ayant accès à des eaux de qualité différente, en adéquation avec les caractéristiques propres à la population mahoraise.
2.2.1 Formulation du modèle

Dans cette section, un modèle mathématique déterministe à compartiments pour dé-
crire la dynamique de transmission de la fièvre typhoïde est formulé de manière à étendre
et à compléter ceux qui existent dans la littérature [4]-[6], [26]-[30]. Le modèle subdi-
vise la population humaine en quatre types de compartiments : S susceptible, I infecté
(personnes infectieuses présentant des signes cliniques), C porteur infecté (porteurs chro-
niques infectieux sans signes cliniques) et R retiré.

Le modèle (Fig.2.1) prend en compte une transmission directe de la typhoïde des in-
dividus infectés aux individus susceptibles. Cependant, la typhoïde est en grande partie
transmise par des bactéries environnementales par l'intermédiaire d'eau, d'aliments, de
boissons contaminées. Pour intégrer ce phénomène biologique, nous considérons un com-
partiment supplémentaire W, qui représente la présence des bactéries dans l'eau.
Nous subdivisons notre population susceptible en deux sous-groupes :

S1 : Population ayant un accès à l'eau de bonne qualité où la contamination se produit à un taux βₚ(I + ωC), c'est-à-dire que la transmission de la maladie se fait d'une manière directe (Personne à Personne avec un taux de transmission spécifique βₚ) entre les individus infectés ou chroniques et les individus susceptibles.

S2 : Population n'ayant pas l'accès à l'eau de bonne qualité où la contamination se pro-
duit à un taux β<sub>p</sub>(I + ωC) + β<sub>w</sub>W où β<sub>w</sub>W est la force d'infection véhiculée par l'eau,
c'est-à-dire que la maladie est transmise par des bactéries environnementales par l'inter-
médiaire d'eau en plus du contact direct de personne à personne.
Nous supposons que les échanges démographiques entre ces deux parties de la popula-
tion sont négligeables mais physiquement elles sont complètement mélangées et donc en
contact direct l'une avec l'autre à Mayotte, ce qui explique aussi le choix de ne pas séparer
les I, C, R en deux catégories.

Les individus de la classe infectée restent contagieux pendant une période de $1/\delta$. Nous supposons que les individus retirés sont immunisés à vie. Une petite fraction ($\theta \in [0, 1]$) deviennent des porteurs chroniques (à vie) et asymptomatiques de la typhoïde (C) qui contribuent à la transmission à un taux relatif $\omega < 1$ par rapport aux infectieux. Le reste des individus connaissent une mortalité induite par la maladie ($\alpha \in [0, 1]$). Les porteurs chroniques et les personnes infectées excrètent tous deux des bactéries dans l'environnement. Toutefois, le taux d'excrétion $\gamma$ par le groupe infectieux est supérieur à celui du groupe de chroniques $\gamma\omega$, cela s'explique par le fait que les porteurs infectieux ne présentent aucun signe d'infection. Les bactéries infectieuses perdent leur viabilité (c'est-à-dire qu'elles se dégradent dans l'eau) à un taux $\xi$.

Nous considérons que les individus naissent (à un taux *b*) susceptibles d'être infectés par la typhoïde avec une probabilité *p* ∈ [0, 1] d'être dans le compartiment *S1* et une probabilité 1 − *p* d'être en *S2*. Le taux de mortalité naturelle est noté par *μ*.
Notons par *N* := *S1* + *S2* + *I* + *C* + *R* la taille de la population totale.
À partir des hypothèses ci-dessus, nous formulons le système d'équations différentielles suivant :

$$
\left\{
\begin{array}{l}
\frac{dS1}{dt} = pbN - \mu S1 - \beta_p (I + \omega C) S1 \\
\frac{dS2}{dt} = (1-p)bN - \mu S2 - (\beta_p (I + \omega C) + \beta_w W) S2 \\
\frac{dI}{dt} = -(\mu + \delta)I + \beta_p (I + \omega C) (S1 + S2) + \beta_w W S2 \\
\frac{dC}{dt} = -\mu C + \delta \theta I \\
\frac{dR}{dt} = -\mu R + \delta (1 - \alpha - \theta) I \\
\frac{dW}{dt} = -\xi W + \gamma (I + \omega C)
\end{array}
\right.
\qquad (2.1)
$$

où les variables et les paramètres sont résumés dans les **tables 2.1** et **2.2**. Le diagramme en compartiments correspondant à ce système est donné par la **figure 2.1**.

<figure><img src="image_19.png" /><figcaption>FIGURE 2.1 – Diagramme de la structure du modèle.</figcaption></figure>
<table><caption>TABLE 2.1 – *Description des variables d'état du modèle épidémique.*</caption><thead><tr><td>Variable d'Etat</td><td>Description</td></tr></thead><tbody><tr><td>S1</td><td>Nombre de personnes susceptibles ayant un accès à l'eau de bonne qualité</td></tr><tr><td>S2</td><td>Nombre de personnes susceptibles n'ayant pas un accès à l'eau de bonne qualité</td></tr><tr><td>I</td><td>Nombre de personnes infectieuses présentant des signes cliniques</td></tr><tr><td>C</td><td>Nombre de personnes infectieuses sans signes cliniques (porteurs chroniques)</td></tr><tr><td>R</td><td>Nombre de personnes ayant une immunité permanente</td></tr><tr><td>N</td><td>Taille de la population totale</td></tr><tr><td>W</td><td>Concentration de bactéries dans l'eau en UFC/ml [31]</td></tr></tbody></table>

### 2.2.2 Calibration des paramètres

La calibration des paramètres d'un modèle est une étape importante pour obtenir des résultats de simulation qui correspondent au mieux aux données réelles. Pour notre modèle décrivant la propagation de la fièvre typhoïde à Mayotte (voir **Fig.2.1**),

* Le taux de natalité *b* et le taux de mortalité naturelle *μ* sont souvent les plus faciles à calibrer car ils sont généralement bien documentés. À Mayotte, ils sont respectivement de 38.9 et 2.9 pour 1000 habitants [9].

Les autres paramètres du modèle présentés dans la table **2.2** ne sont pas identifiés dans la littérature disponible à Mayotte.

* La durée de la contagiosité peut varier en fonction de facteurs tels que l'âge, la souche de la bactérie Salmonella Typhi, la gravité de l'infection et le traitement antibiotique administré. Les valeurs peuvent donc varier en fonction de la population étudiée. Elle a été étudiée dans plusieurs pays d'Afrique et d'Asie, et les valeurs rapportées varient dans une plage de 1 à 5 semaines [30, 32, 33].

* Les taux de transmission directe et indirecte de la fièvre typhoïde varient considérablement selon les contextes épidémiologiques et les populations étudiées. Dans la littérature [34, 35], les valeurs de $\beta_p$ peuvent varier entre 0.1 et 2.3 par personne et par jour, tandis que les valeurs de $\beta_w$ peuvent varier entre 0.01 et 0.6 par personne et par jour. Cependant, il est important de noter que ces valeurs sont souvent spécifiques à chaque contexte et doivent être ajustées en fonction des caractéristiques de la population étudiée.

* Le taux d'excration dans l'approvisionnement en eau ($\gamma$) et le taux de dégradation des particules infectieuses dans l'eau ($\xi$) sont des paramètres qui dépendent de plusieurs facteurs tels que la qualité de l'eau, la présence d'autres microorganismes dans l'eau, etc. Il peut donc être plus difficile de trouver des données précises pour les calibrer. Pour donner un ordre de grandeur, certaines études ont rapporté des valeurs pour $\gamma$ allant de 0.01 à 1.5 par jour [30, 36], tandis que pour $\xi$, les valeurs varient entre 0.01 et 1.5 par jour
<table><caption>TABLE 2.2 – Description des paramètres du modèle.</caption><thead><tr><td>Les Paramètres</td><td>Description</td></tr></thead><tbody><tr><td>b</td><td>Taux de natalité</td></tr><tr><td>p</td><td>Probabilité que les nouveau-nés soient en S1</td></tr><tr><td>μ</td><td>Taux de mortalité naturelle</td></tr><tr><td>1/δ</td><td>Durée de la contagiosité</td></tr><tr><td>β<sub>p</sub></td><td>Taux de transmission directe (de personne à personne)</td></tr><tr><td>β<sub>w</sub></td><td>Taux de transmission indirecte (transmise par l'eau)</td></tr><tr><td>α</td><td>Fraction des infectés qui connaissent une mortalité induite par la maladie</td></tr><tr><td>θ</td><td>Fraction des infectés qui deviennent porteurs chroniques</td></tr><tr><td>ω</td><td>Infectiosité relative des malades chroniques (Contribution des porteurs chroniques dans la propagation)</td></tr><tr><td>γ</td><td>Taux d'excration dans l'approvisionnement en eau</td></tr><tr><td>ξ</td><td>Taux de dégradation des particules infectieuses dans l'eau</td></tr></tbody></table>

[36, 37]. Cependant, ces valeurs doivent être considérées avec prudence car elles peuvent varier considérablement selon les conditions locales.

► les paramètres α, θ et ω sont des paramètres spécifiques au modèle et à la maladie étudiée. Les valeurs de ces paramètres peuvent varier selon les différentes populations étudiées et les caractéristiques de la maladie.

En général, la fraction des infectés qui connaissent une mortalité induite par la maladie (α) peut varier en fonction de la gravité de la maladie et de l'accès aux soins médicaux. Pour la fièvre typhoïde, la mortalité peut varier entre 0.1% et 10% en fonction des conditions sanitaires et de l'accès aux soins [38].

La fraction des infectés qui deviennent porteurs chroniques (θ) peut également varier selon la population étudiée et les caractéristiques de la maladie. Pour la fièvre typhoïde, la proportion de porteurs chroniques peut varier de 1% à 5% [39].

Enfin, l'infectiosité relative des malades chroniques (ω) peut être difficile à quantifier, car elle dépend également de facteurs tels que la dose infectieuse et la durée de l'infection chronique. Les valeurs rapportées dans la littérature pour ce paramètre peuvent varier de 0.5% à 40% [40, 30].

Il convient donc de bien définir ces paramètres et de trouver des données spécifiques à la population étudiée pour pouvoir les calibrer correctement. Ainsi, des études futures seront nécessaires pour estimer ces paramètres afin de calibrer notre modèle général pour la propagation de la fièvre typhoïde à Mayotte. Toutefois, il est possible de réaliser des simulations numériques sous différents scénarios pour caractériser qualitativement le comportement de la transmission de cette maladie hydrique.
### 2.2.3 Discussion et simulations numériques

En nous basant sur la partie précédente portant sur la calibration des paramètres, nous effectuons des simulations numériques pour illustrer le comportement du système (2.1) sous trois scénarios différents. Les simulations ont été effectuées en utilisant le logiciel R en implémentant un schéma d'Euler, dont le pseudo-code (2.2.3) est présenté dans l'algorithme ci-dessous :

<table><caption>Algorithme 1 : Algorithme d'Euler pour le système (2.1)</caption><tr><td></td><td>entrée</td><td>paramètres du modèles : b, μ, p, β<sub>p</sub>, β<sub>w</sub>, ω, α, δ, θ, γ, ξ.<br/>valeurs initiales des variables : s<sub>1</sub><sup>0</sup>, s<sub>2</sub><sup>0</sup>, i<sup>0</sup>, c<sup>0</sup>, r<sup>0</sup>, w<sup>0</sup>.<br/>horizon de simulation H, pas de discrétisation en temps Δt.</td></tr><tr><td></td><td>sortie</td><td>Trajectoire (X<sub>t</sub> = (S1(t), S2(t), I(t), C(t), R(t), W(t)))<sub>0≤t≤H</sub> par un pas de Δt</td></tr><tr><td>1</td><td>début</td><td></td></tr><tr><td>2</td><td>Initialisation</td><td>(S1(0), S2(0), I(0), C(0), R(0), W(0)) ← (s<sub>1</sub><sup>0</sup>, s<sub>2</sub><sup>0</sup>, i<sup>0</sup>, c<sup>0</sup>, r<sup>0</sup>, w<sup>0</sup>)</td></tr><tr><td>3</td><td></td><td>N(0) ← S1(0) + S2(0) + I(0) + C(0) + R(0)</td></tr><tr><td>5</td><td>pour h de 0 à H par pas de Δt faire</td><td></td></tr><tr><td>6</td><td>S1(h + Δt) ← S1(h) + Δt (pbN(h) - μS1(h) - β<sub>p</sub>(I(h) + ωC(h)) S1(h))</td><td></td></tr><tr><td>7</td><td>S2(h + Δt) ← S2(h) + Δt ((1 - p)bN(h) - μS2(h) - (β<sub>p</sub>(I(h) + ωC(h)) + β<sub>w</sub>W(h)) S2(h))</td><td></td></tr><tr><td>8</td><td>I(h + Δt) ← I(h) + Δt (- (μ + δ)I(h) + (β<sub>p</sub>(I(h) + ωC(h))) (S1(h) + S2(h)) + β<sub>w</sub>W(h)S2(h))</td><td></td></tr><tr><td>9</td><td>C(h + Δt) ← C(h) + Δt (- μC(h) + δθI(h))</td><td></td></tr><tr><td>10</td><td>R(h + Δt) ← R(h) + Δt (- μR(h) + δ(1 - α - θ)I(h))</td><td></td></tr><tr><td>11</td><td>W(h + Δt) ← W(h) + Δt (- ξW(h) + γ(I(h) + ωC(h)))</td><td></td></tr><tr><td>12</td><td>N(h + Δt) ← S1(h + Δt) + S2(h + Δt) + I(h + Δt) + C(h + Δt) + R(h + Δt)</td><td></td></tr><tr><td>13</td><td>fin</td><td></td></tr><tr><td>14</td><td>renvoyer</td><td>(X<sub>t</sub> = (S1(t), S2(t), I(t), C(t), R(t), W(t)))<sub>0≤t≤H</sub></td></tr><tr><td>15</td><td>fin</td><td></td></tr></table>

**Premier jeu des paramètres : évolution de l'épidémie sans prise en compte de la démographie.**

Dans un premier temps, nous considérons que la démographie de la population est négligée par rapport à l'évolution de l'épidémie, ce qui signifie que l'impact de l'expansion démographique sur la propagation de la typhoïde est ignoré ($b = \mu = 0$). Cette simplification permet de se concentrer sur les dynamiques de transmission de la maladie sans être influencé par des facteurs externes.

En ce qui concerne les compartiments, la **figure 2.2** représente la population comme étant divisée en plusieurs catégories, chacune étant caractérisée par son état de santé. Les robinets de naissance et de mort sont supposés être fermés, ce qui signifie que la taille de la population reste constante tout au long de la simulation. Initialement, tous les individus sont sains, ce qui est représenté par le premier compartiment en haut $S1 + S2$.

Lorsque les robinets $S1/S2$ et $I$ sont ouverts, l'ensemble de la population s'écoule vers les compartiments $C$ et $R$, tandis qu'une petite fraction des individus infectés s'écoule également vers le compartiment $D$ (Morts), ce qui reflète la mortalité associée à la maladie. Il est important de noter que, tant que $C \neq 0$ ou $W \neq 0$, le robinet $S1/S2$ restera ouvert, même si $I = 0$. Cependant, la propagation de l'épidémie ainsi que sa durée dépendent essentiellement du flux d'entrée/sortie des individus infectés. En effet, plus le flux de sortie
est faible par rapport à l'entrée, plus l'épidémie sera aiguë et sa disparition prendra plus de temps.

<figure><img src="image_20.png" /><figcaption>FIGURE 2.2 – Représentation schématique d'un modèle à compartiment sans démographie.</figcaption></figure>

Dans le cadre de cette simulation, la population initiale est composée de 300000 individus (ordre de grandeur de la population mahoraise) qui sont susceptibles d'être infectés à l'instant $t = 0$. Parmi ces individus, 30% (soit 90000 personnes) n'ont pas accès à de l'eau de bonne qualité, ce qui augmente leur risque d'infection. De plus, nous avons introduit 5 porteurs chroniques dans la population et avons fixé la concentration de bactéries dans l'eau à 20 UFC/ml, et le nombre initial de personnes infectées et retirées à 0. Tout au long des simulations, nous avons utilisé un pas de discrétisation $\Delta t$ de 0.001. Les autres paramètres du modèle sont choisis et fixés selon **2.2.2** et sont répertoriés dans la **table 2.3**.
<table><caption>TABLE 2.3 – *Premier jeu de paramètres*.</caption><thead><tr><th>Les paramètres</th><th>β<sub>p</sub></th><th>β<sub>w</sub></th><th>ω</th><th>δ</th><th>θ</th><th>α</th><th>ξ</th><th>γ</th></tr></thead><tbody><tr><th>Valeurs numériques (1/ semaine)</th><td>0.8</td><td>1.8</td><td>0.07</td><td>0.25</td><td>0.03</td><td>0.001</td><td>0.8</td><td>2</td></tr></tbody></table>

<figure><img src="image_21.png" /><figcaption>FIGURE 2.3 – Évolution des populations S1, S2, I, C, R et la concentration de bactéries dans l'eau W sur un horizon de 100 semaines.</figcaption></figure>

Au fil du temps, le nombre d'individus infectés augmente rapidement jusqu'à atteindre un pic de 130000 personnes (soit 43.3% de la population totale) lors de la 18ème semaine, puis diminue graduellement jusqu'à atteindre zéro en 50 semaines. Au cours de cette période, le nombre de personnes saines tend vers zéro dans les deux sous-ensembles S1 et S2. À la fin de l'épidémie, la population se trouve dans deux états distincts : les porteurs chroniques (représentant 3% de la population) et les retirés (représentant 97% de
la population). En outre, la concentration de bactéries dans l'eau augmente rapidement, atteignant un maximum de 300000 UFC/ml après 18 semaines, avant de diminuer et de se stabiliser à un seuil de 4479 UFC/ml. Malgré la propagation rapide de la maladie, le nombre de décès liés à la typhoïde reste faible, avec seulement 200 décès à la fin de la propagation (Fig.2.3).

Il convient de noter que l'épidémie a pu se propager dans la population même en l'absence d'individus infectés à $t = 0$, en raison de deux sources d'infectiosité : les porteurs chroniques et la concentration de bactéries dans l'eau.

Enfin, il est important de souligner que la population S2 décroit plus rapidement que la population S1 en raison de l'effet supplémentaire de la qualité de l'eau (représentée par le paramètre W). Cette observation met en évidence l'importance de la qualité de l'eau dans la propagation de l'épidémie de typhoïde et souligne la nécessité de mettre en place des mesures efficaces pour prévenir la contamination de l'eau.

## Deuxième jeu de paramètres

Nous maintenons les valeurs des paramètres $\xi$ et $\gamma$, qui décrivent la propagation de la fièvre typhoïde par l'eau, et augmentons les autres paramètres, à l'exception de $\delta$ qui est réduit afin d'obtenir une durée de contagiosité assez longue. Les valeurs modifiées des paramètres sont répertoriées dans la **table 2.4** ci-dessous.

<table><caption>TABLE 2.4 – Deuxième jeu de paramètres.</caption><thead><tr><th>Les paramètres</th><th>β<sub>p</sub></th><th>β<sub>w</sub></th><th>ω</th><th>δ</th><th>θ</th><th>α</th><th>ξ</th><th>γ</th></tr></thead><tbody><tr><th>Valeurs numériques (1/ semaine)</th><td>1.2</td><td>2</td><td>0.4</td><td>0.15</td><td>0.05</td><td>0.01</td><td>0.8</td><td>2</td></tr></tbody></table>
<figure><img src="image_22.png" /><figcaption>FIGURE 2.4 – Évolution des populations S1, S2, I, C, R et la concentration de bactéries dans l'eau W sur un horizon de 300 semaines.</figcaption></figure>

En partant des mêmes conditions initiales, l'analyse des résultats montre que l'épidé-
mie se propage de manière similaire à l'exemple précédent, conformément aux attentes
(Fig.2.4). Cependant, la durée de la maladie est considérablement plus longue, avec un
nombre d'individus infectés qui tend vers zéro après près de 250 semaines. De plus, la
décroissance est beaucoup plus lente, tandis que le pic épidémique est plus élevé et se
produit beaucoup plus rapidement, avec un nombre maximal de 272000 individus infectés
(soit 90.6% de la population). Au fil du temps, la population évolue progressivement vers
les deux états : les porteurs chroniques (5%) et les retirés (93.9%). Toutefois, un petit
pourcentage d'individus subit une mortalité induite par la typhoïde, avec un total de 3000
décès enregistrés au cours de 300 semaines.

Troisième jeu de paramètres : La dimension démographique de la pandémie

Lorsque la démographie de la population est prise en compte, avec un taux de na-
talité élevé, le comportement asymptotique de l'épidémie est considérablement altérée
(Fig.2.5). Tout d'abord, le nombre d'individus infectés et la concentration de bactéries
dans l'eau augmentent rapidement, puis diminuent, ce qui est similaire au modèle sans
démographie. Cependant, l'expansion rapide de la population a un impact significatif sur la propagation de l'épidémie, car les naissances font réaugmenter rapidement les individus infectés, maintenant ainsi les compartiments S1 et S2 remplis d'un côté, et faisant exploser le nombre d'individus infectés (et la concentration de bactéries) de l'autre côté. Les porteurs chroniques, les individus guéris et la taille de la population totale augmentent de manière exponentielle.

Les valeurs des paramètres utilisés dans ce scénario sont indiquées dans le tableau suivant **Table.2.5**. Il convient de noter que nous avons fixé le pourcentage de partage des nouveau-nés entre les deux populations S1 et S2 à 50% ($p = 0.5$).

<table><caption>TABLE 2.5 – *Troisième jeu de paramètres.*</caption><thead><tr><th>Les paramètres</th><th>b</th><th>μ</th><th>β<sub>p</sub></th><th>β<sub>w</sub></th><th>ω</th><th>δ</th><th>θ</th><th>α</th><th>ξ</th><th>γ</th></tr></thead><tbody><tr><td>Valeurs numériques (1/ semaine)</td><td>38.9.10<sup>-3</sup>/52</td><td>2.9.10<sup>-3</sup>/52</td><td>0.8</td><td>1.8</td><td>0.07</td><td>0.15</td><td>0.03</td><td>0.001</td><td>0.8</td><td>2</td></tr></tbody></table>

En d'autres termes, ces modèles mettent en évidence l'interaction complexe entre la dynamique de la population et la propagation de la maladie, soulignant ainsi l'importance de la démographie et de la croissance de la population dans la propagation de la fièvre typhoïde. L'expansion rapide de la population peut accélérer la propagation de l'épidémie, tandis que un taux de natalité plus faible peut atténuer la gravité et la vitesse de transmission.

Cependant, quel que soit le scénario démographique, la population évolue dans une seule direction, car il n'y a pas de possibilité de revenir en arrière et la population restera dans les états infectés, porteurs chroniques ou rétablis.
<figure><img src="image_23.png" /><figcaption>FIGURE 2.5 – Évolution des populations S1, S2, I, C, R, N et W des bactéries sur un horizon de 3000 semaines.</figcaption></figure>

Ces simulations ont permis de mieux comprendre la dynamique complexe de la propagation de la fièvre typhoïde au sein d'une population composée de deux compartiments sains présentant des caractéristiques sanitaires distinctes. Elles ont permis d'explorer différents scénarios et d'analyser l'impact de la démographie, de la transmission directe et indirecte, ainsi que d'autres facteurs sur la propagation de la maladie.

Dans la section suivante, nous allons définir le nombre de reproduction de base, un indicateur clé en épidémiologie permettant d'estimer la dynamique de l'épidémie de la fièvre typhoïde et de projeter son évolution.
### 2.2.4 Le nombre de reproduction de base

Le nombre de reproduction de base, souvent noté $R_0$, est un concept clé en épidémiologie. Il représente le nombre moyen de personnes qu'un individu infecté peut contaminer dans une population susceptible, en l'absence de toute mesure de prévention ou de contrôle de l'épidémie [41]. Ce paramètre est très important dans l'étude de la propagation des épidémies, étant donné qu'il permet de quantifier la virulence initiale d'une maladie infectieuse, et donc de prévoir l'évolution de l'épidémie à court terme. Si $R_0 < 1$, cela signifie que chaque individu infectieux produira en moyenne moins d'une infection secondaire. Et si $R_0 > 1$, chaque individu infectieux produira en moyenne plus d'une infection secondaire et donc que la maladie se répandra dans la population. Ce paramètre $R_0$ est calculé en utilisant l'approche développée par van den Driessche et Watmough [42].

**Théorème 2.1.** Si $b = \mu \neq 0$ (le taux de natalité est égal au taux de mortalité naturelle) alors, le nombre de reproduction de base pour le modèle (2.1) est donné par :

$$R_0 = \frac{1}{(\mu + \delta)p} \left( \beta_p + \frac{\gamma\beta_w(1-p)}{\xi} \right) \left( 1 + \frac{\delta\theta\omega}{\mu} \right).$$

*Preuve.* Notons par $X = (S1, S2, I, C, R, W)$ et $Y = (I, C, W)$. Nous calculons le nombre de reproduction de base en utilisant l'approche *next generation* [42]. Nous définissons le taux de changement des compartiments infectieux à l'équilibre sans maladie (les compartiments infectés sont I, C et W) :

$$\frac{dY}{dt} = \left( \frac{\frac{dI}{dt}}{\frac{dC}{dt}} \right) = \mathcal{F}(X) - \mathcal{V}(X).$$

La partie $\mathcal{F}(X)$ représente les taux d'individus nouvellement infectés rejoignant chaque compartiment. Ces taux dépendent des interactions entre les compartiments et d'autres facteurs influençant la propagation de la maladie.

La partie $\mathcal{V}(X)$ correspond à la différence entre les taux d'individus sortant et entrant dans chaque compartiment par d'autres moyens que l'infection. Cela peut inclure des processus tels que la guérison, la mortalité ou d'autres mécanismes spécifiques à la fièvre typhoïde.

$$\mathcal{F}(X) = \begin{pmatrix} \beta_p(I + \omega C)(S1 + S2) + \beta_w W S2 \\ 0 \\ 0 \end{pmatrix} \quad \text{et} \quad \mathcal{V}(X) = \begin{pmatrix} (\mu + \delta)I \\ \mu C - \delta\theta I \\ \xi W - \gamma(I + \omega C) \end{pmatrix}.$$

Le point d'équilibre sans maladie du système (2.1) $X^* = (S1^*, S2^*, 0, 0, 0, 0)$ est obtenu en résolvant :

$$\frac{dS1}{dt} = \frac{dS2}{dt} = \frac{dI}{dt} = \frac{dC}{dt} = \frac{dR}{dt} = \frac{dW}{dt} = 0,$$

sous la contrainte de l'absence de maladie, c'est-à-dire :

$$I = C = W = 0.$$

Ainsi, $S1^*$ et $S2^*$ sont toute solution strictement positive de :

$$\begin{cases} pb(S1^* + S2^*) = \mu S1^*, \\ (1-p)b(S1^* + S2^*) = \mu S2^*. \end{cases}$$
Ce système admet une solution si et seulement si $b = \mu$ et dans ce cas, il a une infinité de solutions :

$$S2^* = \frac{1-p}{p}S1^*.$$

Ainsi, $X^* = (S1^*, \frac{1-p}{p}S1^*, 0, 0, 0, 0)$.

Ensuite, soit F et V les matrices de $3 \times 3$ (nombre des compartiments infectés) définies par :

$$F = \left[ \frac{\partial \mathcal{F}_i}{\partial Y_j} (X^*) \right] \text{ et } V = \left[ \frac{\partial \mathcal{V}_i}{\partial Y_j} (X^*) \right] \text{ avec } 1 \le i, j \le 3.$$

L'indice de reproduction de base $R_0$ est égal à le rayon spectral de la matrice de *next generation*, $FV^{-1}$.

Dans notre cas, on a :

$$F = \begin{pmatrix} \frac{\beta_p S1^*}{p} & \frac{\beta_p \omega S1^*}{p} & \frac{\beta_w (1-p) S1^*}{p} \\ 0 & 0 & 0 \\ 0 & 0 & 0 \end{pmatrix} \text{ et } V = \begin{pmatrix} \delta + \mu & 0 & 0 \\ -\delta\theta & \mu & 0 \\ -\gamma & -\omega\gamma & \xi \end{pmatrix},$$

par conséquent

$$V^{-1} = \frac{1}{\xi\mu(\delta + \mu)} \begin{pmatrix} \xi\mu & 0 & 0 \\ \xi\delta\theta & \xi(\delta + \mu) & 0 \\ \gamma(\mu + \delta\theta\omega) & \omega\gamma(\delta + \mu) & \mu(\delta + \mu) \end{pmatrix} \text{ et } FV^{-1} = \begin{pmatrix} R_0 & A & B \\ 0 & 0 & 0 \\ 0 & 0 & 0 \end{pmatrix},$$

où

$$R_0 = \frac{S1^*}{(\delta + \mu)} \left( \frac{\beta_p}{p} + \frac{\gamma\beta_w(1-p)}{p\xi} \right) \left( 1 + \frac{\delta\theta\omega}{\mu} \right),$$

$$A = \frac{\omega S1^*}{\mu p} \left( \beta_p + \frac{\beta_w (1-p)\gamma}{\xi} \right),$$

$$B = \frac{\beta_w (1-p) S1^*}{p\xi}.$$

Sans perte de généralité, nous pouvons supposer que $S1^* = 1$, ainsi :

$$R_0 = \frac{1}{(\delta + \mu)} \left( \frac{\beta_p}{p} + \frac{\gamma\beta_w(1-p)}{p\xi} \right) \left( 1 + \frac{\delta\theta\omega}{\mu} \right).$$

□

Cette expression de $R_0$ peut être décomposée pour montrer les contributions de la transmission directe et indirecte. Il est possible de réécrire $R_0 = R_{0,p} + R_{0,w}$, où

$$R_{0,p} = \frac{\beta_p}{p(\delta + \mu)} \left(1 + \frac{\delta\theta\omega}{\mu}\right) \text{ et } R_{0,w} = \frac{\gamma\beta_w(1-p)}{p(\delta + \mu)\xi} \left(1 + \frac{\delta\theta\omega}{\mu}\right).$$

Le terme $R_{0,p}$ représente la contribution de la transmission directe de personne à personne. Il dépend du taux de transmission directe ($\beta_p$), du taux de guérison ($\delta$), du taux de mortalité naturelle ($\mu$), de la fraction des infectés qui deviennent porteurs chroniques ($\theta$) et de l'infectiosité relative des malades chroniques ($\omega$). Plus le taux de transmission directe est élevé et plus la probabilité de devenir porteur chronique est importante, plus
la contribution de la transmission directe sera élevée.

Le terme $R_{0,w}$ représente la contribution de la transmission indirecte par l'eau. Il dépend du taux d'excrétion dans l'approvisionnement en eau ($\gamma$), du taux de transmission indirecte par l'eau ($\beta_w$), de la fraction des nouveaux-nés qui partagent l'eau contaminée ($1 - p$), du taux de dégradation des particules infectieuses dans l'eau ($\xi$), du taux de guérison ($\delta$), du taux de mortalité naturelle ($\mu$), de la fraction des infectés qui deviennent porteurs chroniques ($\theta$) et de l'infectiosité relative des malades chroniques ($\omega$). Si le taux d'excrétion dans l'eau est élevé, le taux de transmission indirecte par l'eau est important et la contribution de la transmission indirecte sera élevée.

À titre d'illustration, la figure 2.6 met en évidence l'évolution du $R_0$ en fonction du taux de transmission par l'eau $\beta_w$ et de l'infectiosité relative des porteurs chroniques $\omega$, tandis que les autres paramètres sont maintenus constants (Table 2.6).

Dans le cas de la fièvre typhoïde, qui est transmise principalement par l'eau contaminée, la valeur de $R_0$ dépend de plusieurs facteurs, tels que le taux de transmission par l'eau, l'infectiosité relative des porteurs chroniques et la qualité de l'eau (Fig.2.6). Plus ces facteurs sont élevés, plus le $R_0$ dépasse la valeur seuil de 1, indiquant ainsi que la maladie se répandra massivement dans la population. En revanche, en l'absence de ces deux facteurs d'infectiosité ($\beta_w = \omega = 0$), le $R_0$ reste inférieur à 1, signifiant que chaque individu infectieux produira en moyenne moins d'une infection secondaire.

Il est donc crucial de prendre en compte ces facteurs dans la mise en place de mesures de prévention et de contrôle de la fièvre typhoïde. Des actions telles que l'amélioration de la qualité de l'eau, l'identification et le traitement des porteurs chroniques, ainsi que la sensibilisation de la population à l'importance de l'hygiène et de la propreté, peuvent contribuer à réduire la valeur de $R_0$ et donc à ralentir la propagation de la maladie.

<table><caption>TABLE 2.6 – *Les valeurs utilisées pour le calcul de R₀.*</caption><thead><tr><th>Les paramètres</th><th>β<sub>p</sub></th><th>μ</th><th>δ</th><th>θ</th><th>ξ</th><th>γ</th></tr></thead><tbody><tr><td>Valeurs numériques</td><td>0.8</td><td>2.9.10<sup>-3</sup>/52</td><td>0.25</td><td>0.03</td><td>0.8</td><td>2</td></tr></tbody></table>
<figure><img src="image_24.png" /><figcaption>FIGURE 2.6 – Évolution du $R_0$ en fonction du taux de transmission par l'eau $\beta_w$ et de l'infectiosité relative des malades chroniques $\omega$.</figcaption></figure>

En conclusion, les simulations réalisées avec différents scénarios de propagation de la maladie de typhoïde ont permis d'identifier les comportements qualitatifs de l'épidémie. Cependant, le taux de reproduction de base $R_0$ est un outil intrinsèque pour évaluer la virulence initiale de la maladie en fonction des différents facteurs épidémiologiques. Tandis que les modèles déterministes permettent de calculer facilement le $R_0$ et de prévoir les conditions de transmission de la maladie. Cependant, le profil déterministe proposé ne colle pas du tout à Mayotte, où la taille de la population est petite, et le nombre des cas déclarés reste faible, fluctue et présente une extinction et une réapparition de la maladie (voir Fig.1.6). Les modèles stochastiques seront donc mieux adaptés à ce type de contexte. Ils permettent de prendre en compte les variations aléatoires et les effets de petite taille de la population, ce qui peut mieux correspondre à la réalité observée.

Ces résultats mettent en évidence l'importance de l'utilisation de simulations stochastiques pour mieux comprendre et prévoir la propagation de maladies infectieuses dans des contextes épidémiologiques complexes comme celui de Mayotte.

## 2.3 Modélisation par un processus markovien déterministe par morceaux

Les modèles stochastiques sont capables de capturer les fluctuations des populations de petite taille contrairement aux modèles déterministes. Ces modèles prennent en compte l'aspect aléatoire des événements impliqués dans la propagation et le développement de maladies infectieuses, et permettent de travailler sur des petits effectifs souvent rencontrés dans les études portant sur des ensembles de populations locales interconnectées d'une
même espèce.

Comme mentionné dans la section **1.2.4** et dans la conclusion précédente, l'évolution de la fièvre typhoïde à Mayotte est caractérisée par une fluctuation entre des valeurs faibles, nécessitant une étude probabiliste de cette dynamique. Dans cette partie, nous introduisons une modélisation hybride de la transmission de la typhoïde en utilisant un processus markovien déterministe par morceaux (**1.4.2**). Cette méthode permet de prendre en compte simultanément des variables continues et des variables discrètes.

### 2.3.1 Formulation du modèle

Nous proposons une approximation hybride originale ($I_t, X_t$) $t \ge 0$ du modèle déterministe (**2.1**) pour formuler notre PDMP. La première composante $I_t$ est discrète et représente le nombre d'individus infectés à l'instant $t$, prenant des valeurs dans $\mathbb{N}$. La deuxième composante $X_t = (S1_t, S2_t, C_t, R_t, W_t)$ est continue à valeurs dans $\mathbb{R}_+^5$. Elle décrit l'évolution des autres variables d'état définies dans la **Table 2.1**.

La **Table 2.7** décrit les sauts possibles du processus $I_t$ partant de l'état $i$, avec $i \in \mathbb{N}$, conditionnellement à $X_t = (s1, s2, c, r, w)$.

<table><caption>TABLE 2.7 – *Sauts possibles du processus I<sub>t</sub>*.</caption><thead><tr><th>Evénements</th><th>Transitions</th><th>Taux</th></tr></thead><tbody><tr><td>Une nouvelle infection</td><td>$i \to i + 1$</td><td>$\beta_p(i + \omega c)(s1 + s2) + \beta_w ws2$</td></tr><tr><td>Une guérison ou<br>passage à l'état<br>chronique ou mort</td><td>$i \to i - 1$</td><td>$(\mu + \delta)i$</td></tr></tbody></table>

Entre deux sauts de la partie discrète ($I_t$)$_{t \ge 0}$, la trajectoire de ($X_t$)$_{t \ge 0}$ déterministe, est solution du système d'équations différentielles (**2.2**) qui dépend de l'état discret $i$. Lors d'un saut de $I_t$, de $i$ vers $i+1$ ou bien $i-1$, la variable continue $X_t$ ne bouge pas.

$$ 
\left\{
\begin{aligned}
\frac{dS1}{dt} &= pbN_t - \mu S1_t - \beta_p(i + \omega C_t) S1_t \\
\frac{dS2}{dt} &= (1-p)bN_t - \mu S2_t - (\beta_p(i + \omega C_t) + \beta_w W_t) S2_t \\
\frac{dC}{dt} &= -\mu C_t + \delta\theta i \\
\frac{dR}{dt} &= -\mu R_t + \delta(1 - \alpha - \theta)i \\
\frac{dW}{dt} &= -\xi W_t + \gamma(i + \omega C_t) \\
N_t &= S1_t + S2_t + i + C_t + R_t.
\end{aligned}
\right.
\qquad (2.2) $$

Le flot du PDMP est déterminé par les équations différentielles du système (**2.2**) reliant les variables continues $S1_t, S2_t, C_t, R_t, W_t$ en fonction du temps $t$, et ce, entre deux sauts du processus $I_t$. Le flot est donné par ces équations différentielles.

L'intensité du PDMP est définie par les taux de transition des sauts du processus $I_t$ décrits dans la **Table 2.7**. Ces taux de transition dépendent de l'état discret $i$ et des valeurs des variables continues $s1, s2, c, r, w$. Lorsque le processus fait un saut de l'état $i$ à $i+1$, le taux de transition est donné par $\beta_p(i + \omega c)(s1 + s2) + \beta_w ws2$. Lorsque le processus fait un saut de l'état $i$ à $i-1$, le taux de transition est $(\mu + \delta)i$.

Le noyau du PDMP est défini comme la probabilité de transition d'un état à un autre
entre les sauts du processus discret, conditionnellement à la trajectoire continue entre ces sauts. Dans le cas présent, le noyau du PDMP correspond à la probabilité conditionnelle de transition de l'état $i$ à $i+1$ ou de l'état $i$ à $i-1$, étant donné les valeurs continues de $X_t$.

En utilisant les taux de transition mentionnés dans la Table 2.7, la probabilité conditionnelle de transition de $i$ à $i+1$ entre les sauts est donnée par :

$$p(i \rightarrow i+1 | s1, s2, c, r, w) = \frac{\beta_p(i + \omega c)(s1 + s2) + \beta_w ws2}{\beta_p(i + \omega c)(s1 + s2) + \beta_w ws2 + (\mu + \delta)i}$$

De même, la probabilité conditionnelle de transition de $i$ à $i-1$ entre les sauts est donnée par :

$$p(i \rightarrow i-1 | s1, s2, c, r, w) = \frac{(\mu + \delta)i}{\beta_p(i + \omega c)(s1 + s2) + \beta_w ws2 + (\mu + \delta)i}$$

En grande population pour $I_0$, le modèle déterministe (2.1) présenté dans la partie 2.2 est une approximation du PDMP. En effet, dans ce cas, les fluctuations dues à la petite taille de la population sont négligeables et le modèle déterministe suffit à décrire l'évolution de la maladie. Cependant, dans le cas d'une population infectée de petite taille, les fluctuations peuvent avoir un impact significatif sur la dynamique de la maladie, ce qui justifie l'utilisation d'un PDMP pour capturer l'effet stochastique de la propagation de la maladie.

### 2.3.2 Discussion et simulations numériques

Nous présentons quelques exemples de trajectoires simulées de notre PDMP en fonction des paramètres décrits dans la Table 2.3, en utilisant les mêmes conditions initiales pour chaque simulation. Le pseudo-code de R est présenté dans l'algorithme 2.3.2.

1er cas : l'évolution de l'épidémie sans démographie (Fig.2.7)
<figure><img src="image_25.png" /><figcaption>FIGURE 2.7 – Évolution des populations S1, S2, I, C, R et la concentration de bactéries dans l'eau W sur un horizon de 100 semaines.</figcaption></figure>

La figure 2.7 illustre l'évolution de la population selon deux mécanismes distincts. Au début, la dynamique de l'épidémie suit un comportement très similaire à celui du modèle déterministe présenté dans la partie 2.2, avec une croissance rapide de l'infec- tion et un pic épidémique important. Cette croissance rapide est rendue possible par le grand nombre de susceptibles au début de l'épidémie. Cependant, en dessous d'un certain seuil de population S1 + S2, la dynamique de la typhoïde ralentit considérablement. Le nombre d'infectés présente des fluctuations avec de petites valeurs, et même s'il atteint zéro et disparaît complètement de la population (I=0), il peut réapparaître en raison de la présence des porteurs chroniques et des bactéries dans l'eau, qui constituent les deux sources d'infectiosité dans notre population. Cette propagation lente est caractérisée par des retours fréquents à l'état 0, ce qui signifie que l'épidémie est s'éteint, mais qu'elle
reprend si les conditions sont à nouveau favorables (à cause de C ici).

<table><caption>Algorithme 2 : Algorithme de simulation pour un PDMP</caption>
  <tbody>
    <tr>
      <td>entrée</td>
      <td>paramètres du modèles : b, μ, p, β<sub>p</sub>, β<sub>w</sub>, ω, α, δ, θ, γ, ξ.<br/>valeurs initiales des variables : s1<sup>0</sup>, s2<sup>0</sup>, i<sup>0</sup>, c<sup>0</sup>, r<sup>0</sup>, 0.<br/>un majorant global M de l'intensité :<br/>horizon de simulation H, pas de discrétisation en temps Δt</td>
    </tr>
    <tr>
      <td>sortie</td>
      <td>Trajectoire (I<sub>t</sub>, X<sub>t</sub> = (S1(t), S2(t), C(t), R(t), W(t)))<sub>0≤t≤H</sub> par pas de Δt</td>
    </tr>
    <tr>
      <td>1</td>
      <td>début</td>
    </tr>
    <tr>
      <td>2</td>
      <td>Initialisation</td>
    </tr>
    <tr>
      <td>3</td>
      <td>(S1(0), S2(0), i, C(0), R(0), W(0)) ← (s1<sup>0</sup>, s2<sup>0</sup>, i<sup>0</sup>, c<sup>0</sup>, r<sup>0</sup>, w<sup>0</sup>)</td>
    </tr>
    <tr>
      <td>4</td>
      <td>t ← 0</td>
    </tr>
    <tr>
      <td>5</td>
      <td>tant que t &lt; H faire</td>
    </tr>
    <tr>
      <td>6</td>
      <td>Tirer le prochain saut virtuel : Tirer T selon la loi Exp(M)</td>
    </tr>
    <tr>
      <td>7</td>
      <td>Ne pas dépasser l'horizon : T ← min{T, H - t}</td>
    </tr>
    <tr>
      <td>8</td>
      <td>Faire avancer le flot :</td>
    </tr>
    <tr>
      <td>9</td>
      <td>pour h de 0 à T, par pas de Δt faire :</td>
    </tr>
    <tr>
      <td>10</td>
      <td>S1(h + Δt) ← S1(h) + Δt (pbN(h) - μS1(h) - β<sub>p</sub>(i + ωC(h)) S1(h))</td>
    </tr>
    <tr>
      <td>11</td>
      <td>S2(h + Δt) ← S2(h) + Δt ((1 - p)bN(h) - μS2(h) - (β<sub>p</sub>(i + ωC(h)) + β<sub>w</sub>W(h)) S2(h))</td>
    </tr>
    <tr>
      <td>12</td>
      <td>C(h + Δt) ← C(h) + Δt (-μC(h) + δθi)</td>
    </tr>
    <tr>
      <td>13</td>
      <td>R(h + Δt) ← R(h) + Δt (-μR(h) + δ(1 - α - θ)i)</td>
    </tr>
    <tr>
      <td>14</td>
      <td>W(h + Δt) ← W(h) + Δt (-ξW(h) + γ(i + ωC(h)))</td>
    </tr>
    <tr>
      <td>15</td>
      <td>N(h + Δt) ← S1(h + Δt) + S2(h + Δt) + i + C(h + Δt) + R(h + Δt)</td>
    </tr>
    <tr>
      <td>16</td>
      <td>fin</td>
    </tr>
    <tr>
      <td>17</td>
      <td>Mettre à jour le temps : t ← t + T</td>
    </tr>
    <tr>
      <td>18</td>
      <td>Déterminer la nature du saut si on n'a pas atteint l'horizon :</td>
    </tr>
    <tr>
      <td>19</td>
      <td>Si t &lt; H alors</td>
    </tr>
    <tr>
      <td>20</td>
      <td>Tirer U selon la loi U[0, 1]</td>
    </tr>
    <tr>
      <td>21</td>
      <td>Si U &lt; (β<sub>p</sub>(i + ωC(T))(S1(T)+S2(T)) + β<sub>w</sub>W(T)S2(T))/M alors</td>
    </tr>
    <tr>
      <td>22</td>
      <td>naissance : i ← i + 1</td>
    </tr>
    <tr>
      <td>23</td>
      <td>fin</td>
    </tr>
    <tr>
      <td>24</td>
      <td>Sinon Si U &lt; (β<sub>p</sub>(i + ωC(T))(S1(T)+S2(T)) + β<sub>w</sub>W(T)S2(T)) + (δ + μ)i)/M alors</td>
    </tr>
    <tr>
      <td>25</td>
      <td>mort : i ← i - 1</td>
    </tr>
    <tr>
      <td>26</td>
      <td>fin</td>
    </tr>
    <tr>
      <td>27</td>
      <td>Sinon ce n'est pas un saut : ne rien faire</td>
    </tr>
    <tr>
      <td>28</td>
      <td>fin</td>
    </tr>
    <tr>
      <td>29</td>
      <td>fin</td>
    </tr>
    <tr>
      <td>30</td>
      <td>renvoyer : (I<sub>t</sub>, X<sub>t</sub>)<sub>0≤t≤H</sub></td>
    </tr>
    <tr>
      <td>30</td>
      <td>fin</td>
    </tr>
  </tbody>
</table>

2<sup>ème</sup> cas : l'évolution de l'épidémie avec démographie (Fig.2.8)

Dans ce scénario avec démographie, nous simulons notre PDMP en utilisant les mêmes valeurs des paramètres que celles utilisées dans le cas déterministe (Tableau 2.5).
<figure><img src="image_26.png" /><figcaption>FIGURE 2.8 – Evolution des populations S1, S2, I, C, R et W des bactéries sur un horizon de 200 semaines.</figcaption></figure>

Dans le cas de la simulation présentée **2.8**, la population suit un mécanisme similaire à celui décrit précédemment, mais avec des variations notables. Tout d'abord, après environ 200 semaines, on observe un changement dans la population : au début de la maladie, nous remarquons une diminution rapide des sous-groupes S1 et S2. Cependant, ceux-ci augmentent à nouveau par la suite en raison de la démographie de la population, notamment les naissances qui contribuent à l'augmentation du nombre d'individus susceptibles d'être infectés. En outre, le nombre d'individus retirés et de porteurs chroniques augmente de manière significative. Les individus infectés présentent également des fluctuations plus marquées, avec plusieurs pics épidémiques plus forts que dans l'exemple précédent **2.7**. Ces résultats suggèrent une transition de la dynamique de la maladie, passant d'un ré-
gime de propagation rapide à un régime de fluctuations modérées et de persistance de la maladie.

Nous avons réussi à simuler notre PDMP pour différents scénarios, mais il est important de noter que nous n'avons pas encore calibré notre modèle sur les données réelles de Mayotte (Fig.1.6) afin de mieux refléter la réalité locale. En effet, notre modèle contient de nombreux paramètres à calibrer spécifiquement pour Mayotte, et nous disposons seulement des cumuls journaliers des nouveaux cas déclarés, ce qui rend impossible l'ajustement précis du modèle à ces observations. Par conséquent, nous proposons de simplifier ce modèle probabiliste dans le but de réaliser des estimations des paramètres à partir des données disponibles.

## 2.4 Conclusion

Dans ce chapitre, nous avons proposé deux modèles pour étudier théoriquement et numériquement la dynamique de la propagation de la fièvre typhoïde à Mayotte. Nous avons commencé par développer un nouveau modèle déterministe qui prend en compte l'effet de la qualité de l'eau sur la transmission de la typhoïde et son impact sur deux populations de caractéristiques sanitaires différentes ayant accès à des eaux de qualité différente. Nous avons pu simuler des trajectoires pour différents scénarios et interpréter leur comportement qualitatif. Ensuite, nous avons proposé une méthode pour calculer le nombre de reproduction de base, qui est crucial pour évaluer la virulence d'une épidémie. Cependant, ce type de modèle déterministe n'est pas adapté à la réalité de l'évolution de la typhoïde à Mayotte, car sa dynamique qualitative ne correspond pas aux données observées. Nous avons donc choisi de passer à une modélisation stochastique par un PDMP, où nous avons considéré une évolution aléatoire des individus infectés et déterministe pour les autres variables d'état. Nous avons pu simuler certaines trajectoires possibles du PDMP pour différentes valeurs des paramètres. Néanmoins, il n'a pas été possible de calibrer le modèle de manière précise pour Mayotte en raison du grand nombre de paramètres impliqués et du manque de données disponibles.

Pour toutes ces raisons, nous proposons dans le chapitre suivant de simplifier notre PDMP en utilisant un modèle de type naissance et mort avec immigration (1.4.3), afin d'essayer de mettre en place une procédure de calibration et d'estimation des paramètres à partir des données fournies par l'ARS de Mayotte. Le but de cette simplification est de créer un modèle avec le minimum de paramètres possible, afin qu'il soit identifiable et puisse être ajusté aux données de manière plus efficace. Nous nous concentrerons sur l'estimation des paramètres clés qui déterminent le comportement du modèle, ce qui nous permettra de mieux quantifier la propagation de la fièvre typhoïde à Mayotte.
# Chapitre 3

Estimation des paramètres d'un processus caché de naissance et mort linéaire avec immigration

## Sommaire

<table><tr><td>3.1</td><td>Introduction</td><td>56</td></tr><tr><td>3.2</td><td>Modèle et schéma d'observations</td><td>57</td></tr><tr><td>3.2.1</td><td>Description du Modèle</td><td>57</td></tr><tr><td>3.2.2</td><td>Schéma d'observations</td><td>59</td></tr><tr><td>3.2.3</td><td>Stratégie d'estimation</td><td>60</td></tr><tr><td>3.3</td><td>Estimation paramétrique</td><td>60</td></tr><tr><td>3.3.1</td><td>Expression analytique des paramètres</td><td>61</td></tr><tr><td>3.3.2</td><td>Modèle de Markov caché</td><td>69</td></tr><tr><td>3.4</td><td>Étude numérique</td><td>74</td></tr><tr><td>3.4.1</td><td>Données synthétiques</td><td>74</td></tr><tr><td>3.4.2</td><td>Données de la typhoïde à Mayotte</td><td>79</td></tr><tr><td>3.4.3</td><td>Application à l'estimation de la probabilité d'atteindre un seuil épidémique dans une période de temps</td><td>80</td></tr><tr><td>3.5</td><td>Conclusion</td><td>84</td></tr></table>
## 3.1 Introduction

Le premier modèle proposé dans le chapitre précédent **2.1**, qui comporte de nombreux paramètres, s'est avéré trop complexe pour être calibré avec précision à l'aide des données fournies par l'agence régionale de santé. En effet, on ne dispose d'aucune information sur les éventuels malades chroniques ni aucune mesure sur les bactéries dans l'eau. Les seules données disponibles sont les cumuls journaliers des nouveaux cas déclarés, qui sont considérés comme isolés selon l'ARS de Mayotte [11]. Afin d'obtenir un modèle ajusté aux données disponibles et d'estimer les paramètres clés de la dynamique de la fièvre typhoïde à Mayotte, nous proposons un modèle stochastique simplifié à un compartiment (**Fig.3.1**), qui ne prend en compte que les effectifs de la population infectée *I*. Dans ce modèle simplifié, nous négligeons les compartiments $S1$ et $S2$ qui sont supposés être infinis par rapport à *I*, et nous n'incluons pas le suivi des porteurs chroniques *C* ni des individus rétablis *R*, puisque nous n'observons que les cumuls de cas isolés. Ce modèle stochastique est un processus de naissance et de mort linéaire avec immigration (LBDPI)(1.4.3). *L'immigration* correspond ici à une contamination exogène, ou contamination de l'environnement, la naissance correspond à une nouvelle contamination de personne à personne, et la mort correspond à l'isolement. Ces processus appartiennent à la famille des processus de Markov de sauts pures et sont utilisés dans un grand nombre de disciplines quantitatives, telles que la biologie évolutive, l'épidémiologie, l'écologie, la génétique des populations, la spéciation et l'extinction, ou la théorie des files d'attente [43, 44, 45, 46, 47].

Ce modèle de LBDPI peut être décrit avec seulement trois paramètres : les taux d'immigration (contamination exogène ou environnementale), de naissance (contamination de personne à personne) et de mort (isolement). Son comportement en temps long est bien conn (voir le Th.1.2), de même que les estimateurs du maximum de vraisemblance basés sur des observations discrètes ou continues dès que les observations disponibles incluent les temps de saut et les positions après saut du processus [48, 49, 50, 51].

Toutefois, dans la pratique épidémiologique, il est rare que l'on dispose de comptages continus ou périodiques des individus infectés. Pour les maladies à déclaration obligatoire, les données les plus couramment collectées sont les cumuls des nouveaux infectés sur des périodes données, généralement des intervalles quotidiens ou hebdomadaires. L'inférence des paramètres pour le processus de LBDPI dans ce cadre d'observation est difficile pour deux raisons principales :

Tout d'abord, lorsque les observations sont discrètes, il est nécessaire d'extraire les probabilités de transition des infectés pour une période de temps donné. La composante immigration invalide la propriété de branchement, de sorte que l'obtention d'une forme analytique pour les probabilités de transition n'est pas simple, et l'obtention d'une forme numérique peut être sujette à des instabilités numériques [52]. Plusieurs chercheurs ont proposé des approches pour estimer les paramètres des processus de naissance et mort généraux observés discrètement. Crawford et Suchard et Crawford et al [47, 53] proposent d'utiliser la théorie des fractions continues pour trouver l'expression explicite de la transformée de Laplace des probabilités de transition afin de considérer le problème d'inférence avec des observations discrètes comme un problème de données manquantes en utilisant l'algorithme d'expectation-maximisation (EM) pour trouver des estimations du maximum de vraisemblance. Xu et al [54] utilisent la même méthode d'algorithme EM pour proposer une approximation de processus de branchement multi-types pour birth-death shift processes. Dans [55], les auteurs proposent d'utiliser deux approches pour estimer les
taux de *naissance* et de *mort* à partir d'un processus de naissance et mort linéaire observé discrètement de type Galton-Watson process et en maximisant la vraisemblance par une méthode d'approximation du point-selle. Malheureusement, aucune de ces approches ne peut être utilisée dans notre cas.

Deuxièmement, dans tous les articles cités précédemment, les observations sont discrètes mais correspondent à des comptages de population infectée. L'extraction du nombre de personnes infectées à partir des cumuls de nouveaux cas déclarés correspond à un modèle de Markov caché (HMM) non standard. En effet, en raison de la dynamique continue sous-jacente, les nouveaux cas déclarés cumulés sur une période donnée dépendent à la fois du nombre de personnes infectées au début et à la fin de la période. Par conséquent, la distribution des émissions du HMM dépend de l'état actuel et de l'état passé de la chaîne cachée. L'algorithme standard de Baum-Welch [25] ne peut donc pas être appliqué directement et doit être adapté à ce cas spécifique.

Ce chapitre est adapté d'un preprint d'article [56], aborde ces deux défis. Tout d'abord, nous décrivons le modèle LBDPI stochastique et notre schéma d'observation spécifique dans la section 3.2. Ensuite, dans la section 3.3, nous présentons la procédure d'estimation des paramètres, où nous établissons une formule analytique pour les paramètres d'un LBDPI en termes de probabilités de transition à pas discret. Nous discutons également les principales propriétés des estimateurs. Nous adaptons ensuite l'algorithme de Baum-Welch à notre forme spéciale de distribution des émissions pour estimer les probabilités de transition de la chaîne cachée à partir des données de cumuls de nouveaux cas déclarés. La section 3.4 est dédiée à l'étude des propriétés numériques des estimateurs, à la fois sur des données synthétiques et réelles. Elle met également en évidence l'application pratique de ces estimations pour évaluer la probabilité d'atteindre un seuil épidémique des infectés à Mayotte. Enfin, la section 3.5 présente quelques remarques finales.

## 3.2 Modèle et schéma d'observations

Dans cette section, nous présentons d'abord le modèle stochastique que nous utilisons pour quantifier le nombre d'individus infectés dans la population, noté *I* (3.2.1). Nous définissons ensuite notre schéma d'observation et énonçons notre problème d'estimation des paramètres (3.2.2). Enfin, nous expliquons notre stratégie pour obtenir des estimateurs (3.2.3).

### 3.2.1 Description du Modèle

Notre nouveau modèle repose sur plusieurs hypothèses simplificatrices par rapport au modèle initial décrit dans 2.1. Tout d'abord, nous supposons que les taux de natalité et de mortalité naturelle sont négligeables, ce qui signifie que l'effet de la démographie de la population est négligeable. De plus, nous considérons que les individus susceptibles (*S*1, *S*2) et la concentration des bactéries dans l'eau (*W*) sont constantes, puisque les individus infectés sont peu nombreux par rapport aux individus susceptibles.

Dans notre modèle simplifié (Fig.3.1), nous nous concentrons principalement sur la dynamique des individus infectés (*I*), en considérant une seule source d'infectiosité (*W*) sans transmission de la typhoïde par les porteurs chroniques (*C*). Les individus retirés (*R*) ne sont pas pris en compte, car les données disponibles ne représentent que des cumuls des cas isolés qui sont considérés guéris après une période de temps inconnue. Cette approche
nous permet de simplifier notre modèle pour être en mesure d'inférer la dynamique. Comme nous l'avons souligné dans la section 1.2.4, la prévalence de la fièvre typhoïde à Mayotte est caractérisée par des niveaux relativement faibles et fluctuants, et des résurgences de la maladie après des plages sans nouvelles infections. Un modèle stochastique de type processus de comptage est donc bien adapté pour modéliser ces données. Parmi ces modèles, l'un des plus simples et pertinent, avec peu de paramètres est le *processus de naissance et mort linéaire avec immigration*. Un tel processus $(I_t)_{t \ge 0}$ comptabilise le nombre d'individus infectés. Chaque individu contamine un nouvel individu à un taux (de *naissance*) $\lambda > 0$ et devient isolé avec un taux (de *mort*) $\mu > 0$ indépendamment des autres individus, ce qui génère des taux de naissance et mort linéaires pour la population globale. La contamination via l'environnement génère également de nouveaux infectés via un processus de Poisson de taux (d'*immigration*) $\nu > 0$.

<figure><img src="image_27.png" /><figcaption>FIGURE 3.1 – Modèle à un seul compartiment : Infecté.</figcaption></figure>

La loi de *I* est caractérisée par son générateur infinitésimal *Q* donné par :

$$Q(i, i + 1) = \lambda i + \nu \qquad \text{pour } i \ge 0,$$

$$Q(i, i - 1) = \mu i \qquad \text{pour } i \ge 1,$$

$$Q(i, i) = -(\lambda + \mu)i - \nu \qquad \text{pour } i \ge 0.$$

Le semi-groupe de transition défini par $(P(t) = (p_{i,j}(t), i, j \in \mathbb{N}))_{t \ge 0}$, où

$$p_{i,j}(t) = \mathbb{P}(I_t = j | I_0 = i), \quad i, j \in \mathbb{N}, t \ge 0,$$

satisfait l'équation de Kolmogorov progressive suivante : pour $t > 0$ et $i, j \in \mathbb{N}$ avec $j \ge 1$

$$\frac{dp_{i,0}(t)}{dt} = \mu p_{i,1}(t) - \nu p_{i,0}(t), \qquad (3.1)$$

$$\frac{dp_{i,j}(t)}{dt} = (\lambda(j - 1) + \nu) p_{i,j-1}(t) + \mu(j + 1)p_{i,j+1}(t) - ((\lambda + \mu)j + \nu) p_{i,j}(t),$$

et avec la condition initiale $p_{i,i}(0) = 1$ et $p_{i,j}(0) = 0$ pour $j \neq i$.

Dans toute la suite, nous supposons que le processus LBDPI est en régime récurrent positif, c'est-à-dire $\nu > 0$ et $\lambda < \mu$. Notre objectif est d'estimer le triple $(\lambda, \mu, \nu)$.
### 3.2.2 Schéma d'observations

La principale difficulté de notre étude provient du schéma d'observation. En effet, le nombre d'individus infectés à un moment donné $I_t$ n'est pas directement observé dans notre jeu de données. Au lieu de cela, nous avons seulement accès aux cumuls de nouveaux cas isolés par jour. Cette formulation est mathématiquement complexe et originale, mais elle est très courante dans les données de santé publique, notamment pour les maladies à déclaration obligatoire.

Notons par $O_{\mid s,t}$ le cumul des nouveaux déclarés correspondant au processus $I$ sur l'intervalle de temps $[s, t]$. Ainsi, $O_{\mid s,t}$ correspond au nombre de sauts d'amplitude $-1$ de $I$ sur l'intervalle de temps $[s, t]$ : $O_{\mid s,t}$ est le cardinal de l'ensemble $\{t \in [s, t]; I_t - I_{t-} = -1\}$, où $I_{t-} = \lim_{r<t, r \to t} I_r$.

Le processus joint $(I_t, O_{\mid 0,t})$ est encore un processus de Markov de saut pur et sa distribution est caractérisée par les probabilités de transition suivantes : pour $t > 0$

$$p_{i,(j,o)}(t) := \mathbb{P}(I_t = j, O_{\mid 0,t} = o \mid I_0 = i), \quad i, j, o \in \mathbb{N}.$$

Elles satisfont également une équation de Kolmogorov (rétrograde) donnée par le lemme suivant.

**Lemme 3.1.** Pour $t > 0$ et $i, j, o \in \mathbb{N}$, on a

$$\frac{dp_{i,(j,o)}(t)}{dt} = -((\lambda + \mu)i + \nu)p_{i,(j,o)}(t) + (\lambda i + \nu)p_{i+1,(j,o)}(t) + (\mu i)p_{i-1,(j,o-1)}(t), \quad (3.2)$$

pour $i \le j + o$, avec $p_{i,(j,o)}(t) = 0$ quand $i > j + o$ et avec les conditions initiales $p_{i,(j,o)}(0) = \delta_{i=j}\delta_{o=0}$.

*Preuve.* Nous utilisons la caractérisation infinitésimale de la distribution du processus joint.

Pour $h > 0$ et $i \le j + o$, on a

$$\begin{align*}
p_{i,(j,o)}(t+h) &= \mathbb{P}(I_{t+h} = j, O_{\mid 0,t+h]} = o \mid I_0 = i) \\
&= \mathbb{P}(I_h = i, O_{\mid 0,h]} = 0, I_{t+h} = j, O_{\mid h,t+h]} = o \mid I_0 = i) \\
&\quad + \mathbb{P}(I_h = i+1, O_{\mid 0,h]} = 0, I_{t+h} = j, O_{\mid h,t+h]} = o \mid I_0 = i) \\
&\quad + \mathbb{P}(I_h = i-1, O_{\mid 0,h]} = 1, I_{t+h} = j, O_{\mid h,t+h]} = o-1 \mid I_0 = i) + 0(h) \\
&= \mathbb{P}(I_h = i, O_{\mid 0,h]} = 0 \mid I_0 = i) \mathbb{P}(I_{t+h} = j, O_{\mid h,t+h]} = o \mid I_h = i) \\
&\quad + \mathbb{P}(I_h = i+1, O_{\mid 0,h]} = 0 \mid I_0 = i) \mathbb{P}(I_{t+h} = j, O_{\mid h,t+h]} = o \mid I_h = i+1) \\
&\quad + \mathbb{P}(I_h = i-1, O_{\mid 0,h]} = 1 \mid I_0 = i) \mathbb{P}(I_{t+h} = j, O_{\mid h,t+h]} = o-1 \mid I_h = i-1) + 0(h)
\end{align*}$$

Puisque $(I_t)_{t \ge 0}$ est un processus de Markov homogène, on obtient :

$$p_{i,(j,o)}(t+h) = (1 - ((\lambda + \mu)i + \nu)h)p_{i,(j,o)}(t) + (\lambda i + \nu)hp_{i+1,(j,o)}(t) \\ + (\mu i)hp_{i-1,(j,o-1)}(t) + 0(h).$$

Par conséquent

$$\frac{dp_{i,(j,o)}(t)}{dt} = -((\lambda + \mu)i + \nu)p_{i,(j,o)}(t) + (\lambda i + \nu)p_{i+1,(j,o)}(t) + (\mu i)p_{i-1,(j,o-1)}(t). \quad \square$$
Nous définissons $\Delta t$ comme étant un pas de temps fixe, généralement d'une durée d'un jour ou d'une semaine. Notre objectif principal est d'estimer les trois paramètres $(\lambda, \mu, \nu)$ à partir des observations disponibles, qui sont représentées par la suite $(O_{[n\Delta t,(n+1)\Delta t]})_{n\in\mathbb{N}}$.

### 3.2.3 Stratégie d'estimation

Notre problème d'estimation fait partie de la classe des problèmes à information causée. En effet, une suite donnée d'observations $(O_{[n\Delta t,(n+1)\Delta t]})_{n\in\mathbb{N}}$ peut correspondre à plusieurs trajectoires différentes possibles de $I$. Toutefois, notre objectif n'est pas de reconstruire la trajectoire la plus probable, mais simplement d'estimer les paramètres $(\lambda, \mu, \nu)$.

Malheureusement, nous n'avons pas pu obtenir une forme analytique des probabilités de transitions $p_{i,(j,o)}(t)$ à partir de l'équation (3.2). Ainsi, il n'est pas possible de mettre en œuvre une méthode d'estimation directe par maximum de vraisemblance avec ce schéma d'observations. Pour pallier ce problème, nous procédons en deux étapes :

La première étape 3.3.1 consiste à étudier la chaîne de Markov à pas discret $(I_{n\Delta t})_{n\in\mathbb{N}}$. Nous obtenons une formule analytique pour sa matrice de transition en fonction des paramètres inconnus $(\lambda, \mu, \nu)$. Cette formule peut être inversée pour obtenir une expression des paramètres $(\lambda, \mu, \nu)$ en termes de probabilités de transition bien choisies. Nous obtenons ensuite des estimateurs par plug-in à partir de n'importe quel estimateur des probabilités de transition, et nous établissons que la consistance et la normalité asymptotique des estimateurs des probabilités de transition se transfèrent aux estimateurs de $(\lambda, \mu, \nu)$.

La deuxième étape 3.3.2 tient compte du cadre à information manquante en interprétant la chaîne de Markov à temps discret $(I_{n\Delta t}, O_{[n\Delta t,(n+1)\Delta t]})_{n\in\mathbb{N}}$ comme un modèle de Markov caché, où $(I_{n\Delta t})$ est la chaîne cachée et $(O_{[n\Delta t,(n+1)\Delta t]})$ les observations. Cependant, dans cette structure, le processus d'émission est non standard car le nombre cumulé de nouveaux cas déclarés $O_{[n\Delta t,(n+1)\Delta t]}$ dépend de l'ensemble de la trajectoire en temps continu de $I$ sur l'intervalle $[n\Delta t, (n+1)\Delta t]$ et pas seulement de sa valeur finale $I_{(n+1)\Delta t}$. En utilisant la propriété de Markov, on peut obtenir les probabilités d'émission de $O_{[n\Delta t,(n+1)\Delta t]}$ conditionnellement à $I_{n\Delta t}$ et $I_{(n+1)\Delta t}$ seulement. Notre première idée a été d'exécuter la procédure HMM standard du package R HMM sur la chaîne bidimensionnelle $Z_n = (I_{n-1}, I_n)$ pour estimer les probabilités de transition de la chaîne cachée ($Z_n$). Malheureusement, la méthode était numériquement instable et entachée d'erreurs dues aux nombreux zéros dans la matrice de transition de ($Z_n$). Au lieu de cela, nous avons complètement réécrit les probabilités forward et backward dans notre contexte particulier et adapté l'algorithme de Baum-Welsh en conséquence. Cette deuxième étape nous permet d'obtenir des estimations des probabilités de transition cachées de la chaîne de Markov $(I_{n\Delta t})_{n\in\mathbb{N}}$ sur la base des observations $(O_{[n\Delta t,(n+1)\Delta t]})$ uniquement. Nous utilisons ensuite les formules explicites de la première étape pour obtenir des estimateurs du triple $(\lambda, \mu, \nu)$.

## 3.3 Estimation paramétrique

Dans cette section, nous détaillons notre méthodologie d'estimation paramétrique. Tout d'abord, nous établissons une formule analytique pour les paramètres d'un LBDPI en termes de probabilités de transition dans la partie 3.3.1. Nous adaptons ensuite l'algorithme de Baum-Welsh à notre forme particulière de loi d'émission pour estimer les
probabilités de transition de la chaîne cachée à partir des cumuls de nouveaux cas dé-
clarés dans **3.3.2**. Enfin, nous obtenons des estimations des paramètres originaux par
plug-in.

**3.3.1 Expression analytique des paramètres**

La première étape de notre procédure d'estimation consiste à établir une formule ana-
lytique pour les paramètres (λ, μ, ν) en tant que fonctions des probabilités de transition
en pas discret de la chaîne de Markov (I<sub>nΔt</sub>)<sub>n∈N</sub>.

Nous présentons d'abord les lemmes suivants qui seront utiles pour établir une solution
analytique de l'équation de Kolmogorov (3.1).

**Lemme 3.2.** Soit $G_i(s,t)$ une fonction différentiable sur $\mathbb{R} \times [0, +\infty[$, pour tout $i \in \mathbb{N}$, la solution de l'équation aux dérivées partielles linéaire

$$
\left\{
\begin{array}{l}
\displaystyle \frac{\partial G_i(s,t)}{\partial t} = (\mu - \lambda s)(1-s) \frac{\partial G_i(s,t)}{\partial s} + \nu(s-1)G_i(s,t), \\
G_i(s,0) = s^i,
\end{array}
\right.
\qquad (3.3)
$$

est donnée par :

$$
G_i(s,t) = \begin{cases}
\frac{(\mu - \lambda)^{\frac{\nu}{\lambda}} \left( s(\mu e^{(\lambda - \mu)t} - \lambda) - \mu (e^{(\lambda - \mu)t} - 1) \right)^i}{\left( \mu - \lambda e^{(\lambda - \mu)t} + \lambda s (e^{(\lambda - \mu)t} - 1) \right)^{\frac{\nu}{\lambda} + i}} & \lambda \neq \mu, \\
\left( \frac{\lambda t + s(1 - \lambda t)}{1 + \lambda t - \lambda ts} \right)^i (1 + \lambda t - \lambda ts)^{-\frac{\nu}{\lambda}} & \lambda = \mu.
\end{cases}
$$

Preuve. Pour résoudre l'équation (3.3), nous utilisons la méthode des caractéristiques
(voir par exemple ([57], chapitre 1)).

Étape 1. Détermination des caractéristiques. Considérons l'équation caractéristique suivante :

$$
\frac{dt}{ds} = - \frac{1}{(\mu - \lambda s)(1 - s)} = \begin{cases} \frac{1}{\lambda - \mu} \left( \frac{1}{1 - s} - \frac{\lambda}{\mu - \lambda s} \right) & \text{if } \lambda \neq \mu \\ -\frac{1}{\lambda(1 - s)^2} & \text{if } \lambda = \mu. \end{cases}
$$

En intégrant par rapport à s, nous avons

$$
t = k - \begin{cases} \frac{1}{\lambda - \mu} \ln \left( \frac{1-s}{\mu - \lambda s} \right) & \text{si } \lambda \neq \mu \\ \frac{1}{\lambda(1-s)} & \text{si } \lambda = \mu. \end{cases} \qquad (3.4)
$$

L'équation (3.4) définit une famille de caractéristiques paramétrées par la constante d'intégration *k*.

Etape 2. Tout au long des caractéristiques de (3.4), $G_i$ satisfait

$$
\frac{dG_i}{ds} = \frac{\nu G_i}{\mu - \lambda s}.
$$
L'intégration par rapport à s donne

$$
G_i = C(k) (\mu - \lambda s)^{-\frac{\nu}{\lambda}},
$$

où $C(k)$ est une autre constante d'intégration, fonction du paramètre $k$. En utilisant
l'équation (3.4), nous obtenons

$$
G_i(s,t) = \begin{cases} C \left( t + \frac{1}{\lambda - \mu} \ln \left( \frac{1-s}{\mu - \lambda s} \right) \right) (\mu - \lambda s)^{-\frac{\nu}{\lambda}} & \text{si } \lambda \neq \mu \\ C \left( t + \frac{1}{\lambda (1-s)} \right) (\lambda (1-s))^{-\frac{\nu}{\lambda}} & \text{si } \lambda = \mu. \end{cases} \tag{3.5}
$$

Nous utilisons maintenant la condition initiale $G_i(s, 0) = s^i$ pour calculer $C(k)$. En fixant $t = 0$ dans l'équation (3.5), nous trouvons

$$
s^i = \begin{cases} C \left( \frac{1}{\lambda - \mu} \ln \left( \frac{1-s}{\mu - \lambda s} \right) \right) (\mu - \lambda s)^{-\frac{\nu}{\lambda}} & \text{si } \lambda \neq \mu \\ C \left( \frac{1}{\lambda (1-s)} \right) (\lambda (1-s))^{-\frac{\nu}{\lambda}} & \text{si } \lambda = \mu. \end{cases} \tag{3.6}
$$

Posons

$$
w = \begin{cases} \frac{1}{\lambda - \mu} \ln \left( \frac{1-s}{\mu - \lambda s} \right) & \text{si } \lambda \neq \mu \\ \frac{1}{\lambda (1-s)} & \text{si } \lambda = \mu. \end{cases}
$$

En mettant

$$
s = \begin{cases} \frac{e^{(\lambda-\mu)w}\mu - 1}{e^{(\lambda-\mu)w}\lambda - 1} & \text{si } \lambda \neq \mu \\ 1 - \frac{1}{\lambda w} & \text{si } \lambda = \mu, \end{cases}
$$

dans (3.6) pour obtenir la fonction $C(w)$ requise,

$$
C(w) = \begin{cases} \left( \frac{e^{(\lambda-\mu)w}\mu - 1}{e^{(\lambda-\mu)w}\lambda - 1} \right)^i \left( \mu - \lambda \frac{e^{(\lambda-\mu)w}\mu - 1}{e^{(\lambda-\mu)w}\lambda - 1} \right)^{\frac{\nu}{\lambda}} & \text{si } \lambda \neq \mu \\ \left(1 - \frac{1}{\lambda w}\right)^i w^{-\frac{\nu}{\lambda}} & \text{si } \lambda = \mu. \end{cases} \quad (3.7)
$$

En combinant les équations (3.5) et (3.7), nous pouvons maintenant écrire la forme com-
plète de $G_i(s,t)$, pour $\lambda \neq \mu$,

$$
G_i(s,t) = \frac{(\mu - \lambda)^{\frac{\nu}{\lambda}} \left( s(\mu e^{(\lambda-\mu)t} - \lambda) - \mu (e^{(\lambda-\mu)t} - 1) \right)^i}{\left( \mu - \lambda e^{(\lambda-\mu)t} + \lambda s (e^{(\lambda-\mu)t} - 1) \right)^{\frac{\nu}{\lambda} + i}},
$$

et pour $\lambda = \mu$

$$
G_i(s,t) = \left( \frac{\lambda t + s (1 - \lambda t)}{1 + \lambda t - \lambda t s} \right)^i (1 + \lambda t - \lambda t s)^{-\frac{\nu}{\lambda}},
$$

□
**Lemme 3.3.** Pour $i, j, n \ge 0$, on a

$$ \sum_{k=0}^{n} \binom{i+k}{i} \binom{j+n-k}{j} = \binom{i+j+n+1}{n}. $$

*Preuve.* Pour tout $x \in ]-1, 1[$ et $\alpha \in \mathbb{R}$, on a

$$ (1-x)^{-\alpha-1} = \sum_{k=0}^{\infty} \binom{\alpha+k}{\alpha} x^k. $$

Ainsi, d'une part, nous avons :

$$ (1-x)^{-(i+j+1)-1} = \sum_{n=0}^{\infty} \binom{i+j+n+1}{n} x^n. \qquad (3.8) $$

et d'autre part :

$$ (1-x)^{-(i+j+1)-1} = (1-x)^{-i-1}(1-x)^{-j-1}, $$

où

$$ (1-x)^{-i-1} = \sum_{k=0}^{\infty} \binom{i+k}{i} x^k, $$

$$ (1-x)^{-j-1} = \sum_{k=0}^{\infty} \binom{j+k}{j} x^k. $$

En utilisant le produit de Cauchy de deux séries infinies, on trouve :

$$ (1-x)^{-(i+j+1)-1} = \sum_{n=0}^{\infty} \left( \sum_{k=0}^{n} \binom{i+k}{i} \binom{j+n-k}{j} \right) x^n. \qquad (3.9) $$

En identifiant les coefficients de (3.8) et (3.9), nous obtenons le résultat souhaité. $\square$

Nous présentons maintenant l'expression explicite des probabilités de transition du processus $I_t$. Bien que nous ayons supposé que $\lambda < \mu$, nous donnons des formules générales qui incluent également le cas où $\lambda \ge \mu$.

**Théorème 3.1.** Pour $t \ge 0$, les probabilités de transition de la chaîne de Markov à temps continu $I_t$ sont données pour $i, j \in \mathbb{N}$ par :

$$ p_{i,j}(t) = q(t)^r \sum_{l=0}^{\min(i,j)} \binom{i}{l} \binom{r+i+j-l-1}{j-l} \left(\frac{\mu}{\lambda}\right)^{i-l} (1-q(t))^{i+j-l} \left(\frac{1-\left(\frac{\mu}{\lambda}+1\right)(1-q(t))}{1-q(t)}\right)^l, \qquad (3.10) $$

où

$$ r = \frac{\nu}{\lambda} \quad \text{et} \quad q(t) = \begin{cases} \frac{\mu - \lambda}{\mu - \lambda e^{(\lambda - \mu)t}} & \lambda \neq \mu, \\ \frac{1}{1 + \lambda t} & \lambda = \mu. \end{cases} $$

*Preuve.* La première étape consiste à transformer l'équation de Kolmogorov **3.1** en une équation différentielle partielle linéaire standard du premier ordre équivalente. Dans ce
but, nous introduisons la fonction génératrice des moments du processus $I_t$: pour tout $(s,t) \in \times[0, +\infty]^2$, et $i \in \mathbb{N}$

$$G_i(s, t) = \sum_{j=0}^{+\infty} s^j p_{i,j}(t).$$

En utilisant l'équation de Kolmogorov progressive et après quelques réarrangements, nous trouvons :

$$\begin{aligned} \frac{\partial G_i(s, t)}{\partial t} &= \sum_{j=0}^{\infty} s^j \frac{\mathrm{d}p_{i,j}(t)}{\mathrm{d}t} \\ &= (\lambda s^2 - (\lambda + \mu)s + \mu) \sum_{j=1}^{\infty} js^{j-1}p_{i,j}(t) + \nu(s-1) \sum_{j=0}^{\infty} s^j p_{i,j}(t). \end{aligned}$$

Ainsi, la fonction génératrice des moments satisfait l'équation différentielle partielle linéaire du premier ordre suivante :

$$\frac{\partial G_i(s, t)}{\partial t} = (\mu - \lambda s)(1 - s) \frac{\partial G_i(s, t)}{\partial s} + \nu(s - 1)G_i(s, t),$$

avec la condition initiale $G_i(s, 0) = s^i$.

Maintenant, en utilisant le lemme 3.2, nous avons : pour $\lambda \neq \mu$,

$$G_i(s,t) = \frac{(\mu - \lambda)^{\frac{\nu}{\lambda}} \left( s(\mu e^{(\lambda - \mu)t} - \lambda) - \mu (e^{(\lambda - \mu)t} - 1) \right)^i}{\left( \mu - \lambda e^{(\lambda - \mu)t} + \lambda s(e^{(\lambda - \mu)t} - 1) \right)^{\frac{\nu}{\lambda} + i}}$$

et pour $\lambda = \mu$

$$G_i(s,t) = \left( \frac{\lambda t + s(1 - \lambda t)}{1 + \lambda t - \lambda ts} \right)^i (1 + \lambda t - \lambda ts)^{-\frac{\nu}{\lambda}}.$$

Le cas $i = 0$ donne la fonction génératrice suivante pour les probabilités de transition $p_{0,j}(t)$, $t \ge 0$ :

$$G_0(s,t) = \begin{cases} \left(\frac{\mu-\lambda}{\mu-\lambda e^{(\lambda-\mu)t}}\right)^{\frac{\nu}{\lambda}} \left(1-\frac{\lambda(1-e^{(\lambda-\mu)t})}{\mu-\lambda e^{(\lambda-\mu)t}}\right)^{-\frac{\nu}{\lambda}} & \lambda \neq \mu, \\ \left(\frac{1}{1+\lambda t}\right)^{\frac{\nu}{\lambda}} \left(1-\frac{\lambda t}{1+\lambda t}\right)^{-\frac{\nu}{\lambda}} & \lambda = \mu. \end{cases}$$

Il s'agit de la fonction génératrice de la loi binomiale négative généralisée avec des paramètres

$$r = \frac{\nu}{\lambda} \quad \text{et} \quad q(t) = \begin{cases} \frac{\mu - \lambda}{\mu - \lambda e^{(\lambda - \mu)t}} & \lambda \neq \mu, \\ \frac{1}{1 + \lambda t} & \lambda = \mu. \end{cases}$$

Par conséquent, pour $j \in \mathbb{N}$ et $t \ge 0$, on obtient

$$p_{0,j}(t) = \binom{j+r-1}{r-1} q(t)^r (1-q(t))^j. \qquad (3.11)$$
Ensuite, la forme générale de la fonction $G_i$ pour $i \neq 0$ peut être écrite comme suit :

$$G_i(s, t) = G_0(s, t)N_i(s, t), \qquad (3.12)$$

où

$$N_i(s, t) = \begin{cases} \left( \frac{\mu - \lambda s - \mu(1-s)e^{(\lambda-\mu)t}}{\mu - \lambda s - \lambda(1-s)e^{(\lambda-\mu)t}} \right)^i & \lambda \neq \mu, \\ \left( \frac{1 + (\lambda t - 1)(1-s)}{1 + \lambda t(1-s)} \right)^i & \lambda = \mu, \end{cases}$$

est la fonction génératrice des moments d'un processus de naissance et mort linéaire sans immigration. En principe, l'équation (3.12) reflète l'indépendance entre le processus de contamination par l'environnement ($G_0$) et l'auto-contamination de la population ($N$). Le processus de naissance et mort linéaire sans immigration possède la propriété de branchement qui se traduit par la forme puissance de sa fonction génératrice des moments, voir par exemple [58] ou [59]. On peut facilement retrouver ses probabilités de transition, notées $p_{i,j}^{(N)}(t)$, $t \ge 0$ à partir du développement en série de $N_i(s, t)$ :

$$p_{i,j}^{(N)}(t) = \sum_{l=0}^{\min(i,j)} \binom{i}{l} \binom{i+j-l-1}{i-1} \alpha(t)^{i-l} \beta(t)^{j-l} (1-\alpha(t)-\beta(t))^l, \quad (3.13)$$

où $\alpha(t) = \frac{\mu}{\lambda}(1 - q(t))$ et $\beta(t) = 1 - q(t)$.

Ainsi, les probabilités $p_{i,j}(t)$ sont données par la convolution suivante de l'équation (3.12) :

$$p_{i,j}(t) = \sum_{k=0}^{j} p_{0,k}(t) p_{i,j-k}^{(N)}(t),$$

et les équations (3.11), (3.13) donnent finalement

$$p_{i,j}(t) = q(t)^r \sum_{l=0}^{\min(i,j)} \binom{i}{l} \left(\frac{\mu}{\lambda}\right)^{i-l} (1-q(t))^{i+j-l} \left(\frac{1-\left(\frac{\mu}{\lambda}+1\right)(1-q(t))}{1-q(t)}\right)^l \\ \times \sum_{k=0}^{j-l} \binom{r-1+k}{r-1} \binom{i-1+j-l-k}{i-1},$$

En utilisant le lemme 3.3, nous obtenons le résultat souhaité

$$p_{i,j}(t) = q(t)^r \sum_{l=0}^{\min(i,j)} \binom{i}{l} \binom{r+i+j-l-1}{j-l} \left(\frac{\mu}{\lambda}\right)^{i-l} (1-q(t))^{i+j-l} \\ \times \left(\frac{1-\left(\frac{\mu}{\lambda}+1\right)(1-q(t))}{1-q(t)}\right)^l.$$

□

**Corollaire 3.1.** En dérivant la fonction génératrice des moments par rapport à $s$ et en prenant $s = 1$, nous pouvons retrouver l'expression bien connue de l'espérance du LBDPI, $m(t) = \mathbb{E}_i[I_t]$ :

$$m(t) = \mathbb{E}[I_t | I_0 = i] = \begin{cases} \frac{\nu}{\lambda-\mu}(e^{(\lambda-\mu)t}-1) + ie^{(\lambda-\mu)t} & \lambda \neq \mu, \\ \nu t + i & \lambda = \mu. \end{cases}$$
**Remarque 3.1.** Le résultat ci-dessus sur l'espérance est cohérent avec celui obtenu dans [60] en utilisant une autre méthode qui consiste à dériver une équation pour $m(t+h)$, en conditionnant sur $I_t$, puis en résolvant directement une équation différentielle sur $m(t)$. De plus, nous constatons que lorsque $\lambda \ge \mu$ et $\nu > 0$, L'espérance tend vers l'infini. Lorsque $\lambda > \mu$, l'espérance augmente exponentiellement vite. Lorsque $\lambda = \mu$, la moyenne conditionnelle de $I_t$ augmente linéairement. Cependant, dans le cas où $\lambda < \mu$ et $\nu > 0$, l'espérance se stabilise et converge vers une constante $\nu/(\mu - \lambda)$.

Il est bien connu que la matrice de transition de la chaîne $(I_{n\Delta t})_{n \in \mathbb{N}}$ est $P(\Delta t)$, nous avons donc une forme analytique pour les probabilités de transition en temps discret. Nous voulons maintenant inverser l'équation (3.10) afin d'obtenir une expression des paramètres d'intérêt $(\lambda, \mu, \nu)$ en fonction des probabilités de transition. Comme il n'y a que trois paramètres inconnus, il devrait suffire de sélectionner trois probabilités de transition spécifiques. Au vu des données réelles décrites dans la section 3.4.2, nous choisissons d'utiliser $p_{0,0}(\Delta t)$, $p_{0,1}(\Delta t)$ et $p_{1,0}(\Delta t)$ car elles correspondent aux transitions les plus fréquemment observées et ont les expressions analytiques les plus simples. Pour simplifier les notations, dans la suite du texte, nous noterons $p_{0,0}$, $p_{0,1}$ et $p_{1,0}$ au lieu de $p_{0,0}(\Delta t)$, $p_{0,1}(\Delta t)$ et $p_{1,0}(\Delta t)$, et de manière générale $p_{i,j}$ au lieu de $p_{i,j}(\Delta t)$.

Nous présentons maintenant le théorème suivant qui donne les estimateurs ponctuels de $\lambda, \mu$ et $\nu$ du processus $I_t$ dans le cas où $\lambda \neq \mu$.

**Théorème 3.2.** Les paramètres $\lambda, \mu$ et $\nu$ sont donnés par les relations suivantes $(\lambda, \mu, \nu) = g(p_{0,0}, p_{0,1}, p_{1,0})$ avec

$$ 
\begin{align*}
\lambda &= g_1(p_{0,0}, p_{0,1}, p_{1,0}) = \frac{\ln\left(\frac{u}{q}\right) (q-1)}{\Delta t (q-u)} \\
\mu &= g_2(p_{0,0}, p_{0,1}, p_{1,0}) = \frac{\ln\left(\frac{u}{q}\right) (u-1)}{\Delta t (q-u)} \\
\nu &= g_3(p_{0,0}, p_{0,1}, p_{1,0}) = \frac{\ln(p_{0,0}) \ln\left(\frac{u}{q}\right) (q-1)}{\ln(q) \Delta t (q-u)},
\end{align*}
$$

où

$$ 
\begin{gather*}
q = \frac{p_{0,1}}{p_{0,0} \ln(p_{0,0})} W \left( \frac{p_{0,0} \binom{\frac{p_{0,0}}{p_{0,1}}+1}{p_{0,1}} \ln(p_{0,0})}{p_{0,1}} \right), \\
u = 1 - \frac{p_{1,0}}{p_{0,0}},
\end{gather*}
$$

et W est la fonction W de Lambert (parfois aussi appelée la fonction Oméga [61, 62]).

*Preuve.* Nous utilisons l'équation (3.10) pour $(i, j) \in \{(0, 0), (0, 1), (1, 0)\}$ et $t = \Delta t$ pour obtenir

$$ 
\begin{align*}
p_{0,0} &= q^r \\
p_{0,1} &= q^r r (1-q) \\
p_{1,0} &= q^r \left(\frac{\mu}{\lambda}\right) (1-q),
\end{align*}
$$
avec $r = \frac{\nu}{\lambda}$ et $q = q(\Delta t) = \frac{\mu - \lambda}{\mu - \lambda e^{(\lambda - \mu)\Delta t}}$.

Les deux premières équations ci-dessus sont équivalentes au système

$$ \left\{ \begin{array}{l} \displaystyle \frac{\ln(p_{0,0})}{p_{0,0}} = r \ln(q) \\[1em] \displaystyle \frac{p_{0,1}}{p_{0,0}} = r(1-q), \end{array} \right. \qquad (3.14) $$

ce qui donne

$$ \frac{p_{0,0}}{p_{0,1}} \ln(p_{0,0}) q e^{\frac{p_{0,0}}{p_{0,1}} \ln(p_{0,0}) q} = \frac{\left(\frac{p_{0,0}}{p_{0,1}} + 1\right)}{p_{0,1}} \ln(p_{0,0}), $$

on a donc,

$$ q = \frac{p_{0,1}}{p_{0,0} \ln(p_{0,0})} W \left( \frac{\left( \frac{p_{0,0}}{p_{0,1}} + 1 \right) \ln(p_{0,0})}{p_{0,1}} \right), \qquad (3.15) $$

où W est la fonction de Lambert. A partir des équations (3.14) et (3.15), nous pouvons facilement trouver que

$$ r = \frac{\ln(p_{0,0})}{\ln(q)}. $$

Définissons maintenant

$$ u = q e^{(\lambda - \mu) \Delta t} = 1 - \frac{p_{1,0}}{p_{0,0}}. $$

L'équation (3.10) pour $(i, j) \in \{(0,0), (0,1), (1,0)\}$ peut maintenant être réécrite comme suit

$$ \left\{ \begin{array}{l} r = \frac{\nu}{\lambda} \\ q = \frac{\mu - \lambda}{\mu - \lambda e^{(\lambda - \mu)\Delta t}} \\ u = \frac{(\mu - \lambda) e^{(\lambda - \mu)\Delta t}}{\mu - \lambda e^{(\lambda - \mu)\Delta t}}. \end{array} \right. \qquad (3.16) $$

Étudions maintenant l'identifiabilité du système **3.16**.

Soit $\lambda, \lambda', \mu, \mu', \nu, \nu' > 0$ tel que $\lambda \neq \mu$ et $\lambda' \neq \mu'$. Supposons que

$$ \left\{ \begin{array}{l} \displaystyle \frac{\nu}{\lambda} \quad = \frac{\nu'}{\lambda'} \\[1em] \displaystyle \frac{\mu - \lambda}{\mu - \lambda e^{(\lambda - \mu)\Delta t}} \quad = \frac{\mu' - \lambda'}{\mu' - \lambda' e^{(\lambda' - \mu')\Delta t}} \\[1em] \displaystyle \frac{(\mu - \lambda) e^{(\lambda - \mu)\Delta t}}{\mu - \lambda e^{(\lambda - \mu)\Delta t}} \quad = \frac{(\mu' - \lambda') e^{(\lambda' - \mu')\Delta t}}{\mu' - \lambda' e^{(\lambda' - \mu')\Delta t}}. \end{array} \right. \qquad (3.17) $$

En prenant le rapport des deux dernières équations, on obtient

$$ \mu = \lambda - \lambda' + \mu'. $$

En substituant l'expression de $\mu$ dans la deuxième équation de (3.17), on obtient

$$ \lambda (1 - e^{(\lambda' - \mu')\Delta t}) = \lambda' (1 - e^{(\lambda' - \mu')\Delta t}), $$
puisque $\lambda' \neq \mu'$, nous déduisons que $\lambda = \lambda'$, donc $\mu = \mu'$.

A partir de la première équation de (3.17), nous obtenons également que $\nu = \nu'$.

Par conséquent, le système d'équations 3.16 est identifiable.

Enfin, la résolution du système (3.16) donne lieu aux expressions de $\lambda, \mu$ et $\nu$ en fonction des probabilités de transition $p_{0,0}, p_{0,1}$ et $p_{1,0}$ données dans le théorème 3.2. $\square$

**Corollaire 3.2.** Soit $\hat{p}_{0,0}^n, \hat{p}_{0,1}^n$ et $\hat{p}_{1,0}^n$ des estimateurs consistants et asymptotiquement normaux de $p_{0,0}, p_{0,1}$ et $p_{1,0}$ avec une matrice de variance-covariance asymptotique $\Sigma'$. Alors, les estimateurs $(\hat{\lambda}^n, \hat{\mu}^n, \hat{\nu}^n) = g(\hat{p}_{0,0}^n, \hat{p}_{0,1}^n, \hat{p}_{1,0}^n)$ obtenus par plug-in à partir des formules du théorème 3.2 sont également consistants et asymptotiquement normaux avec une matrice de variance-covariance asymptotique donnée par $\Sigma = Dg\Sigma' Dg^T$, où $Dg$ est la matrice jacobienne de l'application $g$ et $Dg^T$ est sa transposée.

*Preuve.* Les $\hat{p}_{0,0}^n, \hat{p}_{0,1}^n$ et $\hat{p}_{1,0}^n$ sont des estimateurs consistants de $p_{0,0}, p_{0,1}$ et $p_{1,0}$, c'est-à-dire, quand $n \to \infty$

$$ 
\begin{array}{ccc}
\hat{p}_{0,0}^n & \xrightarrow{p.s.} & p_{0,0} \\
\hat{p}_{0,1}^n & \xrightarrow{p.s.} & p_{0,1} \\
\hat{p}_{1,0}^n & \xrightarrow{p.s.} & p_{1,0}.
\end{array}
 $$

Étant donné que $g$ est une fonction continue (elle est définie par la composition et le produit de fonctions continues), nous déduisons que

$$ 
\begin{array}{ccc}
\hat{\lambda}^n & \xrightarrow{p.s.} & \lambda \\
\hat{\mu}^n & \xrightarrow{p.s.} & \mu \\
\hat{\nu}^n & \xrightarrow{p.s.} & \nu.
\end{array}
 $$

Maintenant, pour la propriété de la normalité asymptotique, les $\hat{p}_{0,0}^n, \hat{p}_{0,1}^n$ et $\hat{p}_{1,0}^n$ sont des estimateurs asymptotiquement normaux de $p_{0,0}, p_{0,1}$ et $p_{1,0}$ avec une matrice de variance-covariance asymptotique $\Sigma'$, c'est-à-dire, quand $n \to \infty$

$$ \sqrt{n} \begin{pmatrix} \hat{p}_{0,0}^n - p_{0,0} \\ \hat{p}_{0,1}^n - p_{0,1} \\ \hat{p}_{1,0}^n - p_{1,0} \end{pmatrix} \xrightarrow{\mathcal{L}} \mathcal{N}_3(0, \Sigma'). \quad (3.18) $$

Puisque $g$ est une fonction différentiable, la $\delta$-méthode à 3.18, donne, quand $n \to \infty$

$$ \sqrt{n} \left( g \left( \hat{p}_{0,0}^n, \hat{p}_{0,1}^n, \hat{p}_{1,0}^n \right) - g \left( p_{0,0}, p_{0,1}, p_{1,0} \right) \right) \xrightarrow{\mathcal{L}} \mathcal{N}_3(0, Dg\Sigma' Dg^T), $$

i.e :

$$ \sqrt{n} \begin{pmatrix} \hat{\lambda}^n - \lambda \\ \hat{\mu}^n - \mu \\ \hat{\nu}^n - \nu \end{pmatrix} \xrightarrow{\mathcal{L}} \mathcal{N}_3(0, \Sigma). $$

Le calcul de la matrice jacobienne de $g$ est donné dans l'annexe B. $\square$

**Application :** Considérons par exemple $\hat{p}_{i,j}^n$, les estimateurs du maximum de vraisemblance de $p_{i,j}$ pour les observations $(I_{k\Delta t}, 0 \le k \le n)$ et pour $(i, j) \in \{(0,0), (0,1), (1,0)\}$. Ces estimateurs correspondent au rapport entre le nombre de transitions de $i$ à $j$ et le nombre total de transitions de $i$ dans la suite $(I_{k\Delta t}, 0 \le k \le n)$, voir 3.4.1. Sous l'hypothèse de récurrence positive ($\lambda < \mu$ et $\nu > 0$) de la chaîne $(I_{k\Delta t})$, ils sont consistants et
asymptotiquement normaux avec une matrice de variance-covariance asymptotique [63, **Théorème 2.1**] donnée par

$$ \Sigma' = \begin{pmatrix} \frac{p_{0,0}(1-p_{0,0})}{\pi_0} & -\frac{p_{0,0}p_{0,1}}{\pi_0} & 0 \\ -\frac{p_{0,0}p_{0,1}}{\pi_0} & \frac{p_{0,1}(1-p_{0,1})}{\pi_0} & 0 \\ 0 & 0 & \frac{p_{1,0}(1-p_{1,0})}{\pi_1} \end{pmatrix}, $$

où $\pi$ est la loi invariante donnée par l'équation **1.3**. Plus précisément, on a, quand $n \to \infty$,

$$ \sqrt{n} \left( \begin{pmatrix} \hat{p}_{0,0}^n - p_{0,0} \\ \hat{p}_{0,1}^n - p_{0,1} \\ \hat{p}_{1,0}^n - p_{1,0} \end{pmatrix} \right) \xrightarrow{\mathcal{L}} \mathcal{N}_3(0, \Sigma'), $$

Ainsi, ces estimateurs satisfont aux hypothèses du corollaire **3.2**. Cependant, ils ne peuvent pas être utilisées dans notre contexte car les $(I_{k\Delta t}, 0 \le k \le n)$ sont cachés. Nous nous tournons donc vers une autre procédure d'estimation pour $p_{0,0}$, $p_{0,1}$ et $p_{1,0}$ via un modèle de Markov caché.

### 3.3.2 Modèle de Markov caché

Nous passons maintenant à l'estimation des probabilités de transition $p_{0,0}$, $p_{0,1}$ et $p_{1,0}$ de la chaîne cachée $(I_{n\Delta t})_{n \in \mathbb{N}}$ à partir des observations disponibles $(O_{\lfloor n\Delta t, (n+1)\Delta t \rfloor})_{n \in \mathbb{N}}$. Pour simplifier les notations, nous notons dans cette section $I_n$ au lieu de $I_{n\Delta t}$ et $O_n$ au lieu de $O_{\lfloor n-1 \rfloor \Delta t, n\Delta t \rfloor}$.

L'objectif de cette section est d'adapter les algorithmes HMM standard pour obtenir des estimations des probabilités de transition de la chaîne cachée $(I_n)$ à partir des observations $(O_n)$. Cependant, comme expliqué dans la section précédente **3.2.3**, $O_n$ dépend à la fois de $I_{n-1}$ et de $I_n$, et conditionnellement à ces deux états $I_{n-1}$ et $I_n$, les observations $O_n$ sont mutuellement indépendantes. Ainsi, le schéma d'émission n'est pas standard (Fig.3.2).

<figure><img src="image_28.png" /><figcaption>FIGURE 3.2 – Schéma d'émission du HMM.</figcaption></figure>

Pour contourner la difficulté découlant du schéma d'émission non standard, nous utilisons la chaîne de Markov bidimensionnelle correspondant à des paires d'états cachés.
Notons par $(Z_n)_{n \in \mathbb{N}^*} = (I_{n-1}, I_n)_{n \in \mathbb{N}^*}$ la chaîne de Markov non observée à deux dimensions consécutifs à valeurs dans l'espace d'état caché $\mathbb{N}^2$. Il est alors facile de prouver que $(Z_n, O_n)$ est un HMM standard avec les caractéristiques suivante.

**Lemme 3.4.** Le processus $(Z_n, O_n)_{n \in \mathbb{N}^*}$ est un modèle de Markov caché dont les caractéristiques sont données par le triple $M = (Q, \psi, \rho)$, où

i) La matrice de probabilité de transition $Q$ de la chaîne cachée $(Z_n)$ est

$$Q_{(i,j),(i',j')} = \mathbb{P}(Z_{n+1} = (i',j') | Z_n = (i,j)) = p_{i',j'} \delta_{i'=j},$$

pour $(i, j), (i', j') \in \mathbb{N}^2$.

ii) La probabilité d'émission $\psi$ du processus $O$ étant donné le processus $Z$ est

$$\psi_{(i,j)}(o) = \mathbb{P}(O_n = o | Z_n = (i,j)) = \frac{p_{i,(j,o)}}{p_{i,j}},$$

pour $i, j, o \in \mathbb{N}$, où $p_{i,(j,o)} = p_{i,(j,o)}(\Delta t)$ est caractérisée par l'équation 3.2.

iii) La loi initiale $\rho$ du processus d'état $Z$ est

$$\rho_{i,j} = \mathbb{P}(Z_1 = (i,j)) = p_{i,j} \pi_i,$$

pour $i, j \in \mathbb{N}$, où $\pi$ est la loi invariante de $I$ donnée par l'équation 1.3.

*Preuve.* i) La probabilité de transition $Q$ : soient $(i, j), (i', j') \in \mathbb{N}^2$, nous avons

$$\begin{aligned} Q_{(i,j),(i',j')} &= \mathbb{P}(Z_{n+1} = (i',j') | Z_n = (i,j)) \\ &= \frac{\mathbb{P}(I_n = i', I_{n+1} = j', I_{n-1} = i, I_n = j)}{\mathbb{P}(I_{n-1} = i, I_n = j)}. \end{aligned}$$

En utilisant la propriété de Markov de la chaîne $(I_n)$, nous obtenons

$$\begin{aligned} Q_{(i,j),(i',j')} &= \frac{\mathbb{P}(I_{n+1} = j' | I_n = i') \mathbb{P}(I_n = i', I_{n-1} = i)}{\mathbb{P}(I_n = j, I_{n-1} = i)} \delta_{i'=j} \\ &= p_{i',j'} \delta_{i'=j}. \end{aligned}$$

ii) La probabilité d'émission $\psi$ : par la définition de la probabilité conditionnelle.

iii) La loi initiale $\rho$ : soit $i, j \in \mathbb{N}$

$$\begin{aligned} \rho_{i,j} &= \mathbb{P}(Z_1 = (i,j)) \\ &= \mathbb{P}(I_1 = j | I_0 = i) \mathbb{P}(I_0 = i) \\ &= p_{i,j} \pi_i. \end{aligned}$$

□

L'algorithme standard pour estimer les paramètres les plus probables $M = (Q, \psi, \rho)$ du HMM compte tenu des observations est l'algorithme forward-backward ou algorithme de Baum-Welch, un cas particulier de l'algorithme d'espérance-maximisation (EM)[25]. Il s'agit d'un algorithme itératif. À l'étape $n$, étant donné les paramètres actuels $M^n = (Q^n, \psi^n, \rho^n)$, la vraisemblance des observations $\mathbb{P}(O_1 = o_1, \dots, O_T = o_T | M^n)$ est maximisée et le modèle $M^n$ est mis à jour en $M^{n+1}$ en conséquence. Le principal avantage de l'algorithme de Baum-Welch est que la vraisemblance peut être maximisée explicitement.
Malheureusement, l'exécution de la procédure HMM standard du package R HMM sur $(Z_n, O_n)$ produit une instabilité numérique et des erreurs dues aux nombreux zéros dans la matrice de transition $Q$ de $Z$. Nous avons donc réécrit les fonctions principales en adaptant les formules récursives à la forme particulière de $Q$ afin d'obtenir des estimateurs explicites des $p_{i,j}$.

Nous commençons par définir la probabilité forward $\alpha_{(i,j)}(t)$ comme la probabilité conjointe d'observer les $t$ premiers vecteurs et d'être dans l'état $(i, j)$ à l'instant $t$. Cette probabilité peut être calculée par la récurrence suivante : pour tout $(i, j) \in \mathbb{N}^2$ et $2 \le t \le T$

$$ \left\{ \begin{aligned} \alpha_{(i,j)}(1) &= p_{i,j}\pi_t\psi_{(i,j)}(o_1), \\ \alpha_{(i,j)}(t) &= \psi_{(i,j)}(o_t)p_{i,j}\sum_{i' \in \mathbb{N}} \alpha_{(i',i)}(t-1). \end{aligned} \right. \qquad (3.19) $$

Ensuite, nous définissons la probabilité backward $\beta_{(i,j)}(t)$ comme la probabilité conditionnelle d'observer les vecteurs après le temps $t$, étant donné que l'état au temps $t$ est $(i, j)$. La probabilité backward peut être évaluée à l'aide de la formule du récurrence : pour tout $(i, j) \in \mathbb{N}^2$ et $2 \le t \le T$

$$ \left\{ \begin{aligned} \beta_{(i,j)}(T) &= 1, \\ \beta_{(i,j)}(t-1) &= \sum_{j' \in \mathbb{N}} p_{j,j'} \psi_{(j,j')} (o_t) \beta_{(j,j')} (t). \end{aligned} \right. \qquad (3.20) $$

En utilisant ces deux probabilités, la vraisemblance est simplement calculée en additionnant tous les produits forward et backward, comme indiqué ci-dessous.

$$ \mathbb{P}_M(O_{1:T} = o_{1:T}) := \mathbb{P}(O_{1:T} = o_{1:T} | M) = \sum_{i \in \mathbb{N}} \sum_{j \in \mathbb{N}} \alpha_{(i,j)}(t) \beta_{(i,j)}(t). $$

Considérons maintenant la probabilité conditionnelle d'une transition de l'état $(i, j)$ au temps $t$ à l'état $(i', j')$ au temps $t+1$, étant donné les observations $O_{1:T} = o_{1:T}$ et le modèle $M$

$$ \begin{aligned} \xi_{(i,j),(i',j')}(t) &= \mathbb{P}(Z_t = (i,j), Z_{t+1} = (i',j') | O_{1:T} = o_{1:T}, M) \\ &= \frac{\alpha_{(i,j)}(t)p_{(i',j')}(\psi_{(i',j')}(o_{t+1})\beta_{(i',j')}(t+1))}{\sum_{i \in \mathbb{N}} \sum_{j \in \mathbb{N}} \alpha_{(i,j)}(t)\beta_{(i,j)}(t)} \delta_{i'=j}. \end{aligned} \qquad (3.21) $$

Soit $\gamma_{(i,j)}(t) = \sum_{j' \in \mathbb{N}} \xi_{(i,j),(j,j')}(t)$, la probabilité conditionnelle d'être dans l'état $(i, j)$ à l'instant $t$, étant donné les observations $O_{1:T} = o_{1:T}$ et le modèle $M$ :

$$ \gamma_{(i,j)}(t) = \frac{\alpha_{(i,j)}(t)\beta_{(i,j)}(t)}{\sum_{i \in \mathbb{N}} \sum_{j \in \mathbb{N}} \alpha_{(i,j)}(t)\beta_{(i,j)}(t)}. \qquad (3.22) $$

En termes d'interprétation, $\sum_{t=1}^{T-1} \xi_{(i,j),(i',j')}(t)$ et $\sum_{t=1}^{T-1} \gamma_{(i,j)}(t)$ représentent respectivement des estimateurs du nombre de transitions de $(i, j)$ vers $(i', j')$ et du nombre de transitions issues de $(i, j)$ dans la suite $Z_{1:T}$, fondés sur l'observation de $O_{1:T} = o_{1:T}$.
Le théorème suivant présente les mises à jour itératives de l'algorithme de Baum-Welch pour estimer les caractéristiques du HMM, étant donné le modèle à l'étape $n$ et la suite des observations $O_{1:T} = o_{1:T}$.

**Théorème 3.3.** Étant donné les paramètres $M^n = (Q^n, \psi^n, \rho^n)$, les estimations du maximum de vraisemblance de $M^{n+1} = (Q^{n+1}, \psi^{n+1}, \rho^{n+1})$ étant donné les observations $(o_1, \dots, o_T)$
sont données, pour $i, j, o \in \mathbb{N}$, par

$$
\begin{align*}
Q_{(i,j),(i',j')}^{n+1} &= \frac{\sum_{t=1}^{T-1} \xi_{(i,j),(i',j')}^{n}(t)}{\sum_{t=1}^{T-1} \gamma_{(i,j)}^{n}(t)} \delta_{i'=j}, \\
\psi_{(i,j)}^{n+1}(y) &= \frac{\sum_{t=1}^{T} 1_{o_t(o)} \gamma_{(i,j)}^{n}(t)}{\sum_{t=1}^{T} \gamma_{(i,j)}^{n}(t)}, \\
\rho_{i,j}^{n+1} &= \gamma_{(i,j)}^{n}(1), \\
p_{i,j}^{n+1} &= \frac{\sum_{i' \in \mathbb{N}} Q_{(i',i),(i,j)}^{n+1}}{\sum_{i' \in \mathbb{N}} 1_{Q_{(i',i),(i,j)}^{n+1} \neq 0}},
\end{align*}
$$

où $\alpha_{(i,j)}^n(t), \beta_{(i,j)}^n(t), \xi_{(i,j),(i',j')}^n(t)$ et $\gamma_{(i,j)}^n(t)$ sont définies dans les équations **3.19**, **3.20**, **3.21** et **3.22**, respectivement.

Preuve. Notons par $\mathcal{Z} = \mathbb{N}^2$ et soit $z = (i,j) \in \mathcal{Z}$, étant donné les paramètres du modèle $M^n$, on cherche à maximiser la probabilité d'observations $\mathbb{P}(O_{1:T} = o_{1:T} | M^{n+1})$ en utilisant l'algorithme EM (Bilmes, 1998). La fonction auxiliaire est donnée par :

$$
\Gamma_{O_{1:T}}(M^n, M^{n+1}) = \sum_{z_{1:T} \in \mathcal{Z}^T} \log \mathbb{P}(O_{1:T} = o_{1:T}, Z_{1:T} = z_{1:T} | M^{n+1}) \mathbb{P}(Z_{1:T} = z_{1:T} | O_{1:T}, M^n), \quad (3.23)
$$

où

$$
\mathbb{P}(O_{1:T} = o_{1:T}, Z_{1:T} = z_{1:T} | M^{n+1}) = \rho_{z_1}^{n+1} \left( \prod_{t=1}^{T-1} Q_{z_t, z_{t+1}}^{n+1} \right) \left( \prod_{t=1}^{T} \psi_{z_t}^{n+1}(o_t) \right).
$$

En remplaçant la probabilité conjointe de $(O_{1:T}, Z_{1:T})$ par son expression dans **3.23**, nous pouvons réécrire la fonction $\Gamma_{O_{1:T}}$ comme suit :

$$
\begin{align*}
\Gamma_{O_{1:T}}(M^n, M^{n+1}) &= \sum_{z_{1:T} \in \mathcal{Z}^T} \log \rho_{z_1}^{n+1} \mathbb{P}(Z_{1:T} = z_{1:T} | O_{1:T}, M^n) \\
&\quad + \sum_{z_{1:T} \in \mathcal{Z}^T} \left( \sum_{t=1}^{T-1} \log Q_{z_t, z_{t+1}}^{n+1} \right) \mathbb{P}(Z_{1:T} = z_{1:T} | O_{1:T}, M^n) \\
&\quad + \sum_{z_{1:T} \in \mathcal{Z}^T} \left( \sum_{t=1}^{T} \log \psi_{z_t}^{n+1}(o_t) \right) \mathbb{P}(Z_{1:T} = z_{1:T} | O_{1:T}, M^n) \\
&= \Gamma_{O_{1:T}}^{\rho}(M^n, M^{n+1}) + \Gamma_{O_{1:T}}^{Q}(M^n, M^{n+1}) + \Gamma_{O_{1:T}}^{\psi}(M^n, M^{n+1}).
\end{align*}
$$

i) Maximisation de $\Gamma_{O_{1:T}}^{\rho}(M^n, M^{n+1})$ par la méthode de multiplicateur de Lagrange pour contrainte $\sum_{(i',j')} \rho_{i',j'}^{n+1} = 1$,

$$
\mathcal{L}(\rho, \theta) = \sum_{(i,j)} \log \rho_{i,j}^{n+1} \mathbb{P}(Z_1 = (i,j) | O_{1:T}, M^n) + \theta \left( \sum_{(i',j')} \rho_{i',j'}^{n+1} - 1 \right).
$$

En dérivant, on obtient,

$$
\begin{align*}
\frac{\partial \mathcal{L}(\rho, \theta)}{\partial \rho_{i,j}^{n+1}} &= \frac{\mathbb{P}(Z_1 = (i,j) | O_{1:T}, M^n)}{\rho_{i,j}^{n+1}} + \theta \\
&= 0,
\end{align*}
$$
et

$$ 
\begin{aligned}
\frac{\partial \mathcal{L}(\rho, \theta)}{\partial \theta} &= \sum_{(i', j')} \rho_{i', j'}^{n+1} - 1 \\
&= 0.
\end{aligned}
$$

Les deux équations donne,

$$ \theta = -1 \text{ et } \rho_{i,j}^{n+1} = \mathbb{P}(Z_1 = (i,j) | O_{1:T}, M^n) = \gamma_{(i,j)}^n(1). $$

ii) La fonction de Lagrange associée à la maximisation de $\Gamma_{O_{1:T}}^Q(M^n, M^{n+1})$ avec la contrainte $\sum_{(i', j')} Q_{(i,j),(i',j')}^{n+1} = 1$ est définie par :

$$ 
\begin{aligned}
\mathcal{L}(Q, \theta) &= \sum_{(i,j)} \sum_{(i', j')} \left( \sum_{t=1}^{T-1} \log Q_{(i,j),(i',j')}^{n+1} \right) \mathbb{P}(Z_t = (i,j), Z_{t+1} = (i',j') | O_{1:T}, M^n) \\
&\quad + \sum_{(i,j)} \theta_{(i,j)} \left( \sum_{(i',j')} Q_{(i,j),(i',j')}^{n+1} - 1 \right).
\end{aligned}
$$

En dérivant et en mettant les dérivées premières des Lagrangiens à zéro, on obtient,

$$ \sum_{t=1}^{T-1} \frac{\mathbb{P}(Z_t = (i,j), Z_{t+1} = (i',j') | O_{1:T}, M^n)}{Q_{(i,j),(i',j')}^{n+1}} + \theta_{(i,j)} = 0, $$

et

$$ \sum_{(i', j')} Q_{(i,j),(i',j')}^{n+1} - 1 = 0. $$

La multiplication par $Q_{(i,j),(i',j')}^{n+1}$ et la sommation sur $(i', j')$ de la première équation permettent d'obtenir,

$$ \theta_{(i,j)} = - \sum_{t=1}^{T-1} \mathbb{P}(Z_t = (i,j) | O_{1:T}, M^n). $$

Ainsi, $Q_{(i,j),(i',j')}^{n+1}$ est donnée par :

$$ Q_{(i,j),(i',j')}^{n+1} = \frac{\sum_{t=1}^{T-1} \mathbb{P}(Z_t = (i,j), Z_{t+1} = (i',j') | O_{1:T}, M^n)}{\sum_{t=1}^{T-1} \mathbb{P}(Z_t = (i,j) | O_{1:T}, M^n)} = \frac{\sum_{t=1}^{T-1} \xi_{(i,j),(i',j')}^n(t)}{\sum_{t=1}^{T-1} \gamma_{(i,j)}^n(t)} \delta_{i'=j}. $$

iii) Le Lagrangien associé à la maximisation de $\Gamma_{O_{1:T}}^\psi(M^n, M^{n+1})$ avec la contrainte $\sum_o \psi_{(i,j)}^{n+1}(o) = 1$ est donné comme suit :

$$ 
\begin{aligned}
\mathcal{L}(\psi, \theta) &= \sum_{(i,j)} \left( \sum_{t=1}^{T} \log \psi_{(i,j)}^{n+1}(o) \right) \mathbb{P}(Z_t = (i,j) | O_{1:T}, M^n) \mathbb{1}_{o_t = o} \\
&\quad + \sum_{(i,j)} \theta_{(i,j)} \left( \sum_{o} \psi_{(i,j)}^{n+1}(o) - 1 \right).
\end{aligned}
$$

Comme précédemment, on dérive et on trouve,

$$ \sum_{t=1}^{T} \frac{\mathbb{P}(Z_t = (i,j) | O_{1:T}, M^n) \mathbb{1}_{o_t = o}}{\psi_{(i,j)}^{n+1}(o)} + \theta_{(i,j)} = 0, $$
et

$$ \sum_{o} \psi_{(i,j)}^{n+1}(o) - 1 = 0. $$

Par conséquence,

$$ \theta_{(i,j)} = - \sum_{t=1}^{T} \mathbb{P}(Z_t = (i,j) | O_{1:T}, M^n), $$

et

$$ \psi_{(i,j)}^{n+1}(o) = \frac{\sum_{t=1}^{T} \mathbb{P}(Z_t = (i,j) | O_{1:T}, M^n) \mathbb{1}_{o_t = o}}{\sum_{t=1}^{T} \mathbb{P}(Z_t = (i,j) | O_{1:T}, M^n)} = \frac{\sum_{t=1}^{T} \mathbb{1}_{o_t = o} \gamma_{(i,j)}^{n}(t)}{\sum_{t=1}^{T-1} \gamma_{(i,j)}^{n}(t)}. $$

iv) Enfin, pour établir un estimateur de $p_{i,j}$, il suffit de remarquer, à partir de la définition de la matrice de transition $Q$ de la chaîne $Z$, que :

$$ \sum_{i'} p_{i,j} = \sum_{i'} Q_{(i',i),(i,j)}. $$

Ainsi,

$$ p_{i,j}^{n+1} = \frac{\sum_{i' \in \mathbb{N}} Q_{(i',i),(i,j)}^{n+1}}{\sum_{i' \in \mathbb{N}} \mathbb{1}_{Q_{(i',i),(i,j)}^{n+1} \neq 0}}. $$

□

Après $m$ itérations, nous obtenons les paramètres optimaux du modèle et les estimations $\hat{\lambda}$, $\hat{\mu}$ et $\hat{\nu}$ sont données par

$$ (\hat{\lambda}, \hat{\mu}, \hat{\nu}) = g(p_{0,0}^m, p_{0,1}^m, p_{1,0}^m). $$

Comme d'habitude pour l'algorithme EM, le résultat dépend fortement de la condition initiale, et il faut également régler le nombre d'itérations et trouver un moyen approprié de tronquer les sommes infinies. Tous ces points sont abordés dans la section suivante.

## 3.4 Étude numérique

Nous étudions maintenant les performances numériques de notre procédure d'estimation sur des données synthétiques et réelles. Tous les codes ont été développés en R. Pour mettre en œuvre la fonction Lambert-W, nous utilisons le package lambertW R [64]. Pour l'algorithme de Baum-Welch, nous avons adapté les fonctions du package HMM [65].

### 3.4.1 Données synthétiques

Nous commençons par étudier les performances des estimateurs $\hat{\lambda}^n$, $\hat{\mu}^n$ et $\hat{\nu}^n$ sur des données simulées. Afin de quantifier la perte d'information due aux états cachées, nous allons utiliser deux cadres différents. Dans un premier temps, nous allons considérer le cadre où le processus $(I_n)$ est entièrement observé. Nous allons alors calculer les estimateurs en utilisant toutes les données disponibles. Dans un second temps, nous allons utiliser le cadre du HMM. En utilisant les données simulées, nous allons également étudier la convergence des estimateurs en fonction de la taille de l'échantillon $n$.
Observations complètes

Lorsque la chaîne de Markov à temps discret ($I_n$) est observée, on peut utiliser les estimations du maximum de vraisemblance pour les probabilités de transition $p_{i,j}$,[66]

$$ \hat{p}_{i,j}^n = \frac{\sum_{k=0}^{n-1} 1_{\{I_k=i, I_{k+1}=j\}}}{\sum_{k=0}^{n-1} 1_{\{I_k=i\}}} $$

Nous avons simulé un LBDPI partant de $I_0 = 0$ avec $\lambda = 0.03, \mu = 0.1$ et $\nu = 0.01$, jusqu'à un horizon temporel (continu) $H$ pour plusieurs valeurs de $\Delta t$. Pour estimer les probabilités de transition entre les états discrets par la méthode du maximum de vraisemblance (EMV), nous utilisons la fonction markovchainFit du package markovchain pour les modèles de chaînes de Markov à temps discret [67].

Le tableau (3.1) présente les valeurs estimées des paramètres avec l'erreur standard asymptotique calculée à partir du corollaire 3.2. Nous avons effectué plusieurs simulations

<table><caption>TABLE 3.1 – Estimation des paramètres et de leur erreur standard lorsque le nombre de personnes infectées est observé sur une période $\Delta t$ et un horizon temporel $H$.</caption><thead><tr><th></th><th colspan="3">H = 1000</th></tr></thead><tbody><tr><td>Pas du temps Δt</td><td>1</td><td>7</td><td>30</td></tr><tr><td>nombre d'observations n</td><td>1001</td><td>143</td><td>34</td></tr><tr><td>λ̂<sup>n</sup></td><td>0.027 (0.100)</td><td>0.081 (0.053)</td><td>0.049 (0.065)</td></tr><tr><td>μ̂<sup>n</sup></td><td>0.084 (0.035)</td><td>0.105 (0.046)</td><td>0.070 (0.128)</td></tr><tr><td>ν̂<sup>n</sup></td><td>0.012 (0.484)</td><td>0.007 (0.106)</td><td>0.007 (0.076)</td></tr><tr><td></td><td colspan="3">H = 5000</td></tr><tr><td>Pas du temps Δt</td><td>1</td><td>7</td><td>30</td></tr><tr><td>nombre d'observations n</td><td>5001</td><td>715</td><td>167</td></tr><tr><td>λ̂<sup>n</sup></td><td>0.041 (0.045)</td><td>0.052 (0.024)</td><td>0.064 (0.029)</td></tr><tr><td>μ̂<sup>n</sup></td><td>0.088 (0.016)</td><td>0.132 (0.020)</td><td>0.117 (0.058)</td></tr><tr><td>ν̂<sup>n</sup></td><td>0.009 (0.216)</td><td>0.012 (0.047)</td><td>0.009 (0.034)</td></tr><tr><td></td><td colspan="3">H = 50000</td></tr><tr><td>Pas du temps Δt</td><td>1</td><td>7</td><td>30</td></tr><tr><td>nombre d'observations n</td><td>50001</td><td>7143</td><td>1667</td></tr><tr><td>λ̂<sup>n</sup></td><td>0.030 (0.014)</td><td>0.023 (0.007)</td><td>0.020 (0.009)</td></tr><tr><td>μ̂<sup>n</sup></td><td>0.098 (0.004)</td><td>0.097 (0.006)</td><td>0.102(0.018)</td></tr><tr><td>ν̂<sup>n</sup></td><td>0.010 (0.068)</td><td>0.010 (0.014)</td><td>0.010 (0.011)</td></tr></tbody></table>

du processus de naissance et de mort avec immigration pour différentes tailles d'échantillon $n$, correspondant à des tailles petites, moyennes et très grandes. Les estimations $\hat{\lambda}^n, \hat{\mu}^n$ et $\hat{\nu}^n$ sont calculées avec ces échantillons et comparées aux valeurs réelles du paramètre.

La table (3.1) montre que les estimations convergent vers les vraies valeurs des paramètres et que les erreurs standards asymptotiques sont faibles. Nous constatons également que plus la taille de l'échantillon $n$ est grande (pour un $\Delta t$ petit et un $H$ fixe, ou pour un $H$ grand et un $\Delta t$ fixe), plus les estimations sont précises, ce qui est cohérent avec la théorie. Cependant, pour de petits échantillons (moins de 500), les estimations oscillent beaucoup, ce qui suggère que les estimateurs ne sont pas très précis. Cela peut être dû à une forte
variabilité des données ou à un échantillonnage insuffisant.

En conclusion, ces résultats confirment l’efficacité de nos estimateurs et leur capacité
à estimer les paramètres du processus de naissance et de mort avec immigration, en
particulier pour des échantillons de taille suffisamment grande.

Observations du cumul des nouveaux cas déclarés

Nous étudions maintenant le cas où seuls les cumuls de nouveaux cas déclarés ($O_n$) sont observés. Pour cela, nous avons simulé un processus de LBDPI partant de $I_0 \sim \pi$ avec des paramètres $\lambda = 0.03$, $\mu = 0.1$ et $\nu = 0.01$ jusqu’à un horizon temporel continu $H$ (Fig.3.3). Ensuite, nous avons extrait les cumuls de nouveaux cas déclarés, ce qui correspond au nombre de sauts d’amplitude $-1$ de $I$ sur l’intervalle $]n\Delta t, (n+1)\Delta t]$ , afin d’estimer les paramètres $\lambda, \mu$ et $\nu$ à l’aide de notre procédure d’estimation (Fig.3.4).

<figure><img src="image_29.png" /><figcaption>FIGURE 3.3 – Trajectoire du processus naissance-mort avec immigration (LBDPI) $I_t$ sur un horizon $H = 5000$ jours.</figcaption></figure>

<figure><img src="image_30.png" /><figcaption>FIGURE 3.4 – Le cumul des nouveaux déclarés $O_n$ correspondant au processus $I$ pour une valeur de $\Delta t = 1$ jour.</figcaption></figure>
Pour le HMM, nous choisissons les paramètres initiaux dans les plages suivantes :
$\lambda^{(0)} \in [0.028, 0.034]$, $\mu^{(0)} \in [0.08, 0.15]$ et $\nu^{(0)} \in [0.008, 0.015]$ avec $length.out = 15$ (en générant 15 valeurs initiales pour chaque paramètre), ce qui nous permet d'initialiser les caractéristiques du HMM : $M^{(0)} = (Q^{(0)}, \psi^{(0)}, \rho^{(0)})$, et nous fixons les paramètres $\Delta t = 1$, $maxIterations = 500$ (nombre maximum d'itérations) et $delta = 10^{-9}$ (critère d'arrêt). Nous exécutons l'algorithme HMM pour chaque valeur de $(\lambda_i^{(0)}, \mu_i^{(0)}, \nu_i^{(0)})$ pour $i = 1, \dots, 15$, et nous choisissons le triplet avec la plus grande vraisemblance.

$$ (\hat{\lambda}, \hat{\mu}, \hat{\nu}) = \max_{(\lambda^n, \mu^n, \nu^n)} \mathbb{P}(O|M^n). $$

Sur le plan numérique, on ne peut pas considérer des espaces d'états infinis. Par conséquent, nous avons tronqué les matrices de transition sur l'espace d'état fini $\{0, 1, \dots, N\}$ comme suit : nous calculons les coefficients $(p_{i,j}, 0 \le i \le N, 0 \le j \le N - 1)$ de la matrice de transition par la formule (3.10), et la dernière colonne est définie par :

$$ p_{i,N} = 1 - \sum_{j=0}^{N-1} p_{i,j}, $$

pour $i = 0, \dots, N$.

Comme nous ne nous intéressons qu'aux probabilités de transition de 0 à 0, de 0 à 1 et de 1 à 0, on peut s'attendre à ce que cette procédure de troncature n'ait pas un grand impact sur l'estimation. (Table.3.2).

<table><caption>TABLE 3.2 – L'impact du choix de N sur l'estimation.
Le framework HMM est implémenté dans R, en utilisant une machine de 32 GO RAM.</caption><thead><tr><th rowspan="2">Troncature de l'espace d'état caché N</th><th colspan="3">H = 5000</th></tr><tr><th>3</th><th>4</th><th>5</th></tr></thead><tbody><tr><td>λ<sup>n</sup></td><td>0.0293</td><td>0.0306</td><td>0.0301</td></tr><tr><td>μ<sup>n</sup></td><td>0.095</td><td>0.11</td><td>0.105</td></tr><tr><td>ν<sup>n</sup></td><td>0.009</td><td>0.011</td><td>0.0105</td></tr><tr><td>Temps d'exécution (h)</td><td>0.98</td><td>1.765</td><td>3.595</td></tr></tbody></table>

Nous pouvons remarquer que l'augmentation de la taille de l'espace tronqué N a un impact sur la complexité de l'algorithme EM, avec un temps d'exécution accru. En effet, à chaque incrément de troncature, le schéma itératif de Baum-Welch prend presque le double de temps. Cependant, en termes de précision, nous avons constaté que pour une troncature de l'espace N égale à 5, les estimations sont relativement bonnes par rapport à une troncature égale à 3 ou 4. En général, nous pouvons conclure que la troncature de l'espace n'a pas un grand impact sur l'estimation, mais plutôt sur le temps d'exécution. Enfin, nous fixons N = 5 et nous illustrons l'impact de l'horizon temporel H sur notre schéma d'estimation. (Table.3.3).
<table><caption>TABLE 3.3 – L’impact de l’horizon temporel *H* sur l’estimation.</caption><thead><tr><th></th><th colspan="4">N = 5</th></tr><tr><th>Horizon temporel <i>H</i></th><th>1000</th><th>2000</th><th>3000</th><th>5000</th></tr><tr><th>nombre d'observations <i>n</i></th><th>1001</th><th>2001</th><th>3001</th><th>5001</th></tr></thead><tbody><tr><th>&lambda;&#x0302;<sup><i>n</i></sup></th><td>0.028</td><td>0.0284</td><td>0.031</td><td>0.0301</td></tr><tr><th>&mu;&#x0302;<sup><i>n</i></sup></th><td>0.085</td><td>0.085</td><td>0.12</td><td>0.105</td></tr><tr><th>&nu;&#x0302;<sup><i>n</i></sup></th><td>0.008</td><td>0.008</td><td>0.012</td><td>0.0105</td></tr></tbody></table>

Il convient de souligner que le choix de l'horizon temporel *H* a une influence substantielle sur les estimations obtenues. Les résultats de nos simulations ont montré que l'utilisation d'un horizon temporel plus long permettait d'obtenir des estimations plus fiables des paramètres. En effet, en augmentant l'horizon temporel, nous obtenons une plus grande quantité d'informations qui peut être utilisée pour estimer les paramètres. Ainsi, l'augmentation du nombre d'observations *n* permet d'améliorer la convergence locale de l'algorithme EM vers les vraies valeurs des paramètres.

Cependant, il est important de noter que l'augmentation du nombre d'observations entraîne également une augmentation du temps de calcul nécessaire pour effectuer les estimations. En résumé, le choix de l'horizon temporel et du nombre d'observations est un compromis entre la qualité des estimations et le temps de calcul nécessaire pour les obtenir.

## Discussion sur l'identifiabilité des paramètres :

Le problème d'identifiabilité dans un HMM se réfère à la capacité de distinguer différentes configurations du modèle qui produisent des distributions d'observation similaires. En d'autres termes, il peut y avoir plusieurs jeux de paramètres différents qui donnent lieu à des modèles HMM avec les mêmes performances d'observation. Cela peut se produire lorsque certaines combinaisons de paramètres du modèle HMM sont équivalentes en termes de générer les mêmes probabilités d'observation. Par conséquent, lors de l'estimation des paramètres du modèle HMM à partir des données observées, il peut y avoir une incertitude quant à la véritable configuration des paramètres, ce qui rend le modèle non identifiable.

Il est important de noter que l'identifiabilité du modèle HMM dépend de l'ensemble de données utilisées pour estimer les paramètres [68, 69]. Dans certains cas, une plus grande quantité de données peut aider à résoudre le problème d'identifiabilité en fournissant des informations supplémentaires permettant de distinguer entre différentes configurations de paramètres. De plus, dans notre cas, nous avons effectué des estimations sur des données simulées qui sont rassurantes quant à ce point d'identifiabilité, car les paramètres désirés ont été retrouvés.

En outre, en prenant en compte la structure spécifique de notre chaîne multidimensionnelle *Z* et en réécrivant les fonctions du HMM en conséquence, nous avons pu restreindre l'espace des solutions possibles et obtenir une estimation plus fiable des paramètres. Cette adaptation de la modélisation nous permet de mieux tenir compte des caractéristiques particulières de nos données et d'améliorer l'identifiabilité du modèle HMM.
### 3.4.2 Données de la typhoïde à Mayotte

Nous étudions maintenant un jeu de données réelles de l'Agence régionale de santé (ARS) française de Mayotte qui représente les cumuls quotidiens de nouveaux cas déclarés de fièvre typhoïde sur l'île entre 2018 et 2022.

Nous disposons ainsi de 1816 observations, avec des valeurs comprises entre 0 et 4 (Table.3.4).

<table><caption>TABLE 3.4 – Le nombre de nouveaux cas déclarés de fièvre typhoïde à Mayotte entre 2018 et 2022.</caption><tr><td>valeur d'observation</td><td>0</td><td>1</td><td>2</td><td>3</td><td>4</td></tr><tr><td>nombre des observations</td><td>1596</td><td>167</td><td>34</td><td>12</td><td>7</td></tr></table>

<figure><img src="image_31.png" /><figcaption>FIGURE 3.5 – Cumul hebdomadaire de nouveaux cas déclarés de fièvre typhoïde entre 2018 et 2022 à Mayotte.</figcaption></figure>

Le modèle stochastique proposé est particulièrement bien adapté aux situations de faible prévalence de la maladie que l'on constate à Mayotte. En outre, des périodes de plusieurs semaines, peuvent s'écouler sans qu'aucun nouveau cas ne soit déclaré, suivies d'une augmentation soudaine du nombre de nouveaux cas. Cette observation suggère que le modèle LBDPI en régime récurrent positif, est approprié pour ce contexte. Il convient de noter que la fièvre typhoïde est une maladie à déclaration obligatoire à Mayotte, ce qui signifie en principe que tous les nouveaux cas déclarés sont observés.

Nous avons opté pour une approche de paramétrage initial basée sur la littérature existante, les experts de la fièvre typhoïde à Mayotte et nos données synthétiques. Plus précisément, nous nous sommes appuyés sur des études telles que [28, 30, 40, 70], qui ont abordé la modélisation de la fièvre typhoïde dans divers contextes géographiques. Nous avons également sollicité l'avis d'experts de la typhoïde à Mayotte pour affiner nos choix. Enfin, nous avons exploré nos données synthétiques pour évaluer la performance du modèle pour différentes combinaisons de paramètres initiaux. Nous avons finalement choisi
des paramètres initiaux dans les intervalles suivants $\lambda^{(0)} = [0.05, 0.08]$, $\mu^{(0)} = [0.11, 0.25]$ et $\nu^{(0)} = [0.015, 0.03]$, en générant une séquence de paramètres de longueur 20 à l'aide de la fonction length.out. Ceci nous permet d'initialiser les caractéristiques de notre modèle, à savoir $P^{(0)}, Q^{(0)}, \psi^{(0)}$ et $\rho^{(0)}$, de la même manière que celle utilisée précédemment pour les données synthétiques.

Pour chaque valeur initiale de $\lambda^{(0)}, \mu^{(0)}$ et $\nu^{(0)}$ dans notre grille, nous avons appliqué un critère d'arrêt basé sur le calcul du résidu, fixé à $\delta\ell = 10^{-9}$, et un nombre maximum d'itérations, défini à $\text{maxIterations} = 500$, pour l'algorithme de Baum-Welch adapté du théorème **3.3**. Nous avons ensuite extrait les valeurs estimées pour les paramètres $\lambda, \mu$ et $\nu$. Enfin, le triplet de valeurs qui maximise la vraisemblance des observations a été sélectionné comme l'estimation finale (voir le tableau des différentes valeurs de log-vraisemblance obtenues en Annexe (C)).

Les estimations obtenues sont présentées ci-dessous.

$$ \hat{\lambda} = 0.054, \quad \hat{\mu} = 0.132 \text{ et } \hat{\nu} = 0.017. $$

Notons que la condition de récurrence positive $\lambda < \mu$ est satisfaite. L'estimation de $\lambda$ indiquant que la transmission de personne à personne joue un rôle significatif dans la propagation de la maladie à Mayotte met en évidence l'importance des interactions directes entre individus infectés et susceptibles dans la transmission de la fièvre typhoïde. Cela souligne l'importance de mesures telles que l'hygiène personnelle, le lavage des mains et la sensibilisation à la maladie pour réduire la transmission de personne à personne.

La durée moyenne entre le moment où une personne devient infectée et le moment où elle déclare sa maladie est de 7.5 jours, ce qui peut être influencé par plusieurs facteurs, tels que la gravité des symptômes, la connaissance et la sensibilisation des individus vis-à-vis de la maladie, l'accès aux soins de santé, et les pratiques de déclaration des cas. Un délai plus court indique une détection et une déclaration plus rapides des cas de typhoïde, ce qui peut contribuer à une réponse plus efficace en termes de prise en charge des individus infectés et de prévention de la propagation de la maladie.

Enfin, en ce qui concerne la transmission par l'eau, l'estimation du taux d'infection dû à la transmission de la bactérie Salmonella Typhi suggère que la qualité de l'eau joue un rôle dans la propagation de la fièvre typhoïde à Mayotte. Cela souligne l'importance de garantir l'accès à une eau potable sûre et de mettre en place des mesures de contrôle de la qualité de l'eau pour réduire les risques de transmission par cette voie. L'amélioration du système d'assainissement, y compris la gestion adéquate des eaux usées et des déchets, ainsi que la promotion de bonnes pratiques d'hygiène, sont des mesures essentielles pour limiter la propagation de la fièvre typhoïde et d'autres maladies hydriques à Mayotte.

Il convient de souligner que ces résultats sont basés sur des estimations spécifiques à Mayotte et doivent être interprétés dans le contexte de cette région. Des études supplémentaires et des mesures de surveillance continues sont nécessaires pour mieux comprendre la dynamique de la propagation de la fièvre typhoïde à Mayotte et pour orienter les interventions de santé publique visant à réduire la prévalence de la maladie.

### 3.4.3 Application à l'estimation de la probabilité d'atteindre un seuil épidémique dans une période de temps

Dans cette section, notre objectif est d'étudier la probabilité qu'un cumul d'individus infectés atteint ou dépasse un seuil épidémique (représenté par l'état j) au cours d'une période de temps t, en fonction du nombre initial d'infectés ($I_0$). Cette probabilité
nous permet d'évaluer la progression de l'infection à partir de l'état initial $i$. En d'autres termes, en utilisant notre modèle LBDPI calibré sur les données de Mayotte, nous cherchons à estimer la probabilité d'observer un nombre d'infectés égal ou supérieur à $j$ dans l'intervalle de temps $[0, t]$, sachant que le processus débute avec $I_0 = i$ individus infectés. Cette mesure est essentielle pour évaluer la dynamique de propagation de l'infection et pour étudier les différents scénarios possibles en fonction des valeurs de $i, j$ et $t$. On désigne par $j > 0$ le seuil épidémique de la fièvre typhoïde à Mayotte, qui correspond en moyenne à un nombre d'infections hebdomadaires de l'ordre de la dizaine. Nous définissons la probabilité d'atteindre l'état $j$ dans un intervalle de temps $[0, t]$, sachant que le processus débute à l'état $i$ par :

$$h_{i,j}(t) = \mathbb{P}_i(T_j \le t) = \mathbb{P}(T_j \le t | I_0 = i),$$

pour $0 \le i \le j$, avec $T_j$ est le temps d'atteinte l'état $j$ par le processus $(I_s)_{s \ge 0}$.

**Théorème 3.4.** La probabilité $h_{i,j}(t)$ satisfait le système d'équations suivant :

$$ 
\left\{ 
\begin{array}{l}
\displaystyle \frac{dh_{0,j}(t)}{dt} = -\nu h_{0,j}(t) + \nu h_{1,j}(t), \\
\\
\displaystyle \frac{dh_{i,j}(t)}{dt} = -((\lambda + \mu)i + \nu)h_{i,j}(t) + (\lambda i + \nu)h_{i+1,j}(t) + \mu i h_{i-1,j}(t) \quad 1 \le i \le j-1, \\
\\
h_{j,j}(t) = 1,
\end{array}
\right.
\qquad (3.24) $$

avec la condition initiale : $h_{i,j}(0) = \delta_{i=j}$.

*Preuve.* Soit $t > 0$ et $0 \le i \le j - 1$, on a :

$$h_{i,j}(t) = \mathbb{P}_i(T_j \le t) = \mathbb{P}_i(T_j \le t, S_1 \ge t) + \mathbb{P}_i(T_j \le t, S_1 < t),$$

où $S_1$ est le premier temps de saut du processus $(I_s)_{s \ge 0}$.

D'une part, nous avons :

$$ 
\begin{aligned}
\mathbb{P}_i(T_j \le t, S_1 \ge t) &= \mathbb{P}_i(T_j \le t | S_1 \ge t)\mathbb{P}_i(S_1 \ge t) \\
&= \mathbb{1}_{i=j}e^{-((\lambda+\mu)i+\nu)t}.
\end{aligned}
\qquad (3.25) $$

D'autre part :

$$ 
\begin{aligned}
\mathbb{P}_i(T_j \le t, S_1 < t) &= \mathbb{P}_i(T_j \le t | S_1 < t)\mathbb{P}_i(S_1 < t) \\
&= \int_0^t \mathbb{P}_i(T_j \le t | S_1 = s)((\lambda + \mu)i + \nu)e^{-((\lambda+\mu)i+\nu)s} ds.
\end{aligned}
$$

Par la propriété de Markov, on obtient :

$$ \mathbb{P}_i(T_j \le t | S_1 = s) = \frac{\lambda i + \nu}{(\lambda + \mu)i + \nu} \mathbb{P}_{i+1}(T_j \le t-s) + \frac{\mu i}{(\lambda + \mu)i + \nu} \mathbb{P}_{i-1}(T_j \le t-s). \quad (3.26) $$

En combinant (3.25) et (3.26), on obtient :

$$ h_{i,j}(t) = \mathbb{1}_{i=j}e^{-((\lambda+\mu)i+\nu)t} + \int_0^t ((\lambda i + \nu)h_{i+1,j}(t-s) + \mu i h_{i-1,j}(t-s))e^{-((\lambda+\mu)i+\nu)s} ds. $$

Ainsi,

$$ \frac{dh_{i,j}(t)}{dt} = -((\lambda + \mu)i + \nu)h_{i,j}(t) + (\lambda i + \nu)h_{i+1,j}(t) + \mu i h_{i-1,j}(t). $$
Pour $i = 0$, nous avons :

$$ \frac{dh_{0,j}(t)}{dt} = -\nu h_{0,j}(t) + \nu h_{1,j}(t). $$

Pour $i = j$, on a :

$$ h_{j,j}(t) = \mathbb{P}_j(T_j \le t) = 1. $$

A l'instant $t = 0$, on a :

$$ h_{i,j}(0) = \mathbb{P}_i(T_j \le 0) = \delta_{i=j}. $$

□

## Implémentation : Schéma de Runge Kutta d'ordre 4

Dans ce qui suit, nous décrivons l'algorithmme Runge-Kutta d'ordre 4 (RK4) pour calculer numériquement la probabilité de saturation $h_{i,j}(t)$.

Dans un premier temps, nous remarquons que nous pouvons mettre le système (3.24) sous la forme suivante :

$$ \begin{cases} \frac{dH_j(t)}{dt} = F(H_j(t)), \\ H_j(0) = e_j, \end{cases} $$

où :

$$ H_j(t) = \begin{pmatrix} h_{0,j}(t) \\ h_{1,j}(t) \\ \vdots \\ h_{i,j}(t) \\ \vdots \\ h_{j-1,j}(t) \\ h_{j,j}(t) \end{pmatrix}, \quad e_j = \begin{pmatrix} 0 \\ \vdots \\ \vdots \\ \vdots \\ 0 \\ 1 \end{pmatrix}, $$

et

$$ F(H_j(t)) = \begin{pmatrix} -\nu h_{0,j}(t) + \nu h_{1,j}(t) & & & & \\ -(\lambda + \mu + \nu)h_{1,j}(t) + (\lambda + \nu)h_{2,j}(t) + \mu h_{0,j}(t) & -\nu h_{0,j}(t) + \nu h_{1,j}(t) & & & \\ & \vdots & \ddots & & \\ -((\lambda + \mu)i + \nu)h_{i,j}(t) + (\lambda i + \nu)h_{i+1,j}(t) + \mu i h_{i-1,j}(t) & & & & \\ & & \vdots & \ddots & \\ -((\lambda + \mu)(j-1) + \nu)h_{j-1,j}(t) + \mu(j-1)h_{j-2,j}(t) + (\lambda(j-1) + \nu)h_{j-3,j}(t) & & & & \\ & & & 0 & \end{pmatrix}. $$

Nous présentons ci-après l'algorithmme RK4.
Algorithme 1 : Algorithme RK4

**entrée** : paramètres du modèles : λ, μ, ν.
Seuil de saturation *j*, période du temps *t*,
pas de discrétisation en temps Δ*t*.
le vecteur initial : *H*<sub>j</sub>(0) = *e*<sub>j</sub>.
**sortie** : La probabilité de saturation *H*<sub>j</sub>(*t*).

1 début
2 **Initialisation**
3 *H*<sub>j</sub>(0) ← *e*<sub>j</sub>
4 **pour k de 0 à t par pas de Δt faire**
5 | *K*<sub>1</sub> = *F*(*H*<sub>j</sub>(*k*))
6 | *K*<sub>2</sub> = *F*(*H*<sub>j</sub>(*k*) + Δ*t*/2*K<sub>1</sub>)
7 | *K*<sub>3</sub> = *F*(*H*<sub>j</sub>(*k*) + Δ*t*/2*K<sub>2</sub>)
8 | *K*<sub>4</sub> = *F*(*H*<sub>j</sub>(*k*) + Δ*t*K<sub>3</sub>)
9 | *H*<sub>j</sub>(*k* + 1) ← *H*<sub>j</sub>(*k*) + Δ*t*/6(*K*<sub>1</sub> + 2*K*<sub>2</sub> + 2*K*<sub>3</sub> + *K*<sub>4</sub>)
10 **fin**
11 **renvoyer** : *H*<sub>j</sub>(*t*).
12 **fin**

**Jeu de paramètres** : Nous fixons la fenêtre de temps *t* à 7 jours (une semaine), le seuil épidémique à *j* = 15 cas et prenons un pas de discrétisation d'ordre *h* = 10⁻³.
Les paramètres *λ*, *μ* et *ν* sont estimés dans la section précédent à partir du modèle de Markov caché **3.2.3** :

$$
\hat{\lambda} = 0.054, \hat{\mu} = 0.132, \hat{\nu} = 0.017.
$$

<figure><img src="image_32.png" /><figcaption>FIGURE 3.6 – Évolution de la probabilité $h_{i,j}(t)$ en fonction d'infectés initiaux $i$.</figcaption></figure>

Ces résultats **Fig.3.6** montrent que la probabilité estimée d'observer un cumul d'infec-
tés égal ou supérieur à 15 cas à Mayotte au cours d'une semaine augmente à mesure que
le nombre initial d'infectés augmente. Cela reflète la tendance naturelle de la maladie à
se propager davantage lorsque le nombre initial d'infectés est plus élevé. Pour $i = 0$, la probabilité est très faible mais n'est pas nulle en raison de la possibilité de contamination par l'eau, ce qui peut entraîner de nouveaux cas même en l'absence initiale d'infectés.

Ces résultats peuvent être utiles pour évaluer les risques associés à différents niveaux d'infection initiale et pour guider les décisions en matière de prévention et de surveillance de la fièvre typhoïde à Mayotte.

## 3.5 Conclusion

Nous avons proposé une méthodologie pour estimer les paramètres inconnus d'un LBDPI dans le cas mathématiquement inhabituel mais épidémiologiquement courant où seuls les cumuls de nouveaux cas déclarés sur des périodes données sont disponibles. Notre approche repose à la fois sur des formules analytiques et des itérations numériques. Nous avons étudié les propriétés des estimateurs sur des données synthétiques et appliqué notre procédure à un ensemble de données réelles sur la fièvre typhoïde à Mayotte qui a motivé ce travail. De plus, notre analyse a permis d'estimer la probabilité d'atteindre un seuil épidémique des infectés dans une période de temps, ce qui est essentiel pour prendre des décisions éclairées et mettre en place des mesures de surveillance et de prévention appropriées pour contenir la propagation de la maladie.

Cette approche est un premier pas vers la compréhension de la propagation de la typhoïde à Mayotte. Nous avons proposé un modèle simplifié avec quelques paramètres interprétables.

Dans les chapitres suivants, nous allons poursuivre notre étude en élargissant notre modèle pour prendre en compte des facteurs supplémentaires qui peuvent affecter la propagation de la fièvre typhoïde à Mayotte. Tout d'abord, nous allons inclure la durée d'incubation de la maladie et les individus exposés, ce qui nous permettra de quantifier plus précisément la transmission de la maladie dans la population. Ensuite, nous allons examiner l'effet de la pluviométrie sur l'évolution de la maladie, car il a été observé que les épidémies de fièvre typhoïde à Mayotte surviennent de préférence en saison humide. Ces analyses supplémentaires nous permettront d'obtenir une vision plus complète de la dynamique de la maladie en vue de développer des stratégies de prévention et de contrôle plus efficaces pour réduire l'impact de la fièvre typhoïde sur la santé de la population.
# Chapitre 4

## Modélisation et estimation des paramètres d'un modèle multi-chaîne cachée

### Sommaire

<table><tr><td>4.1</td><td>Introduction</td><td>86</td></tr><tr><td>4.2</td><td>Modèle et propriété</td><td>87</td></tr><tr><td>4.2.1</td><td>Définition du modèle</td><td>87</td></tr><tr><td>4.2.2</td><td>Propriétés mathématiques du modèle</td><td>88</td></tr><tr><td>4.3</td><td>Estimation des paramètres</td><td>92</td></tr><tr><td>4.3.1</td><td>Observations</td><td>93</td></tr><tr><td>4.3.2</td><td>Stratégie d'estimation</td><td>94</td></tr><tr><td>4.3.3</td><td>Expression analytique des paramètres</td><td>94</td></tr><tr><td>4.3.4</td><td>Modèle de Markov Multi-Chaîne Cachée</td><td>97</td></tr><tr><td>4.4</td><td>Simulation numérique</td><td>102</td></tr><tr><td>4.4.1</td><td>Observations complètes</td><td>102</td></tr><tr><td>4.4.2</td><td>Observations du cumul des nouveaux cas déclarés</td><td>103</td></tr><tr><td>4.5</td><td>Conclusion</td><td>107</td></tr></table>
## 4.1 Introduction

Dans ce chapitre, nous modifions le modèle présenté dans le chapitre précédent pour prendre en compte la durée d'incubation de la maladie, qui est la période entre la contamination par la bactérie Salmonella Typhi et l'apparition des premiers symptômes chez un individu. Cette durée peut varier de quelques jours à plusieurs semaines, et elle peut être influencée par divers facteurs individuels et environnementaux, tels que la charge bactérienne dans l'organisme, l'âge, l'état de santé et d'autres facteurs.

En intégrant cette durée d'incubation dans notre modèle, nous pourrons mieux comprendre la dynamique de transmission de la fièvre typhoïde à Mayotte. Nous serons en mesure de déterminer le délai entre la contamination par la bactérie et le début d'apparition des premiers signes cliniques, tels que la fièvre, les maux de tête, la fatigue, les douleurs abdominales, la constipation ou la diarrhée. Nous nous appuyons sur les caractéristiques médicales de la maladie [71] et sur les données épidémiologiques **1.2.4** fournies par l'ARS de Mayotte pour mener notre étude.

Dans ce nouveau modèle, nous proposons une extension du modèle **3.2** afin de l'affiner et de le rendre plus réaliste. Nous introduisons un nouveau compartiment dans le modèle, correspondant aux individus exposés à la bactérie Salmonella Typhi, noté *E*. Les individus exposés sont définis comme ceux qui ont été en contact avec la bactérie, mais qui ne présentent pas encore de symptômes de la maladie et ils sont non contagieux [71].

<figure><img src="image_33.png" /><figcaption>FIGURE 4.1 – Modèle à deux compartiments : Exposé-Infecté.</figcaption></figure>

Le nouveau paramètre $\alpha$ du modèle exposé-infecté représente le taux d'incubation. Le chapitre présent a pour objectif de modéliser la dynamique des individus exposés et infectés, d'étudier les propriétés théoriques du modèle et de mettre en œuvre une approche paramétrique reposant sur un processus Markovien à temps continu pour estimer les paramètres de ce nouveau modèle. Ces paramètres sont les taux de contamination de personne à personne $\lambda$, de contamination par l'environnement $\nu$, d'incubation $\alpha$ et d'isolement $\mu$. Nous utilisons le même schéma d'observation que précédemment pour estimer ces paramètres **3.2.2**.

Le reste du chapitre est organisé comme suit, dans la section **4.2**, nous décrivons la formulation du modèle stochastique exposé-infecté, ainsi que ses propriétés mathématiques. Dans la section **4.3**, nous présentons la procédure d'estimation des paramètres de notre modèle exposé-infecté en suivant une stratégie d'estimation en deux étapes. Dans la section **4.4**, nous mettons en place des simulations numériques pour illustrer les résultats obtenus.