DOMAINE DES SCIENCES ET TECHNOLOGIES MENTION EN MATHEMATIQUES ET INFORMATIQUE Mémoire en vue de l’obtention du Diplôme de Master en Mathématiques… [614567]
UNIVERSITÉ D’ANTANANARIVO
DOMAINE DES SCIENCES ET TECHNOLOGIES
MENTION EN MATHEMATIQUES ET INFORMATIQUE
Mémoire en vue de l’obtention du Diplôme de
Master en Mathématiques
Option : Mathématiques Appliquées
Spécialité : Calcul Numérique
PRATIQUE DE LA METHODE MCMC
EN
ESTIMATION BAYÉSIENNE
Présenté par
RATISARIJAONA Mickaël David
Le 23 Janvier 2018
Devant la Commission d’examen formée par :
Le Président : Monsieur RABEHERIMANANA Toussaint Joseph
Professeur Titulaire à l’Université d’Antananarivo
Le Rapporteur : Monsieur RAZAFIMAHALEO Martial
Maître de Conférences à l’Université d’Antananarivo
L’Examinateur : Monsieur ANDRIATAHINY Harinaivo
Maître de Conférences à l’Université d’Antananarivo
Remerciements
Toutd’abord,jetiensàremerciernotreDieutoutpuissantdem’avoirdonnélaforce
et le courage durant mes études universitaires, de m’avoir guidé et surtout aidé pour
mener à terme ce mémoire de Master.
Ce travail n’aurait pas été accompli sans la participation et la présence de plusieurs
personnes :
Je tiens premièrement à remercier le Professeur RABEHERIMANANA Toussaint
Josephqui nous a fait l’honneur de présider le jury de la soutenance de ce mémoire.
Je voudrais exprimer mes remerciements les plus sincères à mon encadreur, le Doc-
teurRAZAFIMAHALEOMartial ,Maîtredeconférencesàl’Universitéd’Antananarivo,
qui m’a beaucoup conseillé et aidé dans ce travail.
J’aimeraisaussiremercierchaleureusementl’examinateurdecemémoire,leDocteur
ANDRIATAHINY Harinaivo , Maître de conférences à l’Université d’Antananarivo, et
le premier responsable de la spécialité Calcul Numérique, d’avoir accepté d’offrir ses
connaissances pour l’amélioration de ce mémoire.
JetienségalementàremerciertouslesenseignantsdenotreFacultédesSciencesno-
tamment ceux de la spécialité Calcul Numérique qui m’ont donné une très bonne for-
mation dans ce domaine d’étude.
Je tiens vivement à offrir un très grand remerciement tout particulier à mes parents,
àtousles membresdemafamille, àmesamis,pourleur encouragementetsoutientout
au long de mes études et essentiellement durant le travail de ce mémoire.
I
Table des matières
Introduction 1
I Rappel en inférence statistique 2
1 Concept de base . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2
1.1 Rappel . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2
1.2 Définitions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2
2 Estimation ponctuelle et région de confiance . . . . . . . . . . . . . . . . . 3
2.1 Estimation ponctuelle . . . . . . . . . . . . . . . . . . . . . . . . . . 3
2.2 Qualités d’un estimateur . . . . . . . . . . . . . . . . . . . . . . . . . 4
2.3 Vraisemblance et Information de Fisher . . . . . . . . . . . . . . . . 5
2.4 Estimateur du Maximum de Vraisemblance . . . . . . . . . . . . . . 6
2.5 Région de confiance . . . . . . . . . . . . . . . . . . . . . . . . . . . 8
II Inférence bayésienne 10
1 Motivation du choix Bayésien . . . . . . . . . . . . . . . . . . . . . . . . . . 10
1.1 Principe Classique . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
1.2 Principe Bayésien . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
1.3 Exemple concluant le choix . . . . . . . . . . . . . . . . . . . . . . . 11
2 Inférence bayésienne . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11
2.1 Estimateur Bayésien . . . . . . . . . . . . . . . . . . . . . . . . . . . 12
2.2 Fiabilité de l’estimateur bayésien . . . . . . . . . . . . . . . . . . . . 15
2.3 Région de crédibilité . . . . . . . . . . . . . . . . . . . . . . . . . . . 15
2.4 Loi a priori . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16
2.5 Exemple . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19
III Chaîne de Markov 22
1 Notion de base . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 22
1.1 Définitions et propriétés . . . . . . . . . . . . . . . . . . . . . . . . . 22
1.2 Evolution temporelle . . . . . . . . . . . . . . . . . . . . . . . . . . 23
1.3 Relation de Chapman-Kolmogorov . . . . . . . . . . . . . . . . . . 24
1.4 Propriété de Markov forte . . . . . . . . . . . . . . . . . . . . . . . . 25
2 Classification des états . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
2.1 Irréductibilité . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
2.2 Récurrence . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
2.3 Périodicité . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28
3 Comportement asymptotique . . . . . . . . . . . . . . . . . . . . . . . . . . 29
3.1 Distribution stationnaire . . . . . . . . . . . . . . . . . . . . . . . . . 29
3.2 Convergence . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31
3.3 Réversibilité . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35
II
TABLE DES MATIÈRES III
IV Méthode numérique 36
1 Méthode de Monte-Carlo . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36
1.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 36
1.2 Approche théorique . . . . . . . . . . . . . . . . . . . . . . . . . . . 37
1.3 Quelques méthodes générales . . . . . . . . . . . . . . . . . . . . . . 38
2 La méthode MCMC . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40
2.1 Motivation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40
2.2 Approche théorique . . . . . . . . . . . . . . . . . . . . . . . . . . . 40
2.3 Quelques différents algorithmes . . . . . . . . . . . . . . . . . . . . 41
3 Application . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43
3.1 Mise en œuvre . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43
3.2 Exemple d’application . . . . . . . . . . . . . . . . . . . . . . . . . . 44
Conlusion 48
Annexe A IV
Annexe B V
Bibliographie IX
Introduction
Dans des situations d’expériences aléatoires, parfois il semble cohérent de penser
que nous pourrions avoir une certaine vision du phénomène que nous observions. Or,
le principe classique ne se repose que sur la considération de ce qui est observé seule-
ment.Lesobservationsnefournissantseulementquel’imageduphénomène,ilsepour-
raitquecelle-cisoitmauvaise.Malgrécela,l’imageestcorrigéepardesconditionsetles
théorèmes asymptotiques à partir d’un nombre assez suffisamment grands d’observa-
tions.Maislorsquelesobservationssontinsuffisantes,l’approchefréquentisteposecer-
tainsproblèmes.Icidonc,lepremierproblèmedelaméthodefréquentisteestlenombre
d’observation, qui nécessite donc, dans certain cas, l’apport de certaines connaissances
duphénomèneobservé.L’approchebayésienutilisenotreconnaissanceaprioridupro-
blème, le combine par les observations qui, ensuite nous offre une connaissance, dite à
postériori qui contient toutes informations nécessaire pour combler le manque de don-
née en cas d’insuffisance d’observation ( [ RC06] et [RC13]). Malgré cette méthode éton-
nante, le calcul de l’a postériori impose des calculs d’intégrale très rigoureuse qui pose
souventdoncàsontourdesproblèmes.Ceproblèmedelaméthodebayésienneimpose
donc des méthodes de calcul d’approximation d’intégrale, " Les méthodes numériques
". La méthode Monte-Carlo en fait partie. Cette méthode utilise des variables aléatoires
indépendants pour approximer les valeurs des intégrales. Sa fiabilité est assurée par la
loi desgrands nombre. Laméthode MCMC parcontre, pour sapart, utilise deschaînes
deMarkovpourapproximerlescalculsd’intégrales.LeschaînesdeMarkovétantdesva-
riables aléatoires dépendantes, la loi des grands nombres ne peut pas s’y appliquer. Sa
fiabilitéparcontreestassuréeparunthéorèmed’extensiondelaloidesgrandsnombres
appliqué aux chaînes de Markov, le théorème ergodique.
Danslepremierchapitre,nousrappelonsleconceptdelastatistique.Nousexposonsles
problèmesd’estimationsetcommentlesrésoudres.Nousletermineronsparunexemple
derésolutiond’unproblèmed’estimationenutilisantl’EstimateurduMaximumdeVrai-
semblance et présenterons le cas de l’estimation par région de confiance
Dans le second chapitre, nous exposerons les principales différences entre la statistique
classiqueetlastatistiquebayésienne.Nousprésenteronsensuitelemodèlebayésienpour
comprendre comment fonctionne son inférence. Nous terminerons ce chapître par l’ex-
position d’un exemple comparant les deux approches.
Dans le troisième chapitre, nous nous concentrerons sur les chaînes de Markov. Nous
verronslesdifférentestypesd’étatsconstituantsunechaînedeMarkov.Nousfinironsce
chapitre par l’étude de ses aspects asymptotique pour pouvoir démontrer le théorème
ergodique.
Dansledernierchapitre,nousdébuteronsparlacompréhensiondufonctionnementde
laméthodedeMonteCarlopourpouvoircomprendreceluideMCMCetpourquoiilest
lemieuxadaptépourl’approchebayésienne.Nousfinironscechapitrepardesexemples
d’utilisations de cette méthode dans le calcul d’ estimation de bayésienne.
1
Chapitre I
Rappel en inférence statistique
1 Concept de base
1.1 Rappel
En statistique, tout se repose toujours sur une idée simple : soit un ensemble d’objet
Pappelépopulation ,dontcertainsnombredesescaractéristiquessontinconnues.Lebut
consiste à étudier la variabilité des caractéristiques propre à ces objets.
Lorsque toute la population peut être étudiée, nous pouvons alors décrire les données
observées, après études ,et résumer les informations contenues. Cela à l’aide de repré-
sentation graphique (camembert, histogramme,…) et des indicateurs statistiques (indi-
cateurdedispersion,detendancecentrale,…).Lebutestdoncdedécrirelesdonnéesen
formepluslisible.Touscesméthodesrelèventdecequel’onappellela statistiquedescrip-
tive.
Mais lorsque toute la population est non observable, ce qui est dans la majorité des cas
(faute de grandeur de la population, ou faute de moyen financier, …), l’étude ne pourra
se limiter que sur un sous-ensemble d’objet de la population, appelé échantillon . Nous
chercheronsdoncàdéterminerlescaractéristiquesdelapopulation,viaceuxdel’échan-
tillon,par induction .Connaissantlescaractéristiquesdel’échantillon,onessaied’évaluer
celles de la population après observation de celles de l’échantillon. En d’autres termes,
onestimeralespropriétésinconnuesdeladistributiondelapopulationensebasantsur
les propriétés connues de la distribution d’un échantillon tiré de cette population. La
statistique inférentiel , pour sa part, élabore ces méthodes afin de parvenir à un jugement
ou à une décision au vu des résultats d’un échantillon.
1.2 Définitions
Définition 1. On appelle expérience statistique , une expérience qui consiste à recueillir
une observation xd’un élément aléatoire X, à valeur dans un espace Xet dont on ne
connait pas exactement la loi de probabilité P.
Définition 2. On appelle modèle statistique associé à une expérience statistique le triplet
(X;A;P), où :
X:L’espace des observations
A:La tribu des évènements observable associée.
P:L’espace des observations
On dit que le modèle est discretsiXest dénombrable. C’est le cas lorsque l’élément
2
CHAPITRE I. RAPPEL EN INFÉRENCE STATISTIQUE 3
aléatoire observé X2Xa une loi de probabilité discrète.
Onditquelemodèleest continusiX2Rnetque8P2P,Padmetunedensitédans Rn.
Partant de ces définitions, deux grands modèles y découlent pour résoudre le pro-
blème d’inférence. L’un qui considère que Pappartienne aux familles de loi de proba-
bilité que nous connaissons mais que ses caractéristiques sont inconnus, et que l’on es-
timera. L’autre qui considère que Psoit complètement inconnu, et que l’on estimera sa
forme et ses caractéristiques.
Définition 3. On appelle modèle paramétrique, un modèle où l’on suppose que le type
de loi deXest connu, mais qu’il dépend d’un paramètre inconnu. La famille de lois
de probabilité possible pour Xpeut s’écrire alors comme :
P=fPj2Rng
L’inférenceneseferadoncquesurleparamètre ,enl’estimant(ponctuellementoupar
régions de confiance) et en effectuant des tests d’hypothèses portant sur . On fait alors
de lastatistique paramétrique .
Définition4. Onappelle modèlenonparamétrique ,unmodèleoùPnepeutpassemettre
sous la forme comme ci-dessus. Plus concrètement, ce modèle englobe l’ensemble des
méthodes statistiques qui permettent de tirer de l’information pertinente de données
sansfairel’hypothèsequelaloideprobabilitédecesobservationsappartientàunemême
famille paramétrée connue. Dans ce cas, on fait de la statistique non paramétrique .
Dans tout ce qui suivra, on ne se focalisera que dans le cadre paramétrique.
Définition 5. Soit une population caractérisée par une variable aléatoire X.
Onappelle échantillonaléatoire detaillendeX,lasuitede nvariablesaléatoires (X1;X2;:::;Xn)
iid(indépendant et identiquement dsitribuées) et suivent toutes la même loi de X.
Les réalisations de l’échantillon aléatoire (X1;X2;:::;Xn)qui sont observées sont les n
valeurs (x1;x2;:::;xn).
Définition6. Etantdonnéunéchantillonaléatoire (X1;X2;:::;Xn)iid,onappelle statis-
tiquela fonction hdesnvariables aléatoire. Donc, une variable aléatoire aussi.
SoitSn=h(X1;X2;:::;Xn)une statistique, la valeur de cette statistique correspond à
la réalisation de l’échantillon aléatoire (X1;X2;:::;Xn):
sn=h(x1;x2;:::;xn)
2 Estimation ponctuelle et région de confiance
2.1 Estimation ponctuelle
Définition7. Onappelle estimateur d’unparamètreinconnu lastatistiquedéfiniepar:
^n=g(x1;x2;:::;xn)
L’estimation est une réalisation particulière de l’estimateur.
Leproblèmeenestimationserésumeàchercherunestatistiqueutilisécommeestima-
teur qui soit dotée de bonnes propriétés en échantillon de taille finie, en d’autre terme,
quiconduiseàlameilleureapproximationpossibleduparamètre .L’estimateur ^nde,
étantunevariablealéatoire,sesréalisationspourraientdoncs’écarterplusoumoinsdela
CHAPITRE I. RAPPEL EN INFÉRENCE STATISTIQUE 4
vraievaleurde quel’onchercheàestimer.Lebutestdoncdetrouverunestimateurqui
minimisecetécart.Certes,cesestimationsfluctuentautourde E[^n](l’espérancemathé-
matique de l’estimateur ^n) avec une dispersion caractérisée par Var[^n](la variance de
l’estimateur ^n).Ilesttrèsévidentquelesestimateursdontl’espérancecolleaveclavraie
valeurduparamètre ,etdontlavarianceestlaplusminimale,sontlesplusintéressants.
De ce fait, nous pouvons poser ces propriétés comme critère de choix pour plusieurs
estimateurs. Ces critères permettent alors de comparer différents estimateurs possible
du même paramètre afin de pouvoir choisir le plus optimale.
2.2 Qualités d’un estimateur
Biais :
Définition 8. On appelle biaisd’un estimateur ^ndela valeur :
b
^n
=Eh
^n i
=Eh
^ni
Un estimateur est dit sans biais (ou non biaisé) si :
b(^n) = 0
Un estimateur est dit asymptotiquement sans biais (ou asymptotiquement non biaisé) si :
lim
n!1b(^n) = 0
Efficacité :
Définition 9. Un estimateur ^ndeest ditefficace (ou efficient ) s’il est non biaisé et de
variance minimale parmi tous les estimateurs non biaisés.
Il est dit asymptotiquement efficace (ou asymptotiquement efficient) si :
lim
n!1^n= 0estefficaceouefficient
Convergence et Consistance :
Définition 10. Un estimateur ^ndeest ditconvergent si :
lim
n!1E[^n] = 0
Il est dit consistant s’il converge en probabilité (convergence au sens faible) vers :
8">0 lim
n!1P(j^n j") = 0
Erreur quadratique :
Ilestclairquelespropriétéscitéesci-dessussonttrèsimportantespourameneràun
choix optimal d’un estimateur. Certes, cela n’est pas suffisant puisqu’on pourrait trou-
ver un estimateur, même biaisé mais plus convainquant qu’un autre non biaisé. D’où la
nécessité d’introduire la notion d’erreur.
CHAPITRE I. RAPPEL EN INFÉRENCE STATISTIQUE 5
Définition 11. Soit un estimateur ^nde.
L’erreur quadratique moyenne deest la quantité qui mesure le risque d’utiliser ^npour
estimer, ou encore la précision de l’estimateur ^n, et qui est définie par :
E
^n 2
=Varh
^ni
+b2
^n
(I.1)
Démonstration :Soitd’abordladifférence(erreur)entrel’estimateur ^netlavraivaleur
du paramètre .
^n =^n Eh
^ni
+Eh
^ni
=^n Eh
^ni
+b
^n
L’erreur quadratique s’exprime alors :
^n 2
=
^n Eh
^ni2
+ 2
^n Eh
^ni
b
^n
+b2
^n
L’erreur quadratique moyenne s’exprime alors :
E
^n 2
=E
^n Eh
^ni2
+ 2
^n Eh
^ni
b
^n
+b2
^n
=E
^n Eh
^ni2
+ 2b
^n
Eh
^n Eh
^nii
+Eh
b2
^ni
=E
^n Eh
^ni2
+ 2b
^nh
Eh
^ni
Eh
^nii
+Eh
b2
^ni
AvecE
^n Eh
^ni2
=Varh
^ni
.D’où (I.1).
Nouspouvonstrèsbienremarquerqu’aucasoùplusieursestimateurssontnonbiai-
sés, le meilleur estimateur sera celui qui aura la variance la plus minime, ce qui prouve
la nécessité de prise en compte de l’efficacité.
2.3 Vraisemblance et Information de Fisher
Notion de vraisemblance
Définition12. Lavraisemblance (ouLikelihood enanglais)mesureuneadéquationentre
la distribution observée sur un échantillon aléatoire et une loi de probabilité supposée
décrire une réalité sur la population dont l’échantillon est issue.
Partant de cette définition, la vraisemblance est donc un outil fondamental qui utilise
les observations sur les échantillons afin de résoudre le problème de l’inférence.
Considérons donc un échantillon aléatoire (X1;X2;:::;Xn)où lesXisontiid.
Définition 13. Lafonction vraisemblance est la statistique définie par :
Dans le cas discret :
`(jx1;x2;:::;xn) =P(X1=x1;X2=x2;:::;Xn=xn;) =nY
i=1P(Xi=xi;)
Dans le cas continu :
`(jx1;x2;:::;xn) =f(X1;X2;:::;X n)(x1;x2;:::;xn;) =nY
i=1fXi(xi;)
CHAPITRE I. RAPPEL EN INFÉRENCE STATISTIQUE 6
Principedevraisemblance :L’informationapportésparuneobservationde xsurest
entièrementcontenuedanslafonctionvraisemblance `(jx).Deplus,si x1etx2sontdeux
observations qui dépendent du même paramètre , et telles qu’il existe une constante c
satisfaisant :
`1(jx1) =c`2(jx2)
Alors pour tout , elles apportent la même information sur et doivent conduire à la
même inférence.
Information de Fisher
Soit un échantillon aléatoire X= (X1;X2;:::;Xn)où lesXisontiid.
Définition 14. On appelle quantité d’information de Fisher apportée par les observations
de l’échantillon X, sur le paramètre , la quantité, s’il existe, définie par :
Dans le cas unidimensionnel :
I() =E@2ln (f(xj))
@2
Dans le cas multidimensionnel :
I() =E@2ln (f(xj))
@i@j
Où(i;j= 1;2;:::;n )
Il est moins évident d’interpréter cette notion d’information de Fisher en ce contexte.
Mais par son lien avec l’inégalité de Cramer-Rao , on comprend mieux son utilité.
Soittnun estimateur de g(). On a alors :
E
(tn g())2
(g0() b0
())2
nI()
Donc, sig() =et quetnsoit non biaisé, alors :
E
(tn g())2
1
nI()
Var[tn]1
nI()
Ce résultat illustre le concept que l’information de Fischer se compare à une variance.
Ilestévidentque,pluslesinformationsquenousapportelesobservationssontenabon-
dance, il sera plus plausible que notre estimation sera plus proche du vrai valeur à esti-
mer et fera moins d’erreur que le cas contraire.
2.4 Estimateur du Maximum de Vraisemblance
Comme nous l’avions bien pu le constater, dans la notion de vraisemblance, la fonc-
tion de vraisemblance est une fonction du paramètre , qui est inconnue et propre à la
population,maisquiestégaleauxproduitsdesprobabilités(oudensitédanslecasconti-
nue)dechaqueréalisationdel’échantillonaléatoire.Puisquec’estlerésultatduproduit
deplusieursprobabilités,ilestnécessairementcomprisentre 0et1.Decefaitdonc,plus
la valeur de cette fonction est proche de 0, plus la valeur de , qui donne le résultat est
moins plausible. Partant de ce postulat donc, la valeur de qui maximisera la fonction
de vraisemblance, sera la plus crédible. C’est le principe de l’estimation du maximum
devraisemblance .
CHAPITRE I. RAPPEL EN INFÉRENCE STATISTIQUE 7
Définition15. Soitunéchantillonaléatoire X= (X1;X2;:::;Xn)tiréd’unepopulation,
où lesXisontiid.
On appelle estimateur du maximum de vraisemblance (EMV), l’estimateur définie par :
^= max
2`(jX)
Danslapratique,onpréfèremaximiserlelogarithmedelafonctiondevraisemblance( ou
le log-vraisemblance ) afin de transformer les produits en sommes. Puisque maximiser la
fonction de vraisemblance équivaut à maximiser le log-vraisemblance, car la fonction
logarithme est strictement croissante. On obtient donc :
^= max
ln`(jX)
Sousdeshypothèsesderégularité,l’estimateurdumaximumdevraisemblancepos-
sède de très bonnes propriétés. En effet, il est :
Convergent(au sens faible)
Pluslenombred’échantillonaugmente,plusl’estimateurserapprocheradelavraieva-
leur du paramètre.
Asymptotiquement non biaisé
Pluslenombred’échantillonaugmente,plusl’écartmoyenentrel’estimateuretlavraie
valeur du paramètre tend vers 0. C’est-à-dire que son biais tend vers 0. Un biais nul
signifie qu’en moyenne, l’estimation est exacte.
Asymptotiquement efficace
Plus le nombre d’échantillon augmente, plus sa variance tends vers sa borne inférieur.
Asymptotiquement normalement distribué
^! N
;1
I()
n!1
Exemple 1. Soit l’expérience statistique suivante qui consiste à estimer la proportion
d’avoir pile d’une pièce de monnaie en lançant nfois la pièce. Les données sont des va-
riables aléatoires X1;X2;:::;Xniidde même loi de Bernoulli B().
La fonction de vraisemblance s’exprime comme :
`(jX) =P(X1=x1;X2=x2;:::;Xn=xn;) =Pn
i=1xi(1 )n Pn
i=1xi
Par la méthode du maximum de vraisemblance on a :
^= max
2`(jX) = max
ln`(jX)
Avec
ln [`(jX)] = lnh
Pn
i=1xi(1 )n Pn
i=1xii
=nX
i=1xiln() +
n nX
i=1xi!
ln(1 )
CHAPITRE I. RAPPEL EN INFÉRENCE STATISTIQUE 8
Or
@lnh
`(^jX)i
@= 0
Déterminons donc la valeur de ^en résolvant cette équation.
@lnh
`(^jX)i
@= 0,Pn
i=1xi
n Pn
i=1xi
1 = 0
)(1 )nX
i=1xi
n nX
i=1xi!
= 0
)nX
i=1xi nX
i=1xi n+nX
i=1xi= 0
D’où
^=1
nnX
i=1xi
Parlaméthodedumaximumdevraisemblancedonc,oninduitquecetteproportionest
sensiblement égale à la moyenne empirique des échantillons recueillis.
2.5 Région de confiance
Les estimations ponctuelles ne tiennent pas en compte les erreurs dues aux fluctua-
tions d’échantillonnage. Elles n’apportent donc pas d’informations sur la précision des
résultats. Dans ce cas, il est nécessaire de nous offrir un intervalle de valeur contenant,
avec une certaine probabilité, la vraie valeur du paramètre, afin d’évaluer la confiance
que l’on peut avoir sur cette valeur : c’est l’ estimation par intervalle de confiance . Dans le
cas général, c’est-à-dire lorsque la dimension du paramètre est supérieure ou égale à 2,
on parle alors de région de confiance .
Définition16. Etantdonnés X1;:::;Xnunéchantillontiréd’unepopulationet 2]0; 1[;
on appelle région de confiance pour le paramètre , de niveau de confiance 1 , la
famille non vide de parties de ,C(x1;:::;xn)telle que :
82;P
2C(x1;:::;x n)
= 1
Ceci s’interprète comme : il y a un risque de probabilité que le paramètre ne se
trouverait dans la région de valeur C(x1;:::;x n).
Pour la détermination de cette région de confiance, le recourt à des théorèmes li-
mites1et des fonctions asymptotiquement pivotales2sont indispensables ([ FM07]).
Exemple2. Etantdonnéesunéchantillon X1;:::;Xni.i.ddeloiinconnuedontlamoyenne
est inconnu et sa variance connue, alors grâce à la LFGN et le TCL, on peut dire que
l’intervalle :
Xn 1;96:pn;Xn+ 1;96:pn
est un intervalle de confiance asymptotique à 95 pourcent. C’est-à-dire :
1. Loi forte des grands nombres, théorème centrale limite.
2. La fonction converge vers une loi indépendante du paramètre lorsque la taille de l’échantillon tend
vers l’infini.
CHAPITRE I. RAPPEL EN INFÉRENCE STATISTIQUE 9
82;P
2
Xn 1;96:pn;Xn+ 1;96:pn
! 0;95
n!1
Remarque 1. La région de confiance, ne peut être déterminée qu’en grande échantillon
(siladistributiond’échantillonnagen’estpasgaussienne).Puisquel’onnesaitvraiment
pas à quelle vitesse a lieu cette convergence, en pratique, on accepte l’hypothèse que
cette convergence est atteinte lorsque la taille de l’échantillon soit supérieure à 30. Plus
la région de confiance est serrée, plus l’estimation est intéressante.
Chapitre II
Inférence bayésienne
1 Motivation du choix Bayésien
Rappelonsd’abordquelebutprincipaledelastatistiqueestdetiréuneinduction(par
estimation),aprèsobservationd’unphénomènealéatoire,lesprincipalescaractéristique
de ce qui génère ce phénomène afin de donner une analyse d’un phénomène passé ou
bien de faire une prédiction d’un phénomène (de même nature) à venir.Il y a eu tou-
jours une théorie concurrente entre l’approche classique (ou fréquentiste) et l’approche
bayésienne puisque c’est deux approches ont leur propre démarche même d’aborder la
résolution du problème d’inférence.
1.1 Principe Classique
Dans l’approche classique (aussi dite fréquentiste), on fait l’inférence directe sur les
paramètres ayant générés le phénomène via l’observation du phénomène. En d’autre
terme,l’uniquesourced’informationestl’informationprovenantdesobservations.Dans
la notion de vraisemblance, on a :
`(jx) =f(xj)
Où`(jx)est la fonction de vraisemblance et f(xj)la densité de données observées.
Lafonctiondevraisemblanceétantunefonctionde .Onessaidelanormaliserpouren
faire une fonction de densité sur puis l’utiliser pour estimer .Par exemple, par l’es-
timation de type maximum de vraisemblance, on cherche la valeur ^qui maximise la
fonction de vraisemblance. C’est-à-dire :
^= max
2`(jX)
Icidonc,l’inversionde (xj)vers (jx)estpurementformellealorsque,dansl’approche
bayésienne, l’inversion est plus cohérente par l’utilisation de la formule de Bayes1.
1.2 Principe Bayésien
Avecl’approchefréquentiste,leparamètreinconnu estconsidérénonaléatoiremais
juste comme étant une simple variable, alors que dans cette nouvelle approche, l’idée
même de base est de considéré le paramètre inconnu comme aléatoire et admettant
une densité que nous noterons ()appelédensité à priori .
Lebutestdoncd’utilisercettedensitéàpriorietdelecombinerauxinformationsappor-
téesparlesobservations(ladensitédesobservations f(xj))afind’obtenirunenouvelle
1. Voir Annexe A 1.
10
CHAPITRE II. INFÉRENCE BAYÉSIENNE 11
densité dite à postériori, que l’on note (jx)), par l’usage de la formule de Bayes. Nous
avons donc :
(jx) =f(xj)()
f(x)=f(xj)()R
f(xj)()d
Comme nous pouvons le constater, ici l’inversion de cause (paramètre) à effet (observa-
tion) se fait de manière plus satisfaisante car notre connaissance à priori du problème
poséestactualiséeaprèslesobservationspourainsiavoiruneconnaissanceàpostériori.
La densité à postériori est une sorte de mise à jour de notre densité à priori après les
observations.
1.3 Exemple concluant le choix
Soitunefemmevoulantsavoirsiouiounonelleestenceinte.Elleprocèdedoncàun
test.Supposonsmaintenantquedesétudesontpudémontrerqueletestindiquepositif
9=10surdesfemmesréellementenceinteetnégatif 95=100surdesfemmesnonenceintes.
Notons les évènements : E(elle est enceinte), Tp(le test est positif), Tn(le test est néga-
tif). En terme probabiliste, nous avons donc :
P(TpjE) = 0;9etP
TnjE
= 0;95
Maintenant,aprèsavoirfaitletest,ellechercheàinduiresiouiounonelleestréellement
enceinte.ParlaformuledeBayes,ilestpossibledecalculerlaprobabilitéqu’ellelesoit:
P(EjTp) =P(TpjE)P(E)
P(Tp)
Ici,P(E)estlaprobabilitéd’êtreenceinte(indépendammentdurésultattest).Lafemme,
en connaissance de sa santé, estimera donc cette probabilité à priori. Si elle utilisait des
méthodes contraceptives, elle choisira un a priori faible puisqu’elle a une fine raison de
s’inquiéter.Maisparcontre,sielleaeudesrapportsnon-protégésdanslesvoisinagesde
la période d’ovulation, et ayant maintenant des vomissements fréquents utiliserait un a
priori plus élevé. Le résultat du test sera donc accentué, ou atténué par cette estimation
à priori.
L’apportdecetteinformationàprioripourraitdoncrenforcerounuancer(corriger)l’image
que nous apporte les observations. Et c’est cette estimation à priori que les méthodes
classiques nient systématiquement.
2 Inférence bayésienne
Comme nous l’avions bien vu, l’approche bayésienne se différencie de l’approche
fréquentiste du fait que le paramètre n’est plus considéré comme tout simplement in-
connu mais un phénomène aléatoire dont son comportement est plus ou moins connu.
Soit donc la définition suivante :
Définition 17. On appelle modèle bayésien la donnée paramétrique où, f(xj)est la den-
sité des observations et (x)la loi du paramètre (que l’on estimera) indépendant des
observations.
Etant donné ces deux densités, nous pouvons construire :
CHAPITRE II. INFÉRENCE BAYÉSIENNE 12
La densité jointe
'(;x) =f(xj)()
La loi marginale de x
f(x) =Z
f(xj)()d
La densité à postérieur de
(jx) ='(;x)
f(x)=f(xj)()R
f(xj)()d
On peut remarquer (jx)/f(xj)()
2.1 Estimateur Bayésien
Notion de coût, de décision et de risque
Leproblèmeauquelons’intéresseiciestceluid’unindividuplongédansunenviron-
nementdonné(nature)etqui,surlabased’observations,estconduitàmenerdes actions
et à prendre des décisions qui auront un coût. Considérons donc les espaces suivantes :
X: L’espace des observation.
: L’espace des états de nature(l’espace des paramètres en statistique).
A: L’espace des actions ou décisions. On note aune action.
D: L’ensemble des règles de décision. On note une décision.
L’inférence consiste à choisir une règle de décision 2Dconcernant 2sur la base
d’uneobservation x2X,xetétantliésparlaloi f(xj).Ici,larèglededécisionestun
estimateur,l’actionestuneestimation(valeurdel’estimateuraupointd’observation x).
Pour choisir une décision, on construit une relation de préférence en considérant une
mesure coût ou perte encourue lorsqu’on prend la décision (x)et que l’état de nature
est.
Définition 18. On appelle fonction coût , toute fonction LdeAdansR.
Lafonction L(;a)évaluelecoûtd’unedécisionaquandleparamètrevaut .Ellepermet
doncdemesurerlaperteencourueparunemauvaisedécision,c’est-à-direunemauvaise
évaluation de .Donc, il s’agit d’une fonction de .
On dira qu’une décision est bonnesi elle conduit à un coût nul. Autrement dit, une
bonne décision est solution de l’équation :
L(;(x)) = 0
étant inconnu, on ne peut évidemment pas résoudre cette équation. Classer les déci-
sionsparlaconsidérationducoûtestdoncimpossible.Celui-cineprendpasencompte
l’information apportée par le modèle f(xj).
Cesremarquesconduisentàconsidérerlamoyennedelaperte,c’estle risquefréquentiste .
Définition 19. On appelle risque fréquentiste , lecoût moyendu coût d’unerègle de déci-
sion :
R(;) =E[L(;(x))] =Z
XL(;(x))dP(x)
CHAPITRE II. INFÉRENCE BAYÉSIENNE 13
Etantdonnéquel’approcheBayésiennemetàdispositiondustatisticienuneloiapriori
, nous pouvons considérer la moyenne du risque fréquentiste suivant la loi a priori,
c’est lerisque bayésien .
Définition 20. On appelle risque bayésien (ourisque de Bayes ), l’espérance du risque fré-
quentiste suivant la loi à priori et que l’on note r(;):
r(;) = E[R(;)]
=Z
XR(;)()d
=Z
Z
XL(;(x))f(xj)dx()d
=Z
Z
XL(;(x))(jx)f(x)dxd
On définit ainsi alors le coût a posteriori (;(x))comme étant la moyenne du coût
par rapport à la loi a posteriori :
(;(x)) =E(:jx)[L(;(x))] =Z
L(;(x))(jx)d
Il s’agit d’une fonction de x. D’où la proposition suivante :
Proposition 1. Le risque de Bayes r(;)n’est autre que la moyenne du coût a posteriori
(;(x))suivant la loi marginal f(x).
Démonstration : On sait que :
r(;) =Z
Z
XL(;(x))(jx)f(x)dxd
Orf(xj)() =(jx)f(x)
On a donc :
r(;) =Z
XZ
L(;(x))(jx)df(x)dx=Z
X(;(x))f(x)dx
L’estimateur de Bayes
Vu précédemment, la prise d’une décision (le choix d’un estimateur) va engendrer
un coût que l’on va quantifier à l’aide de la fonction perte. Logiquement, on cherchera
toujours la décision qui aura le coût le plus minimal. Autrement dit, on cherche une
décision qui minimise en moyenne la fonction de coût.
Définition 21. On appelle Estimateur bayésien (ouEstimateur de Bayes ) associé à un coût
Let à une distribution a priori , toute décision (x)qui minimise le risque de Bayes
r(;).
(x) = min
2Dfr(;)g
Remarque 2. L’estimateur de Bayes vari selon le coût Lassocié à celui-ci.
Onpeuttrouverplusieurssortesdefonctiondecoûtetnouspouvionsaussienconstruire
avecplusieursmanières.Maisilexistedesfonctionsdecoûtsusuelsutiliséspourl’infé-
rence bayésienne : le coût absolu, le coût quadratique et le coût 0-1. Parmi eux, le coût
quadratique (appelé aussi le coût classique) reste le plus utilisé du fait de sa convexité.
CHAPITRE II. INFÉRENCE BAYÉSIENNE 14
Dans le cas unidimensionnel
Définition 22. On appelle coût quadratique la fonction définie par :
L(;(x)) = ( (x))2
ThéorèmeII.1.1. L’estimateurdeBayes deassociéàlaloiapriori etaucoûtquadratique
est la moyenne a posteriori de :
(x) =E(:jx)[] =Z
(jx)d
Démonstration : Par définition, l’estimateur bayésien (x)minimise le risque de Bayes
r(;),quin’estautrequelamoyenneducoûtaposteriori (;(x))suivantlaloimar-
ginalf(x)(d’après la proposition II.1 ). Donc minimise le coût a posteriori (;(x)).
On a donc :
(;(x)) = E(:jx)[L(;)]
=E(:jx)
( (x))2
=E(:jx)
2
+ 2(x)E(:jx)[] +2(x)
=E2
(:jx)[] + 2(x)E(:jx)[] +2(x) +E(:jx)
2
E2
(:jx)[]
=
(x) E(:jx)[]2+Var(:jx)[]
Ceci n’est minimal que si (x) =E(:jx)[]
Dans le cas multidimensionnel
Dans le contexte multidimensionnel, c’est-à-dire où = (1;2;:::;n), l’estimateur de
Bayes s’exprime comme un vecteur :
(x) = (1(x);2(x);:::;n(x))
Où
i(x) =E(:jx)[i] =Z
i(ijx)di
Les(ijx)sontobtenusenintégrant (jx)surtouteslescomposantesde autrequei.
Remarque3. Lecalculdel’estimateurbayésiendevientdeplusenplusdifficiledèsque
le nombre de dimension augmente.
Avant les avancées des méthodes numériques2, c’était seulement les densités appar-
tenant aux familles conjugués (que nous verrons plus loin) qu’on pouvait trouver les
estimateurs bayésiens.
2. Que nous développerons dans le chapitre IV
CHAPITRE II. INFÉRENCE BAYÉSIENNE 15
2.2 Fiabilité de l’estimateur bayésien
Admissibilité
Définition 23. Un estimateur 1domineun autre estimateur 2si :
8;R(;2)R(;1)et90;R(;2)>R(;1)
Définition 24. Un estimateur est diteadmissible s’il domine tout estimateur autre que
lui.
ThéorèmeII.2.2. Siladistributionàpriori estpositivesur ,strictementpositivesurunsous
ensemble de , de risque de Bayes fini, et la fonction de risque R(;)est une fonction continu
depour tout, alors l’estimateur de Bayes est admissible.
Démonstration : Supposons que cet estimateur de Bayes ne soit pas admissible. Cela
implique qu’il existe un autre estimateur 0qui domine . Donc :
8;R(;)R(;0)et90;R(;)>R(;0)
Cela implique que :
Z
XR(;0)()d<Z
XR(;)()d,r(;0)<r(;0)
Cequiestimpossiblepuisquel’estimateurdeBayesdevraitminimiserlerisquedeBayes.
Doncdevrait être nécessairement admissible.
Remarque 4. Par ce théorème, on peut en déduire que l’estimateur de Bayes associé à
une loi de probabilité à priori et au coût quadratique est admissible.
Propriétés asymptotique
Comme de l’approche bayésien est le plus utilisées dans le cas des tailles d’échan-
tillon très petit ([ RC06]), il est rare qu’on se soucie du comportement asymptotique de
cetestimateur.Nousn’exposeronspasassezcepoint,puisquecetteapprochenedépend
pas vraiment des propriétés asymptotiques des échantillons. Mais en possède de très
bonne propriétés asymptotiques.
La distribution des estimateurs bayésiens (la loi à postériori) basées sur n variables
aléatoires i.i.d tend à devenir indépendant de la distribution a priori quand n tend vers
l’infini et c’est qui lui procure de très bonnes propriétés. En effet, il est consistant et
asymptotiquement normalement distribuées ( [ JR10], [RC06] et [IK85] )
2.3 Région de crédibilité
Pareillement à l’approche classique, l’approche bayésien peut nous proposer la no-
tion de région de confiance qui est plus naturelle3que son voisin fréquentiste. En effet,
danslecadrebayésien,larégiondeconfiance,représentel’imageconcrètedelaconfiance
que l’on peut avoir sur la valeur du paramètre concernée. Puisqu’au sens fréquentiste,
unerégiondeconfiance C,d’aprèsladéfintion16etlaremarque1,signifiequelaproba-
bilité pour que soit dans Csoit égale à 1 si l’on répétait l’expérience statistique un
grand nombre de fois. Une région de confiance classique n’aura donc de sens que pour
ungrandnombred’expériencealorsqu’avecl’approchebayésienneproposeunerégion
de confiance conditionnellement qu’aux observations.
3. RCO6
CHAPITRE II. INFÉRENCE BAYÉSIENNE 16
Définition 25. Une région Cdeest dite -crédible si et seulement si :
P(2CjX)1
Làdonc,larégiondeconfiancebayésienne(dite régiondecrédibilité pourledifférencié
de la région de confiance classique) s’interprète comme : la probabilité pour que soit
dansC,étantdonné(conditionnellement)l’expériencedéjàréalisé,estsupérieurà 1
. Ici donc, on n’a plus besoin de grand nombre d’expérience pour définir une région de
crédibilité. Le recourt aux méthodes spécifiques, comme la construction des fonctions
asymptotiquement pivotales ne sont donc plus nécessaires.
Partant de cette définition donc, on peut constater qu’il existe plusieur région -
crédible. Mais identiquement à celui de son voisin fréquentiste, vu dans la Remarque
1, celui dont le volume est le plus serrée sera la plus intéressante. Ce qui nous amène à
la définition suivante :
Définition 26. C
est dite une région HPD (Highest Postériori Density) ou PFDP (Plus
Forte Densité à Postériori) si et seulement si :
C
=; (jX)h
Avec
h = supfh;P(; (jX)hjX)1 g
C’est la région -crédible qui est la plus intéressante ([ JR10]).
Remarque 5. La détermination de la région de crédibilité ne peut se faire sans la loi à
postériori en main.
Le calcul de la région de crédibilité peut se faire soit par des méthodes analytiques (qui
nepeuvents’appliquerquedansdetrèsrarecas)soitpardesméthodesdesimulations(le
plus utilisé actuellement)
2.4 Loi a priori
Comme nous l’avions bien vu dans la section précédente, une fois la loi a priori
connue, la détermination de la loi a posteriori et de l’estimateur de Bayes se fait quasi
automatique. La détermination de la loi a priori est donc une étape très cruciale dans
l’inférence bayésienne. Mais dans certaine mesure, c’est aussi la plus difficile.
Les paramètres de la loi a priori sont appelés hyperparamètres . Ils sont déterminés par la
connaissanceapriorisurleparamètre,donccesontnormalementdesquantitésconnues.
Evidemment, dans la pratique, il est rare que l’information a priori soit suffisamment
précisepourconduireàunedéterminationexactedelaloiapriori,ausensoùplusieurs
lois de probabilité peuvent être compatibles avec cette information.
Il est donc nécessaire le plus souvent de faire un choix (partiellement) arbitraire4de loi
a priori, ce qui peut avoir un impact considérable sur l’inférence qui en découle.
Approche informative
Famille conjuguée
Définition 27. Une familleFde distribution de probabilité sur est diteconjuguée (ou
fermée par échantillonnage ) par une fonction de vraisemblance f(xj)si,
Pour toute distribution a priori 2F, la distribution a posteriori (:jx)2F.
Une des difficultés de l’approche bayésienne est le calcul de la loi a posteriori. Le
calcul est alors plus facile quand la loi a priori et a posteriori sont de la même forme.
L’induction de la forme de l’estimateur depuis la loi à priori se fait donc directement.
4. Mais restant toujours dans le rationnel ([ RC06] et [RC13])
CHAPITRE II. INFÉRENCE BAYÉSIENNE 17
Famille exponentielle
Définition 28. On appelle Famille exponentielle à s-paramètre, les familles de loi de dis-
tribution de la forme :
f(xj) = exp"sX
i=0i()Ti(x) B()#
h(x)
Où lesTisont des statistiques. ietBdes fonctions de .
Les familles exponentielles peuvent être décrites sous une autre forme dite Canonique
en reparamétrisant ieni:
f(xj) = exp"sX
i=0iTi(x) B()#
h(x)
Proposition 2. Soitf(xj)appartenant à une famille exponentielle. Alors une famille à priori
conjuguée pour f(xj)est de la forme :
;() =K(;) exp[ B()]où>0et
2
AvecK(;)une constante de normalisation.
La loi a posteriori sera de la forme :
(jx)/exp[(+T(x)) (+ 1)B()]
Démonstration : On sait que (jx)/f(xj)(), on a donc :
(jx)/exp [T(x) B()]h(x)K(;) exp[ B()]
/exp [T(x) + B() +B()]h(x)
/exp [(T(x) +) (1 +)B()] =T(x)+;1+()
Voiciuntableauillustrantunbrefexemplemontrantlesrelationsentrefamilleexponen-
tielle et famille conjuguée ( dite aussi naturelle).
Tableau 1. Exemples de lois conjuguées
f(xj)() (jx) E(:jx)[]
N(;2)N(; 2)N
x
2+
2;1
2+1
22!
x
2+
2
(n;)
( ; )
( +n; +x) +n
+x
P()
( ; )
( +x; + 1) +x
+ 1
B(n;)B^eta( ; )B^eta( +x; +n x) +x
+ +n
CHAPITRE II. INFÉRENCE BAYÉSIENNE 18
Approche non informative
Lorsqu’on ne dispose pas d’information assez fiable pour pouvoir choisir une loi a
priori, le recourt à des lois a priori dite non informative sont possible. Laplace fut le pre-
mier à utiliser cette technique.
Invariance par translation
L’invariance par translation est un peu similaire à la règle d’ équiprobabilité des évène-
ments élémentaires. Elle consiste à amener une loi a priori invariante par translation,
c’est-à-dire :
() =( 0)80
Pour quevérifie cette égalité, il faut qu’il soit une fonction constante. Donc, suit une
loi uniforme sur .
Lois a priori de Jeffreys
La règle de Jeffreys consiste construire une loi a priori non informative en utilisant l’in-
formation de Fisher : I(). En effet, l’information de Fisher représente la quantité d’in-
formation sur dans les observations. De ce fait, il parait évident que les valeurs de
pour lesquelles I()est plus grande doivent être plus probables a priori. Cela équivaut
à minimiser l’influence de la loi à priori, donc choisir le non informatif que possible.
La méthode de Jeffreys consiste donc à considérer des lois a priori de la forme :
Dans le cas unidimensionnel
()/p
I(), AvecI() =Eh
@2
@2ln(f(x))i
Dans le cas multidimensionnel
()/p
I(), AvecIij() =Eh
@2
@i@jln(f(x))i
eti;j= 1;2;:::;n
Loi a priori impropre
La loi a priori étant par définition une loi de probabilité de , donc son intégrale sur
toutsonensembledevraitêtreégaleà1.Ordanscertainscas,commedel’approchenon
informative, la loi a priori peut ne pas respecter cette définition. C’est-à-dire :
Z
()d=1
Dans ce cas, on dit que est une distribution a priori impropre.
Même si la loi a priori est impropre, elle peut être amenée à des réponses bayésiennes.
En d’autre terme, il se pourrait que ()ne soit pas intégrable mais que f(xj)()le
soit. Et que, de ce fait, la loi a posteriori peut être trouvée.
CHAPITRE II. INFÉRENCE BAYÉSIENNE 19
Modèle hiérarchique
Rappelons que les hyperparamètres (paramètres de la distribution à priori ) de-
vraient être des valeurs connu à partir de l’information à priori. Mais on peut encore
pousserl’approchebayésienneencoreplusloin,enconsidérantleshyperparamètresin-
connu et suivent pour leur part aussi, leur propre distribution à priori. Donnant ainsi
donc un modèle hiérarchique . La loi à priori se décompose donc et prend la forme sui-
vante :
() =Z
1(j1)1(1)d1
Où1est appelé hyperparamètre de niveau 1.
La distribution à postériori s’exprime alors comme suit :
(jx) =f(xj)()
f(x)=f(xj)R
1(j1)1(1)d1R
R
1f(xj)(j1)1(1)d1
En généralisant l’approche, on a la définition suivante.
Définition29. Onappelle modèlebayésienhiérarchique unmodèlestatistiquebayésienoù
la loi à priori est décomposé en un ou plusieurs distribution conditionnelles :
(j1);1(1j2);:::;k(kjk+1);k+1(k+1)
Oùiest appelé hyperparamètre de niveau i.
Cette approche tend à modéliser le manque d’information sur les hyperparamètres
de la distribution à priori par une nouvelle loi sur ces paramètres.
Remarque 6. En utilisant ce modèle, le calcul analytique de la distribution à postériori
devient incontestablement très difficile en fonction que le niveau augmente. Le recourt
à des méthodes de calcul numérique approximative est inévitable.
2.5 Exemple
Exemple 3. Reprenons l’exemple 2. Mais supposons maintenant qu’on serait d’obtenir
des informations (avis d’expert, connaissance obtenu dans le passé, …) pouvant aider
à donner une approximation (moyenne et variance 2) de la proportion à estimer.
On aimerait donc combiner l’information apportée à priori et celle des observations.
Il nous faut donc trouver une loi à priori adapté à. Cette loi a priori devrait avoir
comme support [0; 1]puisqu’il s’agit d’une probabilité. Nous allons prendre5comme
distribution à priori la loi beta de première espèce 1(a;b).
Les hyperparamètres aetbsont des valeurs qui doivent être connus. Ils sont trouvés
grâce aux valeurs et2données par les experts. Puisque 1(a;b), on a :
E[] =a
a+b= et Var [] =ab
(a+b)2(a+b+ 1)=2
Commenouslemontrelatableau1,ladistributionàpostérioriestaussiunedistribution
betadepremièreespèce 1(a+x;b+n x)oùx=Pn
i=1xi.Etl’estimateurbayésienassocié
s’exprime :
n=a+Pn
i=1xi
a+b+n
5. Voir [MHM15] pour la justification de ce choix
CHAPITRE II. INFÉRENCE BAYÉSIENNE 20
Rappelons que l’EMV associé à cette expérience s’exprime :
^EMV =1
nnX
i=1xi
On remarque que les deux estimateurs se rapprochent lorsqu’on a beaucoup de don-
nées. C’est lorsque les données sont en petites quantités que leur différence peut être
très importante. Et c’est là qui est intéressant avec l’estimateur de Bayes. Même avec
peud’observation,cettemanquepeutêtrecomplétéeparl’informationàpriori,puison
peutendéduireunestimateuradmissible.Etdanscecas,l’estimateurdeBayesestplus
convainquant que l’EMV.
Maintenant, nous allons supposer qu’on n’a aucune information sur cette proportion.
Danscecas,onestdansunesituationnoninformative,cequinousamèneàl’ignorance
complète sur cette proportion . L’incapacité d’évaluer à priori nous permet de consi-
dérer que notre distribution à priori suit une loi uniforme sur [0; 1], et qui n’est autre
que la loi 1(1;1). Et l’estimateur bayésien associé s’exprime donc :
n=1 +Pn
i=1xi
2 +n
Qui est vraiment très proche de l’EMV.La loi à priori étant non-informative, les infor-
mations dont nous disposons ne viennent que des observations. La loi à postériori est
très proche de la distribution d’échantillonnage qui ramène donc à la même inférence,
d’après le principe de vraisemblance.
Exemple4. Ceciestuneapplicationdel’exemple3.Danscetteexpérience,onlanceune
pièce de monnaie 10 fois etdès qu’on observe une pile on note 1 et face onnote 0. Nous
avonsdoncl’échantillon X=X1;:::;X 10oùchaquevariablealéatoire Xi20; 1.Lavrai-
semblance s’exprime de la façon suivante :
P(X=k) =Ck
npk(1 p)n k
Supposonsquelapièceestnontruquée(c’est-à-direque = 0;5)maisqu’onnelesait
pas.Etqu’après10lancées,onaeul’échantillonsuivant: X=1;0;1;0;1;1;0;1;1;1soit7
piles et 3 faces. Ici on a donc :
P(X= 7) = 120p7(1 p)3
.Approche classique :
Par la méthode du maximum de vraisemblance on a :
^p=1
nnX
i=1xi= 0;7
Sur R (pour le code, voir Annexe B 1.), on a les résultats suivants :
-^p= 0;7
-IC0;95= [0;347 ; 0;933]
.Approche bayésienne :
Ondoitd’abordpasserparladéterminationdelaloiàpriori (p).Puisqu’ils’agitd’éva-
luerpindépendamment des observations. Naturellement, on peut :
CHAPITRE II. INFÉRENCE BAYÉSIENNE 21
– Soit ignorer complètement quelle valeur cette proportion pourrait avoir et donner
les mêmes chances à toutes les valeurs possibles de p(à priori non-informatif). On aura
donc comme loi à priori 1(p) = (1; 1).
-Soitapporternotrecroyanceàpriorisur p(àprioriinformatif).Puisqu’ils’agitd’une
pièce de monnaie, on peut supposer qu’à priori, en moyenne, cette proportion vaut 0,5.
Et on prendra une marge d’erreur de 0,01. On aura donc comme loi à priori 2(p) =
(12; 12).
Maintenant on passe à la détermination de la loi à postériori. Puisque ce sont tous
deux des lois conjuguées, on a :
-1(pjX) = (8; 4)comme loi à postériori de 1(p)
-2(pjX) = (19; 15)comme loi à postériori de 2(p)
Le tableau suivant résume l’inférence faite :
Tableau 2. Comparaison ente estimation bayésienne et classique
Situation Non informative informative
Loi à priori 1(p) = (1; 1) 2(p) = (12; 12)
Loi à postériori 1(pjX) = (8; 4) 2(pjX) = (19; 15)
Estimateur de Bayes 1=8
8 + 40;6672=19
19 + 150;559
Intervalle de crédibilité à 0,95 [m1;M1] = [0;39; 0;89] [m2;M2] = [0;392; 0;719]
EMV EMV ^= 0;7 EMV ^= 0;7
Intervalle de confiance à 0,95 [0;347; 0;933] [0;347; 0;933]
La détermination de l’intervalle de crédibilité a été faite sous R (voir Annexe B 2.).
Remarque 7. Remarquons bien qu’en situation non-informative, l’estimateur de Bayes
est vraiment proche de l’EMV. Mais lorsqu’on a ajouté notre connaissance à priori que
pourraitavoirlaproportion,onvoitquel’estimationobtenuparlaméthodebayésienne
est plus proche de la valeur théorique de p.
En termes d’estimation par intervalle de confiance, l’intervalle de crédibilité est plus
dérrée que celle de l’intervalle de confiance classique.
On voit bien que même avec peu de donner, cette manque peut-être compléter par les
informations6à priori
6. Evidemment, nécessite la vérification de la fiabilité de l’information apportée avant son utilisation.
Laraisonpourlaquelledelademanded’avisd’expertdudomaineouutilisationderésultatdesanciennes
études.
Chapitre III
Chaîne de Markov
1 Notion de base
1.1 Définitions et propriétés
Avant de présenter le formalisme des chaînes de Markov, on rappellera brièvement
quelques notions sur les processus.
Définition 30. Unprocessus aléatoire est un phénomène qui se déroule dans le temps et
qu’on peut repérer son évolution à chaque instants tpar une variable aléatoire Xt.Xt
désigne l’état du phénomène à l’instant t.
Définition31. Onappelle processusstochastique unprocessusaléatoiredonttouslesétats
du phénomène appartiennent à un même ensemble.
De manière plus formelle, un processus stochastique est une suite de variable aléatoire
indicée par un ensemble T(le temps) en général infini (dénombrable ou non), à valeurs
dans un espace mesurable, et toutes définies sur le même espace de probabilité.
PropriétéIII.1.3. (PropriétédeMarkov"faible" ):Pourunprocessusquis’estdéroulédans
le passéXt(tt0), l’état du phénomène futur Xt(t > t 0)ne dépend que de l’état du
phénomène présent Xt0). Là, on parle aussi de processus à mémoire limité .
Définition 32. On appelle processus markovien un processus stochastique possédant la
propriété de Markov. Si de plus, l’ensemble des états E du phénomène est discret, alors
c’est une chaîne de Markov .
Unesuite (Xn)(n0)devariablealéatoireàvaleursdansunensembleauplusdénombrable
Eestunechaîne de Markov d’espaced’étatsEsietseulementsipourtout n2N;in2Etelsque
P(Xn+1=in+1jXn=in;:::;X 0=i0) =P(Xn+1=in+1jXn=in)
Si de plus le terme de droite est indépendant de n, alors la chaîne de Markov est dite
homogène . C’est-à-dire ,8n2Neti;j2E:
P(Xn+1=jjXn=i) =P(X1=jjX0=i)
Définition 33. On appelle probabilité de transition de i vers j la probabilité de passer de
l’étatià l’instantnvers l’étatjà l’instantn+ 1, et on la dénote pij.
pij=P(Xn+1=jjXn=i)
On note par 0laloi( oudistribution ) deX0tel que :
8i2E 0(i) =P(X0=i)
22
CHAPITRE III. CHAÎNE DE MARKOV 23
Proposition 3. Soit(Xn)(n0)une chaîne de Markov. Pour (i0;:::;in)2En+1
On a :
P(Xn=in;:::;X0=i0) =0(i0)n 1Y
k=0pikik+1
Démonstration :
P(Xn=in;:::;X0=i0) = P(Xn 1=in 1;:::;X0=i0):P(Xn=injXn 1=in 1;:::;X0=i0)
D0apreslad efinition 32 = P(Xn 1=in 1;:::;X0=i0):P(Xn=injXn 1=in 1)
=P(Xn 1=in 1;:::;X0=i0):pinin 1
:::
=P(X0=i0)pi0i1pi1i2:::pinin 1
=0(i0)n 1Y
k=0pikik+1
Définition 34. On appelle Matrice de transition la matriceP= (pij)(i;j2E):
0
BBB@pi0i0pi0i1pi0in
pi1i0pi1i1pi1in…………
pini0pini1pinin1
CCCA
UnechaînedeMarkovestcaractériséeparlaloi 0deX0etparsamatricedetransition.
Proposition 4. Toute matrice de transition d’une chaîne de Markov vérifie les propriétés sui-
vantes :
(i)Pour tout couple (i;j)deE,0pij1
(ii)Pour touti2E, on aP
j2Epij= 1
Toute matrice vérifiant ces deux propriétés sont appelés matrice stochastique .
Démonstration :
(i)Très évident puisque 8(i;j)2E2, lespijsont des probabilités.
(ii)Comme étant la somme des probabilités sur toutes les valeurs possibles d’une
variable aléatoire, il est aussi évident que cette somme vaut 1.
1.2 Evolution temporelle
Une fois la matrice de transition et la distribution initiale des états spécifiés, il est
possible de calculer l’évolution de cette distribution de probabilité avec le temps.
CHAPITRE III. CHAÎNE DE MARKOV 24
Proposition5. Soit(Xn)(n0)unechaînedeMarkovdematricetransition P,etsoit0laloide
X0etncelle deXn. Alors8n2N, on a :
n=0Pn
Démonstration : Nous allons procéder par récurrence :
Notons parR(n)la proposition n=0Pn. Soitj2E. On sait que :
1(j) =P(X1=j) =X
i2EP(X1=jjX0=i):P(X0=i)
=X
i2Epij0(i)
= (0P)j
D’où1=0P, doncR(1)vrai.
Maintenant Supposons que R(n)est vrai et montrons que R(n+1)l’est aussi. Soit j2E.
On sait que :
n+1(j) =P(Xn+1=j) =X
i2EP(Xn+1=j)jXn=i):P(Xn=i)
=X
i2Epijn(i)
=X
i2Epij(0Pn)i
= (0Pn+1)j
D’oùn+1=0Pn+1, doncR(n+1)vrai.
Ainsi8n2N, on an=0Pn
Soit(Xn)n0est une chaîne de Markov, dont l’ensemble des états est Eet la matrice
de transition P= (pij)(i;j)22E. Pourn0eti;j2E, on désigne par pij(n)la probabilité,
partant de l’état ià l’instant 0 et d’être à l’état jà l’instantn. En d’autres termes :
8i;j2E,8n>0,p(n)
ij=P(Xn=jjX0=i)
1.3 Relation de Chapman-Kolmogorov
La relation de Chapman-Kolmogorov nous dit que : pour qu’une chaîne de Markov
a pu passer de l’état ià l’étatjenm+nétapes, il a bien fallu en métapes d’aller de i
à un certain état kpuis ennétapes d’aller de cet état kàj. Autrement dit : 8i;j2E,
8n;m2N
P(Xm+n=jjX0=i) =X
k2EP(Xm=kjX0=i):P(Xn=jjX0=k)
Ou encore
p(m+n)
ij =jjX0=i) =X
k2Ep(m)
ikp(n)
kj
Démonstration : C’est immédiat grâce à l’associativité du produit matricielle :
P(m+n)=Pm+n=PmPn=P(m)P(n)
CHAPITRE III. CHAÎNE DE MARKOV 25
Proposition 6. Soientn0;r1deux entiers. Alors :
P(Xn+r=in+r;:::;Xn+1=in+1jXn=in;:::;X0=i0) =n+rY
k=n+1pik 1ik
Démonstration :
P(Xn+r=in+r;:::;Xn+1=in+1jXn=in;:::;X0=i0) =
=P(Xn+r=in+r;:::;Xn+1=in+1jXn=in)
=Qn+r
k=n+1pik 1ik
NouspouvonsformaliserlapropriétédeMarkovfaiblesousuneversionéquivalente
etplusgénérale.EnNotantpar A unévènementappartenantàlatribudupassé Fn 1=
(Xn;:::;X0)etparA+unévènementappartenantàlatribudufutur (Xn+1;Xn+2;:::).
Nous obtenons :
SiP(A;Xn=i)>0, alors :
P(A+jA ;Xn=i) =P(A+jXn=i)
On voit que cette dernière est plus générale que celle de la propriété de Markov faible.
1.4 Propriété de Markov forte
SoientAunévènementet Zunevariablealéatoire,onnoterapar Pi(A) =P(AjX0=i)
etEi[Z] =E(ZjX0=i). Ici, on s’intéressera plus sur la chaîne de même matrice de
transition que (Xn)et de loi initiale i. On l’appelle aussi « chaîne issue de i ».
SoitTun temps d’arrêt adapté à la suite (Xn)n0. Rappelons qu’un événement A2FT
(tribu engendrée par la réunion de Fn) s’il a la propriété suivant :
8n2N;A\T=n2Fn
Théorème III.4.4. ( Propriété de Markov Forte )
SoitTun temps d’arrêt à valeur dans [0;1[adapté à la chaîne de Markov. Soient iun état,
A2FTetP(T <1;A;Xn=i)>0. Alors :
P(XT+1=j1;:::;XT+r=jrjT <1;A;XT=i) = P(X1=j1;:::;Xr=jrjX0=i)
Démonstration : Posons par A0=fT <1;A;XT=ig,B=fXT+1=j1;:::;XT+r=jrg
et pourn0,Bn=fXn+1=j1;:::;Xn+r=jrg
P(A0\B)=P
n0P(T=n;A;XT=i;B)=P
n0P(T=n;A;Xn=i;Bn)
=P
n0P(BnjT=n;A;Xn=i)P(T=n;A;Xn=i)
PuisquefT=n;Ag2FT, alors donc
P(BnjT=n;A;Xn=i) =P(BnjXn=i)et on a :
CHAPITRE III. CHAÎNE DE MARKOV 26
P(A0\B)=P
n0P(BnjXn=i)P(T=n;A;Xn=i)
=P
n0P(Xn+1=j1;:::;Xn+r=jrjXn=i)P(T=n;A;Xn=i)
=P
n0pij1pj1j2:::pjr 1jrP(T=n;A;Xn=i)
=pij1pj1j2:::pjr 1jrP
n0P(T=n;A;Xn=i)
=pij1pj1j2:::pjr 1jrP(T <1;A;XT=i)
=pij1pj1j2:::pjr 1jrP(A0)
D’après la formule de Bayes :
pij1pj1j2:::pjr 1jr=P(A0\B)
P(A0)=P(BjA0)
2 Classification des états
Soit une chaîne de Markov de matrice de transition Pet de distribution initiale pi0.
2.1 Irréductibilité
Définition 35. On dit que l’état jestaccessible à partir de l’état i, et on note i j, s’il
existen0tel quen(j) =0(i)Pnc’est-à-dire p(n)
ij>0.
Définition 36. On dit que les états ietjcommuniquent si on a à la fois i jetj i, et
on la notei!j
Définition 37. On dit qu’une chaîne de Markov est irreductible si tous ses états commu-
niquent entre eux. C’est-à-dire 8i;j2E,i!j.
On peut formuler l’irréductibilité comme tel :
8i;j2E;8n0;9n02N=P(Xn+n0=jjXn=i)>O
Remarque :n0dépend des états ietj
2.2 Récurrence
8i2EOn désigne par Tiletemps d’atteinte de l’étatià partir de l’instant 1 :
Ti= inf
n1fXn=ig
On notera par Nile nombre de fois que la chaîne visite l’état ien comptant le point de
départ :
Ni=X
n2N1Xn=i
Définition 38. On dit qu’un état i2Eestrécurrent(ou persistent) si, partant de cet état,
la chaîne y reviendras presque sûrement (p.s). En d’autre terme, la probabilité que la
chaîne de Markov retourne à son état initiale en un temps fini est égale à 1 :
Pi(Ti<1) =P(Ti<1jX0=i) = 1
CHAPITRE III. CHAÎNE DE MARKOV 27
Si de plus, l’espérance du temps de retour est fini, alors l’état iest ditrécurrent positive .
C’est-à-dire :
Ei(Ti) =E(TijX0=i)<1
Aucascontraire,s’iln’estpasrécurrentalorsonditqu’ilest transient(outransitoire) .C’est-
à-dire :
Pi(Ti<1) =P(Ti<1jX0=i)<1
Remarque : Nous pouvons tout de suite remarquer que dès que la chaîne parvient à
revenir à son état initial, Ti<1etNi>1. En d’autre terme :
Pi(Ti<1) =P(Ti<1jX0=i) =P(Ni>1jX0=i) =Pi(Ni>1)
UnechaînedeMarkovestdite transiente (respectivement récurrente )sitoussesétatssont
transients (respectivement récurrents ).
Proposition 7. Les propriétés suivantes permettent de caractériser les états récurrents et tran-
sients :
(i)Siiest transient alors Pi(Ni<1) = 1,P
n2Np(n)
ii<1et conditionnellement à fX0=ig,
Niest de loi géométrique de paramètre Pi(Ti=1).
(i)Siiest transient alors Pi(Ni=1) = 1etP
n2Np(n)
ii=1
Démonstration : L’intuition est comme ceci : si on revient presque sûrement au moins
une fois, on peut répéter l’argument une fois qu’on est revenu, et on revient donc au
moins deux fois et ainsi de suite.
Soitk2Net supposons que Pi(Ti<1)6= 0:
Pi(Ni=k) = Pi(Ni=kjTi<1)Pi(Ti<1)
=Pi" Ti 1X
n=01Xn=i+1X
n=Ti1Xn=i=kjTi<1!#
Pi(Ti<1)
=Pi"
1 +1X
n=01XTi+n=i=kjTi<1!#
Pi(Ti<1)
=Pi"1X
n=01XTi+n=i=k 1jTi<1#
Pi(Ti<1)
On peut maintenant appliquer la propriété de Markov forte : conditionnellement à
(Ti<1), la chaîne (XTi+n)à même loi que la chaîne partant de i. Donc :
Pi(Ni=k) = Pi"1X
n=01XTi+n=i=k 1#
Pi(Ti<1)
=Pi(Ni=k 1)Pi(Ti<1)
Ainsi ,8k2N,Pi(Ni=k) =Pi(Ni=k 1)Pi(Ti<1). D’où
La suite (Pi(Ni=k))k2Nest une suite géométrique de raison q=Pi(Ti<1)
(i)Siiest transient, on a q=Pi(Ti<1)<1et on en déduit que
Pi(Ni=k) =Pi(Ni= 1)ak 1,avecPi(Ni= 1) = Pi(Ti=1) = 1 q
CHAPITRE III. CHAÎNE DE MARKOV 28
Donc, conditionnellement à X0=i,Nisuit bien une loi géométrique de paramètre
Pi(Ti=1). De ce fait, elle possède une espérance finie. Donc :
Pi(Ni<k) = 1)X
n2Np(n)
ii<1
Celanousditaufinalquesiiesttransientalorslachaîneypassepresquesûrementqu’un
nombre fini de fois.
(ii)Siiest récurrent, on a q=Pi(Ti<1) = 0et donc Pi(Ni2N) = 0, d’où
Pi(Ni=1) = 1)X
n2Np(n)
ii=1
Etcelanousramèneàdirequesi iestrécurrentalorslachaîneypassepresquesûrement
en un nombre infini de fois.
Proposition 8. Sii!jalors ils ont la même classe de récurrence.
Démonstration :
Supposons que l’état iest récurrent.
Commei!jalors9n1; n22Ntels quep(n1)
ij>0etp(n2)
ji>0. Alors :
X
np(n1+n+n2)
jjX
np(n1)
ijp(n)
iip(n2)
ji=p(n1)
ijp(n2)
jiX
np(n)
ii=1
Ainsi l’état jest aussi récurrent.
Maintenant supposons que l’état jest transient.
Commei!jalors9n1; n22Ntels quep(n1)
ij>0etp(n2)
ji>0. Alors :
1>X
np(n1+n+n2)
jjX
np(n1)
ijp(n)
iip(n2)
ji=p(n1)
ijp(n2)
jiX
np(n)
ii
Ainsi l’état iest aussi transient.
2.3 Périodicité
L’étude ici se focalise sur les conditions de temps qui sépare deux retours au même
état est ou n’est pas multiple d’un temps minimum. Pour cela, on introduit la notion de
période.
Définition 39. Soiti2E. On appelle période de i, et on le note d(i), le P.G.C.D de tous
les entiersn1pour lesquels p(n)
ii>0.
d(i) =pgcd
n1;p(n)
ii>0
Sid(i)>1alors on dit que iestpériodique de période d(i).
Sid(i) = 1alorsonditque iestapériodique .Onditqu’unechaînedeMarkovest apériodique
si tous ses états sont apériodiques.
Proposition 9.
(i) Sip(1)
ii>0alors l’étatiest apériodique.
(ii) Si deux états i et j communiquent alors ils ont la même période.
CHAPITRE III. CHAÎNE DE MARKOV 29
Démonstration :
(i)Soit l’étatitel quep(1)
ii>0. On a :
d(i) =pgcd
n1;p(n)
ii>0
=pgcd
1;pgcd
n>1;p(n)
ii>0
= 1D’où l’apériodicité.
(ii)Comme les états ietjcommuniquent, alors :
9M;N2Ntels quep(M)
ij>0etp(N)
ji>0. Et8k1on a :p(M+nk+N)
iip(M)
ijh
p(k)
jjin
p(N)
ji
Donc8k1tel quep(k)
jj>0, on a :p(M+nk+N)
ii>0,8n1. On a donc d(i)divise
M+nk+N,8n1, ce qui implique que d(i)divisek, qui est multiple de d(j), donc
d(i)divised(j). De la même façon on montre que d(j)divised(i).
Ainsid(i) =d(j).
3 Comportement asymptotique
Soit(Xn)n0estunechaînedeMarkovirréductiblesurunensembledénombrable E
dematricedetransition P.Danscequisuit,ontravailleraquesurleschaînedeMarkov
homogène.
3.1 Distribution stationnaire
Définition 40. Une distribution de probabilité surEest dite stationnaire si :
8j2E j=X
i2Eipij
De façons plus général, une mesure surEsatisfaisant j=P
i2Eipij8j2E, est
appeléeune mesure invariante de la chaîne.
On désignera l’expression ci-dessous comme étant le nombre moyen de passage de
la chaîne sur un état ientre deux passages d’un autre état k:
'(k)
i=Ek TkX
n=11Xn=i!
(III.1)
Proposition 10. Si de plus la chaîne est récurrente, alors on a 8k2E:
(i)'(k)
k= 1
(ii)'(k)est une mesure invariante
(iii)8i2E, on0<'(k)
i<1
(iv)est l’unique mesure invariante telle que '(k)
k= 1
Démonstration :
(i)Pour 1n<Tk,Xn6=ket queXTk=k, d’où l’évidence.
(ii)Soienti;j2E, on a :
'(k)
i=Ek 1X
n=11Xn=i;nTk!
=1X
n=1Pk(Xn=i; nTk)
=X
j2E1X
n=1Pk(Xn 1=j; nTk)pji=X
j2Epji1X
m=0Pk(Xm=j; mTk 1)
=X
j2EpjiEk Tk 1X
m=01Xm=j!
=X
j2EpjiEk TkX
m=11Xm=j!
=X
j2Epji'(k)
j
CHAPITRE III. CHAÎNE DE MARKOV 30
D’où l’invariance de la mesure.
(iii)Comme'(k)est une mesure invariante, alors 8n0:
'(k)
i=X
j2E'(k)
jPj(Xn=i) (III.2)
(i)et (III.2) impliquent que : 8j2Eon a 1 ='(k)
k'(k)
jPj(Xn=i)Et comme la chaîne
est irréductible, alors 9n=Pj(Xn=i)>0
D’où8j2E; '(k)
j<1
D’une autre part : 9nj'(k)
iPk(Xn=i)>0, par l’irréductibilité de la chaîne.
(iv)Soit une autre mesure invariante telle que k= 1. Alors :
k=X
i6=k kpij+pkjpkj;8j2E
En minorant iparpkion en déduit que :
kX
i6=kpkipij+pkj=Pk(X2=j;2Tk) +Pk(X1=j;1Tk)
Et par récurrence nous obtenons :
kn+1X
m=1Pk(Xm=j; mTk) =E0
@min(n+1;Tk)X
m=11Xm=j1
A
Quandn!1, k'(k)
j,8j2E. Posons par la mesure défini par : = '(k).
La mesure est une mesure invariante et on a k= k '(k)
k= 0.
k=P
j2EjPj(Xn=k) = 0;8n0. La chaîne étant irréductible, alors nécessaire-
ment
8j2E; j= 0) ='(k)
Théorème III.1.5. Les propriétés suivantes sont équivalentes :
(i)Il existe une distribution stationnaire
(ii)9k2Etels que Ek(Tk)<1
(iii)8k2E,(ii)est vrai.
Démonstration :
(ii))(i): Supposons que (ii)est vérifié, donc kest récurrent positive donc la chaîne
étant irréductible, est aussi récurrente positive. Donc d’après la proposition précédente
'(k)est l’unique mesure invariante qui est égal à 1 en k. Or :
X
j2E'(k)
j=Ek TkX
n=1X
j2E1Xn=j!
=Ek TkX
n=11!
=Ek(Tk)<1
On peut donc définir la mesure j='(k)
jP
j2E'(k)
jqui est une mesure de probabilité inva-
riante donc une distribution stationnaire.
(i))(iii): Soitune distribution stationnaire et k2E. Alors
définit par
j=j
k
CHAPITRE III. CHAÎNE DE MARKOV 31
est une mesure invariante et que
j= 1.D’après la proposition précédente
j='(k). Et
donc :
Ek(Tk) =X
j2E
j=X
j2Ej
k=P
j2Ej
k=1
k<1
(iii))(ii): Évident.
Et donc, si la distribution stationnaire existe, alors son unicité découle de l’unicité de
'(k)
jet on a :
8i2E ; k=1
Ek(Tk)
Nouspouvonsainsiconclurequel’irréductibilitéetlarécurrencepositivesontlescondi-
tions nécessaires et suffisantes de l’existence et l’unicité de la distribution stationnaire.
3.2 Convergence
Définition 41. Une chaîne de Markov est dite ergodique si elle est à la fois irréductible,
apériodique et récurrente positive.
ThéorèmeIII.2.6. Soit(Xn)n0unechaînedeMarkovergodiqueetsoit sonuniquedistribu-
tion stationnaire. Alors pour toute distribution initiale , on a :
lim
n!1P(Xn=i) =i
Démonstration :
Soit(Xn;Yn)n0unechaînedeMarkovsur EE,dedistributioninitiale =
etde
probabilité de transition p
(i;j);(k;l)=pikpjl.
Supposons que XnetYnsont deux chaînes indépendantes de distributions initiales res-
pectivesetet toutes les deux apériodiques, irréductibles, récurrentes positives et de
même matrice de transition P. Mais qu’en est-il de (Xn;Yn)n0?
Soit un état k. Notons par :
k=fn2NjPk(Xn=k)>0g
La relation de Chapman-Kolmogorov implique que si n;m2 k)n+m2 k. La
chaîne (Xn)n0étant apériodique : pgcd( k) = 1
Nous allons avancer que 9n0j8tn0; t2 k(1)
Pour le montrer, supposons d’abord que 9n; m2 k)qui sont premiers entre eux. Par
le théorème de Bézout, 9p; q2Njnp mq=1(On peut supposer que cela vaut 1).
Posonsn0=qnmet soit l’entier rtel que 0rn, et nous avons :
n0+r=n0+r(np mq) =mq(n r) +rpn2 k
Alors pour t=n0+r+ns;0rnon at2 k.Ce qui démontre (1).
Fixons maintenant des états i;j;k;l2E. Par l’irréductibilité de P,9r2Ntel que
Pi(Xr=k)>0.
D’après ce que nous venons de voir, nous avons donc 8nn0:
Pi(Xr+n=k)Pi(Xr=k)Pk(Xn=k)>0
)Pi(Xr=k);8nn0+r
CHAPITRE III. CHAÎNE DE MARKOV 32
En suivant le même raisonnement, on peut affirmer aussi que :
9m0r2N)Pi(Xr=k)8mm0+s
Donc, on peut déduire que
9t0= max(m0+s; n 0+r);P(i;j)[(Xt;Yt) = (k;l)]8tt0
Et cela implique l’irréductibilité de la chaîne (Xn;Yn)n0. Ceci implique aussi que
l’état (i;j)estapériodiquepuisque pgcd(ftjtt0g) = 1.Parsonl’irréductibilité,elleest
aussi donc apériodique.
Commelachaîneadmetaussiévidementladistribution
commedistributioninva-
riante, alors le théorème III.1.5 nous dis qu’elle est récurrente positive.
Considérons l’ensemble D=f(i;i); i2EgEEEet notons par TDle temps de
premier passage sur cet ensemble. Où TD= inf(n>0jXn=Yn).
Nous allons montrer d’abord que XnetYnont même loi pour nTD. Pour ce faire,
nous allons considérer le processus (Zn)n0défini par :
Zn(
Xn; pournTD
Yn; pourn>T D
On a :
P(Z0=i0;:::;Zn=in) =i0n 1Y
k=0pikik+1;8(i0;:::;in)2En+1; n0
Donc d’après la Proposition 3, (Zn)n0est une chaîne de Markov de distribution ini-
tialeet matrice de transition P, et est donc égale à la loi à (Xn)n0. AinsiXnetZnont
même loi pour n2N. C’est-à-dire, pour n2N;P(Zn=i) =P(Xn=i).
Commei=P(Yn=j), donc :
P(Xn=i) j=P(Zn=i) P(Yn=i)
=P(Zn=i;nTD) +P(Zn=i;n>TD)
P(Yn=i;nTD) P(Yn=i;n>TD)
=P(Zn=i;nTD) P(Yn=i;nTD)
P(Zn=i;nTD)
P(nTD)
Comme (Xn;Yn)n0est récurrent positif, alors TDest presque sûrement fini, donc :
lim
n!1P(nTD) = 0
d’où :
lim
n!1jP(Xn=i) ij= 0)lim
n!1jP(Xn=i)j=i
Proposition 11. Soit(Xn)n0une chaîne de Markov ergodique, '(k)la mesure vu dans III.1,
f:E!Rune fonction et N(n)
ila somme telles que :
(i)N(n)
i=Pn
l=01Xl=i
(ii)P
i2Ejf(i)j'(k)
i<1
CHAPITRE III. CHAÎNE DE MARKOV 33
Alors pour n’importe quelle distribution initiale on a :
1
N(n)
inX
l=1f(Xl) !X
i2Ef(i)'(k)
i
n!1
Démonstration :SoitTk= 1; 2; 3;:::lestempsderetoursuccessifsenl’état ketquela
chaîne parte de cette état :
k+1= inf
n kfXn=kg
Posons par :
Up=l= p+1X
l= p+1f(Xl)
Supposons que fsoit positive sur E. On a donc :
E[Up] = Ek2
4l= p+1X
l= p+1f(Xl)3
5
=Ek"TkX
l=1f(Xl)#
; par l0homog eneitedelacha ^ine
=Ek"TkX
l=1X
i2Ef(i)1Xl=i#
=X
i2Ef(i)E0"TkX
l=11Xl=i#
=X
i2Ef(i)'(k)
i
Ceci est fini d’après l’hypothèse (ii).
On peut remarquer que fUpgp1sont des variables aléatoires indépendantes. Et par la
loi des grands nombre, on a :
1
nnX
p=1Up !X
i2Ef(i)'(k)
i
n!1
Remarquons que :
N(n)
kn N(n)
k+1
Car si la chaîne de longueur a pu visité N(n)
ifois l’étatk, donc il est nécessaire que N(n)
k
( la longueur minimale de visite de l’état kenN(n)
kfois) soit inférieur à la longueur de
la chaîne(qui est égale à n). Et si la chaîne n’a pas pu visiter N(n)
k+ 1fois l’étatk, donc
ellen’apasencoreatteintlalongueurminimalepourpouvoirvisiterl’état kenN(n)
k+ 1
fois(qui est égale à N(n)
k+1).
On a donc :P
N(n)
k
l=1f(Xl)
N(n)
kPn
l=1f(Xl)
N(n)
kP
N(n)
k+1
l=1f(Xl)
N(n)
k
CHAPITRE III. CHAÎNE DE MARKOV 34
Orlachaîneestergodique,doncrécurrentepositive,donc N(n)
k!1quandn!1.Or:
lim
n!1P
N(n)
k
l=1f(Xl)
N(n)
k= lim
n!1P
N(n)
k+1
l=1f(Xl)
N(n)
k=X
i2Ef(i)'(k)
i
d’où
1
N(n)
inX
l=1f(Xl) !X
i2Ef(i)'(k)
i
n!1
Enfinsifestdesignearbitraireonpose f+= max(0;f)etf = max(0; f)etonobtient
quef=f+ f , et par l’hypothèse (ii)cette différence est fini. On arrive donc à :
1
N(n)
inX
l=1f(Xl) !X
i2Ef(i)'(k)
i
n!1
Théorème III.2.7. (Théorème ergodique)
Soit(Xn)n0unechaînedeMarkovergodiquededistributionstationnaire ,etsoitfunefonction
telle que :
E[jf(Xk)j] =X
i2Ejf(i)j'(k)
i<1
Alors pour n’importe quelle distribution initiale on a :
1
nnX
l=1f(Xl) !X
i2Ef(i)(k)
i
n!1
Démonstration : En appliquant la proposition 11 pour la fonction constante 1, et parce
que la chaîne est ergodique, on a :
1
N(n)
inX
l=11(Xl) =n
N(n)
i !X
i2E1(i)'(k)
i=X
i2E'(k)
i<1
n!1
RappelonsquedanslethéorèmeIII.1.5, i='(k)
iP
i2E'(k)
iestl’uniquedistributionstation-
nairedecettechaîne.Donccequinousramèneà i= i.Onvoitque i/'i.Etpuisque
E[jf(Xk)j]<1, on a :
1
N(n)
inX
l=11(Xl) = !X
i2E1(i)'(k)
i<1
n!1
On en déduit que :
lim
n!11
nnX
l=1f(Xl) = lim
n!1N(n)
i
n1
N(n)
inX
l=1f(Xl) =P
i2Ef(i)'(k)
iP
i2E'(k)
i
Commei='(k)
iP
i2E'(k)
i, on obtient :
lim
n!11
nnX
l=1f(Xl) =X
i2Ef(i)i
CHAPITRE III. CHAÎNE DE MARKOV 35
3.3 Réversibilité
Définition 42. Une chaîne de Markov est dite réversible si sa matrice de transition Pad-
met un vecteur = (i)i2E2[0;1[Evérifiant la propriété suivante :
ipij=jpji8i;j2E
Tout vecteur 6= 0vérifiant cette propriété est appelé vecteur réversible par rapport à la
matrice de transition P.
Théorème III.3.8. SoitPune matrice stochastique et 2[0;1[Eun vecteur non nul. Alors :
(i) Siest réversible par rapport à P, alors est une mesure invariante.
(ii) Siest réversible par rapport à P etP
j2Ej<1, alors la mesure définie pari=
iP
j2E
jest une distribution stationnaire.
(iii) Siest une distribution stationnaire, alors 8n2N:
P(X0=i0;:::;Xn=in) =P(X0=in;:::;Xn=i0);8(i0;:::;in)2En+1
Démonstration :
(i)Supposons que est réversible par rapport à P. Alors :
X
i2Eipij=X
i2Ejpji=jX
i2Epji=j
D’où l’invariance de .
(ii)Supposons que est réversible par rapport à PetP
j2Ej<1. D’après (i),est
une mesure invariante. Donc la mesure définie par i=iP
j2E
jest une mesure de
probabilité, donc une distribution stationnaire.
(iii)Soitest une distribution stationnaire. On a :
P(X0=i0;:::;Xn=in) =i0pi0i1pi1i2:::pin 1in
=pi0i1i1pi1i2:::pin 1in
=pi0i1pi1i2i2:::pin 1in
=:::
=pi0i1pi1i2:::pin 1inin
=P(X0=in;:::;Xn=i0)
Tout cela implique que si une chaîne de Markov est réversible par rapport à une proba-
bilité, alors cette probabilité est nécessairement invariante.
Chapitre IV
Méthode numérique
Nousavonsvudansleprécédentchapitrelesbasesfondamentalessurleschaînesde
Markov. A présent donc, nous avons en main les éléments essentiels pour comprendre
lapartieChaînedeMarkovdel’abréviationMCMC.Maintenant,allonsvoircequel’on
entend par Monte Carlo.
Maisrappelonstoutd’abordcequel’onavudanslechapitreII.Pourfairedel’estimation
bayésienne, nous avons d’abord besoin de calculer la densité à postériori (jx)associé
à sa loi à priori , qui est le rapport entre la densité jointe f(xj)()et la loi marginale
dex(qui est l’intégrale de la densité jointe sur tout l’espace du paramètre ). Mais sou-
vent, la loi marginale de xest difficile à calculer analytiquement.De plus, si le nombre
de dimension du paramètre à estimer augmente, le calcul de l’estimation bayésienne
devient de plus en plus difficile, voire même impossible analytiquement, surtout avec
des modèles hiérarchiques de plusieurs niveaux. Tout cela implique que faire de l’esti-
mationbayésienneesttrèsrigoureuxencalcul,quecesoitparestimationponctuelleou
parrégiondeconfianceCertes,aveclesfamillesconjuguées,ceproblèmeestcontourné.
Mais quand n’est-il des densités n’appartenant pas à cette famille? C’est par ce grand
problème qu’on a nécessairement besoin de méthode de calcul d’approximation d’inté-
grale, les " Méthodes Numériques ".
1 Méthode de Monte-Carlo
1.1 Introduction
Exemple 5. Calcul de la valeur par la méthode de Monte-Carlo.
Soit un disque de rayon rdans un carré de côté 2r. Le principe est de générer aléatoire-
ment plusieurs points dans le carré.
Figure 1. Simulation Monte Carlo pour approximer
36
CHAPITRE IV. MÉTHODE NUMÉRIQUE 37
Dupointdevufréquentiste,laprobabilité pqu’unpointappartienneaucercleéquivaut
au nombre de point dans le cercle divisé par le nombre totale de point générer dans
le carré, avec un nombre de point suffisamment grand pour recouvrir tout le carré. On
remarqueaussiquecetteprobabilitéestégaleàlasurfaceducerclediviséparlasurface
du carré.
Notons donc par Scerla surface du cercle, Scarla surface du carré, Ncer(n)le nombre
des points générés aléatoirement tombés dans le cercle et enfin nle nombre total des
points générés aléatoirement dans le carré. On obtient :
lim
n!1Ncer(n)
n=p=Scer
Scar=r:r
4r2=
4
De ce résultat, on peut obtenir la valeur approximée de .
Maintenant, supposons que le cercle dans ce même carré ne soit plus un cercle mais
une figure quelconque et que le calcul de sa surface est impossible à résoudre analyti-
quement.
Figure 2. Simulation Monte Carlo pour calculer une surface quelconque
Notonspar Sfiglasurfacedenotrenouvellefigureet Nfig(n)lenombredespointsgéné-
rés aléatoirement tombés à cette figure. On a donc :
lim
n!1Nfig(n)
n=p=Sfig
Scar)Sfig= lim
n!1Nfig(n)
nScar
Decerésultat,onpeutobtenirlavaleurapproximéedelasurfacedenotrenouvellefigure
et quel que soit sa forme.
Cesdeuxexemplesnousillustrentquenouspouvonsutiliserl’aléatoirepourapproxi-
mer des nombres, ou pour calculer des surfaces dans un plan, ou des volumes dans un
espace, donc aussi pour les calculs des intégrales.
1.2 Approche théorique
Définition 43. Lesméthodes de Monte-Carlo (puisqu’il y en a plusieurs) regroupe une
familledetechniquesd’échantillonnagealéatoirenumériqueayantpourbutdecalculer
des intégrales. Le problème se résume à évaluer l’intégrale :
I=Z
Xh(x)f(x)dx=Ef[h(X)]
Oùfest une fonction de densité de support XRm,hune fonction quelconque.
CHAPITRE IV. MÉTHODE NUMÉRIQUE 38
Pourapproximerlavaleurde I,ongénèreungrandnombrede nvariablespseudos-
aléatoiresXi,i2[n]i:i:ddedensitéfetdeproposercommeapproximationlamoyenne
empirique :
hn=1
nnX
i=1h(xi)
Proposition 12. Avec un très grand nombre de variable aléatoire généré, l’estimation par cette
méthode possède de très bonnes propriétés :
Grâce à la loi forte des grands nombre :
hn !n!1Ef[h(X)]p:s
SiEf[h(X)]estfini,lavitessedeconvergencede hnestdeO(pn)etsavarianceasympto-
tique est :
Var[h(X)] =1
nZ
X(h(x) Ef[h(X)])2f(x)dx
Avec la même condition, grâce au théorème central limite :
pn(h(x) Ef[h(X)])p
Var[h(X)] ! N (1;0)
n!1
Cela conduit à pouvoir construire un intervalle de confiance sur l’approximation de Ef[h(X)].
1.3 Quelques méthodes générales
Nous allons voir deux des méthodes les plus connus.
Acceptation-Rejet
Soitlafonctioncible f(x),ladensitéassociéeàlavariablealéatoiredontnousvoulons
simuler une suite de réalisation.
Supposons qu’ils existent1une loiq(x)et une constante ktelles que :8×2; kq (x)
f(x).qest appelé loi instrumental ou proposal distribution (en anglais ).
Figure 3. Exemple d’une loi instrumentale pour une fonction f
1. Que l’on connait et que l’on ait accès facilement.
CHAPITRE IV. MÉTHODE NUMÉRIQUE 39
La méthode d’Acceptation-Rejet consiste à :
(1) Echantillonner xideq(x)
(2) Echantillonner uiuniformément dans [0;kq(x1)]
(3) Siuif(xi), accepterxi
(4) Revenir à (1)
(5) Répéter (1)-(4) jusqu’à un certain nombre d’itération fixé au préalable
(6) Retourner tous les xiacceptés
Remarque 8. Plusqsera proche de f, moins souvent l’algorithme va rejeter.
Cette méthode doit vérifier certaines conditions qui la donnent de ce fait des inconvé-
nients. En effet :
q(x) doit être non nulle pour chaque valeur de xtelle quef(x)est non nulle
Siq(x)ne donne pas une borne serrée, cette méthode peut être très inefficace
Trouver une bonne loi instrumentale n’est pas toujours facile
Comme nous pouvons bien le voir, l’une des inconvénients de la méthode d’Accep-
tation Rejet réside dans le cas d’un mauvais choix de loi instrumentale, qui implique
beaucoupderejetd’échantillonetquientraîneralaconvergencelentedelaméthode.La
méthode suivante nous apporte une approche pour résoudre cet inconvénient.
Échantillonnage Préférentielle
L’échantillonnagePréférentiel ouImportanceSampling (enanglais ),contrairementàl’Acceptation-
Rejet,cetteméthodenerejettepasleséchantillonsmaislesutilisestousmaisenlesassi-
milantdespoidsquidétermineraleurimportancedanslecalculd’espérancedelafonc-
tion cible.
Cette nouvelle méthode consiste à générer des échantillons distribués suivant une loi
instrumentale quelconque, puis de le transformer pour qu’il possède les propriétés ca-
ractéristiques de notre fonction cible. On a donc :
Ef[h(X)] =Z
h(x)f(x)dx
=Z
h(x)f(x)
q(x)q(x)dx
=Eq
h(X)f(X)
q(X)
'1
nnX
i=1h(xi)f(xi)
q(xi)
L’Echantillonnage Préférentiel consiste à :
(1) Echantillonner nfois deq(x)pour avoir x1;x2;:::;xn.
(2) Calculer les n poids w1;w2;:::;wnavec :wi=f(xi)
q(xi)
(3) Retourner1
nPn
i=1h(xi)f(xi)
q(xi)
Maintenant, supposons que f soit de la forme suivante :
f(xi) =~f(xi)
Favec F =Z
~f()d
CHAPITRE IV. MÉTHODE NUMÉRIQUE 40
Quiestidentiqueàlaformedeladistributionàpostériorietquetrouveranalytiquement
Fest difficile. Le poids widevient :
wi=nh
~f(xi)=q(xi)i
Pn
i=1h
~f(xi)=q(xi)i
Remarque 9. Contrairement à l’Acceptation-Rejet, l’Echantillonnage Préférentiel ne re-
jette aucun échantillon. Elle a une variance plus réduite que celle de la méthode d’Ac-
ceptationRejet.Mais identiquementàla méthoded’Acceptation-Rejet,plus laloiinstru-
mentale est proche de la loi cible, plus la méthode est efficace.
2 La méthode MCMC
2.1 Motivation
Vuprécédemment,laméthodedeMonte-Carloconsisteutilisedesvariablesi.i.dgé-
nérés dans l’espace des paramètres (le support de la distribution cible). Mais comme
dansl’Acceptation-Rejetoul’Echantillonnagepréférentielle,laméthodeestplusefficace
lorsque la loi instrumentale soit plus proche de la loi cible. Or choisir une telle loi ins-
trumentale n’est loin d’être évident. Et plus la dimension de notre loi cible augmente,
échantillonner des variables i.i.d sur une espace à grande dimension devient très cou-
teux en calcul et ce qui rend la méthode moins efficace. Ceci est dû à la malédiction de la
dimension (oucurseofdimension ).Donc,pourpouvoirrésoudrecesinconvénients,ilfaut
fairedeséchantillonnagesintelligents(voir[ MHM16])afindeminimiserleséchantillon-
nages inutiles (dans les parties ou notre loi est nul).
2.2 Approche théorique
Définition 44. Les méthodes de Monte-Carlo par Chaîne de Markov (Markov Chain Monte
Carloen anglais et abrévié MCMC) regroupent une famille de méthode produisant une
chaînedeMarkovergodiquedontladistributionstationnaireestladistributiond’intérêt.
Puisque la chaîne est distribuée selon , on peut alors utiliser l’estimation Monte-
Carlo afin d’approximer E[h(X)]. Pour que cela marche, il faut que la chaîne produite
remplisse les conditions2pour converger vers sa loi stationnaire (la loi cible ), et puis
par le théorème ergodique, l’estimation converge (voir [ RC96]).
Toutes les méthodes MCMC peuvent être généralisées par cet algorithme :
(1) Commencer l’algorithme à la position "actuelle" dans l’espace des paramètres actuel
(2) Proposer un "saut" une nouvelle position dans l’espace des paramètres suivant
(3) Accepter ou Rejeter le saut proposé
(4) Si le saut est accepté, revenir à (1)
(5) Si le saut est rejeté, rester sur la même position et revenir à (1)
(6) Après un nombre défini de sauts, retourner toutes les positions acceptées
Les principales différences entre les algorithmes MCMC se trouver aux instructions
(2) et (3) de cet algorithme. C’est-à-dire, leur "façon de sauter" et comment "ils décident
de sauter".
2. Irréductible, apériodique et récurrente positive
CHAPITRE IV. MÉTHODE NUMÉRIQUE 41
2.3 Quelques différents algorithmes
Plusieurs sont les algorithmes MCMC, mais nous allons voir deux des plus connus.
Nous allons voir en premier l’algorithme Metropolis (dont le fonctionnement est relati-
vementtrèssimpleàcomprendre),ensuitel’algorithme Metropolis-Hasting (quiestl’une
des versions générale des algorithmes MCMC).
L’algorithme de Metropolis
L’algorithme Metropolis est un algorithme MCMC qui utilise une distribution nor-
malepourproposerunsaut.Cettedistributionnormaleaunevaleurmoyenne quiest
égale à la position actuelle et prend une "largeur de proposition" comme écart-type.
Cette largeur de proposition est un paramètre de l’algorithme de Metropolis et a un
impact significatif sur la convergence. Une largeur de proposition plus élevée va sau-
ter plus loin et couvrir plus d’espace dans la distribution postérieure, mais vas souvent
rejeté des échantillons,donc vas mettre un certain moment pour converger. Cependant,
une largeur de proposition plus petite ne couvrira pas la plus grande partie de l’espace
aussi rapidement et pourrait donc prendre plus de temps à converger.
Soit doncnotre loi cible. L’algorithme de Metropolis est la suivante :
Algorithme 1. Algorithme de Metropolis
Entrées: Le point de départ 0de l’algorithme; La longueur N de la chaîne.
Sortie: La chaîne (0;1;:::N)
1. i = 1
2.pouriNfaire
3. Sur une position i, proposerpropselon la loiN(i;2)
4. Générer le critère d’acceptation
p= min
1;(prop)
(i)
5. Générer uuniformément dans [0; 1]
6.Siupalorsi+1 prop
7.Sinoni+1 i
.fin
fin
Avant d’entamer la validation cet algorithme, vérifions tout d’abord que la chaîne
qu’il produit possède les bonnes propriétés pour converger :
(i)Ladistributionquigénère propétantuneloinormale,dontlesupportrecouvrecelui
de,lachaînepeutdoncvisitertouslesétatsdusupportde .Donc8i;j2N;9n2N
tel quep(n)
ji>0. Donc la chaîne produite est irréductible.
(ii)On peut remarquer que i+1 propest possible à tout moment. Donc, il existe au
moinsunétat itelquep(n)
ii>0.Laproposition9et (i)impliquequelachaîneproduite
est apériodique.
(iii)On peut remarquer que i+1 propComme la distribution stationnaire est la loi
cible, le théorèmeIII.1.5 et (i)impliquent que tous les états sont récurrents positives.
(i);(ii);(iii)impliquent que la chaîne est ergodique, et par le théorème III.2.6, elle
possède belle et bien les conditions nécessaires pour converger vers .
CHAPITRE IV. MÉTHODE NUMÉRIQUE 42
L’algorithme de Metropolis-Hasting
Celui de Metropolis utilise une distribution normale pour proposer un saut, celui
de Metropolis-Hasting, pour sa part, utilise une distribution arbitraire3q(jactuel),
de façons qu’elle soit non nulle sur tout le support de la distribution cible et qui est
conditionnelle par l’état présent de la chaîne.
Soitdoncnotreloicibleet qladistributioninstrumentale.L’algorithmedeMetropolis-
Hasting est la suivante :
Algorithme 2. Algorithme de Metropolis-Hasting
Entrées: Le point de départ 0de l’algorithme; La longueur N de la chaîne.
Sortie: La chaîne (0;1;:::N)
1. i = 1
2.pouriNfaire
3. Sur une position i, proposerpropselon la loi q(propji)
4. Générer le critère d’acceptation
p= min
1;(prop)q(ijprop)
(i)q(propji)
5. Générer uuniformément dans [0; 1]
6.Siupalorsi+1 prop
7.Sinoni+1 i
.fin
fin
Mais tout d’abord, vérifions que la chaîne qu’il produit possède aussi les propriétés
pour être vers la distribution cible :
(i)8i; jappartenant au support de ,q(jji)>0, sinon la loi instrumentale a été
mal choisit. Donc 8i;j2N;9n2Ntel quep(n)
ji>0. Donc la chaîne produite est irré-
ductible.
(ii)Identiquement à celui de Metropolis, i+1 propest possible à tout moment. La
chaîne produite est apériodique.
(iii)La probabilité d’acceptation a été construite pour respecter la condition de réversi-
bilité. En effet, on a
()P0=()q(j0)p0
Si
p0= 1)p0=()q(j0)
(0)q(0j)<1
On aura donc
(0)P0=(0)q(0j)p0
=(0)q(0j)()q(j0)
(0)q(0j)
=()P0
Dans le cas contraire, on aura le résultat similaire. Donc est réversible par rapport la
chaîne, donc par le théorème III.3.8, est une distribution stationnaire par rapport à la
chaîne
3. Maisdoitêtrefacileàsimuler,disponibleanalytiquementetqu’elleaitunedispersionassezimpor-
tante pour que la chaine peut ainsi explorer tout le support de la distribution cible.
CHAPITRE IV. MÉTHODE NUMÉRIQUE 43
(iv)D’après (iii), la distribution stationnaire est , le théorème III.1.5 et (i)impliquent
que tous les états sont récurrents positives.
(i);(ii);(iii);(iv)impliquentquelachaîneproduiteseraasymptotiquementdistribuése-
lon.
Bienquecesoientdesalgorithmesrelativementsimples,iln’estpasaussitôtévident
commentilsnousaidentàéviterleproblèmedecalculerd’intégraledelaloimarginale.
En fait, en divisant l’a postériori du paramètre proposé par l‘à postériori du paramètre
courant. Par la formule de Bayes, la loi marginale f(x)du rapport se simplifie :
(propjx)
(ijx)=f(xjprop)(prop)
f(x)
f(xji)(i)
f(x)=f(xjprop)(prop)
f(xji)(i)
Lecôtédroitdecettedernièreégaliténecontientquelesdensitésjointes,quenouspou-
vons obtenir grâce aux observations et aux informations à priori. Ainsi, en divisant l’a
postériori d’un point par l’a postériori d’un autre, nous échantillonnons des régions de
probabilité à postériori plus souvent.
3 Application
Notonsbienquelaprécisiondesestimationssoutiréesparcesméthodesdépendsur-
toutsurlecomportementasymptotiquedeschaînesqu’ilsontgénéré.Maisriennenous
apprend, ni sur comment avoir une très bonne précision, ni sur comment savoir si la
chaîne ait effectivement atteint son régime stationnaire
3.1 Mise en œuvre
Période de chauffe
Toutes les simulations obtenues par les méthodes MCMC peuvent être divisées en
deux phases : la première, la phase où les états visités par la chaîne sont dépendants de
l’état initial de la chaîne, et la seconde, la phase de simulation où la chaîne a atteint sa
distributionstationnaire.Lescandidatsgénéréslorsdelapremièreétape,appelée période
de chauffe (ouBurn-inen anglais), ne peuvent pas être utilisés lors de l’échantillonnage,
puisqu’ils ne sont pas considérés comme distribués selon la distribution cible. Au-delà
decettepériodedechauffe,lesvaleursgénérésdelachaînesontconsidéréscommedis-
tribué selon la loi cible, et donc, peuvent être échantillonnés. La durée de la période de
chauffevarieselonlaloicible,cequiveutdirequelenombred’itérationn0pourcettepé-
rioden’estpasconstantetc’estàl’utilisateurdel’algorithmedevérifierempiriquement
silapériodedechauffeestterminéeets’ilestdoncpossibledecommencerl’échantillon-
nage.
Remarque10. Uneerreurdansl’estimationden0entraineuneaugmentationdel’erreur
d’estimation, si l’échantillonnage commence durant la période de chauffe, ou bien une
augmentation du temps de calcul, la période de chauffe est prolongée inutilement.
Diagnostic de convergence
Bien que la convergence théorique soit assurée par le théorème ergodique, nous ne
savonspasencorecombiend’itérationsnousfaudrait-ilpouravoirlaconvergenceenloi.
Il existe plusieurs tests concernant la convergence des MCMC, mais nous avons choisis
CHAPITRE IV. MÉTHODE NUMÉRIQUE 44
celui de Gelman-Rubin qui se base sur la simulation de plusieurs chaînes en parallèle
afin de déterminer numériquement s’ils ont atteint leur régime stationnaire.
Ensimulant mchaînesdemêmelongueur n,nousdéfinissonspar (j)
ilavaleurdela jme
simulation à l’instant i, et les statistiques suivantes :
G=m
n 1mX
j=1
j 2etS=1
mmX
j=1S2
j
Avec
j=1
nnX
i=1(j)
i; =1
mmX
j=1 j; S2
j=1
n 1nX
i=1
(j)
i j2
Lastatistique G(lavarianceinter-chaîne)estuneestimationdelavarianceglobaleentre
les simulations et S(la variance intra-chaînes) est une estimation de la variance dans
chacunedessimulations.L’estimationdelaRéductiond’échellepotentielle,quel’onnote
^R, correspond au rapport suivant :
^R=vuutn 1
nS+m+ 1
nmG
S=r
n 1
n+m+ 1
nmG
S
Lorsque les simulations sont toutes en régime stationnaire, ce rapport est égal à 1 (voir
[G.R92]). Dans le cas contraire, on retarde l’échantillonnage. En pratique, si ce rapport
est en dessous de 1;2, on peut accepter l’hypothèse de stationnarité en loi.
3.2 Exemple d’application
Reprenonsl’exemple4,maismaintenantsupposonsquel’onn’apasaccèsàladistri-
bution à postériori (que l’on ne peut pas approcher par les méthodes analytique). Mais
qu’on voudrait quand même l’approximer.
Nous allons approximer les 2 distributions à postériori avec l’algorithme de Metropolis
dont la loi instrumentale est une normale d’écart-type = 0;01, afin d’aboutir à l’esti-
mation bayésienne.
Nous avions eu 2 lois à priori, l’une non-informative (la loi 1= (1;1)) et l’autre in-
formative(laloi 1= (12;12)).Pourcalculerleurloiàpostérioricorrespondante,nous
n’aurons plus besoin d’évaluer le dénominateur de la loi à postériori, c’est-à-dire leur
prédictive).
Nousavonsutilisétroistypesdesimulationspourapproximercesloiscibles.Lescodes
R pour ces simulations sont dans l’Annexe B 3.
Les histogrammes sont les approximations de la simulation et la courbe en rouge sont
les distributions cible. Voici donc les graphiques obtenu :
Avec 10 000 itérations de MCMC avec 500 itérations de période de chauffe.
CHAPITRE IV. MÉTHODE NUMÉRIQUE 45
Figure 4. Approximation des 2 à postériori avec 10 000 échantillons
Là on voit que dans les deux cas, l’approximation n’est pas encore très bien fiable. Ce
qui nécessite le prolongement de l’échantillonnage.
Avec 100 000 itérations de MCMC avec 5 000 itérations de période de chauffe.
Figure 5. Approximation des 2 à postériori avec 100 000 échantillons
Avec100000échantillons,onremarquequel’approximationcommenceàêtretrèsproche
des distributions cible. L’erreur peut-être vraiment négligeable.
Avec 1 000 000 itérations de MCMC avec 10 000 itérations de période de chauffe.
CHAPITRE IV. MÉTHODE NUMÉRIQUE 46
Figure 6. Approximation des 2 à postériori avec 1 000 000 échantillons
Avec 1 millions d’échantillons, on voit que l’approximation suit parfaitement la distri-
bution cible.
On remarque tout de suite que plus la taille de l’échantillon est important, plus l’ap-
proximation est précise.
Maintenantqu’onaunetrèsbonneapproximation(nousallonsutiliserlatroisième)des
distributions cible, on peut passer à l’inférence. Le tableau suivant compare le résultat
analytique et le résultat obtenu par les simulations.
Tableau 3. Comparaison ente estimation par méthode analytique et par simulation
Situation Non informative Informative
Loi à priori 1(p) = (1;1) 2(p) = (12;12)
Loi à postériori 1(pjX) = (8;4)2(pjX) = (19;15)
Estimateur de Bayes 1=8
8 + 40;6672=19
19 + 150;559
analytique
Intervalle de crédibilité [m1;M1] = [0;39; 0;89] [m2;M2] = [0;392; 0;719]
à 0,95
Estimateur de Bayes MCMC10;667MCMC20;559
par MCMC
Intervalle de crédibilité [0;393; 0;892] [0;392; 0;721]
à 0,95 par MCMC
Voyons par ce tableau à quel point l’approximation par la méthode de MCMC est très
prochedurésultatdonnéeparlaméthodeanalytique.Nouspouvonstrèsbienconstater
que la détermination de la loi à postériori n’est plus nécessaire, seule la vraisemblance
CHAPITRE IV. MÉTHODE NUMÉRIQUE 47
et l’à priori suffisent pour aboutir à faire de l’inférence bayésienne grâce aux méthodes
MCMC.
Onsaitmaintenantqu’enfonctiondel’augmentationdunombred’échantillonquel’ap-
proximation devienne à son tour plus précise. Mais on ne sait pas quelle est la quantité
d’échantillon nécessaire que l’on devrait avoir afin d’arrêter l’augmentation du nombre
d’échantillon.C’estlàqu’intervientletestdeGelmanetRubin.SurR,onabesoind’ins-
tallerunpackagenommer"coda"afindepouvoirl’utiliserfacilementsanspourautant
construire un algorithme du calcul de la réduction d’échelle potentiel (voir [ R.C11]).
Maispourfairedel’inférencepurementbayésienne,ilexistedeslogicielsadaptéàcela:
WinBUGS,OpenBUGS,….Cesonttoutesdeslogicielsspécialementconçuspourréaliser
et faciliter des inférences en utilisant les approches bayésiennes. Ces logiciels utilisent
des algorithmes MCMC pour approximer les calculs. Nous allons utiliser OpenBUGS
mais sous R. Pour cela, il faudra installer4le package " R2OpenBUGS " afin de pouvoir
l’utiliser. Pour pouvoir réaliser notre inférence, nous n’auront besoin que des données,
la loi à priori et la vraisemblance. En sortie nous auront : une estimation ponctuelle, un
intervalle de crédibilité à 0,95, la valeur de la réduction d’échelle potentielle (ici, nous
allonslancer3chaînesenmêmetemps).C’est-à-direlesrésultatsattenduspourfairede
l’estimation bayésienne.
Voici donc les résultats obtenu par l’utilisation d’OpenBUGS sous R. Les codes R sont
dans l’Annexe B 6.
Tableau 4. Résultat des simulations sous OpenBUGS via R
Simulations Situation Estimation Intervalle de ^R
Ponctuelle crédibilité à 0,95
10 000 Non informatif 0,668 [0,390;0,890] 1,001
échantillon
et 3 chaines
en parallèles informatif 0,559 [0,392;0,719] 1,001
50 000 Non informatif 0,667 [0,389;0,891] 1,001
échantillon
et 3 chaines
en parallèles informatif 0,559 [0,392;0,719] 1,001
Dans les 2 cas, dans les deux situations, on les même valeurs de ^R, ça veut dire que la
chaîneadéjàatteintsastationnarité.Donc,avec10000échantillons,onpeutdéjàaccepter
l’hypothèsequelachaîneestrépartitselonladistributioncible,cequisignifiequ’onpeut
passer à l’inférence.
4. IlsepourraitqueRdonnedesmessagesd’erreuraprèsl’appelledupackage.Sic’estlecas,installer
d’abord le package " MASS " et ensuite réinstaller le package.
Conclusion
Nousavonsvuquel’approchebayésienprésenteunavantageparrapportàsonvoi-
sin,lefréquentiste,dufaitqu’ilpeucontournerlenécessitédesapproximationsasymp-
totiques, qui est la base même de l’approche fréquentiste. Grâce à la loi à priori, toutes
les informations sont réunies dans la distribution à postériori. Et grâce à elle, on peut
passer direct à l’inférence (estimation ponctuelle et région de crédibilité).
Par contre, nous avons vu que ce privilège présente un grand problème. Le calcul de la
loiàpostériori,àcausedesondénominateur f(x) =R
f(xj)(),présentetoujoursde
problème ce qui rend les méthodes bayésienne mises de côtés durant plusieurs années.
Elle requiert des calculs d’intégrales, souvent impossibles à la main.
Le but de cet étude est de présenté l’algorithme MCMC qui est un algorithme stochas-
tiquecommeceluideMonteCarlo,maislemieuxadaptéauxméthodebayésiennepuis-
qu’on peut approximer directement la loi à postériori sans avoir à calculer son dénomi-
nateur. Cet algorithme fait le rapport entre deux à postériori de deux points distincts
dans l’espace des paramètres à échantilloner. De ce fait, leur dénominateur commune
sesimplifientautomatiquement.Laformecomplètedeladistritubtionàpostériorin’est
doncplusnécessairepourpouvoirl’approximer.Seulementqu’àunmultiplicativeprès,
grâce à la proportionnalité de l’à postériori et sa densité jointe : (jx)/f(xj)().
48
Annexe A
1. Théorème de Bayes
Dans le cas discret :
Soit deux évènements A et B. On a :
P(BjA) =P(AjB)P(A)
P(B)
Dans le cas continu :
Soit deux variable etxtelles que2etx2 . On a :
(jx) =f(xj)()
f(x)=f(xj)()R
f(xj)()
2. Théorème Centrale Limite
Soit(Xk)1kn)une suite de variable aléatoir i.i.d, de mêmes moyenne et variance 2.
On a alors :
Sn=X1+:::+Xn ! N
n; n2
n!1
3. Loi Forte des Grands Nombre
Soit(Xk)1kn)unesuitedevariablealéatoiri.i.d.Sideplus,ilssontdecarréintégrable,
alors on a :
lim
n!11
nnX
i=1Xi=E[X0]; p:s
IV
Annexe B
1.
# Entrée des données
X<- c(1,0,1,0,1,1,0,1,1,1)
# Inférence
binom.test(table(X)[2:1])
2.
# Intervalle de crédibilité à 95% obtenu par la loi à posteriori beta(8,4)
m1 = qbeta(0.025,8,4)
M1 = qbeta(0.975,8,4)
# Intervalle de crédibilité à 95% obtenu par la loi à posteriori beta(19,15)
m2 = qbeta(0.025,19,15)
M2 = qbeta(0.975,19,15)
3.
# Entrée des données
X<- c(1,0,1,0,1,1,0,1,1,1)
# Les lois à priori
pi1 <- function(x) dbeta(x,1,1 # A priori non informatif beta(1,1)
pi2 <- function(x) dbeta(x,12,12) # A priori informatif beta(12,12)
# Distribution d’échantillonnage ou la vraisemblane
vraisemblance<-function(x) (x^sum(X))*(1-x)^(length(X)-sum(X))
# Les densités jointes associées
dens.jointe1 <- function(x) vraisemblance(x)*pi1(x)
dens.jointe2 <- function(x) vraisemblance(x)*pi2(x)
###### Construction de l’algorithme MCMC
# La loi instrumentale
qM <- function (x) rnorm (1 , x , 0.01)
# Marche de l’algorithme
pasM <- function ( x , f , q ) {
xp <- q( x ) # Proposition du nouveau point
alpha <- min (1,f(xp)/f(x)) # Probabilité d’acceptation
if (runif ( 1 ) < alpha ) # Accepte le nouveau point
{x <- xp} # avec une probabilité alpha
drop(x)
}
V
CHAPITRE IV. MÉTHODE NUMÉRIQUE VI
# L’algorithme de Metropolis
M <- function ( x , f , q , nbpas ) {
res <- matrix ( NA , nbpas , length(x))
for ( i in seq_len ( nbpas ))
res [ i ,] <- x <- pasM ( x , f , q )
drop ( res )
}
# Période de Chauffe
BurninM<-function(x,f,q,nbiter){
for (i in length(nbiter)){
x<-pasM(x,f,q)
}
drop(x)
}
# Première approximation des deux lois à postériori
res1.1 <- M( BurninM(runif(1) , dens.jointe1 , qM , 500) , dens.jointe1 , qM ,
10000)
res1.2 <- M( BurninM(runif(1) , dens.jointe2 , qM , 500) , dens.jointe2 , qM ,
10000)
# Deuxième approximation des deux lois à postériori
res2.1 <- M( BurninM(runif(1) , dens.jointe1 , qM , 5000) , dens.jointe1 , qM ,
100000)
res2.2 <- M( BurninM(runif(1) , dens.jointe2 , qM , 5000) , dens.jointe2 , qM ,
100000)
# Troisième approximation des deux lois à postériori
res3.1 <- M( BurninM(runif(1) , dens.jointe1 , qM , 10000) , dens.jointe1 , qM ,
1000000)
res3.2 <- M( BurninM(runif(1) , dens.jointe2 , qM , 10000) , dens.jointe2 , qM ,
1000000)
4.
# Graphique de la première approximation
layout(matrix(1:2,1,2))
hist (res1.1, 100 , freq = FALSE , ylim = c( 0 , 4 ), las = 1 ,
xlim = c ( 0 , 1 ) , main = "A postériori 1" , xlab="x" , ylab="" )
curve ( dbeta(x,8,4), add=TRUE,col = "red" , 0 , 1 , n = 1000 , las=1 )
hist (res1.2, 100 , freq = FALSE , ylim = c( 0 , 5 ), las = 1 ,
xlim = c ( 0 , 1 ) , main = "A postériori 2" , xlab="x" , ylab="" )
curve ( dbeta(x,19,15), add=TRUE,col = "red" , 0 , 1 , n = 1000 , las=1 )
# Graphique de la deuxième approximation
layout(matrix(1:2,1,2))
hist (res2.1, 100 , freq = FALSE , ylim = c( 0 , 4 ), las = 1 ,
xlim = c ( 0 , 1 ) , main = "A postériori 1" , xlab="x" , ylab="" )
curve ( dbeta(x,8,4), add=TRUE,col = "red" , 0 , 1 , n = 1000 , las=1 )
hist (res2.2, 100 , freq = FALSE , ylim = c( 0 , 5 ), las = 1 ,
xlim = c ( 0 , 1 ) , main = "A postériori 2" , xlab="x" , ylab="" )
curve ( dbeta(x,19,15), add=TRUE,col = "red" , 0 , 1 , n = 1000 , las=1 )
CHAPITRE IV. MÉTHODE NUMÉRIQUE VII
# Graphique de la troisième approximation
layout(matrix(1:2,1,2))
hist (res3.1, 100 , freq = FALSE , ylim = c( 0 , 4 ), las = 1 ,
xlim = c ( 0 , 1 ) , main = "A postériori 1" , xlab="x" , ylab="" )
curve ( dbeta(x,8,4), add=TRUE,col = "red" , 0 , 1 , n = 1000 , las=1 )
hist (res3.2, 100 , freq = FALSE , ylim = c( 0 , 5 ), las = 1 ,
xlim = c ( 0 , 1 ) , main = "A postériori 2" , xlab="x" , ylab="" )
curve ( dbeta(x,19,15), add=TRUE,col = "red" , 0 , 1 , n = 1000 , las=1 )
5.
# Inférence après approximation des lois à postériori
# Estimation ponctuelle
mean(res3.1) # A posteriori 1
mean(res3.2) # A posteriori 2
# Intervalle de crédibilité à 95%
quantile(res3.1 , probs=c(0.025 , 0.975)) # A posteriori 1
quantile(res3.2 , probs=c(0.025 , 0.975)) # A posteriori 2
6.
# Entrée des données
X<- c(1,0,1,0,1,1,0,1,1,1)
# Spécification des modèles
# Situation Non informative
modele1 <- function() {
nbface <- sum(X[]) # Nombre de réalisations du côté pile
nbface ~ dbin(p,n) # Vraisemblance
p ~ dbeta(1,1) # A priori non-informatif
}
write.model(modele1,"modele1.txt")
# Situation Informative
modele2 <- function() {
nbface <- sum(X[]) # Nombre de réalisations du côté pile
nbface ~ dbin(p,n) # Vraisemblance
p ~ dbeta(12,12) # A priori informatif
}
write.model(modele2,"modele2.txt")
# Spécification du paramètre à estimer
parametres <- "p"
# Entrée des données dans le modèle
n <- length(X)
data <- list("X","n")
# Spécification du début de l’agorithme
inits1 <- function() {p=rbeta(1,1,1)} # Pour modèle 1
inits2 <- function() {p=rbeta(1,12,12)} # Pour modèle 2
CHAPITRE IV. MÉTHODE NUMÉRIQUE VIII
# Simulation MCMC dont 10 000 échantillons et 500 pour le période de chauffe
res1.1 <- bugs(data, inits1, parametres, "modele1.txt", n.chains=3, n.iter=10000
, n.burnin=500)
res1.2 <- bugs(data, inits2, parametres, "modele2.txt", n.chains=3, n.iter=10000
, n.burnin=500)
# Simulation MCMC dont 50 000 échantillons et 1000 pour le période de chauffe
res2.1 <- bugs(data, inits1, parametres, "modele1.txt", n.chains=3, n.iter=50000
, n.burnin=1000)
res2.2 <- bugs(data, inits2, parametres, "modele2.txt", n.chains=3, n.iter=50000
, n.burnin=1000)
# Résultats du premier échantillonnage
print(res1.1 , digits.summary=3)
print(res1.2 , digits.summary=3)
# Résultats du second échantillonnage
print(res2.1 , digits.summary=3)
print(res2.2 , digits.summary=3)
Bibliographie
[AM97] A.Monfort, CoursdeStatistiqueMathématique ,EcoleNationaledelaStatistique
et de l’Administration Economique et Centre d’Etudes des Programmes Econo-
miques ,Economica, 3ème édition, 1997
[B.R98] Brooks, S. and Roberts, Assessing convergence of Markov chain Monte Carlo algo-
rithms, 1998.
[EC14] Ertur Cem, Rappels : inférence statistique , notes de cours, 2014.
[DB16] Delyon Bernard, Estimation paramétrique , Cours de Master 2, 2016.
[FM07] F. Malrieu, Intervalles de confiance , Préparation à l’épreuve de modélisation –
Agrégation Externe de Mathématiques, Université de Rennes I, 2007.
[G.R92] Gelman, A. et Rubin, D. B. , Inference from iterative simulation using multiple se-
quences, Statistical Science, 7, 473–483, 1992.
[IK85] Ibragimov, A. and Khas’minskii, R.Z., On non-parametric estimation of the value of
a linear functional in gaussian white noise ,Theory of Probability and its Applications,
29(1), 1985.
[JL02] J.Lacroix, ChaînesdeMarkovetProcessusdePoisson,notesdecours ,UniversitéPierre
et Marie Curie, 2002.
[JR10] J. Rousseau, Statistique Bayésienne , notes de cours, ENSAE ParisTech, 2010.
[KW96] KrauthWener, IntroductionToMonteCarloAlgorithms ,notesdecours,EcoleNor-
male Supérieure CNRS-Laboratoire de Physique Statistique, 1996.
[LB05] Lecoutre,B., Etsivousétiezunbayésienquis’ignore? ,RevueMODULAD,92-105,
2005.
[MHM15] MichaelHalls-Moore, BayesianInferenceofaBinomialProportion-TheAnalyti-
cal Approach , Webographie, 2015.
[MHM16] Michael Halls-Moore, Markov Chain Monte Carlo for Bayesian Inference – The
Metropolis Algorithm , Webographie, 2016
[NB07] Nils Bergkund, Chaînes de Markov , notes de cours, Universités d’Orléans, 2007.
[P.B07] Parent E. et Bernier J., Le Raisonnement Bayésien – Modélisation et inférence ,
Springer, Paris, 2007.
[P.B.B07] Parent E., Bernier J. et Boureux J., Pratique du calcul Bayésien , Springer, Paris,
2010.
[RC96] Robert, C. , Methodes de Monte Carlo par chaines des Markov , Economica, 1996.
[RC06] Robert, C. , Le Choix Bayésien – Principes et pratique , Springer, Paris, 2006.
[R.C04] Robert, C. et Casella, G. , Monte Carlo Statistical Methods , Springer, New York,
second édition, 2004.
[R.C11] Robert, C. et Casella, G. , Méthode de Monte Carlo avec R , Springer, 2011.
[RC13] Robert,C., Desspécificitésdel’approchebayésienneetdesesjustificationsenstatistique
inférentielle , Webographie, 2013.
IX
BIBLIOGRAPHIE X
[SJ16] Sprenger,J., Bayésianismeversusfréquentismeeninférencestatistique ,TilburgCenter
for Logic, Ethics and Philosophy of Science, Le Bayésianisme aujourd’hui, p.167-
192, 2016.
[YI12] YildromI., BayesianInference:Metropolis-HastingsSampling ,Rochester,NY14627,
2012.
Imprétrant : RATISARIJAONA Mickaël David
Tel: (+261) 34 76 307 30
e-mail: r.jaonamickael@gmail.com
Titre : Pratique de la méthode MCMC en estimation bayésienne
Résumé : Dans ce mémoire, nous avons pu exposer l’avantage de l’approche bayésien
par rapport à l’approche fréquentiste en termes d’estimation ponctuelle et intervalle de
confiance, dans le cadre paramétrique de la statistique inférentielle. Mais ce privilège
présente un inconvénient en termes de calcul. Notre but est de présenter l’algorithme
MCMC, qui est une méthode d’approximation de calcul comme Monte Carlo, mais le
mieux adapté pour l’inférence bayésienne.
Mots clés : Inférence Statistique, Inférence bayésienne, Chaîne de Markov, Méthode de
Monte Carlo, Méthode de Monte Carlo par Chaîne de Markov.
Abstract :Inthismemoir,wehavecouldexposetheadvantageofthebayesianapproach
in comparison with the frequentist approach, concerning the ponctual estimation and
intervalofconfidence,withintheparametriccontextofthestatisticalinference.Butthat
privilegeshowsaninconvenienceconcerningcalculation.OurgoalistoshowtheMCMC
algorithm,thatisanapproximatemethodofnumericalcalculationlikeMonteCalro,but
the better adapted for the bayesian inference.
Keywords : Statistical inference, Bayesian inference, Markov Chain, Monte Carlo Me-
thod, Markov Chain Monte Carlo Method.
Encadreur : Dr RAZAFIMAHALEO Martial
Maître de Conférences
Département de Mathématiques et Informatique
Faculté des Sciences
Université d’Antananarivo
Tel : (+261) 34 09 178 95
E-mail : martialeo@hotmail.fr
Copyright Notice
© Licențiada.org respectă drepturile de proprietate intelectuală și așteaptă ca toți utilizatorii să facă același lucru. Dacă consideri că un conținut de pe site încalcă drepturile tale de autor, te rugăm să trimiți o notificare DMCA.
Acest articol: DOMAINE DES SCIENCES ET TECHNOLOGIES MENTION EN MATHEMATIQUES ET INFORMATIQUE Mémoire en vue de l’obtention du Diplôme de Master en Mathématiques… [614567] (ID: 614567)
Dacă considerați că acest conținut vă încalcă drepturile de autor, vă rugăm să depuneți o cerere pe pagina noastră Copyright Takedown.