Corrigé du TP des 19 et 20 septembre 2002
 


Exemple d'utilisation

Les fonctions suivantes sont les fonctions de test utilisées dans les différents exercices : essai, essaiprim, test, bord.



% Corrige de l'exercice 1
% -----------------------

% Cet exercice permet de garder un echantillon x.

clear;
close all;
sortie=0;
while (sortie==0),
    r=input('Type de variable recherchee : ''bin'', ''exp'', ''beta'', ''poi'', ''FIN'' ?');
    switch r,
    case 'bin'
        disp('    Variable aleatoire binomiale de parametres n et p,');
        n=input('    n ?');
        p=input('    p ?');
        m=input('    Taille de l''echantillon ?');
        u=dbinom([0:n],n,p);
        [x,y]=escalier([0:n],u,0,0,1,'k',1);
        titre=sprintf('Fonction de repartition de la loi binomiale de parametres %d et %1.3f\n',n,p);
        v=rbinom(m,n,p);
        [z,t]=repartition(v,0,1,'r',0);
        titre=strcat(titre,sprintf('\n et fonction de repartition empirique avec %d tirages\n',m));
        title(titre);
        clear x;
        x=u;
    case 'exp'
        disp('    Variable aleatoire exponentielle de parametre lambda,');
        lambda=input('    lambda ?');
        m=input('    Taille de l''echantillon ?');
        u=ones(size([0:0.02:5]))-exp(-[0:0.02:5]);
        x=courbe((1/lambda)*[0:0.02:5],u,2,'k',1);
        titre=sprintf('Fonction de repartition de la loi exponentielle de parametre %1.3f\n',lambda);
        v=rexpweib(m,lambda);
        [z,t]=repartition(v,0,2,'r',0);
        titre=strcat(titre,sprintf('\n et fonction de repartition empirique avec %d tirages\n',m));
        title(titre);
        clear x;
        x=v;
    case 'beta'
        disp('    Variable aleatoire beta de parametres i et j,');
        i=input('    i ?');
        j=input('    j ?');
        m=input('    Taille de l''echantillon ?');
        u=pbeta([0:0.01:1],i,j);
        x=courbe([0:0.01:1],u,3,'k',1);
        titre=sprintf('Fonction de repartition de la loi beta de parametres %d et %d\n',i,j);
        v=rbeta(m,i,j);
        [z,t]=repartition(v,0,3,'r',0);
        titre=strcat(titre,sprintf('\n et fonction de repartition empirique avec %d tirages\n',m));
        title(titre);
        clear x;
        x=v;
    case 'poi'
        disp('    Variable aleatoire de Poisson de parametre lambda,');
        lambda=input('    lambda ?');
        m=input('    Taille de l''echantillon ?');
        u=[exp(-lambda)];
        for j=1:floor(5*lambda),
            u(j+1)=u(j)*lambda/j;
        end
        [x,y]=escalier([0:floor(5*lambda)],u,0,0,4,'k',1);
        titre=sprintf('Fonction de repartition de la loi de Poisson de parametre %1.3f\n',lambda);
        v=rpoiss(m,lambda);
        [z,t]=repartition(v,0,4,'r',0);
        titre=strcat(titre,sprintf('\n et fonction de repartition empirique avec %d tirages\n',m));
        title(titre);
        clear x;
        x=v;
    case 'FIN'
        sortie=1;
    otherwise
        disp('Entrez une des variables proposees !');
    end
end



% Corrige de l'exercice 2
% -----------------------

% Cet exercice permet de garder un echantillon x

