Category Archives: Randomness

De la créativité des gangsters

Pendant mon séjour récent en Nouvelle Angleterre, j’ai survolé le livre de Leonard Mlodinow, the drunkard’s walk. Et au hasard de mes lectures, je suis tombé sur la petite histoire suivante

Autrement dit, en utilisant les cinq derniers chiffres d’une quantité économique comme le us treasury balance, on aurait un générateur de nombre aléatoire…
Par contre la suite est un peu plus surprenante,


http://freakonometrics.hypotheses.org/files/2015/12/sammy.gifLa loi de Benford (que j’avais pu évoquer ici ou ) parle des premiers chiffres, mais cette fois on parle deslast five digits. Donc visiblement l’évoquation n’est pas pertinente ici. Mais qui sait ? Ca reste malgré tout une histoire intéressante. Considérons – histoire de tester cette légende – les deux sources de données suivantes,http://treasurydirect.gov/ et http://economagic.com/. Sur ce dernier site, un petit travail de mise en forme des données est nécessaire.

b1=read.table("http://freakonometrics.free.fr/debtus1.txt",
   header=TRUE,sep="\t")
b2=read.table("http://freakonometrics.free.fr/debtus2.txt",
   header=TRUE)
X1=as.character(b1$Dollar.Amount)
n1=nchar(X1)
Y1=substr(X1,n1-8,n1-3)
X1=as.numeric(substr(Y1,1,2))*1000+as.numeric(substr(Y1,4,6))
x=X1/100000
X2=b2$DEBT
Y2=trunc(as.numeric(X2))
X2=as.character(Y2)
n2=nchar(X2)
Y2=substr(X2,n2-4,n2)
y=as.numeric(Y2)/100000
y=y[y<1]

Pour rappel, un générateur aléatoire (standard) vérifie deux propriétés importantes

  • les nombres doivent être tirés suivant une loi uniforme sur [0,1], i.e. ici, si on divise les nombres à 5 chiffres par 10000,
  • les tirages doivent être indépendants entre eux.

La première propriété semble assez naturelle, et correspond à l’histoire racontée dans un commentaire posté ici (expliquant comment un casino avait été au bord de la faillite car une roulette faisait sortir certains chiffres trop souvent, et j’essayais de comprendre comment utiliser l’information qu’un chiffre sort plus souvent). La seconde est probablement encore plus importante.

  • Visualisation de la distribution

La première idée est de visualiser la densité de nos séries de chiffres. Pour éviter les problèmes de bord (et comme c’est juste en introduction) on va utiliser un histogramme, et pas une estimation à noyau.

hist(x,col="red")
hist(y,col="blue")

On obtient pour la première série la courbe rouge, et pour la seconde la courbebleue,

On note qu’a priori, pour la première série, l’hypothèse d’uniformité n’est peut être pas la plus réaliste…

  • Test de Kolmogorov-Smirnov

On peut aussi mettre en œuvre un test de Kolmogorov-Smirnov afin de tester si la loi uniforme est adaptée.:

> ks.test(x,"punif")
 
	One-sample Kolmogorov-Smirnov test
 
data:  x 
D = 0.1047, p-value = 0.01645
alternative hypothesis: two-sided 
 
> ks.test(y,"punif")
 
	One-sample Kolmogorov-Smirnov test
 
data:  y 
D = 0.0456, p-value = 0.3581
alternative hypothesis: two-sided

On retrouve ici confirmée l’intuition précédante: la loi uniforme est pertinente pour la seconde série, pas la première.

  • Les autocorrélations de la série

Travaillons uniquement sur la seconde série. On peut étudier l’autocorrélation de notre série de nombres, ou peut-être un peu plus malin, sur les quantiles gaussiens associés (les autocorrélations étant intéressantes pour les séries gaussiennes),:

plot(acf(y))
plot(acf(qnorm(y)))

ie. pour la série brute

et pour la série normalisée,

Bref, on pourrait être tenté de valider l’hypothèse d’indépendance entre les tirages.
  • Run test (de Wald–Wolfowitz)

L’idée est de comparer une série de chiffres à la médiane, s’ils sont plus grands, on note + (ou A) et sinon – (ou B). On crée alors une série du genre “+++−−++−−++++++−−−” et on compte les séries de + et les séries de -, les runs,

library(lawstat)
runs.test(y,plot=TRUE)

Runs Test - Two sided

data:  y
Standardized Runs Statistic = -0.2462, p-value = 0.8055

Bref, la légende me semble à prendre avec des pincettes (car fonction de la source considérée), même si l’idée est intéressante (si l’on met de côté les aspects d’aléa moral). Et l’analyse sur la loi de Benford ne semble pas valide: les derniers chiffres sur les grands nombres ne se comportent pas du tout comme les premiers.

Cursed numbers ?

In Lost, Hugo “Hurley” Reyes played the numbers 4, 8, 15, 16, 23 and 42 at the lottery, and ended up winning the $114-million jackpot. And over the ensuing weeks, everyone around him seems to suffer increasingly bad luck: Hurley’s grandfather dies of a heart attack, his brother’s wife walks out, his mother breaks her ankle while the house Hurley bought her goes up in flames, and Hurley himself is falsely arrested.
Anyway, last week (here) 4 numbers (out of 6) appeared at the lottery in LA. As pointed out by Xian (here), the odds were not that small, i.e. it is a 1‰ chance,

http://freakonometrics.blog.free.fr/public/perso/probaloto.png
Hence, with one lottery per week, the return period is 16 years. Note this percentage is very close to what we did observe on the French lottery (below the statistics in ‰, from here, in a zip file),

> loto=read.table("loto.csv",dec=",",header=TRUE,sep=";")
> ntirage=nrow(loto)
> loto=loto[51:ntirage,]
> ntirage=nrow(loto)
> N=as.matrix(loto[,c("boule_1","boule_2","boule_3",
   "boule_4","boule_5","boule_6")])
