La méthode de
Monte Carlo


La méthode de la propagation de l’incertitude (GUM) présentée dans la section précédente a permis de montrer qu’elle s’appliquait dans le cas où la non-linéarité du modèle du processus de mesure est négligeable afin d’identifier correctement les dérivées partielles (sensibilités) et lorsque les écarts d’erreur sont relativement petits. Cependant, dans le cas où ces hypothèses ne sont pas vérifiées, la méthode de la propagation de distributions (GUM - S1) dite Méthode de Monte Carlo (MCM) [1] va être appliquée. Elle va permettre de traiter les cas où : la probabilité de la grandeur de sortie n’est pas connue (le niveau de l’intervalle de confiance ne peut donc être spécifié correctement), les grandeurs d’entrée sont décrites par des lois asymétriques, les différentes grandeurs d’entrée contribuent de manière inhomogène à l’incertitude-type composée (composantes majoritaires), l’incertitude-type d’une des grandeurs d’entrée est de l’ordre de grandeur de la mesurande. Dans certains cas particuliers, la méthode de Monte Carlo peut également être utilisée en complément, dans le cas où nous souhaitons faire un test de validation de l’utilisation de la loi de propagation de l’incertitude (GUM), où il n’est pas évident que les conditions d’application soient vérifiées.

L’objectif de la méthode de Monte Carlo est de reconstituer artificiellement un phénomène aléatoire. Elle va simuler un échantillon fictif de réalisation de ce phénomène à partir d’hypothèses sur les variables aléatoires, ainsi la densité de probabilité sera construite à partir d’un échantillonnage qui comprend toutes ces variables aléatoires simulées. Par conséquent, ce n’est plus les incertitudes-types que nous propageons, mais directement les variables aléatoires générées à partir de fonctions de densité de probabilité (PDF) prédéfinies. C’est à dire que la MCM attribue des densités de probabilité aux grandeurs d’entrée du modèle, de façon à en déduire une densité de probabilité pour la grandeur de sortie (mesurande). Ainsi, les analyses statistiques effectuées sur un très grand échantillon (par exemple un million de variables aléatoires pour chaque grandeur d’entrée) seront alors très proches de la population. Par conséquent l’espérance et l’écart-type de la population avoisineront la moyenne et l’incertitude-type (écart-type) de la grandeur de sortie [2, 3].

En prenant un modèle mathématique f du mesurande Y et contenant n grandeurs d’entrée Xi introduisant des erreurs complexes, l’écriture du modèle mathématique est donnée par :


(1)

Ces grandeurs d’entrée Xi sont représentées cette fois-ci par des fonctions de densité de probabilité. C’est à dire que leur densité de probabilité va être construite par une multitude de variables qui seront piochées aléatoirement grâce à un générateur de nombres aléatoires (RNG : random number generator). Les nombres générés sont qualifiés de pseudo-aléatoires car c’est un processus stochastique1 qui se sert d’une valeur de référence (graine) issue de propriétés physiques pour choisir la suite de nombres aléatoires (par exemple à partir d’un bruit électronique, de l’heure, d’un nombre, etc.). Le RNG utilisé dans le cas de nos travaux est l’algorithme de Mersenne Twister [4] implémenté sous Matlab ayant comme valeur de graine 0. Cette graine permet de travailler ainsi sur le même espace de probabilité pour différentes simulations de Monte Carlo, permettant ainsi de répéter les simulations sur le même modèle mathématique avec des rééchantillonnages indépendants. L’algorithme de Mersenne Twister est conçu de telle sorte qu’il soit très peu probable d’obtenir la même suite de nombres aléatoires dans différentes simulations indépendantes de Monte Carlo (auto-corrélation). Par ailleurs, en plus d’être un RNG rapide, lorsqu’il produit une suite de nombres aléatoires, la suite ne se répètera qu’à partir du 219937 − 1 ème tirage. Sa périodicité est donc suffisante pour l’utilisation de Monte Carlo dans nos simulations, puisque le nombre de tirages M ne dépasse pas le million de nombres aléatoires (nombre de tirages préconisé par le GUM-S1).

Ce RNG va produire en sortie une distribution uniforme (loi rectangulaire) de nombres aléatoires bornée sur [0, 1] notée R[0, 1]. Ainsi, à partir de cette distribution, nous pouvons établir n’importe quelle fonction de densité de probabilité pour une grandeurs d’entrée Xi construite à partir de M variables générées aléatoirement qri (où r est le numéro du run avec en tout M tirages et i la grandeur d’entrée considérée avec en tout n grandeurs d’entrée Xi). Le tableau 1 présente différents exemples de fonctions de densité de probabilité utilisées en métrologie.