clear;
close all;
sortie=0;
while (sortie==0),
    r=input('Type de generation de variable recherchee : ''discrete'', ''rejet'', ''inversion'', ''FIN'' ?');
    switch r,
    case 'discrete'
        disp('    Generation de v.a. discrete de vecteur de probabilite p');
        p=input('    p ?');
        m=input('   Taille de l''echantillon ?');
        x=[];
        for i=1:m,
            x=[x,rdiscrete(p)];
        end
        [z,t]=escalier([1:length(p)],p,0,0,1,'k',1);
        titre=sprintf('Fonction de repartition de la loi dicrete de parametres %s\n',num2str(p));
        [z,t]=repartition(x,0,1,'r',0);
        titre=strcat(titre,sprintf('\n et fonction de repartition empirique avec %d tirages\n',m));
        title(titre);
    case 'rejet'
        disp('    Generation de v.a. continue de densite f<=c par methode de rejet sur l''intervalle [a,b]');
        nom=input('    f ?');
        a=input('    Borne inferieure ?');
        b=input('    Borne superieure ?');
        c=input('    Borne de la densite ?');
        m=input('    Taille de l''echantillon ?');
        x=[];
        for i=1:m,
            x=[x,rejet(nom,a,b,c)];
        end
        y=courbeintf(nom,a,b,2,'k',1);
        titre=sprintf('Fonction de repartition de la loi continue de densite %s\n',nom);
        [z,t]=repartition(x,0,2,'r',0);
        titre=strcat(titre,sprintf('\n et fonction de repartition empirique avec %d tirages\n',m));
        title(titre);
    case 'inversion'
        disp('    Generation de v.a. continue de fonction de repartition F sur l''intervalle [a,b]');
        nom=input('    F ?');
        a=input('    Borne inferieure ?');
        b=input('    Borne superieure ?');
        m=input('    Taille de l''echantillon ?');
        x=[];
        for i=1:m,
            x=[x,inversion(nom,a,b)];
        end
        y=courbeF(nom,a,b,3,'k',1);
        titre=sprintf('Fonction de repartition %s\n',nom);
        [z,t]=repartition(x,0,3,'r',0);
        titre=strcat(titre,sprintf('\n et fonction de repartition empirique avec %d tirages\n',m));
        title(titre);
    case 'FIN'
        sortie=1;
    otherwise
        disp('Entrez une des techniques proposees !');
    end
end



% Corrige de l'exercice 3
%------------------------

% On suppose que l'on dispose d'un vecteur de tirages aleatoires x.

f=input('Numero de la figure ?');
t=input('Nom de la variable aleatoire ?');
et=input('Esperance theorique (''?'' si inconnue) ?');
n=length(x);
sn=cumsum(x);
d=diag(1:n);
e=inv(d);
snsn=sn*e;
figure(f);
clf;
plot(snsn);
hold on;
if et~='?'
    plot([0,n],[et,et],'r-');
end

titre=sprintf('Loi des grands nombres pour une variable aleatoire %s',t);
title(titre);



% Corrige de l'exercice 4
%-------------------