> P=rep(NA,nrow(N))
> for(s in 1:nrow(N)){
+ P[s]=sum(N[s,1]%in%c(4, 8, 15, 16, 23, 42)+
+          N[s,2]%in%c(4, 8, 15, 16, 23, 42)+
+          N[s,3]%in%c(4, 8, 15, 16, 23, 42)+
+          N[s,4]%in%c(4, 8, 15, 16, 23, 42)+
+          N[s,5]%in%c(4, 8, 15, 16, 23, 42)+
+          N[s,6]%in%c(4, 8, 15, 16, 23, 42))
+ }
> table(P)/nrow(N)*1000
P
         0          1          2          3          4
435.732113 405.366057 137.271215  19.966722   1.663894

But what about the full sequence…? Imagine that in France, at the official lottery, the exact sequence played by Hugo appears. What a coincidence.The probability that the sequence appears, assuming that there are 48 possible numbers in the lottery, is

http://freakonometrics.blog.free.fr/public/perso/lotto2.png
i.e. the expected number of draws we need before seeing that sequence for the first time is almost a billion.
Now if we look at all official lotteries around the world, say 100 per week, what is the probability to see Hurley’s sequence shows up – at least once – in 25 years (assuming that after 25 years, no one will remember Lost and those cursed numbers) ? It looks like it is a 1% chance…
http://freakonometrics.blog.free.fr/public/perso/proba-lost.png
So let us wait and see…

Mandelbrot, fractals and counterexamples in applied probability

Benoît Mandelbrot died yesterday. Like most of the blogs dealing with applied mathematics, it looks like I have to mention this event. Unfortunately, I don’t know much about fractals…

The first time I heard about Mandelbrot and chaos was when I have been working on fractional time series (see eghere). Murad Taqqu gave a very interesting short course in Paris, and I have been using it in two papers (actually one more should appear soon in Climate Change).

The second time was in Québec (city), five years ago, when Roger Nelsen gave a talk on copulas with fractal suppport (here). By that time, we were finishing our paper with Alessandro (in mathematical finance, and limiting theorems). I remember adding a Remark in the paper (that can be found here), since using that kind of copulas was a nice way to show that, without sufficient regularity conditions, the limit we were looking for had no sense.

A few months after, with Johan, in a paper on pitfalls on lower tail dependence for Archimedean copulas (here), we used again this fractal construction (here applied to Archimedean copulas) to find a nice counterexample to a (false) theorem on regular variation. Again, the goal was to understand how the dependence structure of https://perso.univ-rennes1.fr/arthur.charpentier/latex/mandel01.png given https://perso.univ-rennes1.fr/arthur.charpentier/latex/mandel02.png and https://perso.univ-rennes1.fr/arthur.charpentier/latex/mandel03.png changed when https://perso.univ-rennes1.fr/arthur.charpentier/latex/mandel04.png goes to 0. In the case of Archimedean copulas, the copula of the conditional pair is still Archimedean (with another generator, except for Clayton copula). The graph below show how https://perso.univ-rennes1.fr/arthur.charpentier/latex/mandel05.png changes, as https://perso.univ-rennes1.fr/arthur.charpentier/latex/mandel04.png decreases… Actually, I draw https://perso.univ-rennes1.fr/arthur.charpentier/latex/mandel06.png since Archimedean generators are not unique.

I have also decided to plot https://perso.univ-rennes1.fr/arthur.charpentier/latex/mandel07.png.

Here, we see that there is no way of talk about a possible limit for the conditional copula because of a fractal behavior in 0 of the generator (even if my fractal are not as nice as the one you can find on the internet….). So thanks Benoît for giving us a nice toy to build interesting counterexamples !

Marches aléatoires en dimension 1, 2 ou 3, ou histoire éthylique

Revenons à nos histoires d’alcoolique un instant (les précédentes sont en ligne iciet ). L’avant dernière fois (ici), nous avions parlé de la marche aléatoire, et du temps mis avant de trouver la sortie de la rue de la soif.

Et j’avais évoqué la probabilité la probabilité qu’il ne récupère pas les clés (i.e. qu’il ne repasse pas par l’origine) s’il peut sortir un jour de la rue. En probabilité, on appelle cela la probabilité de repasser par l’origine et cela vérifie une propriété étrange dès lors que l’on change de dimension… On va suivre notre alcoolique et noter https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy01.png sa position à la date https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy02.png, puis poser

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy03.png

autrement dit la probabilité de se retrouver à sa position initiale à la date n (et de retrouver ainsi ses clés), le https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy04.png pouvait être en n’importe quelle dimension, on posera

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy05.png

correspondant à la probabilité de premier retour à la position initiale, et enfin

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy06.png

qui est la probabilité de revenir un jour au point de départ, c’est à dire https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy04.png… C’est cette dernière grandeur qui nous intéresse….

  • De la probabilité de récupérer ses clés dans la rue

La rue est ici infinie (contrairement à l’autre fois, où on pouvait en sortir, soit vers la place de Lices, soit vers la place Sainte Anne). Autrement dit, on s’intéresse à la probabilité qu’une marche aléatoire sur https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy08.png issue de https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy04.png repasse un jour par https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy04.png. S’il va à droite ou à gauche avec probabilités https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethyl99.png et https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethyl98.png respectivement, notons que https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy11.png, car il faut un nombre de pas pairs pour revenir à l’origine, et

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy12.png

Il est possible de montrer que

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy13.png

et en particulier

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy14.png

On peut relier les https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy15.png et les https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy16.png en notant que

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy17.png

Or https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy16.png constitue une distribution de probabilité, dont on peut définir la fonction génératrice

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethyl18.png

On notera que que l’on peut écrire

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy19.png

soit, en utilisant l’expression donnée auparavant

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethyl19.png

Soit finalement

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethyl20.png

ce qui se traduit par https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy21.png si la série https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy22.png converge. Aussi,

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethyl25.png

donc si https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy26.png, notre alcoolique repassera presque sûrement par le point d’origine, alors que si https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethyl97.png, il ne repasse pas par le bar où il a oublié ses clés avec probabilité https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy27.png|.

  • De la probabilité de récupérer ses clés sur une place

