Comparaison et sélection de modèles intégro-différentiels (IDE): les principes.

Dans deux articles précédents (Modélisation spatio-temporelle avec une équation intégro-différentielle (IDE) avec QGis et R et Modélisation spatio-temporelle avec une équation intégro-différentielle (IDE) avec QGis et R(suite))nous avons vu comment calculer et ajuster un modèle intégro-différentiel à un ensemble de données spatio-temporelles, en utilisant R et QGis.

Pour le même ensemble d’observations nous avons vu que l’on pouvait jouer sur différents paramètres, principalement en calculant un noyau du modèle invariant dans l’espace ou un noyau variant dans l’espace.

Les différents résultats de modèles apparaissent, à priori, aussi plausibles les uns comme les autres. Nous verrons dans cette article comment évaluer la qualité du modèle pour pouvoir décider lequel retenir.

Comparaison des résultats du modèle, mais à quoi?

Avant de pouvoir parler d’évaluation des modèles, nous devons décider à quoi nous comparerons notre modèle. Les procédures de modélisation spatio-temporelle hiérarchique nous donnent la distribution prédictive du processus spatio-temporel latent, Y , avec un ensemble d’observations (mesures), Z. Comme nous sommes le plus souvent intéressés par ce processus latent, nous aimerions évaluer notre modèle en fonction de sa capacité à fournir des représentations raisonnables de Y .
Mais par définition, ce processus est caché ou latent – ce qui signifie qu’il n’est pas observé directement – et nous ne pouvons donc pas évaluer directement notre processus modélisé par rapport au processus réel à moins de le faire par simulation. Par ailleurs, nous pouvons évaluer notre modèle à l’aide de distributions prédictives de données, où nous comparons les prédictions de données (et non de Y), fondées sur notre modèle, aux données réellement observées.

Le plus souvent, nous devrons comparer les observations de validation du monde réel, disons Zv, aux observations prédites à partir de notre modèle, disons Zp. La question maintenant est, à quelles observations comparons-nous Zp? La capacité de généralisation d’un modèle est une propriété qui indique dans quelle mesure il peut prédire un ensemble de données d’essai (aussi appelé ensemble de données de validation) qui est différent des données utilisées pour ajuster le modèle. (Veuillez noter que les mots « test » et « validation » sont souvent utilisés de façon interchangeable dans ce contexte; nous préférons utiliser « validation ».) Donc, supposons que nous avons utilisé un échantillon de données, Z, pour ajuster notre modèle. Avant de décrire les différentes possibilités de sélection des données de validation Zv, notez que les processus spatio-temporels ont certaines propriétés qui doivent être prises en compte lors de la comparaison des prédictions du modèle avec les observations du monde réel. En particulier, comme pour les séries temporelles, les processus spatio-temporels ont une dépendance unidirectionnelle du temps et, comme les processus spatiaux, ils ont divers degrés de dépendance spatiale. Ces dépendances doivent être prises en compte dans la mesure du possible lors de l’évaluation d’un modèle spatio-temporel. En général, le choix des observations de validation Zv peut alors être l’un des suivants.

Validation à partir des données d’observation

Il peut être instructif d’utiliser les observations de l’ensemble de données de calcul (Zv = Z) pour évaluer notre modèle, en particulier lors de l’évaluation de la capacité du modèle à s’adapter aux données et pour vérifier les hypothèses du modèle au moyen de diagnostics.
Toutefois, dans le contexte de la prévision, il n’est généralement pas recommandé d’utiliser les données de calcul pour valider la capacité prédictive du modèle, car l’erreur du modèle est généralement optimiste en ce sens qu’elle sous-estime l’erreur prédictive qui serait observée dans un échantillon indépendant. Il n’est peut-être pas surprenant que l’ampleur de cet optimisme soit liée à la mesure dans laquelle une valeur prévue de l’ensemble de données de calcul influe sur sa propre prédiction.

Validation par échantillonnage dans les données d’observation.

Il est souvent utile de construire des échantillons de validation à partir des données d’observations. Dans le calcul du modèle on omet ces échantillons qui peuvent correspondre soit à des pas temporels complets, soit à un échantillonnage aléatoire parmi tous les pas temporels. Bien que cette deuxième solution puisse être appropriée, une évaluation plus appropriée des modèles spatio-temporels résulte de l’omission de « blocs» de données. Ceci est dû au fait que la structure de dépendance spatio-temporelle doit être très bien caractérisée pour combler adéquatement de larges lacunes dans les processus spatio-temporels (processus particulièrement dynamiques).

Validation des prévisions.

L’une des méthodes de validation les plus utilisées pour les données dépendantes du temps consiste à omettre les données du dernier pas de temps dans les observations et d’utiliser le modèle pour prévoir cette période future. Pour prédire l’évolution des caractéristiques spatiales dans le temps, le modèle spatio-temporel doit tenir compte adéquatement de la dépendance spatio-temporelle . Par conséquent, la validation des prévisions fournit un étalon pour de telles évaluations.

