Est-il possible de retrouver les équations qui gouvernent la dynamique d’un système environnemental exclusivement à partir de séries de mesures ?

Dans le passé, la démarche pour obtenir les équations qui décrivent des processus dynamiques se faisait par un va-et-vient entre observation et théorie, les observations servant de base pour construire notre représentation théorique du monde, puis d’élément de validation ou d’invalidation.
C’est dans les domaines de l’ingénierie électrique et des statistiques, autour des années 1980, qu’a émergé l’idée d’obtenir des équations directement à partir de séries temporelles d’observations. Les premiers développements étaient alors essentiellement linéaires et donc peu adaptés aux phénomènes réels. Ce n’est qu’au cours des années 1990 que les premiers modèles ont pu être obtenus pour des dynamiques non linéaires, et ce pour des cas théoriques et expérimentaux. Ces premiers modèles, reconstruits automatiquement à partir de séries temporelles, permettaient de reproduire la dynamique originale de façon très satisfaisante, mais pour des raisons d’équifinalité1 , les équations obtenues n’étaient pas nécessairement celles des systèmes originaux.
En 2015, un modèle à 3 équations a pu être obtenu pour décrire la dynamique de l’épidémie de peste qui a sévi à Bombay (aujourd’hui Mumbai) au début du XX^e siècle. Ce modèle était inattendu, car permettant de formaliser mathématiquement, directement à partir d’observations, le couplage dynamique entre le nombre de décès humains et le nombre de cas d’infection de deux groupes de rongeurs (rats noirs et rats bruns), et dans une formulation très différente de celles des modèles communément utilisés en épidémiologie. De plus, une interprétation éco-épidémiologique de chacun des termes de ce modèle a pu être proposée, laissant entrevoir l’idée que les équations responsables de la dynamique d’un système peuvent être directement extraites de données observationnelles.

Deux chercheurs du CESBIO ont mis en place un ensemble d’expérimentations numériques pour tester la possibilité de remonter aux équations originales en partant de séries temporelles issues de l’intégration numérique d’équations aux dérivées ordinaires. L’algorithme GPoM2 (Generalized polynomial modelling) qui avait été utilisé pour obtenir le modèle éco-épidémiologique de peste évoqué ci-dessus a été réutilisé à cet effet en suivant la même procédure de recherche de modèle. Cet algorithme s’appuie sur la technique de modélisation globale initiée au cours des années 1990. Il vise à obtenir des équations aux dérivées ordinaires polynomiales directement à partir de séries observationnelles. Initialement, cette technique était destinée à être appliquée à une variable unique. L’algorithme GPoM utilisé dans cette étude s’appuie sur le même formalisme théorique3 , mais son fonctionnement a été modifié et généralisé pour travailler avec plusieurs variables4 .

Une première série d’expérimentations numériques a été conduite, visant à tester le potentiel de l’approche sur un système particulier, le système chaotique de Rössler-1976. Ce système dynamique a été choisi pour sa capacité à générer une certaine complexité dynamique à partir d’une formulation très simple (trois variables et une seule non-linéarité pouvant donner lieu à une trajectoire imprévisible à long terme). L’approche a également été testée en situations dégradées, en modifiant la longueur des séries temporelles, leur échantillonnage, les conditions initiales, le régime dynamique (notamment périodique ou chaotique), ou encore en bruitant les observations et en perturbant le système. Ces expérimentations ont prouvé la possibilité de retrouver les équations originales de ce système particulier.

Une deuxième série d’expérimentations a alors été menée pour explorer la généralité du résultat en appliquant le même outil sur de nombreux systèmes dynamiques. Vingt-sept autres systèmes ont été testés, toujours de petite dimension (jusqu’à cinq variables), mais tous non triviaux et très diversifiés dans leurs propriétés dynamiques, géométriques, algébriques, statistiques et topologiques. Ces systèmes incluaient des modèles de convection, de climat, de particules dans une boîte, de croissance tumorale, de dynamo terrestre, d’oscillateurs couplés, de dynamique de population, ainsi que de nombreux cas purement mathématiques. L’algorithme de modélisation étant polynomial, quatre systèmes non polynomiaux ont également été inclus dans cet ensemble afin d’identifier les risques de mésinterprétation.
Les résultats de cette seconde série d’analyses ont permis de montrer la puissance de l’outil, la plupart des équations originales étant retrouvées, souvent complètes, parfois partielles, selon la concision des systèmes considérés (nombre de termes), et généralement sans détections erronées. Les résultats des tests appliqués aux systèmes non polynomiaux les plus complexes n’ont pas été faussement associés à des modèles polynomiaux (toutes les équations ayant été rejetées) tandis que les formulations obtenues pour les systèmes non polynomiaux moins complexes correspondaient elles à une approximation formelle en séries de Taylor. En outre, les résultats se sont avérés robustes au bruit, le niveau de tolérance pouvant toutefois varier fortement d’un système à l’autre.

• 1La qualité d’un modèle peut être estimée en analysant son ajustement à un jeu de données. La qualité de cet ajustement est généralement basée sur des critère statistiques. Si plusieurs modèles différents conduisent à la même qualité d’ajustement, il y a équifinalité. Le critère ne permet donc pas de discriminer le modèle le plus réaliste. Il faut alors soit renoncer à trouver le meilleur modèle, soit trouver un nouveau critère mieux adapté.
• 2Plateforme de modélisation GPoM développée au CESBIO au cours des projets GloMo, AMoGlo, SpatioGloMo et MoMu dans le cadre des programmes LEFE-INSU et Défi InFiNiTi-CNRS.
• 3Celui de la théorie des systèmes dynamiques non linéaires (ou théorie du chaos).
• 4La principale difficulté pour retrouver les équations originales d’un système tient à la sélection de la structure algébrique qui doit être à la fois absolument drastique et extrêmement rigoureuse. Or le nombre de structures possibles est absolument gigantesque. Pour une formulation quadratique avec trois variables, le nombre de structures possibles est de 230 (un peu plus d’un milliard). Tout l’enjeu de la procédure de cet algorithme est donc de réduire ce nombre de modèles possibles à un tout petit nombre. Dans l’exemple choisi ici, l’algorithme va permettre de réduire les 230 modèles possibles à environ mille modèles potentiels, un tel nombre devenant testable par intégration numérique. La dernière étape consistera à comparer le comportement dynamique des modèles dans l’espace des phases.

L’ensemble de ces résultats renforcent l’idée qu’il est possible de remonter aux équations originales d’un système lorsque celui-ci est polynomial et suffisamment concis (dans le meilleur des cas, jusqu’à 9 termes pour une reformulation complète), et qu’une formulation approchée peut également être obtenue lorsque les équations originales sont proches d’une formulation polynomiale. Pour cette raison, tout en gardant en tête le potentiel de l’approche à obtenir une formulation concise des équations, il apparaît envisageable de proposer une interprétation contextuelle (biologique, physique, chimique, épidémiologique, etc.) pour les termes des modèles obtenus avec cet outil.

Ce travail a été financé par le programme LEFE du CNRS-INSU dans le cadre des projets AMoGlo (Assimilation de données dans un modèle chaotique obtenu par modélisation globale), SpatioGloMo (Modélisation globale spatialisée) et MoMu (Modélisation globale multivariée) et par le programme Défi InFiNiTi du CNRS dans le cadre du projet Musc & SlowFast (Modélisation multiéchelle des systèmes lents-rapides).