La place est ici infinie de part et d’autre. Autrement dit, on s’intéresse à la probabilité qu’une marche aléatoire sur https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy09.png issue de https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy04.png repasse un jour par https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy04.png. Il a alors 1 chance sur 4 d’aller devant lui, autant d’aller en arrière, d’aller à gauche ou d’aller à droite. On peut noter que là encore
https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy11.png, car il faut un nombre de pas pairs pour revenir à l’origine, et

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethyl26.png

Or on peut réécrire cette somme en notant que

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy28.png

et que

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethyl29.png

aussi,

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethyl30.png

que l’on peut approcher (via la formule de Stirling) sous le forme

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy31.png

La série https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy22.png est alors divergente, et donc https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethyl31.png. Aussi, même sur une place, notre alcoolique reviendra presque sûrement à sa place initiale (et ce une infinité de fois, car une fois qu’il y est revenu et qu’il en repart, il repassera presque sûrement à cette position, etc).

  • De la probabilité de récupérer ses clés dans l’espace

Bon, cette fois notre alcoolique passe en dimension 3…. j’admets, il faut de l’imagination. On s’intéresse à la probabilité qu’une marche aléatoire sur https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy10.png issue de https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy04.png repasse un jour par https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethy04.png. Il a alors 1 chance sur 6 d’aller devant lui, autant d’aller en arrière, d’aller à gauche ou d’aller à droite, d’aller en haut ou d’aller en bas.

Il y a désormais 6 voisins, et en raisonnant comme avant (le maillage est un peu plus complexe), on obtient que

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethyl32.png

Je laisse les amateurs de combinatoires faire des calculs ou des approximations, mais on peut montrer que

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethyl33.png

La série \sum u_n est alors convergente, et

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ethyl34.png

qui vaut un peu plus de 1/3. Autrement dit, on a seulement 1 chance sur 3 de revenir à l’origine et de retrouver ses clés…
Bref, en dimension 1 ou 2, on repasse presque sûrement une infinité de fois à l’origine, mais en dimension 3, on a 2 chances sur 3 de ne jamais y repasser…. J’ai toujours trouvé ça étonnant….

Histoire éthylique, et arrêt optimal