Hélas,pas de recette

La validation et sélection d’un modèle s’apparente plus à un travail de détective qu’à l’application d’un protocole. Nous allons suivre ici une démarche « type », en sachant que le type de données, le nombre d’observations, le nombre de pas de temps, et tant d’autres choses, vous obligeront chaque fois à la modifier.

Mais, même si nous n’avons pas de recette, il y a quelques étapes incontournables:

1-Constituer les groupes de données qu’on utilisera pour la validation et celles qu’on va utiliser pour calculer et ajuster le modèle.

2-calculer et ajuster le modèle sur une partie des observations

3-Effectuer les prédictions pour le ou les groupes de données de validation

4-Calculer les scores du modèle en analysant les erreurs de prédiction sur les données de validation

Les étapes 2 et 3 nous les avons vues dans les deux articles précédents. Nous allons donc voir la première et dernière étape. D’abord d’un point de vue général, puis on verra l’implémentation en R pour QGis.

Sous-échantillonnage pour la validation

Cette étape dépend des données sur lesquelles on travaille. Nous continuons à travailler sur les données d’observations de l’oiseau Troglodyte aux États-Unis. Nous avons des observations annuelles sur 21 années avec une quarantaine de points chaque année.

Nous allons extraire une partie des observations que nous laisserons de côté lors du calcul du modèle. Une fois le modèle ajusté on calculera les prédictions pour ces observations et nous aurons alors des prédictions et des observations que nous pourrons comparer.

Selon le but recherché par la modélisation, on choisira une méthode d’échantillonnage ou une autre.

Nous allons utiliser ici deux méthodes:

Nous allons extraire un pas de temps complet que nous mettrons de côté et que nous appellerons par la suite « Bloc ».

Sur les 20 pas de temps restants, on va extraire de manière aléatoire, un certain pourcentage d’observations, que nous mettrons aussi de côté et que nous appellerons par la suite « Random ».

Avec toutes les observations restantes (le volume bleu de l’image), nous allons calculer et ajuster le modèle.

La validation sur le pas de temps complet permet d’avoir une idée de l’utilité du modèle à prévoir les campagnes d’observations futures.

La validation sur l’échantillonnage aléatoire permet d’avoir une idée de l’utilité du modèle à « boucher les trous » des observations manquantes.

Mais, pour l’instant, nous ne faisons que valider notre modèle. Comme nous avons vu dans l’article précédent avec les noyaux variant spatialement, nous pouvons aussi avoir différents paramètres et donc différents modèles. On aura alors à choisir le « meilleur ». Petite digression ici: n’oublions jamais l’aphorisme de Box qui disait que nous savons déjà que tout modèle est faux, mais que ce que nous ignorons c’est à quel point il l’est.

La démarche est d’utiliser les mêmes groupes de données (observations,bloc et random) pour calculer chaque modèle et d’utiliser les scores de validation de chacun pour comparer et choisir le modèle final.

Analyse des erreurs de prédiction

Nous allons voir le graphique des prévisions en fonction des observations qui permet de visualiser l’ajustement du modèle. Malheureusement, notre vue ne suffit pas toujours pour comparer deux nuages de points et choisir celui qui est mieux ajusté.

Nous allons donc calculer une série de « scores » mathématiques qui vont nous aider à décider.

Corrélation observations-prédictions

Ce n’est pas tant la valeur de la corrélation qui nous intéresse, même si c’est déjà une indication, mais disons que la visualisation du nuage de point s’apparente plus à l’analyse exploratoire qu’à la validation à proprement parler.

Elle peut nous permettre de percevoir des problèmes ou des artefacts qui seraient difficilement mis en évidence par des méthodes plus « mathématiques ».

Ces méthodes sont le calcul de « scores » que nous décrivons à la suite.

Le biais et la variance

Le biais est l’erreur résultant de la différence entre la ou les valeurs attendues d’un modèle et la ou les valeurs réelles (ou « correctes ») pour lesquelles nous voulons faire des prédictions sur plusieurs itérations. Dans les concepts scientifiques d’exactitude et de précision, le biais est très semblable à l’exactitude.

La variance est définie comme l’erreur résultant de la variabilité entre différentes prédictions de données dans un modèle. En ce qui concerne la variance, la ou les valeurs correctes n’ont pas autant d’importance que l’étendue des différences de valeur entre les prédictions. La variance joue également un plus grand rôle lorsque nous effectuons plusieurs essais de création de modèles.

L’image suivante est très parlante. La cible est l’objectif de notre modèle. Les impacts sont les prédictions de notre modèle.