Tableau 1 - Exemple d’attribution d’une loi de distribution d’une grandeur d’entrée X à partir du générateur de variable aléatoire R[0, 1]. Les fonctions de densité de probabilité (PDF) sont représentées sous forme d’histogramme représentant la répartition des M = 105 variables aléatoires qr générées.


Par conséquent, à travers le modèle de mesure f du mesurande Y, les grandeurs d’entrée Xi ne vont plus directement propager une incertitude-type u(Xi), mais elles vont propager des variables qri piochées une à une aléatoirement à partir de leurs fonctions de densité de probabilité (PDF) :


(2)

C’est à dire que pour un premier scénario (run r=1), une première variable aléatoire est piochée aléatoirement pour chaque grandeur d’entrée (q11, q12, ..., q1n), ainsi toutes ces variables aléatoires vont façonner un premier scénario du mesurande Y1. Ensuite, une fois les M scénarios effectués, le mesurande Yr=1,...,M va décrire une densité de probabilité de sortie à partir des M valeurs de sortie calculées sur chaque scénario. Le nombre de combinaisons et d’interactions possibles entre les paramètres est élevé. C’est pourquoi l’outil statistique de Monte Carlo est nécessaire pour identifier, parmi les nombreuses configurations (scénarios) possibles, celles qui sont les plus critiques dans le résultat du mesurande du fait des interactions entre les grandeurs d’entrée. Le modèle f n’a pas obligatoirement d’expression mathématique explicite. Dans nos travaux, ce modèle va être considéré comme un modèle déterministe, c’est à dire qu’un tirage sur chacune des grandeurs d’entrée (géométrie, mesures environnementales, etc.) va façonner une configuration unique et indépendante sur chacun des M scénarios.

Les différentes étapes de la procédure de Monte Carlo sont résumées dans la figure 1. Tout d’abord, une PDF est attribuée à chaque grandeur d’entrée Xi. Le choix de la PDF repose sur l’ajustement (« fit ») d’un ensemble d’observations répétées ou sur les connaissances disponibles (tolérances du constructeur, spécifications techniques, etc.). Ensuite, un générateur de nombres pseudo-aléatoires est utilisé pour échantillonner à partir de ces PDF afin d’obtenir toutes les variables aléatoires qri. De cette façon, le comportement physique des grandeurs d’entrée est pris en compte grâce à une bonne représentation de différentes valeurs possibles pour chaque grandeur d’entrée. Dans nos travaux, le nombre de tirages M = 105 est choisi afin d’avoir un bon compromis entre la fiabilité des résultats et le temps de calcul. Le modèle f est alors appliqué à ces M-échantillons afin d’obtenir M valeurs possibles pour la grandeur de sortie Yr. Dans la dernière étape, le mesurande Y est caractérisé par sa valeur moyenne Yr et son incertitude-type u(Yr) qui est tout simplement l’écart-type de sa densité de probabilité construite grâce à ces M valeurs de sortie.



Figure 1 - Ce schéma présente le processus de MCM, où (a) représente une illustration de la propagation de la loi de distribution des n grandeurs d’entrée indépendantes Xi, afin d’obtenir M valeurs possibles pour chaque grandeur d’entrée qri. Puis, en (b), ces M échantillons (qr1, qr2, ..., qrn) r=[1;M] sont transférés dans le modèle de mesure f, afin d’obtenir M valeurs du mesurande Yr, où la distribution de probabilité (c) représente la région de convergence basée sur ces M valeurs de sortie. Enfin, en (d), l’incertitude-type u(Yr) et la moyenne Yr du mesurande sont récupérées.

Par ailleurs, la représentation du PDF associée à la grandeur de sortie Y est très bien définie, l’incertitude élargie U(Y ) sera déterminée par les quantiles de cette distribution de sortie [Ylow; Yhigh]. C’est à dire que pour un niveau de confiance p, nous aurons p valeurs parmi les M valeurs de sortie Yr comprises entre Ylow et Yhigh. Par exemple pour p = 95%, l’incertitude élargie va être caractérisée par la région de couverture U = [Ylow; Yhigh] qui va contenir au moins 95% des valeurs possibles pour le mesurande Y. Cette région de couverture peut être définie soit par une région symétrique autour de la moyenne ou par une région la plus compacte. Dans notre cas, nous prendrons essentiellement une région symétrique autour de la moyenne.

Afin d’illustrer la méthode de Monte Carlo, nous reprenons l’exemple du cylindre où le modèle de mesure est donné par :


(3)