suite aux pressions générales, je vais reprendre mes discussions d’alcoolique…. ou plutôt reprendre des classiques de finance, en expliquant que ce sont simplement des problèmes que se posent les amateurs de boissons fortes (de là à conseiller plutôt de recruter dans les bars qu’à la sortie des grandes écoles…).
Bref, avant d’avoir entamé sa marche aléatoire dans la rue de la soif (ici, correspondant aux problèmes d’options à barrière traduit en termes financier), puis d’avoir un soucis avec ses clés (), puis la maréchaussée (ici), notre héros (car on peut maintenant l’appeler un héros après 4 billets qui lui sont consacrés) avait du choisir son bar… Le problème est loin d’être simple. Il y a 20 bars dans la rue (disons https://perso.univ-rennes1.fr/arthur.charpentier/latex/opt-am-50.png pour faire quelque chose de plus formel). Il arrive de la place sainte Anne, et là, il souhaite choisir le bar le moins cher. Le soucis est qu’il n’a pas le droit de faire demi-tour1 et il ne connaît pas les prix pratiqués dans les différents bars. Il part avec un a priori qui est que le prix d’une pinte est compris entre 3 et 6 euros, que le prix est uniformément réparti entre ces deux prix, et que les prix sont indépendants d’un bar à l’autre. Pour les financiers, il a une option (de commander une bière), et peut l’exercer quand il le souhaite. Une option américaine en quelque sorte. Supposons qu’on soit arrivé au bar https://perso.univ-rennes1.fr/arthur.charpentier/latex/opt-am-01.png. On peut soit payer https://perso.univ-rennes1.fr/arthur.charpentier/latex/opt-am-02.png (qui est supposé aléatoire, uniformément distribué et indépendant des autres bars), soit espérer que l’on puisse payer moins cher plus loin,
Soit https://perso.univ-rennes1.fr/arthur.charpentier/latex/opt-am-03.png la valeur de cette option, alors

https://perso.univ-rennes1.fr/arthur.charpentier/latex/opt-am-04.png

i.e.

https://perso.univ-rennes1.fr/arthur.charpentier/latex/opt-am-05.png

soit

https://perso.univ-rennes1.fr/arthur.charpentier/latex/opt-am-06.png

où https://perso.univ-rennes1.fr/arthur.charpentier/latex/opt-am-07.png désigne la loi du prix de la bière (soit ici une loi uniforme) avec une condition terminale de la forme

https://perso.univ-rennes1.fr/arthur.charpentier/latex/opt-am-09.png

car il a soif, et ne quittera pas la rue sans avoir bu un verre !
Classiquement, par backward induction, on peut résoudre ce programme, à partir de la loi de https://perso.univ-rennes1.fr/arthur.charpentier/latex/opt-am-10.png. Posons https://perso.univ-rennes1.fr/arthur.charpentier/latex/opt-am-12.png. Alors

https://perso.univ-rennes1.fr/arthur.charpentier/latex/opt-am-13.png

et

https://perso.univ-rennes1.fr/arthur.charpentier/latex/opt-am-14.png

soit simplement

https://perso.univ-rennes1.fr/arthur.charpentier/latex/opt-am-16.png

soit enfin

https://perso.univ-rennes1.fr/arthur.charpentier/latex/opt-am-18.png

Je laisse les plus courageux simplifier les calculs. La “frontière d’exercice” est alors obtenue par récurrence. Numériquement, le code est alors simplement

> n=20
> u=rep(NA,n)
> b=6;a=3
> u[n]=(b+a)/2
> for(k in (n-1):1){
+ u[k]=1/(b-a)*(u[k+1]*(b-u[k+1])+(u[k+1]^2-a^2)/2)
+ }

Dès qu’on atteint la barrière, on s’assoit au bar. On note que plus on avance dans la rue, moins on est exigent: au tout début, on ne s’assoit pas à moins de 3 euros 30… mais plus on avance, plus on relève le seuil d’exigence. Le calcul sous forme intégrale donne ici

> u=rep(NA,n)
> b=6;a=3
> u[n]=(b+a)/2
> for(k in (n-1):1){
+ g=function(x){pmin(x,rep(u[k+1],length(x)))/(b-a)}
+ u[k]=integrate(g,lower=a,upper=b)$value
+ }

J’avais déjà abordé ce problème dans un précédant billet, sur les options américaines, mais on peut maintenant aller un peu plus loin… que se passe-t-il si on suppose que les prix sont discret (par exemple par tranches de 50 centimes ou 1 euro) ? L’avantage avec ces méthodes numériques est que l’on peut très facilement enlever des hypothèses, par exemple ici on aurait

> h=2
> K=(b-a)*h+1
> PRIX=seq(a,b,by=1/h)
> u2=rep(NA,n)
> b=6;a=3
> u2[n]=(b+a)/2
> for(k in (n-1):1){
+ g=function(x){pmin(x,rep(u[k+1],length(x)))}
+ u2[k]=sum(g(PRIX)*1/K)}

pour des seuils à 1 euros (les seuls prix possibles étant 3,4,5 ou 6 euros).

Ou la frontière suivante si les prix varient par tranche de 50 centimes.

Compte tenu de la discrétisation, notons que la vraie frontière devient alors ici

Bref, comme toujours, les problèmes d’alcooliques rejoignent les problèmes d’exercice optimal d’options américaines, problème classique en finance de marché…
1 pour rendre cette histoire crédible, à chaque bar rencontré il demande le prix d’une pinte. S’il estime que c’est trop cher, il s’exclame “mais c’est bien trop cher ici !” et s’en va. Sinon il commande et s’installe. Cette exclamation rend improbable – à ses yeux – l’idée de revenir finalement s’installer au bar….

Histoire éthylique, et loi de l’arcsinus

Avant de revenir sur la marche aléatoire et le problème de la ruine d’un assureur comme problème de temps d’arrêt (je n’ai pas insisté en cours car les temps d’arrêt nécessitent un peu de recul sur les processus), je voulais prendre 5 minutes pour faire un exemple qui devrait intéresser les étudiants. On va se placer un peu avant les billets sur l’étude de l’ébriété rennaise, ici et .

  • De la probabilité de sortir un jour de la rue de la soif

Supposons que notre homme soit encore rue de la soif, un jeudi soit à 23 heures, dans un état d’ébriété assez avancé. Souhaitant quitter cette rue, il fait un pas (disons 1 mètre pour faire simple), puis hésite sur sa direction: il a 1 chance sur 2 de continuer dans la même direction, et 1 chance sur 2 de faire demi-tour.  Il se trouve à une distance L'image “https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-05.png” ne peut être affichée car elle contient des erreurs. de la place Sainte Anne et à une distance L'image “https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-06.png” ne peut être affichée car elle contient des erreurs.de la place des Lices. On se demande s’il arrivera un jour à quitter la rue de la soif, et surtout au bout de combien de temps (ou après avoir parcouru quel distance)… Au delà du temps moyen passé à errer dans la rue, on pourra se demander la probabilité pour qu’il soit sorti de la rue pour 8 heures, histoire d’avoir une chance d’arriver à l’heure en cours.

Il s’agit du problème de la marche aléatoire (qui, ma foi, est assez approprié au problème de l’ébriété) de type pile ou face. On note tau l’instant (exprimé en nombre de pas) de la sortie de la rue de la soif. Formellement, on note https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-01.png la position à la date https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-02.png. On suppose qu’à la date https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-02.png=0, il est à la position 0 (disons sa chaise à la terrasse du 1929). tau est l’instant où https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-01.png quitte l’intervalle https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-03.png. S’il sort par la partie inférieure (en https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-05.png) il sera place Sainte Anne, alors que s’il sort en https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-06.png, il sort par la place de Lices (la longueur de la rue étant https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-04.png).
Notons que https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-07.png où https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-08.png prend les valeurs https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-09.png avec probabilités https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-10.pnghttps://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-50.png s’il se dirige vers Sainte Anne, https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-11.png s’il part place des Lices.
On cherche

https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-12.png

C’est ce qu’on appelle un temps d’arrêt. Il s’agit d’une variable aléatoire, positive, qui prend valeurs dans https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-555.png.
Le problème est simple à formuler, mais pour le résoudre, c’est une autre histoire…
On peut d’abord montrer que

https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-13.png

ce qui est un peu rassurant (c’est le lemme de Fatou). On peut alors montrer que https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-14.png appartient forcément à https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-15.png, et donc

https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-16.png

aussi,

https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-17.png

ce qui est assez intuitif compte tenu de la symétrie naturel du problème. On a alors d’autant plus de sortir d’un coté que l’on part d’un endroit proche de ce côté.
On peut aussi monter que https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-18.pngest une martingale, nulle à la date initiale. Aussi,

https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-19.png

(je passe les difficultés techniques, mais on peut montrer que ce résultat est rigoureusement juste en travaillant sur la martingale arrêtée). Bon, on a pas mal d’information…. mais on peut aussi regarder la loi. On peut faire quelques simulations, dans un premier temps, afin de visualiser la densité (sous forme d’histogramme), avec ci-dessous la position de notre individu en ordonné, et en abscisse le temps depuis qu’il a quitté sa table.

https://perso.univ-rennes1.fr/arthur.charpentier/image-st-michel-v2.gif

L’histogramme ci-dessous
est le comptage des moments où l’individu a quitté la rue.

https://perso.univ-rennes1.fr/arthur.charpentier/compte-image-st-michel-v2.gif

Je me suis arrêté assez tôt sinon le fichier image serait trop lourd.
Pour étudier la distribution, mieux vaut éviter de travailler sur la fonction de répartition ou la densité du temps d’arrêt. En fait, le plus simple est de récupérer la fonction génératrice de tau. Pour cela, on va utiliser la fonction génératrice des incréments,
https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-20.png

en utilisant la notation du cosinus hyperbolique. On peut alors montrer

https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-21.png

Pour tout https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-22.png,

https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-23.png

En fait, on retrouve ici quelque chose qui est très proche de la probabilité qu’a le mouvement brownien de quitter un corridor (qui ont été beaucoup étudiée de part les applications en finance, non pas sur l’étude des gens bourrés, via les options parisiennes, ici ou , mais qui avait été étudié dans les revues de physique dans les années 60). Bon, comme on quitte le monde discret pour le monde continu, on ne travaille plus sur la fonction génératrice des moments https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-24.png, mais sur la transformée de Laplace https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-25.png dont on peut montrer que la transformée de Laplace du temps de sortie d’un corridor de la forme https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-521.pngest

https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-30.png

qui se rapproche de la formule précédente. Pour faire le lien, notons que

https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-31.png

et que pour https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-32.png,

https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-33.png

Bref, en bricolant, on peut relier assez clairement les deux notions (ce qui n’a rien d’une démonstration, j’en conviens). Mais on peut aussi se poser d’autres questions fondamentales……

  • De la probabilité de récupérer ses clés

J’avais évoqué ici les problèmes de clés de notre indigène. Mais j’avais oublié de dire qu’il avait oublié ses clés sur sa table en quittant le bar. On peut se demander quelle est la chance pour qu’il ne repasse jamais à sa table…. Autrement dit, on souhaite calculer

https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-522.png

ou éventuellement, si on sait que le bar ferme à l’heure https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-523.png,

https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-524.png

Dans le cas du mouvement Brownien, on peut montrer que

https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-525.png

On utilise la propriété de Markov, qui garantie que https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-526.png est indépendant de https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-527.png. On peut alors montrer que si https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-528.png est le temps d’arrêt défini comme le premier instant où le processus atteint https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-529.png, alors

https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-530.png

Dans le cas du mouvement Brownien, cette propriété est une des lois de l’arcsinus, dont on trouvera des démonstrations ici ou . L’arcsinus apparaît ici dans un calcul d’intégral, en notant que

https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcool-531.png

Compte tenu des similitudes que nous avions observées dans la section précédente entre le cas continu et le cas discret, notre étudiant devrait pouvoir calculer assez simplement la probabilité qu’il a de récupérer ses clés…. enfin, dès qu’il aura dessaoulé un peu afin de faire les calculs….

Histoire éthylique, et problèmes de clés

Pour poursuivre le billet commencé ici, avant d’expliquer comment notre personne imbibée d’alcool est allée du bar à son appartement, je voulais raconter une petite mésaventure qu’il a eu.

Lors de ses multiples essais de clés, des agents de la maréchaussée ont trouvé étrange qu’il ait besoin d’autant d’essais pour trouver la bonne clé.

On part d’une distribution a priori d’une chance sur deux,

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ivre-01.png

On sait que

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ivre-02.png

alors que

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ivre-03.png

et 0 sinon.
On peut alors utiliser la formule de Bayes pour déterminer la probabilité qu’une personne soit ivre sachant qu’il lui a fallu https://perso.univ-rennes1.fr/arthur.charpentier/latex/ivre-04.pngessais,

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ivre-06.png

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ivre-07.png

Je laisse les plus courageux écrire la forme de cette probabilité, mais ça n’est pas sorcier….
Si on regarde ce qui se passe avec https://perso.univ-rennes1.fr/arthur.charpentier/latex/ivre-05.png=15 clés, par exemple, on obtient la probabilité suivante, en fonction du nombre d’essais nécessaire pour ouvrir la porte,

Pour la partie de droite, on a un résultat assez intuitif: s’il faut plus d’essais que le nombre de clés dans le trousseau, c’est qu’on est ivre. En revanche, tant qu’on n’atteint pas la valeur fatidique de https://perso.univ-rennes1.fr/arthur.charpentier/latex/ivre-08.png, alors plus le nombre de clés essayé est grand, plus cela “prouve” que l’on a de forte chances d’être à jeun. On retrouve ça sur le graphique ci-dessous,

Moralité, si on se fait arrêté par la police après avoir eu besoin de 12 essais pour trouver la bonne clé, ce petit laïus permet de prouver que l’on a peu de chances d’être ivre (malgré les apparences, j’en convient).
Mais pour l’instant, nous avons supposé que la maréchaussée intervenait lors de l’ouverture. Que se passe-t-il si elle intervient avantl’ouverture de la porte ? Autrement dit, on cherche ici

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ivre-09.png

dont l’évolution pour notre trousseau de https://perso.univ-rennes1.fr/arthur.charpentier/latex/ivre-05.png=15 clés donne

ou, en fonction de https://perso.univ-rennes1.fr/arthur.charpentier/latex/ivre-05.png et https://perso.univ-rennes1.fr/arthur.charpentier/latex/ivre-04.png,

https://blogperso.univ-rennes1.fr/arthur.charpentier/public/perso/desproges.jpgOn retrouve alors quelque chose qui pourrait paraître plus intuitif… Quoi que. Si au 12ème essai, la porte n’était pas ouverte, la police a tout intérêt à arrêter notre jeune héros, car la probabilité qu’il soit ivre est supérieure à la probabilité a priori. En revanche, s’il réussit à l’ouvrir au 13ème essai, elle a intérêt à ne pas l’arrêter, car il y a de très fortes chances pour qu’il soit à jeun. Étonnant non ?

Histoire éthylique, et de clés

En tant que statisticien, j’aime beaucoup les alcooliques… ils nous permettent de voir certains problèmes de modélisation sous un autre angle.
Considérons une personne ayant un trousseau de 10 clés, ou disons n pour généraliser. Il souhaite ouvrir la porte de chez lui, mais ne sait plus laquelle est la bonne. Combien d’essais – en moyenne – faudra-t-il avant de tomber sur la bonne clé ?
Pour une personne à jeun (et saine d’esprit), le plus simple est de tenter une clé au hasard, puis, si ça ne marche pas, de prendre celle à sa droite (ou à sa gauche, peut importe, mais il faut une stratégie et s’y tenir, à moins de pouvoir noter quelle clé a été utilisée). Bref, au bout de n essais, on est certain de tomber sur la bonne clé. Si le fait de manière un peu formelle (histoire de justifier une intuition), et que l’on note https://perso.univ-rennes1.fr/arthur.charpentier/latex/etyl-01.png l’évènement “la bonne clé est tirée au bout de https://perso.univ-rennes1.fr/arthur.charpentier/latex/etyl-02.png essais“, alors

https://perso.univ-rennes1.fr/arthur.charpentier/latex/etyl-03.png

Raisonnons alors simplement,

https://perso.univ-rennes1.fr/arthur.charpentier/latex/etyl-04.png

puis

https://perso.univ-rennes1.fr/arthur.charpentier/latex/etyl-05.png

puis

https://perso.univ-rennes1.fr/arthur.charpentier/latex/etyl-06.png

etc. Aussi, on peut légitimement penser que

https://perso.univ-rennes1.fr/arthur.charpentier/latex/etyl-07.png

Moralité (il s’agit de la moyenne d’une loi uniforme)

https://perso.univ-rennes1.fr/arthur.charpentier/latex/etyl-08.png

Pour une personne ivre (ou un peu moins à jeun), à chaque essai, elle prend une clé au hasard dans le trousseau (on pourrait appeler ça la stratégie du poisson rouge). Le nombre d’essais peut être potentiellement infini. Et en reprenant les notations précédentes, on note que si

https://perso.univ-rennes1.fr/arthur.charpentier/latex/etyl-09.png

en revanche

https://perso.univ-rennes1.fr/arthur.charpentier/latex/etyl-10.png

puisqu’à l’étape 2, on n’a pas écarté la clé qui ne marche pas… Aussi,

https://perso.univ-rennes1.fr/arthur.charpentier/latex/etyl-11.png

qui correspond à une loi géométrique, ou loi de Pascal, de paramètre https://perso.univ-rennes1.fr/arthur.charpentier/latex/etyl-12.png. Pour ceux qui ne savent plus comment calculer l’espérance d’une telle loi, on peut prendre trois lignes pour refaire le calcul

https://perso.univ-rennes1.fr/arthur.charpentier/latex/etyl-13.png
soit tout simplement

https://perso.univ-rennes1.fr/arthur.charpentier/latex/etyl-14.png

i.e.
https://perso.univ-rennes1.fr/arthur.charpentier/latex/etyl-15.png

Aussi https://perso.univ-rennes1.fr/arthur.charpentier/latex/etyl-16.png.
Grosso modo, il faut deux fois plus de temps pour trouver une clé quand on est rond comme une queue d’pelle, que quand on est à jeun (sans compter le temps où l’on parle à son trousseau de clé, en espérant qu’il trouvera un moyen de désigner la bonne clé). La prochaine partie du billet parlera de “marche aléatoire“, et de difficultés pour quelqu’un d’ivre mort à quitter a rue Saint Michel….

De l’uniformité des tirages du loto

J’ai abordé à plusieurs reprises les tirages du loto (ici ou ), en mentionnant par exemple que Slutsky avait proposé les utiliser comme générateurs de nombres aléatoires. Or je suis tombé un peu par hasard l’autre jour sur un rapport – qui date un peu, certes – expliquant que la loterie britannique était biaisée, la boule 37 tombant (selon ledit rapport) un peu trop souvent.
On retrouve le même genre de remarque sur le loto français, où certains trouvent surprenant que certains numéros soient sortis plus souvent que d’autres…. Certains notent – sans rire – que le fait que la boule 7 soit sorti 722 fois sur 4858 tirages confirme le pouvoir magique de ce chiffre (en moyenne, une boule n’est sortie que 694 fois). Je laisse les amateurs de pratiques ésotériques s’amuser avec la base qui est téléchargeable ici (on pourra noter que les multiples de 5 apparaissent aussi plus souvent que la moyenne).

  • De la normalité des tirages

Mais je m’égare. Le but était de voir si, sur les tirages du loto en France, certaines boules sortaient éventuellement tropsouvent. Si on prend la base et que l’on regarde les statistiques d’ordre, on obtient le graphique suivant

En l’occurrence, la boule qui est le plus sortie est sortie 745 fois (il s’agit du 34), et celle qui est le moins sortie est sortie 652 fois (le 29). Cette répartition me rappelant la  fonction quantile de la loi normale, j’ai décidé de la représenter aussi sur le graphique. Intuitivement, quand on tire des boules dans une urne (avec remise), on s’attend à trouver un loi binomiale, qui s’approche d’une loi normale quand le nombre de tirages est grand. Avec près de 5000 tirages, on pourrait penser que l’approximation asymptotique ne serait pas absurde.
Mais si l’ajustement n’est pas mauvais pour les boules les plus sorties, on voit qu’il est assez mauvais pour les boules les moins sorties….La première explication à laquelle j’ai pensée est que mon approximation est valide si l’on fait des tirages avec remise. Or formellement, entre chaque tirage, on remet bien toutes les boules dans l’urne, mais à l’intérieur d’un tirage, les tirages se font sans remise (ce qui fait qu’une boule ne peut sortir qu’un fois par tirage). Bref, ceci pourrait expliquer la non-normalité…

  • Le latin cube

Ceci m’a fait pensé à l’idée du latin hypercube, qui est souvent utilisé pour rendre “plus uniforme” les tirages avec des lois de Monte Carlo.
Les figures ci-dessous montrent des simulations de Monte Carlo classiques à gauche, et avec une partition avec 10 tirages uniformes dans 100 carrés, à droite.

C’est l’idée du latin hypercube… Les figures ci-dessous correspondent à des estimations nonparamétriques des densités par des noyaux betas, avec le même paramètre de lissage (pouvant être interprété comme une fenêtre),

Autrement dit, l’utilisation de la méthode latin hyperbube permet d’avoir des choses beaucoup plus uniformes…
C’est finalement la même chose que l’on a ici. Mais ceci fait qu’on devrait être de variance plus faible, mais cela n’explique pas l’asymétrie que l’on observe. Car effectivement, si l’on calcul la skewness (moment centré réduit d’ordre 3), on obtient 0,1928.

  • Faisons un peu de simulations

Afin de voir s’il y avait un soucis, la méthode la plus simple consiste à faire du Monte Carlo. Et afin de valider cette histoire de latin cube, deux techniques ont été considérée ici,

  • tirer 34000 boules (avec remise) parmi {1,2,….,49}, ce qui doit se rapprocher de la normalité,
  • tirer 7 boules (sans remise) parmi {1,2,….,49}, et répéter 4858 fois cette opération, ce qui correspond à un vrai tirage du loto.

En faisant 10000 tirage de ces 4858 tirages du loto, on obtient les résultats suivants, résumés par les valeurs moyenne en trait plein, et les quantiles à 5% et 95% en trait pointillé. Les traits verts sont la simulation de vrais tirages, et les traits mauves correspondent à la version avec remise. On peut noter que la courbe verte présente moins de variance que la courbe mauve, ce qui semble confirmer ce que je racontais dans le paragraphe précédant,

Bref, on note que les boules qui sont peu sorties sont finalement trop sorties par rapport à ce qu’aurait du donné des tirages indépendants. Mais ceci ne prouve rien, si ce n’est que c’est statistiquement surprenant….

"Générer" du hasard, partie 5

La dernière méthode que  je présenterais est est celle des générateurs congruentiels linéaires. On retrouve l’idée de linéarité que l’on avait pas exemple dans le premier générateur proposé (ici), la notion de congruence est liée à des problèmes d’arithmétique et d’écriture de nombres dans les bases. On retiendra une fonction de récurrence de la forme suivante

X_{n+1} = (a cdot X_n + c) mod m

Pour que l’algorithme marche bien, on essaye de partir d’une valeur initiale correspondant à un nombre premier.
Prenons un premier exemple simple (inspiré par wikipedia). Si on prend des valeurs un peu au hasard, par exemple on obtient

  • avec X0 = 10, la suite : 10, 10, 10, 10, 10, …
  • avec X0 = 11, la suite : 11, 35, 123, 19, 235, 3, 91, 243, 203, 227, 59, 211, 171, 195, 27, 179, 139, 163, 251, 147, 107, 131, …
  • avec X0 = 12, la suite : 12, 60, 236, 28, 204, 252, 172, 220, 140, 188, 108, 156, 76, 124, 44, 92, 12, 60, 236, 28, 204, 252, 172, …

Bref, on notera déjà qu’il est effectivement intéressant de partir d’un nombre premier, car on boucle alors beaucoup moins vite. Ah oui, et pour se ramener sur l’intervalle [0,1], on divise par m. On a alors intérêt à avoir m assez grand.
Prenons un second exemple simple (inspiré du générateur proposé par Robert Sedgewick). On prend , et on part de 0.

[1]         0        1 31415822 40519862 31536702 84962343
[7]  59428603 54128063 38284723  4802584 19281465 53057765
[13] 47900065 67928365 55662665 20022965 84329266 25717386
[19] 95163907 67972648 42464008 74270568 69856328 96165288
[25] 74221448 24728808  5671369 13328949 75902129 98182909

La série a l’allure suivante,
xxxxEt si on regarde l’histogramme et la fonction d’autocorrélation sur les 10 000 premières valeurs, on a gardé l’uniformité (ou tout du moins on ne doit pas en être trop loin), et on n’a plus d’autocorrélation.

Parmi les générateurs ultra-classiques, mais mauvais, on retiendra le RANDU, datant des années 60, implémenté alors sur les IBM, X_{n+1} equiv (65539 times X_n) mod 2^{31} Mais cet algorithme présente de très gros défauts numériques. Je n’aborderais pas ici les problèmes de ce générateur (je l’avais fait dans le chapitre sur les simulations dans le tome 2 du livre avec Michel Denuit).
Bref, voilà ce qu’on peut raconter rapidement sur les méthodes de construction de nombres aléatoires, ou comme on dit généralement pseudo-aléatoires. Pour quelques réflexions sur l’utilisation de méthodes déterministes pour générer du hasard, je renvoie au très beau petit livre d’Ivar Ekeland intitulé Chaos. Et histoire d’insister un peu lourdement, ces billets (comme la plupart de ceux présent sur ce blog) n’ont pas grande prétention, si ce n’est présenter des concepts parfois compliqués “avec les mains” (et des dessins). Pour aller plus loin, on trouvera de très bons résumés de ces techniques sur le net, ici ou . Je laisse les plus curieux continuer. Sinon le chapitre 7 du Monte Carlo de Fishman propose une analyse poussée de ces méthodes. Et je voulais aussi noter que ce que j’ai proposé dans ces billets, c’est simplement de tirer des nombres aléatoires uniformément répartis sur [0,1]. Pour générer autre chose que de l’uniformité, je renvoie à la bible de Luc Devroye, non-uniform random variate generation, épuisé chez Springer, mais en ligne ici.

"Générer" du hasard, partie 4

Bon l’idée de récurrence les dernières fois (ici et ) n’était pas forcément mauvaise, mais telle qu’elle était présentée, ça ne marchait pas très bien… Une autre méthode peut être celle proposée par John von Neumann, parfois appelée middle square. On prend un nombre à 4 chiffres, on l’élève au carré, et on prend les 4 chiffres du “milieu“.  Par exemple à partir de 1234, dont le carré est 1522756, puis on prend le carré de 2275, c’est à dire 5175625… etc. On prend alors comme nombre cette suite de nombres divisés par 10000 (pour revenir entre 0 et 1). Sur cette exemple on note que l’on peut assez rapidement boucler, et retomber sur la même suite (ceci reste le problème central de tous les algorithmes récursifs).
XXX

Bref, on n’y est toujours pas… Pourtant, de loi, la suite des nombres peut  sembler aléatoires (en oubliant le fait que l’on s’est restreint à 4 chiffres),

[1]  0.2756 0.5536 0.7296 0.1616 0.1456 0.9936 0.4096 0.7216
[9]  0.7065 0.4225 0.5062 0.3844 0.6336 0.4896 0.7081 0.4056
[17] 0.1136 0.0496 0.6016 0.2256 0.9536 0.5296 0.7616 0.3456
[25] 0.3936 0.2096 0.3216 0.2656 0.4336 0.8008 0.8064 0.8096
[33] 0.5216 0.6656 0.2336 0.6896 0.4816 0.3856 0.8736 0.7696

Mais si on regarde cet exemple de plus prêt, on note que le chiffre 6 apparaît très fréquemment à la fin, par exemple…. Même si cette méthode est souvent présentée (par exemple dans la page de wikipedia consacrée aux méthodes de simulations, ici), elle marche beaucoup moins bien que les méthodes que j’avais présentées auparavant (ici et ) et que l’on cite assez peu.

"Générer" du hasard, partie 3

Pour poursuivre l’idée développée dans le billet précédant (ici), prenons une autre fonction de récurrence.

Toujours dans ces notes de cours (donné il y a un peu plus de 10 ans par Christian Gouriéroux), il y avait cette tent function que l’on retrouve sur le graphique de gauche. Si je fais tourner l’algorithme à partir de π-3, on obtient, avec une simple fonction R

[1]  0.14159265 0.28318531 0.56637061 0.86725877 0.26548246 0.53096491
[7]  0.93807017 0.12385966 0.24771932 0.49543864 0.99087728 0.01824545
[13] 0.03649090 0.07298179 0.14596358 0.29192717 0.58385434 0.83229132
[19] 0.33541736 0.67083471 0.65833057 0.68333886 0.63332228 0.73335543
[25] 0.53328913 0.93342173 0.13315654 0.26631308 0.53262615 0.93474770

Bref, ça marche bien là aussi (à condition là aussi d’avoir une écriture très longue de π). On peut là aussi montrer simplement que l’on génère une suite de nombre uniformément distribués sur [0,1], car cette loi est invariante par la fonction de récurrence (tout comme celle proposée par Whittle). Mais on n’a toujours pas d’indépendance….

“Générer” du hasard, partie 2

Une première piste (qui n’est pas de moi mais qui a été proposé par Peter Whittle en 1963, et que j’ai retrouvé dans des notes prises lors d’un cours de Christian Gouriéroux) est une idée simple de récurrence. On prend un nombre u entre 0 et 1 (si possible irrationnel, sinon ça n’a pas grand intérêt), et on considère la suite. Par exemple si on prend comme valeur initiale π-3, on obtient, avec une simple fonction R

[1]  0.14159265 0.28318531 0.56637061 0.13274123 0.26548246 0.53096491
[7]  0.06192983 0.12385966 0.24771932 0.49543864 0.99087728 0.98175455
[13] 0.96350910 0.92701821 0.85403642 0.70807283 0.41614566 0.83229132
[19] 0.66458264 0.32916529 0.65833057 0.31666114 0.63332228 0.26664457

Bref, ça marche bien (à condition numériquement d’avoir une écriture avec un nombre infini de décimales).  Le plus simple est de passer en base 2 (on peut trouver π en base 2 sur le net, par exemple ici). Bref, on peut obtenir la série suivante,
Plus généralement sur le caractère aléatoire des décimales de π, je renvoie aux publications de David H. Bailey (ici).
Bref, cette méthode pourrait marcher car on génère effectivement des nombres uniformément sur l’intervalle [0,1]… sauf que les tirages ne sont pas du tout indépendants. On observe un beau processus autorégressif AR(1), avec une autocorrélation  en 2-h.
XXBref, en l’état cette méthode n’est pas vraiment géniale pour générer du “hasard“… à suivre donc.

“Générer” du hasard, partie 1

J’ai déjà mis plusieurs billets sur mon blog afin d’évoquer les méthodes de Monte Carlo (ici ou là), car de mon point de vue, ce sont des techniques fondamentales en économétrie par exemple, si on travaille sur des petits échantillons. Mais pour l’instant, j’avais toujours éludé le problème délicat du générateur… J’avais toujours admis que l’on avait à notre disposition un générateur de nombres aléatoires sur [0,1], qui génère à chaque appel des nombres indépendants, supposés uniformément distribués sur [0,1].

Je voulais prendre un peu de temps pour présenter quelques pistes de méthodes simples de méthodes permettant de générer des suites de nombre donnant l’illusion du hasard….
Une méthode plus ancienne (que j’ai connue, et oui) est celle des tables de nombres aléatoires, que l’on retrouvait par exemple dans les anciennes éditions du livre de Philippe Tassi ou de Gilbert Saporta. Je laisse les plus curieux traîner dans les vieux livres pour voir comment lire ce genre de table…. Une autre solution est celle que j’évoquais en introduction de mon (vieux) cours de séries temporelles à Dauphine (ici), à savoir l’idée de Slutsky d’utiliser les tirages du loto pour modéliser des prix (et donc inventer ainsi les processus MA, moving average), dans The Summation of Random Causes as the Source of Cyclical Processes. On peut trouver sur  le net les tirages du loto en France depuis des dizaines d’années (ici par exemple) et je suppose qu’en cherchant un peu, on peut trouver des choses similaires dans beaucoup d’autre pays.

Pour ma part, pour revenir sur une loi uniforme sur [0,1], je vais dans un premier temps diviser par 50.

Dans un second temps, on va concaténer deux boules consécutives (en base 50), puis se ramener sur [0,1].
Aussi, les boules 13 et 37 donnent comme nombre 0.2548.

On a ainsi un peu plus de possibilités sur les valeurs atteintes. Graphiquement, on obtient les résultats suivants,
XXXBref, on a à notre disposition des gens qui ont tiré des boules dans une urne. Alors rigoureusement, il faut faire attention car les tirages sont fait sans remise pour les 7 premières boules, puis on remet les boules pour le tirage suivant… Bref, les tirages ne sont pas rigoureusement indépendant les uns des autres (on ne peut pas avoir 3 fois consécutives le nombre “13”, alors qu’en théorie cela devrait être possible).
Bref, on est un peu limité avec ce genre d’outils, car ce qui m’intéresse c’est que mon ordinateur me génère autant de nombres aléatoires que je le souhaite.

An Open Lab-Notebook Experiment