Cross Validation: A Beginner’s Guide de Caleb Neale, Demetri Workman, Abhinay Dommalapati
https://towardsdatascience.com/cross-validation-a-beginners-guide-5b8ca04962cd

Il est clair que ce que nous recherchons c’est le modèle avec le biais et la variance la plus petite.

Nous allons calculer trois scores:

  • le BIAIS avec la fonction moyenne(observations-prédictions) . Plus petit sera le biais, mieux c’est. La moyenne de nos prédictions sera plus proche de la moyenne de nos observations.
  • la PCV (predictive cross-validation), la validation croisée prédictive, calculée avec la fonction moyenne(observations-prédictions)². Malgré son nom pompeux, ce n’est que les écarts au carré. Plus petit sera la PCV, mieux c’est.
  • la SCV (standardized cross-validation), la validation croisée standardisée, calculée avec la fonction moyenne( (observations-prédictions)²/variance des prédictions). Dans l’idéal, si le modèle est parfait, le résultat de la SCV est de un. Plus la SCV s’approche de 1 est mieux c’est.

Le CRPS (continuous ranked probability score)

Le score de probabilité classé en continu est devenu l’une des règles de notation (scoring) les plus populaires dans les statistiques spatio-temporelles.

Pour comprendre le CRPS sans rentrer dans des notations mathématiques, les figures suivantes feront l’affaire. La figure de gauche montre la fonction de distribution cumulative, F(x), d’une prédiction avec une moyenne de 6,5 et une écart type de prédiction 1. L’observation est Z = 6, et la zone ombrée indique la différence entre la fonction de distribution cumulative de l’observation et la distribution prédictive. La figure de droite montre la valeur du CRPS qui correspond à la surface de la zone sous et sur la courbe.

Wikle, C. K., Zammit-Mangion, A., and Cressie, N. (2019), Spatio-Temporal Statistics with R

Pour comparer deux modèles nous aurons trois distributions cumulatives:

La courbe noire représente des observations. Les courbes rouges et bleue, le résultat de deux modèles. On voit que le modèle rouge s’approche mieux de nos observations. Et c’est celui qui a un CRPS plus petit. La valeur absolue du CRPS dépend de nos mesures. Il est donc inutile de l’analyser par rapport à une valeur absolue, il nous sert pour comparer deux modèles calculés sur le même lot de données.

Règles de classement (scoring) multivariées

Les règles de classement ci-dessus sont des quantités uni-variées qui peuvent être moyennées ou plus généralement appliquées à travers le temps et l’espace dans nos observations. Pour calculer les différents scores nous utilisons les valeurs des observations ainsi que les valeurs des prédictions. Bien que moins fréquentes, certaines règles de notation tiennent explicitement compte de la nature multivariée d’une prédiction multivariée, ce qui peut être important lorsqu’il y a des dépendances dans le modèle de processus (entre des variables dans l’espace ou le temps).

Un petit point sur la notion de score. Supposons p(z) une distribution de probabilité des prédictions de Z, que nous souhaitons comparer à une observation Z avec une certaine fonction F (c’est-à-dire que nous souhaitons valider notre modèle prédictif). Un score sera une mesure de comparaison entre la distribution prédictive et la valeur observée, notée S(p,Z). En général, les scores sont définis de façon à ce que de plus petits scores indiquent de meilleures prédictions. Le score S est « adéquat »si le résultat de la fonction F(p) >=F(q) pour deux distributions prédictives, p(z) et q(z), où q(z) est la distribution prédictive « vraie ». Autrement dit, le score attendu est minimisé lorsque la distribution prédictive coïncide exactement avec la distribution prédictive réelle. La règle de notation est strictement correcte si ce minimum dans le score attendu se produit seulement quand p(z) = q(z) pour tout z, c’est-à-dire quand la distribution prédictive est la même que la distribution vraie. Le concept « d’adéquation » est très intuitif en ce sens qu’il formalise la notion que si notre distribution prédictive coïncidait avec la vraie distribution, q(z), alors elle devrait être au moins aussi bonne, mais jamais moins bonne, qu’un autre distribution prévisionnelle, p(z), pas égale à q(z).

Nous verrons que dans les scripts R de validation on a implémenté deux classements multivariés: l’energy score et le variogram score d’ordre p.

L’energy score et Variogram score