Quatre PDF vont être mises en place pour définir les grandeurs d’entrée R, ε1, h et ε2. Nous attribuons comme PDF, par l’expérience, des lois gaussiennes pour les grandeurs R et h avec respectivement une moyenne de 13.53 cm et 15.32 puis un écart-type de 0.05 cm et 0.06 cm. Concernant les grandeurs ε1 et ε2, nous prendrons pour chacune une PDF de loi uniforme bornée à ±0.01 cm. La procédure de mesure est présentée dans la figure 2. Dans un premier temps M = 100000 valeurs sont récupérées aléatoirement pour chaque grandeur d’entrée à partir de leur PDF (échantillonnage). Ensuite chaque jeu de valeurs (run r = 1 à M) va façonner un scénario indépendant à travers le modèle de mesure f, ainsi nous obtenons M valeurs de sortie du mesurande Vol. A partir de ses valeurs de sortie, nous obtenons une densité de probabilité à laquelle nous pouvons extraire la moyenne Vol = 8810.72 cm3 et l’incertitude-type (écart-type) u(Vol) = 83.65 cm3.



Figure 2 - Processus de MCM pour l’exemple du cylindre.

Pour finir, l’incertitude élargie U(Vol) va être déterminée à partir d’un intervalle élargi [Vollow; Volhigh] en ayant pris comme niveau de confiance à 95%. A partir de la densité de probabilité de sortie, nous obtenons comme intervalle élargi du mesurande [Vollow; Volhigh] = [8647.45, 8975.48] cm3. Nous avons donc 95% de chance que le volume du cylindre se trouve entre 8647.45 cm3 et 8975.48 cm3. L’intervalle déterminé par la méthode de propagation de distributions est cohérente par rapport aux valeurs déterminées avec la méthode de propagation des incertitudes. Nous pouvons également conclure sur les hypothèses de linéarité du modèle mathématique et de la non corrélation des grandeurs d’entrée.

Nous pouvons ainsi constater que l’utilisation de la Méthode de Monte Carlo est plus aisée pour déterminer l’incertitude-type ou élargie d’un mesurande. Cependant, selon l’ampleur du modèle mathématique, les temps de calculs peuvent devenir relativement longs (plusieurs jours, semaines, etc.). Mais il reste un outil adapté lorsque le modèle mathématique devient complexe ou lorsqu’il ne peut pas être directement défini.

Cependant, la méthode de propagation des incertitudes (GUM) a l’avantage de fournir un indice de sensibilité à partir des dérivées partielles (section indice de sensibilité). Cet indice reste néanmoins un indice de sensibilité local autour de la meilleure estimation du mesurande (estimation du premier ordre ou second ordre, ensuite il devient compliqué d’estimer les ordres supérieurs). L’investissement supplémentaire pour fournir un indice de sensibilité avec la méthode Monte Carlo (GUM - S1) apportera une information plus riche, à savoir un indice de sensibilité valable sur tout le domaine de variation du mesurande, qu’on appelle indice de sensibilité totale. Les prochains outils statistiques détaillés dans les sections suivantes vont permettre de récupérer les sensibilités de chaque grandeur d’entrée. Cette analyse de sensibilité va étudier la façon dont l’incertitude du mesurande u(Y) à travers le modèle peut être répartie entre les différentes sources d’incertitude des grandeurs d’entrée u(Xi). Ainsi, nous utiliserons dans un premier temps la méthode des plans de Morris [6] qui va permettre d’identifier un ensemble de facteurs influents et de fixer constants les facteurs jugés non influents. C’est un outil peu coûteux, il va permettre de restreindre des grandeurs d’entrée pour réduire le temps de calcul avec l’outil suivant qui concerne les indices de Sobol [7]. Ce dernier outil va établir directement une analyse quantitative de sensibilité pour estimer l’influence des facteurs d’incertitude des grandeurs d’entrée sur le mesurande.


[1] : 101 J 2008 Evaluation of measurement data – Supplement 1 to the “Guide to the expression of uncertainty in measurement” — Propagation of distributions using a Monte Carlo method BIPM 1st ed
[2] : Cox M G and Siebert B R L 2006 The use of a monte carlo method for evaluating uncertainty and expanded uncertainty Metrologia 43 S178–S188
[3] : Cox M, Harris P and Siebert B R L 2003 Evaluation of measurement uncertainty based on the propagation of distributions using monte carlo simulation Measurement Techniques 46 824–833
[4] : Saito M and Matsumoto M 2008 SIMD-oriented fast mersenne twister : a 128-bit pseudorandom number generator Monte Carlo and Quasi-Monte Carlo Methods 2006(Springer Nature) pp 607–622
[5] : Box G E P and Muller M E 1958 A note on the generation of random normal deviatesThe Annals of Mathematical Statistics 29 610–611
[6] : Morris M D 1991 Factorial sampling plans for preliminary computational experimentsTechnometrics 33 161–174
[7] : Sobol' I 1993 Sensitivity estimates for nonlinear mathematical modelsMathematical Modelling and Computational Experiments 1 407–414

Autres documents
Partagés