disp('Marches au hasard dans le plan R^2');
clear;
close all;
sortie=0;
couleur='rgbky';
while (sortie==0),
    ter=input('  Pas i.i.d. sur le reseau (1) ou pas d''orientation aleatoire (2) ou ''FIN'' ?');
    switch ter,
    case 1
        n=input(' Nombre de realisations pour le dessin (<5) ?');
        n=max(1,min(n,4));
        figure(1);
        clf;
        hold on;
        axis equal;
        axis([-5 5 -5 5]);
        temps='';
        x=[];
        y=[];
        t=[];
        for i=1:n,
            clear x;
            clear y;
            clear t;
            x=[0];
            y=[0];
            t=0;
            while max(abs(x),abs(y))<5,
                t=t+1;
                [pa,pas]=pasreseau(x(t),y(t));
                x=[x,pa];
                y=[y,pas];
            end
            temps=strcat(temps,num2str(t));
            if i~=n,
                temps=strcat(temps,' et   ');
            end
            dessinmarche(x,y,num2str(couleur(i)));
            disp('       <pause>');
            pause;
        end
        titre=sprintf('Marche au hasard classique, temps de sortie des %d trajectoires egaux a %s',n,temps);
        title(titre);
        n=input('Nombre de realisations pour l''estimation du temps de sortie issu de (0,0) ?');
        t=[];
        tt=[];
        for i=1:n,
            clear x;
            clear y;
            clear tt;
            x=[0];
            y=[0];
            tt=0;
            while max(abs(x),abs(y))<5,
               tt=tt+1;
               [pa,pas]=pasreseau(x(tt),y(tt));
               x=[x,pa];
               y=[y,pas];
           end
           t=[t,tt];
       end
       et=sum(t)/n;
       s=sprintf('L''estimateur du temps de sortie issu de (0,0) pour %d tirages est egal a %f.',n,et);
       disp(s);
       disp('Resolution par marche aleatoire du probleme de Dirichlet pour la fonction g');
       g=input('g ?');
       for i=-5:5,
           v(i+6,1)=feval(g,i,-5);
           v(i+6,11)=feval(g,i,5);
           v(1,i+6)=feval(g,-5,i);
           v(11,i+6)=feval(g,5,i);
       end
       t=[];
       x=[];
       y=[];
       tt=[];
       for i=-4:4,
           for j=-4:4,
               clear t;
               t=[];
               for k=1:n,
                 clear x;
                 clear y;
                 clear tt;
                 x=[i];
                 y=[j];
                 tt=0;
                 while max(abs(x),abs(y))<5,
                     tt=tt+1;
                    [pa,pas]=pasreseau(x(tt),y(tt));
                    x=[x,pa];
                    y=[y,pas];
                 end
                 t=[t,feval(g,x(tt),y(tt))];
               end
               v(i+6,j+6)=sum(t)/n;
           end
       end
       figure(2);
       surf(v);
       title(sprintf('Resolution approchee par marche au hasard du probleme de Dirichlet pour la fonction %s',g));
       % Pour une comparaison avec la solution approchee obtenue par resolution effective par methode
       % de Jacobi ou Gauss-Seidel ou relaxation, se referer aux tps de Stephane.
    case 2
        n=input(' Nombre de realisations pour le dessin (<5) ?');
        n=max(1,min(n,4));
        figure(1);
        clf;
        hold on;
        axis equal;
        axis([-5 5 -5 5]);
        temps='';
        x=[];
        y=[];
        t=[];
        for i=1:n,
            clear x;
            clear y;
            clear t;
            x=[0];
            y=[0];
            t=0;
            while max(x.^2+y.^2)<25,
                t=t+1;
                [pa,pas]=pasincline(x(t),y(t));
                x=[x,pa];
                y=[y,pas];
            end
            temps=strcat(temps,num2str(t));
            if i~=n,
                temps=strcat(temps,' et   ');
            end
            dessinmarche(x,y,num2str(couleur(i)));
            disp('       <pause>');
            pause;
        end
        cercle(0,0,5);
        titre=sprintf('Marche au hasard avec pas d''angle aleatoire, temps de sortie des %d trajectoires egaux a %s',n,temps);
        title(titre);
        n=input('Nombre de realisations pour l''estimation du temps de sortie issu de (0,0) ?');
        t=[];
        tt=[];
        for i=1:n,
            clear x;
            clear y;
            clear tt;
            x=[0];
            y=[0];
            tt=0;
            while max(x.^2+y.^2)<25,
               tt=tt+1;
               [pa,pas]=pasincline(x(tt),y(tt));
               x=[x,pa];
               y=[y,pas];
           end
           t=[t,tt];
       end
       et=sum(t)/n;
       s=sprintf('L''estimateur du temps de sortie issu de (0,0) pour %d tirages est egal a %f.',n,et);
       disp(s);
    case 'FIN'
        sortie=1;
    end
end



% Corrige de l'exercice 5
%------------------------

% On utilisera dans cet exercice une variable aleatoire obtenue par inversion
% de la fonction de repartition F. On utilisera un estimateur sans biais de la variance.

clear;
close all;
F=input('Nom de la fonction de repartition a inverser ?');
a=input('Borne gauche de son support ?');
b=input('et borne droite ?');
n=input('Taille de l''echantillon ?');
m=input('Nombre de series de tirages pour la loi empirique ?');
for i=1:n,
    for j=1:m,
        x(i,j)=inversion(F,a,b);
    end
end
sn=cumsum(x)';
d=diag(1:n);
e=inv(d);
snsn=sn*e;
ex=input('Esperance theorique (''?'' si inconnue) ?');
if ex=='?'
    ex=approxesp(F,a,b);
end
va=input('Variance theorique (''?'' si inconnue) ?');
if va=='?'
    va=approxmdeux(F,a,b)-ex*ex;
end
snsn=snsn-ex*ones(size(snsn));
normalise=(1/sqrt(va))*snsn*sqrt(d);
figure(1);
axis([-3 3 0 1]);
M=moviein(length(1:10:n));
cpt=1;
for k=1:10:n,
    [x,y]=repartition(normalise(:,k),0,1,'r',1);
    hold on;
    u=[-3:0.02:3];
    v=pnorm(u,0,1);
    t=courbe(u,v,1,'b',0);
    M(:,cpt)=getframe;
    cpt=cpt+1;
end
movie(M)