Si on a des méthodes de scoring adéquates on peut quantifier la performance prédictive des modèles probabilistes. Bien qu’il existe une vaste sélection de ces règles de scoring pour les quantités univariées, il n’existe que peu de règles de scoring pour les quantités multivariées, et bon nombre d’entre elles exigent que les prévisions soient données sous la forme d’une fonction de densité de probabilité. Le score énergétique (energy score), une généralisation multivariée du score de probabilité classé en continu (CRPS), est le seul score couramment utilisé applicable dans le cas des prévisions d’ensembles importants, où la distribution prédictive multivariée est représentée par un échantillon fini. C’est le cas avec les modèles IDE. Malheureusement, sa capacité à détecter des mauvais calculs de corrélations entre les composants du modèle IDE est quelque peu limitée. Une alternative de règles de scoring est basée sur le concept géostatistique des variogrammes. La sensibilité de ces règles de scoring fondées sur les variogrammes et meilleure pour les moyennes, variances et corrélations mal prédites.

Le score d’énergie multivariée (ES) et le score variogramme d’ordre p (V Sp) sont disponibles en R dans le paquet scoringRules. Les deux fonctions que nous utiliserons sont es_sample et vs_sample. Cependant, pour calculer ces scores, nous devons d’abord simuler les prévisions à partir de la distribution prédictive. Pour ce faire, nous avons besoin non seulement des écarts marginaux de prédiction, mais aussi de toutes les covariances de prédiction. En raison de la taille des matrices de covariance de prédiction, la notation multivariée ne peut être effectuée que sur au plus quelques milliers de prédictions à la fois.

Pour les scores « univariés » (biais,PCV, SCV, CRPS)nous avons travaillé seulement avec les valeurs des observation et de prédictions. Pour pouvoir calculer des scores multivariés, il nous faut travailler sur les matrices de covariances et aussi sur une fonction de probabilités qu’on doit construire avec des simulations de notre modèle. Les scripts proposés effectuent, par défaut, 100 simulations pour construire la fonction de densité des probabilités.

Comme pour le reste des scores, lors de la comparaison de deux modèles, les plus faibles valeurs de ES et de VS iniquement une meilleure adéquation entre les prédictions et les observations.

Les critères d’information

Lorsque l’on estime un modèle statistique, il est possible d’augmenter la vraisemblance du modèle en ajoutant un paramètre. Les critères d’information (d’Akaike, tout comme le critère d’information bayésien, BIC), permettent de pénaliser les modèles en fonction du nombre de paramètres afin de satisfaire le « critère de parcimonie « (« les multiples ne doivent pas être utilisés sans nécessité »). On choisit alors le modèle avec le critère d’information le plus faible.

En statistique, l’estimateur du maximum de vraisemblance est un estimateur statistique utilisé pour inférer les paramètres de la loi de probabilité d’un échantillon donné en recherchant les valeurs des paramètres maximisant cette fonction de vraisemblance.

En termes plus simples, les critères d’information privilégient les modèles qui utilisent moins de paramètres pour arriver à prédire les observations.

Critère AIC

Le critère d’information le plus célèbre est peut-être celui d’Akaike (AIC).

Le calcul de ce critère est basé sur une valeur objective, calculée à partir du modèle, la log-vraisemblance.

Cette valeur est pénalisée par le nombre de paramètres dans le modèle. En comparant deux modèles, le modèle avec l’AIC inférieure est meilleur, toutes choses étant égales par ailleurs. L’AIC a le désavantage de dépendre du nombre de paramètres que l’utilisateur juge pertinents. Mais quand il y a des dépendances non connues entre certains paramètres, ce nombre s’avère faux.

Bien qu’il y ait des corrections dans le calcul de l’AIC qui tentent de régler certains de ces problèmes, il faut faire attention lors de son interprétation. L’AIC n’est pas un critère approprié pour la sélection de modèles entre différents types de modèles hiérarchiques Bayesiens. Par contre, dans notre cas, on comparera non seulement le même type de modèle (IDE) mais des réalisations différentes des mêmes paramètres.

On peut donc avoir une bonne confiance sur la valeur de ce critère d’information et retenir le modèle avec l’AIC le plus petit.

Critère BIC

Un autre critère d’information d’usage courant est le critère bayésien d’information (BIC).
Comme pour l’AIC, nous préférerons les modèles avec des valeurs BIC plus petites.
Bien qu’il s’agisse d’un critère d’information « bayésien », il n’est pas non plus approprié pour la sélection de modèles entre différents BHM (modèles hiérarchiques Bayesiens). Encore une fois, le BIC s’appuie sur les estimations des paramètres et ne fournit aucun moyen d’ajuster le terme de pénalité pour tenir compte du nombre effectif de paramètres dans les modèles ayant des effets aléatoires et une dépendance.

A la différence de l’AIC, la log-vraisemblance est pénalisée par le nombre de paramètres, mais aussi parce qu’on prend en compte le nombre d’observations (taille de l’échantillon).

Ayant fait le tour de l’ensemble des critères de validation implémentées dans les scripts R proposés, il est temps de passer à la pratique.

Chose que nous ferons dans le prochain article.

Laisser un commentaire

Votre adresse e-mail ne sera pas publiée. Les champs obligatoires sont indiqués avec *