TP 4 et 5 : Corrigé - Mathématiques en ECE 2

Download Report

Transcript TP 4 et 5 : Corrigé - Mathématiques en ECE 2

ECE 2 - Mathématiques
Mr Dunstetter - ENC-Bessières 2014\2015
TP 4 et 5 : Chaînes de Markov
TP 4 et 5 : Corrigé
I
Simulation de l'exemple épidémiologique du cours d'introduction
Exercice 1.
En informatique, on représente l'état S par le nombre 1, l'état I par le nombre 2 et l'état M par le nombre
3.
1. n=input('Entrer la durée de l''expérience : ')
X=zeros(1,n+1)
X(1)=1
for k=1:n do
if X(k)==1 then
if rand() <1/2 then
X(k+1)=3
else
X(k+1)=1 \% (inutile)
end
elseif X(k)==2 then
if rand()<1/10 then
X(k+1)=1
else
X(k+1)=2 \% (inutile)
end
else
if rand()<4/5 then
X(k+1)=2
else
X(k+1)=3 \% (inutile)
end
end
end
disp(X)
2. La simulation des chaînes de Markov étant rapidement très lourde (cf ci-dessus pour seulement trois
états), une fonction a été prédénie en Scilab pour la simuler : Soient
• X0 un état xé ;
• A la matrice des probas de transition de la chaîne de Markov ;
L'instruction
X=grand(n,'markov',A',X0)
donne un vecteur ligne X représentant la simulation des variables (X1 , X2 , ..., Xn ) d'une chaîne de
Markov de matrice de transition A et de position initiale x0.
Il faut faire très attention aux points suivants :
• On doit donner A0 la transposée de A car la dénition d'une matrice de transition en Scilab
est la transposée de notre dénition.
• Si on souhaite acher l'état initial de la chaîne de Markov, il faut rajouter la valeur initiale
x0 au vecteur X !
• Cette instruction permet de simuler la variable aléatoire Xn , pas d'obtenir sa loi Un (on verra
plus loin deux méthodes pour le faire).
1
ECE 2 - Mathématiques
Mr Dunstetter - ENC-Bessières 2014\2015
TP 4 et 5 : Chaînes de Markov
A=[ 1/2 , 1/10 , 0 ; 0 , 9/10 , 4/5 ; 1/2 , 0 , 1/5]
X0=1
X=grand(10,'markov',A',X0)
X=[X0,X]
3. Loi empirique de X10 :
(a) p=zeros(1,100)
for k=1:100 do
A=[ 1/2 , 1/10 , 0 ; 0 , 9/10 , 4/5 ; 1/2 , 0 , 1/5]
X0=1
X=grand(10,'markov',A',X0)
p(1,k)=X(10)
end
(b) L'instruction tabul donne une matrice avec en première colonne les valeurs du caractère de la
série statistique et en deuxième colonne l'eectif de la valeur dans la série statistique.
Pour changer les eectifs en fréquence, on récupère la colonne des eectifs, on la divise par sa
somme pour obtenir les fréquences, puis on remplace dans M :
N=M(:,2)
N=N/sum(N)
M(:,2)=N
ou en une seule fois :
M(:,2)=M(:,2)/sum(M(:,2))
puis on fait :
p=zeros(1,1000)
for k=1:1000 do
A=[ 1/2 , 1/10 , 0 ; 0 , 9/10 , 4/5 ; 1/2 , 0 , 1/5]
X0=1
X=grand(10,'markov',A',X0)
p(1,k)=X(10)
end
M=tabul(p,"i")
M(:,2)=M(:,2)/sum(M(:,2))
4. Loi théorique de X10 : On sait calculer Uk le vecteur de la loi de Xk . Calculons les lois de X0 ,
X1 , ..., X10 :
// Loi initiale
U=[1;0;0];
// On crée une matrice [U0 U1 ...Un]
E= zeros(3,11);
E(:,1)=U; // range dans la première colonne de E la loi U0.
for i=2:11 do
U=A*U; // U prend les valeurs de U(i-1)
E(:,i)=U; // range dans la i-ième colonne de E la loi U(i-1)
end;
disp(E)
La dernière colonne doit être parfaitement cohérente avec la loi empirique, surtout celle avec 1000
simulations. On voit apparaître une convergence de la loi de Un vers une loi limite (chaque probabilité
semble converger vers une certaine limite).
5. Représentation graphique des lois : Taper :
scf(1)
n=0:10
plot2d(n', [E(1,:)', E(2,:)',E(3,:)'])
title('les suites (sn), (in) et (mn)')
On remarque alors que la loi de (Xn ) converge vers une loi limite π .
2
ECE 2 - Mathématiques
Mr Dunstetter - ENC-Bessières 2014\2015
II
TP 4 et 5 : Chaînes de Markov
Le Google Page Rank
Nous allons illustrer la méthode utilisée par Google pour classer les pages web par "indice de popularité".
L'idée est d'étudier l'évolution des déplacements d'un individu sur le net, et de regarder la probabilité pi
d'être sur un site i au bout d'une très grand nombre de déplacements (limite de la chaîne de Markov). pi
est alors appelée "indice de popularité de la page i".
Nous commencerons par regarder un cas très simple (et très édulcoré) où il n'existerait sur le Web que 3
pages internet, puis on étudiera le cas général.
1
Le cas particulier
A l'issue de sa requête, un internaute est susceptible d'aller sur trois sites A, B , et C . On suppose que
l'internaute se déplace forcément vers un lien de la page Web où il se trouve.
• le site A contient un lien vers lui-même, un lien vers B , et un lien vers C ;
• le site B contient 5 liens vers lui-même, un vers A et un vers C ;
• le site C comporte un lien vers A, 7 liens vers B et 4 vers C .
Au départ, l'internaute choisit au hasard l'un des trois sites. Par la suite, la probabilité de passer d'un
site (à l'instant n) vers un autre (à l'instant n + 1) est proportionnelle au nombre de liens au premier site
vers le deuxième.
Pour n ∈ N∗ , on désigne par an , bn et cn les probabilités des évènements An , Bn et Cn dénis par :
• An : "à l'instant n l'internaute est sur le site A",
• Bn : "à l'instant n l'internaute est sur le site B ",
• Cn : "à l'instant n l'internaute est sur le site C ".
1. Représenter le modèle par un graphe (diagramme de transition).
2. Pour obtenir les probabilités de l'instant n + 1 en fonction de celles de l'instant 1, on utilise les
probas totales avec le sce de l'instant n :
(An , Bn , Cn ) est un sce, d'après les probabilités totales puis composées :
an+1
De même
3. On en déduit que
= P (An+1 ) = P (An ∩ An+1 ) + P (Bn ∩ An+1 ) + P (Cn ∩ An+1 )
=
P (An )PAn (An+1 ) + P (Bn )PBn (An+1 ) + P (Cn )PCn (An+1 )
=
1
1
1
an + bn + cn
3
7
12

1
cn
an+1 = 31 an + 17 bn + 12





7
bn+1 = 13 an + 57 bn + 12
cn





4
cn+1 = 13 an + 17 bn + 12
cn

 
an+1
1/3
 bn+1  = 1/3
cn+1
1/3

1/3
donc la matrice de transition est S = 1/3
1/3
1/7
5/7
1/7
1/7
5/7
1/7
 
1/12
an
7/12  bn 
4/12
cn
 

1/12
1/3 1/7 1/12
7/12 = 1/3 5/7 7/12.
4/12
1/3 1/7 1/3
3
ECE 2 - Mathématiques
Mr Dunstetter - ENC-Bessières 2014\2015
TP 4 et 5 : Chaînes de Markov
4. function X=deplacements(n)
X0=grand(1,1,'uin',1,3)
A=[ 1/3 , 1/7 , 1/12 ; 1/3 , 5/7 , 7/12 ; 1/3 , 1/7 , 1/3 ]
X=grand(n,'markov',A',X0)
X=[X0,X]
endfunction
puis (dans la commande scilab) :
--> deplacements(100)
Le site no 2 (donc le site B) a l'air beaucoup plus visité que les deux autres, qui semblent eux se valoir.
Cela s'explique par le fait que le site B est celui qui a le plus de lien sur chaque page, il est donc
normal qu'il soit le plus visité.
5. A=[ 1/3 , 1/7 , 1/12 ; 1/3 , 5/7 , 7/12 ; 1/3 , 1/7 , 1/3 ]
U=[1/3 ; 1/3 ; 1/3 ]
E= zeros(3,11)
E(:,1)=U
for i=2:11 do
U=A*U
E(:,i)=U
end
disp(E)
6. n=0:10
plot2d(n', [E(1,:)', E(2,:)',E(3,:)'])
On peut conjecturer que :
a ' 0, 1607143 ,
b ' 0, 625 ,
c ' 0, 2142857
7. A=[ 1/3 , 1/7 , 1/12 ; 1/3 , 5/7 , 7/12 ; 1/3 , 1/7 , 1/3 ]
[a,b]=spec(A)
X=a(:,1)
X=X/sum(X)
permet de visualiser dans X un vecteur propre associé à la valeur propre 1, puis en divisant par sa
somme de trouver le vecteur propre associé à la valeur propre 1 et dont la somme des composantes
vaut 1 : c'est la loi stable que l'on cherchait ; on obtient :
a = 0, 1607143 ,
b = 0, 625 ,
c = 0, 2142857
(les chires après la virgule disponibles ne permettent pas de faire la diérence avec la valeur de
U1 0, l'approximation de la question précédente était donc très très bonne).
2
Le cas général
On modélise l'ensemble des pages du WEB par un ensemble ni d'états E = {1, 2, ..., N }.
(le nombre N est bien très grand, de l'ordre de 1013 , mais reste ni).
On modélise les déplacements de l'internaute par une chaînes de Markov (Xn ) où la variable Xn est le
numéro de la page où se situe l'internaute après n déplacements.
Après avoir fait des études statistiques des comportements sur le Web, Google part du principe que :
• A l'instant 0, l'internaute choisit une page internet au hasard parmi les N pages.
• lorsqu'un internaute est sur une page i :
dans 85% des cas, il se déplacera au hasard vers l'une des li pages pointées par i (en pratique,
une page a au plus un lien vers une autre).
4
ECE 2 - Mathématiques
Mr Dunstetter - ENC-Bessières 2014\2015
TP 4 et 5 : Chaînes de Markov
dans 15% des cas, il se déplacera au hasard vers n'importe laquelle des pages du Web.
La matrice de transition est alors dénie par :
P (j, i) =
0, 85 ×
1
li
+ 0, 15 ×
0, 15 ×
1
N
1
N
si i pointe vers j
sinon
Nous allons illustrer l'algorithme de Page Rank sur un groupe autonome de N = 5 pages diérentes,
ayant des liens pointant les pages les unes vers les autres, selon le diagramme suivant :
1
2
4
3
5
1. A l'aide du graphique ci-dessous, on peut calculer, pour chaque état i, le nombre de pages vers
lesquelles il pointe et obtenir :
0,85
2

0, 03
 0,85 + 0, 03
 3
P =
 0,850, 03

3 + 0, 03
0,85
3 + 0, 03
donc
+ 0, 03
0, 03
0, 03
0,85
2 + 0, 03
0, 03
0, 03
0,85
2 + 0, 03
0, 03
0, 03
0,85
2 + 0, 03
0, 03
0,85
2 + 0, 03
, 0, 03
0, 03
0,85
2 + 0, 03

0, 03
0, 85 + 0, 03


0, 03


0, 03
0, 03

0, 03 0, 455 0, 03 0, 03 0, 03
0, 96/3 0, 03 0, 455 0, 455 0, 88


0, 03 0, 03 0, 03 0, 03
P =
.
 0, 03
0, 96/3 0, 455 0, 03 0, 03 0, 03
0, 96/3 0, 03 0, 455 0, 455 0, 03

2. P=[ 0.03 , 0.455 , 0.03 , 0.03 , 0.03 ; 0.96/3 , 0.03 , 0.455 , 0.455 , 0.88 ; 0.03 , 0.03 , 0.03
[a,b]=spec(P)
X=a(:,1)
X=X/sum(X)
permet d'obtenir pour loi stationnaire :


0, 1819705
0, 3591360



X=
0, 0298912
0, 2345506
0, 1944516
3. P=[ 0.03 , 0.455 , 0.03 , 0.03 , 0.03 ; 0.96/3 , 0.03 , 0.455 , 0.455 , 0.88 ; 0.03 , 0.03 , 0.03
U=[1/5 ; 1/5 ; 1/5 ; 1/5 ; 1/5 ]
E= zeros(5,11)
E(:,1)=U
for i=2:11 do
U=P*U
E(:,i)=U
end
disp(E)
On remarque que la loi de Xn s'approche de la loi stationnaire au fur et à mesure : la loi de la chaîne
de Markov converge donc vers la loi stable.
4. n=0:10
plot2d(n', [E(1,:)', E(2,:)',E(3,:)'])
5
ECE 2 - Mathématiques
Mr Dunstetter - ENC-Bessières 2014\2015
TP 4 et 5 : Chaînes de Markov
5. Les pages les plus populaires sont, dans l'ordre : 2, 4, 5, 1 , 3.
Cela s'explique très bien : la ligne 2 de P (ou bien le graphique) permet de voir que la page 2 est la
seule vers laquelle pointe toutes les autres : elle est donc clairement la plus populaire.
Si on se contente de regarder P , on l'impression qu'on a ensuite un groupe de 2 pages à peut près
égales (4 et 5) avec une préférence pour 5), puis la page 1, puis la page 3.
Ce n'est pas le cas dans la pratique : c'est l'eet de la page 2, très populaire donc très visitée :
comme 5 n'a pas de lien sur 2, malgré un grand nombre de lien, elle est "rattrapée" en popularité
par la page 4 (qui la dépasse même) et la page 1 qui ont un lien sur la page la plus populaire.
III
Autres chaînes de Markov
Exercice 2.
Un joueur joue à un jeu au casino. On suppose qu'il dispose d'une fortune de a euros et que celle du casino
est de b euros. A chaque tour, il gagne un euro avec probabilité p ∈]0, 1[ ou perd un euro (donné au casino)
avec probabilité q = 1 − p. Le jeu s'arrête lorsque le joueur ou le casino se retrouve ruiné.
On peut alors modéliser la situation par une chaîne de Markov homogène (Xn ) où Xn représente la fortune
du joueur à l'issue du n-ième tour.
1. On donne la matrice P pour a + b = 5, avec des état de 0 à 5, il faut ensuite agrandir selon le même
principe :


1
0

0
P =
0

0
0
q
0
p
0
0
0
0
q
0
p
0
0
0
0
q
0
p
0
0
0
0
q
0
p
0


0

0

0
1
2. function P=transition(a,b,p)
q=1-p
P=zeros(a+b+1,a+b+1)
P(1,1)=1
P(a+b+1,a+b+1)=1
for i=1:a+b-1 do
P(i,i+1)=q
end
for i=3:a+b+1 do
P(i,i-1)=p
end
endfunction
3. function []=ruinejoueur(p,a,b,n)
X0=a+1
P=transition(a,b,p)
X=grand(n,'markov',P',X0)
X=[X0,X]
plot2d(0:n,X)
endfunction
4. Avec les paramètres donnés, le joueur est systématiquement ruiné après un certain nombre de parties.
Avec d'autres valeurs on peut remarquer que si p < 0, 5 (assez largement), le joueur est systématiquement ruiné ; si p > 0, 5 (assez largement), le casino est systématiquement ruiné. Pour p = 0, 5
(notamment diérent de 0,5 mais assez proche), cela dépend des fortunes de départ et des fois mais
le jeu se termine à chaque fois.
5. On peut montrer que si b est très grand et si p ≤ 1/2 , le joueur nit ruiné. Vérier ces résultats sur
des simulations.
6
ECE 2 - Mathématiques
Mr Dunstetter - ENC-Bessières 2014\2015
TP 4 et 5 : Chaînes de Markov
De même vérier que si p > 1/2, il est possible que ce soit le casino qui nisse ruiné, même si b est
très grand.
6. function f=ruinejoueurfrequence(p,a,b,n)
P=transition(a,b,p)
X=[]
for i=1:1000 do
Y0=a+1
Y=grand(n,'markov',P',Y0)
Y=[Y0,Y]
X=[X;Y]
end
Z=X(:,n+1)
f=length(find(Z==1))/1000
endfunction
7. La fonction qui suit crée la matrice P , puis eectue 1000 simulations de X1 , ..., Xi , ..., Xn (l'état
de la chaîne de Markov à l'instant i) qu'elle range dans un tableau X (une simulation par ligne, la
colonne i donne donc les 1000 simulations de Xi ).
Ensuite pour chaque i entre 1 et n, elle cherche, parmi les 1000 simulations de Xi , combien valent
0 (qui correspond à un joueur ruiné à l'instant i) et divise ce nombre par 1000 pour obtenir la
fréquence de ruine du joueur à l'instant i dans les 1000 simulations (dans la variable pest).
Enn elle trace le graphe de ces fréquences pour i allant de 1 à n + 1, pour observer si elle converge
vers une valeur (qui serait alors la probabilité de la ruine du joueur pour les valeurs de a, b et p
données).
function []=ruinejoueurproba(p,a,b,n)
P=transition(a,b,p)
for i=1:1000 do
simu=grand(n,'markov',P',a+1)-1
X(i,:)=[a,simu]
end
for i=1:n+1 do
pest(i)=length(find(X(:,i)==0))/1000
end
clf(),plot(1:(n+1),pest)
endfunction
On obtient à chaque fois, à condition d'aller à un assez grand nombre de simulations, une courbe
qui se stabilise autour d'une valeur.
Exercice 3. On veut modéliser la mobilité sociale d'une population. Pour cela, on s'intéresse à la profession des individus d'une génération, comparée avec celle de leurs parents.
L'INSEE fournit par exemple le tableau suivant, à partir d'un sondage eectué sur 18373 individus de
sexe masculin :
Agr.
Indép.
Cadres
Interm.
Emp.
Ouvr.
Total
Agr.
621
237
220
405
214
1126
2823
Indép.
24
550
479
439
183
499
2174
Cadres
10
127
749
325
100
116
1427
Interm
9
138
552
568
182
352
1801
Emp.
7
210
490
643
393
770
2513
Ouvr.
43
621
700
1520
855
3896
7635
où on a utilisé les abréviations : Agr.=Agriculteur, Indép.=Indépendants,Interm.=Professions intermédiaires, Emp.=Employés et Ouvr.=Ouvriers.
7
ECE 2 - Mathématiques
Mr Dunstetter - ENC-Bessières 2014\2015
TP 4 et 5 : Chaînes de Markov
L'eectif placé ligne i colonne j représente le nombre de sondés qui exerce la profession i, sachant que leur
père exerçait la profession j .
Une colonne donne donc les eectifs des catégories professionnelles des sondés, pour une catégorie socioprofessionnelle xée pour le père : on visualise ainsi sur chaque colonne les "transitions intergénérationnelles".
On modélise la situation par une chaîne de Markov (Xn ) à 6 états : Xn donne la catégorie socioprofessionnelle de l'individu de la n-ième génération.
1. La matrice de transition est la matrice :

621/2863
 237/2863

 220/2823
P =
 405/2823

 214/1126
1126/2823
24/2174
550/2174
479/2174
439/2174
183/2174
499/2174
10/1427
127/1427
749/1427
325/1427
100/1427
116/1427

9/1801
7/2513
43/7635
138/1801 210/2513 621/7635 

552/1801 490/253 700/7635 

568/1801 643/2513 1520/7635

182/1801 393/2513 855/7635 
352/1801 770/2513 3896/7635
et on la crée sous scilab avec la fonction suivante :
function A=transition()
A=zeros(6,6)
A(1,:)=[621 , 24 , 10 , 9 , 7 , 43]
A(2,:)=[237 , 550 , 127 , 138 , 210 , 621]
A(3,:)=[220 , 479 , 749 , 552 , 490 , 700]
A(4,:)=[405 , 439 , 325 , 568 , 643 , 1520]
A(5,:)=[214 , 183 , 100 , 182 , 393 , 855]
A(6,:)=[1126 , 499 , 116 , 352 , 770 , 3896]
for i=1:6 do
A(:,i)=A(:,i)/sum(A(:,i))
end
endfunction
2. function Y=socio(x0,n,r)
P=transition()
Y=zeros(1,r)
for k=1:r do
X=grand(n,'markov',P',x0)
Y(1,k)=X(1,n)
end
endfunction
3. n=1000;r=5000;
for x0=1:6 do
Y=socio(x0,n,r);
F=zeros(1:6)
for k=1:6 do
F(k)=length(find(Y==k))/r
end;
scf(x0-1), clf(), bar(1:6,F)
end
Ce programme trace, pour chaque catégorie socio-professionnelle de départ, les fréquences d'appartenance à chaque catégorie socio-professionnelle au bout de 1000 générations, sur 1000 simulations.
On obtient, quel que soit l'était de départ, grosso modo le même graphique : comme les fréquences,
sur un grand nombre de simulation, approximent la loi, on a vérié que quel que soit l'état de départ,
la loi de la chaîne de Markov convergeait vers la même loi limite : ici c'est donc la seule loi stable
par la matrice de transition, contrairement au cas de la ruine du joueur par exemple.
4. On obtient, quel que soit l'était de départ, grosso modo le même graphique : comme les fréquences,
sur un grand nombre de simulation, approximent la loi, on a vérié que quel que soit l'état de départ,
8
ECE 2 - Mathématiques
Mr Dunstetter - ENC-Bessières 2014\2015
TP 4 et 5 : Chaînes de Markov
la loi de la chaîne de Markov convergeait vers la même loi limite : ici c'est donc la seule loi stable
par la matrice de transition, contrairement au cas de la ruine du joueur par exemple. On vérie que
cette loi stable est unique en faisant :
P=transition()
[a,b]=spec(P)
X=a(:,1)
X=X/sum(X)
et on retrouve dans X les fréquences des graphiques précédents.
Un sociologue pourrait en déduire que le brassage socio-professionnel dans notre société est bon,
puisque quel que soit la catégorie socio-professionnelle d'une personne, au bout d'un certain nombre
de générations, il a les mêmes probabilités d'appartenir à chaque autre catégorie socio-professionnelle.
Il en déduirait également que la profession d'agriculteur va de venir de plus en plus minoritaire dans
la société.
Par contre ce modèle est un peut trop simple car il manque de précisions sur le devenir des catégories : il n'est pas sûr par exemple qu'un ls d'ouvrier qui est devenu cadre n'ait pas plus de chances
de voir son ls redevenir ouvrier que la moyenne des autres cadres (car même si son père ne l'est
pas, d'autre personnes comme son grand-père qui sont ouvriers gravitent autour de lui et peuvent
lui donner "envie" de faire ce type de métier, etc....).
Mais le plus gros défaut reste le fait que les probabilités de passer d'une catégorie à une autre sont
en constant changement : elles ne risquent pas de rester identiques pendant 1000 générations (ni
même pendant une seule !)
9