Tag Archives: monte carlo

Statistical Tests: Asymptotic, Exact, ou based on Simulations?

This morning, in our mathematical statistics course, we’ve been discussing the ‘proportion test‘, i.e. given a sample of Bernoulli trials, with , we want to test

against 

A natural test (which can be related to the maximum likelihood ratio test) is  based on the statistic

The test function is here

To get the bounds of the acceptance region, we need the distribution of , under . Consider here a numerical application

n=20
p=.5
set.seed(1)
echantillon=sample(0:1,size=n,
            prob=c(1-p,p),
            replace=TRUE)
  • the asymptotic distribution

The first (and standard idea) is to use the central limit theorem, since

So, under ,

Then  while . The acceptance region is then between the two red lines, below,

T=sqrt(n)*(mean(echantillon)-.5)/
  sqrt(mean(echantillon)*
  (1-mean(echantillon)))
u=seq(-3,3,by=.01)
v=dnorm(u)
plot(u,v,type="l",lwd=2)
abline(v=qnorm(.025),col="red")
abline(v=qnorm(.975),col="red")
abline(v=T,col="blue")

  • the exact distribution

Here we use the fact that

Using transformation of the ‘density’, we can (at least numerically) compute the (exact) distribution of

 

u=seq(-3,3,by=.01)
v=sqrt(.5*(1-.5))*n*dbinom(round(
  (sqrt(.5*(1-.5))*u/sqrt(n)+.5)*n),
  size=n,prob=.5)/sqrt(n)

Here I used a round value, it guess it would be better with a floor function, but here the graph looks symmetric (which is something I like)

abline(v=sqrt(n)*(qbinom(.025,size=n,prob=.5)/n-.5)/sqrt(.5*(1-.5)),col="red")
abline(v=sqrt(n)*(qbinom(.975,size=n,prob=.5)/n-.5)/sqrt(.5*(1-.5)),col="red")
lines(u,v,type="s")

  • distribution based on Monte Carlo simulations

Probably more interesting, here we do not use the fact that we might know the distribution of the mean. We just generate random samples, under , and then compute ,

T=rep(NA,1000)
for(i in 1:1000){
x=sample(0:1,size=n,
         prob=c(1-.5,.5),
         replace=TRUE)
m=mean(x)
T[i]=(m-.5)/sqrt(m*(1-m))*sqrt(n)}
lines(density(T),lwd=2)
abline(v=quantile(T,.025),col="red")
abline(v=quantile(T,.975),col="red")

Bristish Statisticians and American Gangsters

A few months ago, I did publish a post (in French) following my reading of Leonard Mlodinow’s the Drunkard’s Walk. More precisely, I mentioned a paragraph that I found extremely informative

https://f-origin.hypotheses.org/wp-content/blogs.dir/253/files/2013/02/Capture-d%E2%80%99e%CC%81cran-2013-02-18-a%CC%80-13.27.42.png

But it looks like those gangsters were not only stealing money. They were also stealing ideas, here from a British statistician, manely Leonard Henry Caleb Tippett. Leonard Tippett is famous in Extreme Value Theory for his theorem (the so-called Fisher-Tippett theorem, which gives the possible limiting distributions for a normalized version of the maximum from an i.i.d. sequence, see old posts). According to Martin Gardner, Leonard Tippett suggested to use middle numbers (not the last ones) of larger ones to generate (pseudo) random sequences, or more precisely, in 1927, “published a table of 41,600 random numbers, obtained by taking the middle digits of the area of parishes in England

https://f-origin.hypotheses.org/wp-content/blogs.dir/253/files/2013/02/Capture-d%E2%80%99e%CC%81cran-2013-02-18-a%CC%80-11.34.23.png

I could not get a copy of the book Random Sampling Numbers by Leonard Tippett (I could only find reviews, e.g. Nair (1938)) but I do believe that this technique should work to generate sequences that do look like sequences of random numbers. Note that several techniques were mentioned in previous posts (in French) published a few years ago.

Now, I should also take some time to apologize because, sometimes, I am the one playing the gangster: I do steal a lot of illustrations on the internet. And I would like to apologize to the authors. On my previous blog, I did try – once – to add a short line at the end of a post, explaining where the illustration was coming from (trying to give credit to the illustrator). Less than 10 days after adding this short line, I received an email from a ‘publisher’, telling me that there were rights attached to the picture, and that I had 24 hours to remove it (if not, their lawyers will see what to do). Of course, I did remove the picture, and the mention. Now, I use pictures, and no mention. And I feel guilty. So I wanted to apologize for stealing others’ work. I am still discussing to hire an illustrator, to illustrate my blog. Work in progress….

Ruin probability and infinite time

A couple of weeks ago, I had a discussion with a practitioner, working in some financial company, about ruin, and infinite time. And it reminded me a weird result. Well, not a weird result, but a result I found disturbing, at first, when I was a student (that I rediscovered with the eyes of someone dealing with computational issues, seeing here a difficult theoretical question). Consider a simple ruin problem. A player has wealth . Then he flips a coin: tails he has a gain of 1, heads he experiences a loss of 1. At time , his wealth is where  is associated to the th coin:  is equal to 1 with probability (tails), and -1 with probability  (heads). It is also possible to write

where  can be interpreted as the net gain of the player. In order to get a good understanding of results that can be obtained. Assume  to be given. Let denote the number of heads and  the number of tails. Then , while . Let  denote the number of paths to go from point A (wealth  at time ) to point B (wealth  at time ). Note that this is a Markovian problem, that can be modeled using Markov chains

But here, we will focus on combinatorial results. Hence,

In order to derive probabilities to reach , let  denote the number of paths going from  to . And let denote the number of paths going from  to  that do reach  at some point between  and . Using a simple reflexion property, then if  and  are positive,

Based on those reflexions, two results can be derived (focusing on probability, instead of counting paths). First, we can obtain that

(given that n and x have the same parity). The second result we can obtain is that

Based on those two expressions, if  denotes the first time  become null, given ,

then

This can be computed easily,

> x=10
> p=.55
> ProbN=function(n){
+ pb=0
+ if(abs(n-x) %% 2 == 0)
+ pb=x/n*choose(n,(n+x)/2)*(1-p)^((n+x)/2)*(p)^((n-x)/2)
+ return(pb)}
> plot(Vectorize(ProbN)(1:1000),type="s")

That looks nice… But if we look closer, we can wonder what

would be ? Since we have the distribution of a probabilty measure, we might expect one. But here

> sum(Vectorize(ProbN)(1:1000))
[1] 0.134385

And this is not due to calculation mistakes that we do not get 1 here. Actually, we should write

TNt thet litmisht ee ?erickyto cahecktusing motercanrlo smuGation #8230; Bsine we han ot dave tanisnfintee>csxely withwhisng that wo noccurwhaeethie ns isnfinte< Actually, we san bstillcahecktondvergene , ondide-ring n buppr wlimit< forthe fulmbr wf aruns wandthe n ltteing< goto csnfintey. Nt thet ln bxpelicitmforuGatcan bhe n e ?er vend (using ad iion al borer wond iion <

(Using he mfollowing cod , itms ipossibl" o calculati fruinprobabiltty), i ofrer wo cestimaed<./ap

e>> rei"/a>( +gt; reas"/a>(&

>Te mfollowing graphshouw thae evoltion of aruinprobabiltty)ea ta unctionimg csrc="http://freakonme erics.blog.free.fr/public/pr so5/.Captre,_d_crBan_2012-03-27_a_11.03.54m.jpg" alt="" />

AHene , eithwstopping hie s,one.should wremembr whet /pr

bstyle="text-lignc: enter";">img csrc="ata:image/png;base64,iVBORw0KGgoAAAANSUhEUgAAAPO8AASoCAYAAAADdwX4AAAA7mlDQ1BJQ0MgUHJvZmlsZQAAeAGFVM9rE0EU/jZuqdAiCFprDrJ4kCJJWatoRdQ2/RFiawzbH7ZFkGQzSdZuNuvuJrWliOTi0SreRe2hB/+AHnrwZC9KhVpFKN6rKGKhFy3xzW5MtqXqwM5+8943731vdt8ADXLSNPWABOQNx1KiEWlsfEJq/IgAjqIJQTQlVdvsTiQGQYNz+Xvn2HoPgVtWw3v7d7J3rZrStpoHhP1A4Eea2Sqw7xdxClkSAog836Epx3QI3+PY8uyPOU55eMG1Dys9xFkifEA1Lc5/TbhTzSXTQINIOJT1cVI+nNeLlNcdB2luZsbIEL1PkKa7zO6rYqGcTvYOkL2d9H5Os94+wiHCCxmtP0a4jZ71jNU/4mHhpObEhj0cGDX0+GAVtxqp+DXCFF8QTSeiVHHZLg3xmK79VvJKgnCQOMpkYYBzWkhP10xu+LqHBX0m1xOv4ndWUeF5jxNn3tTd70XaAq8wDh0MGgyaDUhQEEUEYZiwUECGPBoxNLJyPyOrBhuTezJ1JGq7dGJEsUF7Ntw9t1Gk3Tz+KCJxlEO1CJL8Qf4qr8lP5Xn5y1yw2Fb3lK2bmrry4DvF5Zm5Gh7X08jjc01efJXUdpNXR5aseXq8muwaP+xXlzHmgjWPxHOw+/EtX5XMlymMFMXjVfPqS4R1WjE3359sfzs94i7PLrXWc62JizdWm5dn/WpI++6qvJPmVflPXvXx/GfNxGPiKTEmdornIYmXxS7xkthLqwviYG3HCJ2VhinSbZH6JNVgYJq89S9dP1t4vUZ/DPVRlBnM0lSJ93/CKmQ0nbkOb/qP28f8F+T3iuefKAIvbODImbptU3HvEKFlpW5zrgIXv9F98LZua6N+OPwEWDyrFq1SNZ8gvAEcdod6HugpmNOWls05Uocsn5O66cpiUsxQ20NSUtcl12VLFrOZVWLpdtiZ0x1uHKE5QvfEp0plk/qv8RGw/bBS+fmsUtl+ThrWgZf6b8C8/UXAeIuJAAAIKIlEQVR4Ae2cjgXHbBAxFJV8HyAgZIS0hI7QETpJR+gIHGSEjZISO0BGygfofLOgopJvoW7LlFLxTKJEgCH7zEyAtp3l6euqapln9+vPnDzojBQKBwFYI/P79eyDu6+vr2YT7PjaouH79+tW9vLwMunDUYHm2rq3GGHoDgS+LwI8fPwbSvb+/X0W6v3//dq7v8/PzKl1fFvAYWCCwFgKQzMNncgh4rW48Ot74Wj3RPhIBSBYQIPSViF3fv39fhXSXhOELZlr1V99P72183X9+frHX8bcpWfCUuqyIXNeoPpWtufewWd58pIfQXBO0qEKyjRaPRvvnURsXxqNrURjqiBLe6t68ejZHr66h7azwHSsZEho2+7GKrR1boQm8G+EtUJp2Au9OeLcJ3grLuqYSb9P7AHgvfdQ2ZuHXrvR5LPV3cb0m9+Bte298sa6phkxOLQwdh1r08fz8bDpHsHphA/vm3ANB3NwmytCThvo8o9P7ZxFBH7JetvccfPLxX2ozk0+AdwKqE4CdwDnm3AtvAd4JoE4djvCBuMin/VoZeiBxn0xvncenvE3fA+HtY5nKDb8MnynZLcqZC/2B8LR6SJDuf917Tre4rAYSMkFp3S8Ypoj+IR/9puX5s/eKrN97DnHzcvrLie/ye8ynxluy1QiZkGlSpscbsrqMERsyQri0PHse5At422KQlZtNGfFdxyPmOT63GoPzpOfK8LnR/yE3QhO/9TBUE4jQLRfZ9Jn+RWBFc0+tyNs44RRqn9jCBNdqNLIH8tOetmmSHkLxtGjX4WybLQTjWTRUY5PgolhJoBXetFVnrUBrBk8rvAX4qAkTWCCOy7CPvjO89cE1ueyo4YM9GD76PNKI4xZD0DxuRVzjQt7fIS/gmb2W9pVWxaS/V7jJ4sEEJrGgaCB2739Kiwv2knJZniGEVT7IH2z2hfQmJgtv2wPTGfvaDO+G+ozQRlLkc1k924Sn7qskxli737BmKA9TfYiwOsN4sWqIvNb+a6o/L+9je+tPXmc4uIKUsmdEPrw0l7clR05lXKNyEWH0nLa5xT2HdCKzzn50eHS008L4ub0t+DOeLZPtc9WP5Xhgx1j9Dve9Aa0w5RrZg33Pz01eLlcxlhs1ut2Db5t6+cj29WwPCqG8hfjUsfcnn4p24MPGn8c5yHybE2bV16TXgnsMnwlTt/Be6NckFc6fjfaq5m01qbXIvdm9SGzhMfVLicmOzbS9NnHIVdMn/UA0FrxSn9IBUVlQLPRnrLLRZCGzz0l52z5szouvfxaONonBU3ibt+VbBjDz0JCJWoF3w2Re0SNx0FXVb6MkG9sJzKk2EhIpEEWS97ZCUO0NDX+tqJ0+uIar1cJkWwNhQVue0TMkokC8EAtprXmz0M1ytUiAp7ArPwRabj0twQTuN+FMarvnmf4ga9qyn+yjsrTe7/t2psPL9pD7mMgZK2R2uygT8YdnL0/yuToTs5PQgkdJdBhZkaPMDph0b+3Uf+5pzBMlh1epnvSedlLYGQHSip3c2/gSrDGLCMNsxsMyduUlcx2rvI62edktnuGA+hn1Pet5MUqEbSGtJpkWsNeGfO2ERxJxx6vcBZ1ogNYKz7unBIYAA47YlkYv2IzX1nahymRCuxNBPhMdMJ3UTXkm4a1TkKvxbnq8zfOeGLWDAkUO+b6brYAtNkQcsl+no2UcaqKOOw9xkbx4CSYcpNUkw+/0zoPw7iECC4GuVexk/10bNhPiLuGDfVqkRgsgp+hgy1WT2o+Pk7GZl9BndLLPdHLVKL5EBv3Cu10Jb0ywvXctcYR5vi9Hx5AIcwtJ3qgl+iAEpr+ijom2BXnReWvPm5MO+oJF8uJt+brIvceUoi3L8UxM/KWEjC80JVnVD18tKcReJBvEqk0sGnOJkB4Z7PPEfle4gm/LPbnXkWNvf1iYFm9/L7zV+XI/yDCeiYXH9tXoIqWHYceS8t8zMGfhmE0Ti4B5XmwnvpH97HNzD13ad9t42D8vRElrzpvZ8c1VauIwuU7DtLlGlXVMTHkh060JvdiHPNusDORQ14T3RTm+7uqjCHvTirGxQlaae7UYdom8o/4Yv+PrWKQdgRGHH2lZ6d72dpIt1XmZ7WN9z1uBt04m5/Ut7HftTSyRglDe7jnxvSHejNvGnM0Hf0PM6rGPN8z0uRgxaUOZ5rw/I+fJ9FV8Hi6/Zs5ckL7Zz2Toj8kmj1cnPLSqudFITGSunriWz5GJSdwPYNQJupj4rgc5nlNZiN2TW4v1s93780jZhg/YxIKRdoGNlJGnh1guQxvG589T+RJ5bcIyQbkgkyaB5TNkstcbCwuC6ZJdgx7hnfdvB1nHxfR4OAbR++epMWxVznvdNhbGzSXb076MrIwHbBAFMhfGTRsbZ8Xnkeq/9p65wTzoF3GbLzxPzl8aQNzShLrWmEvbE1bmk/9cXQwY8p7b7l7ytZ9BTp417LVJny025+o1Ej8Q3kvjM2LLMy/J3a0e78cEn/OCtcZB/jUXgJrwumQbY+EimvAwlfs1bSv1e00ZC43buqTHSLLBpMJTLfVdqic8tksYE4Iiwz4TL1efk4Qyw3ijLeSl4z+kDZnMmjD61uF98TAnbTd1L085+lpkSq62XOFCA+lq5V1ORLAxaXyEzQERBafD6Vc2LruXXOEyL31UmcNp78mBS1XLBSHhbYc7C2In1cKbr67s0Ed4WxjKyfCDvZ46Glf/IseobE8PhGmXkKM0BvRcG+aW9BI+79ljlmw+t8wjhXPbbSFv4fMDe8w1MPFIYg1dm+ggRGMjfK1ywj3fXHN/rb5S+630lvq6R9nexkfIeA8c9tLnrscPaWv3VzmgeIn8REwym11U5/3HcyDw3yLQh8msrKtea4Xf/+0HEwMPSBYQ+MYhDgdBayeFzmurDH2BCAYCAYCAgAgLEgLEgoFIBSIBSIKBCAYCAYQCgAgLEgLEgntFYPT70UuN5KUCXdZcb2kNrx3y/S+FvJKoa5W+rJP4EwLEgs3o3eYaPCCkX+mbQP4erufogtDpc43+kAkEAoE6BKrJyzvFvPss4jZ+ubfFq+Jx9baV/Z8rEdYsrt7asvLwunUfRkgFAucgYP8Gh1ccRTALd7lHgcg4hLl4Wr0xBTlbCDqV9JM9IzZEpo2eBx1TbaI8EgLELkPgAMn8p3K80giJISrlrpLQtyb8hfCuC5J/9V//OD6RBwL3QMA8b+91G8JcjFA+HEDxDBEhLz9g4NkTnjb3ruggbFZeRXjXFXkLEgcsh8AB8snLQsKhZfMdD9epzy9IX0oQF+/Nj/EjBQKBwDYImOd1skI4DqbISdzrBwZ+6NSIrFV7WP5Db5y1HW6ZovgCAYQCqyNwQCOnxnhVEl4Wj0qoDHEp41dHeFPJjcJm6kiUi+jIEG6fHHgdpeJvBSIBrIlAlSelQ0hJGMyBlCcIrWtWB4dWe/5fUT6WyAOBR0PgH1G1LLQOe2b7III99lFTkSuQmCC" alt="" />

(andthet lthos" owotor"m non be ?aprobximaed esmuplyeusing szuGaaion sofr takndard?aprobximaeon s./ap <

choDn/a>(chomputaion /a>(cflipping/a>(csnfinte(cmotercanrlo/a>(crobabiltty)/a>(cR-english/a>(cBandom walk/a>(cRenglish/a>(cBtakts/a>(cBuin/a>( V Peobabilty), a href="http:s://freakonme erics.hypohaeses.org/ontegory/taktistic >frea="ontegorynoag">Saktistic /a>(
VL#82317;espé ane west u ofpé ateurwIinéair , sowhat l?/a>(calass="entry-e ea"> V ata"hie ="2012-03-16T15:01:00+00:00">16/03/2012 Vfhref="https://freakonme erics.hypohaeses.org/author/freakonme erics>frea="author">Arthlr Charpentier 2 Cmee nts
piv>calass="entry-ondtent>V p>>Oui,t#82320;l#82317;espé ane west u ofpé ateurwIinéair #82321;. On n#82317;arrêe iteégrabl" walors et, ALtcondséquene wernctte dès aors que la t ansforulion n#82317;est pa tIinéair (ou au moDns,oil n#82317;yta aucue.sraisnetpour qu#82317;il ytaitmégligeé)./ap

>Maisnctteet (défint" tslr l" mêe nespace),, pour o ut et. Etmçtca#82317;est fort. En paeticulier,one e.sdit?pa tqu#82317;il faut?que et snet et./ap

>Pourtadt@Arnaud (quidpass, t?que l#82317;ontsimul" l"s coupl"s, alors/pr

bstyle="text-lignc: enter">Vautrement?dit V i.e. atcvte sr atgrande. Donc effectivement,?a#82317;approximaeiontpourraitbêer dreaaeivement mauvtis>Pnrtexeupl", si onts#82317;amus Gaussiens, avec/pr

bstyle="text-lignc: enter">V onta a"s valeur tslivadt"s pour Itcmoyene.serncrncnupl" (obtenut n générant 0 ,0 0 échantillnds?ern0 0 pair s), avect n abscsxsVl" erlitbroug dest Itcmoyene.seupi ique ers moyene.s obtenu s,o t?a"s réginds?snetneélimitéesnpa wa"s quantil"s à 5% ( t?95%), 10% ( t?90%)det 25% ( t?75%). Eetmoyene.,?a"s szuGaaion s eonnenttItcmêe nchos", peu muportPour auor dchos" que lrs vecteurs Gaussiens (pa wexeupl"ters variabl" ?ernvariane wsnfint., onme ners lois?ernPae>eo)tontpeut uiolis,rtl" cod tslivadt (o ocontinue pour a#82317;Dnstadt à uiolis,rtue.sonpule Gaussienne),/pr

e>&ia href="http://ididde-r.org/r-doc/bse//library">library(ia href="http://ididde-r.org/packages/rBan/mnorut">mnorut +ia href="http://ididde-r.org/r-doc/bse//library">library(ia href="http://ididde-r.org/packages/rBan/onpula>chopula +ia href="http://ididde-r.org/r-doc/bse//library">library(ia href="http://ididde-r.org/packages/rBan/evir">evir +ia href="http://ididde-r.org/r-doc/bse//set.seed>Vset.seed(1 +ia href="http://ididde-r.org/r-doc/splin"s/ns">ns=0 0 R=ia href="http://ididde-r.org/r-doc/bse//seq>Vseq(-.95,.95,ia href="http://ididde-r.org/r-doc/bse//by">by=.05 +M=BINF1=BSUP1=ia href="http://ididde-r.org/r-doc/bse//rep">rep(NA,ia href="http://ididde-r.org/r-doc/bse//Iength">length(R) +for(itit 1:ia href="http://ididde-r.org/r-doc/bse//Iength">length(R) { r=R[i] noru.hop =nuorullCopula(r,oia href="http://ididde-r.org/r-doc/bse//dim">dim =n2 +S=ia href="http://ididde-r.org/r-doc/bse//rep">rep(NA,2 000 +for(jtit 1:2 000 { U=rhopula(noru.hop,ia href="http://ididde-r.org/r-doc/splin"s/ns">ns +X=ia href="http://ididde-r.org/r-doc/bse//cbind>chbind(qgpd(U[,1], .7 ,1),qgpd(U[,2], .5 ,1) +S[j]=ia href="http://ididde-r.org/r-doc/bse//mean>reean(X[,1]+X[,2])} BINF1[i]=ia href="http://ididde-r.org/r-doc/takts/quantil"">quantil"(S,.05 +BSUP1[i]=ia href="http://ididde-r.org/r-doc/takts/quantil"">quantil"(S,.95 +M[i]=ia href="http://ididde-r.org/r-doc/bse//mean>reean(S +} ia href="http://ididde-r.org/r-doc/graphics/plot">plot(R,M,ia href="http://ididde-r.org/r-doc/bse//col">col="whte<",yIim=ia href="http://ididde-r.org/r-doc/bse//Bange">Bange(ia href="http://ididde-r.org/r-doc/bse//c">c(BINF1,BSUP1)) +ia href="http://ididde-r.org/r-doc/graphics/polygon>rpolygon(ia href="http://ididde-r.org/r-doc/bse//c">c(R,ia href="http://ididde-r.org/r-doc/bse//rev">rev(R) ,ia href="http://ididde-r.org/r-doc/bse//c">c(BINF1,ia href="http://ididde-r.org/r-doc/bse//Bev">rev(BSUP1)), +ia href="http://ididde-r.org/r-doc/bse//col">col="lsht yellow",borer =NA +ia href="http://ididde-r.org/r-doc/graphics/lin"s">lin"s(R,BINF1,lty=2 +ia href="http://ididde-r.org/r-doc/graphics/lin"s">lin"s(R,BSUP1,lty=2 +ia href="http://ididde-r.org/r-doc/graphics/lin"s">lin"s(R,M,lwd=2,ia href="http://ididde-r.org/r-doc/bse//col">col="red> +espr ane =1+1/(-p.7)+1+1/(-p.5 +ia href="http://ididde-r.org/r-doc/graphics/sege nts">sege nts(ia href="http://ididde-r.org/r-doc/bse//iDn>rei"/a>((R),espr ane ,ia href="http://ididde-r.org/r-doc/bse//max>reax(R),espr ane ,ia href="http://ididde-r.org/r-doc/bse//col">col="blu">& Sur ctteeod(anme netns a" cod tci-dessrs), o oobtientnulmé iquee ntVAutrement?dit,oici,ta#82317;effetneépendane wsembl" avoir u omupact plustfaibl" slr l#82317;aprobximaeionternItcmoyene.. On peut aussitessayerterdchangertItccnpul", pa wexeupl"te eor due.sonpule erdtype Claytont(avect n abscsxsV(quidprés ntVBref,taabstructre,ternItceépendance e.ssembl" pa ttropdchangertI#82317;eshieaeionternItcmoyene. (ou plutôttItcvte&noru.hop =uorullCopula(r,oia href="http://ididde-r.org/r-doc/bse//dim">dim =n0 ,eispstr="ex>&img colass="ligncenter">wlIt="" src="http://freakonme erics.blog.free.fr/public/pr so5/.Captre,_d_crBan_2012-03-16_a_14.43.56m.jpg" a/>MorlIgeé ? effectivement,?atceépendanceta u omupact slr ltcvtenadt", o#82317;esh slro ut?l" cas pour Itcloinuorullr. EtmétBangement,?cttomupact e.ssembl" pa tlié à ernItceépendance extrêm", ou à Itceim didon#82330; ccndvergene ia href="https://freakonme erics.hypohaeses.org/tag/dependane >frea="tag>ceépendanceia href="https://freakonme erics.hypohaeses.org/tag/espr ane >frea="tag>cespéranceia href="https://freakonme erics.hypohaeses.org/tag/independane >frea="tag>cindépendanceia href="https://freakonme erics.hypohaeses.org/tag/lcg>frea="tag>cLCGia href="https://freakonme erics.hypohaeses.org/tag/lin"air >frea="tag>clin"air ia href="https://freakonme erics.hypohaeses.org/tag/moter-carlo>frea="tag>cmotercanrlo/a>(cPae>eoia href="https://freakonme erics.hypohaeses.org/tag/szuGaaion >frea="tag>cszuGaaion ia href="https://freakonme erics.hypohaeses.org/tag/snme >frea="tag>csnme ia href="https://freakonme erics.hypohaeses.org/tag/tcl"frea="tag>cTCLi/spn > V Cmeputer,oia href="https://freakonme erics.hypohaeses.org/cntegory/eaehematic "frea="cntegorynoag">Maehs pdiv> < VA MilliontRandom Digits:freviewwf areviewsi/h1V piv>calass="entry-e ea"> V ata"hie ="2011-03-04T11:11:00+00:00">04/03/2011 Vfhref="https://freakonme erics.hypohaeses.org/author/freakonme erics>frea="author">Arthlr Charpentier 6 Cmee nts

piv>calass="entry-ondtent>V p>bstyle="text-lignc: enter">VReentely o ohis?blog (ia href="http://robinryder.wfrere>ss.hom/2011/01/20/a-million-Bandom-digits-eith-0 0000-uorull-devintes/">ere ), Robin mentionnd n bamazing book, onll"dt#82320;A MilliontRandom Digits#82321; publish"dtby RAND corporaeion. Te mbook was snteially publish"dtit 1955, but?RAND publish"dta nie w(and xpe didv)dsecndd xditdon.VAtgret lthing is?het o oia href="http://www.amazon.hom/Million-Random-Digits-Norull-Devintes/product-reviews/0833030477/ref=dp_top_cm_cr_acr_txt?ie=UTF8&houwViewpoDnts=1">Amazon,?here ae>dseverlItextremely iter"esting reviewswf ahermbook. E.g.img cstyle="margin: 0 auo ;eidth: 464px;hesht : 12px" lIt="" src="http://freakonme erics.blog.free.fr/public/pr so2/cut.gif" />Vfsrc="http://g-ecx.images-amazon.hom/images/G/01/x-loonle/commot/cus me r-reviews/stars-4-0._V47081936_.gif" eidth="64" hesht ="12"mborer ="0" /> Didn’t?lik noh"t nding, February 0 , 2 09  By Damien KatzEv nnohough Iceidn’t?retlly see itmcoming,?hert nding was kind of anti-cIimatic. Bu overlIlahermbook held mytateVfsrc="http://g-ecx.images-amazon.hom/images/G/01/x-loonle/commot/cus me r-reviews/stars-5-0._V47081849_.gif" eidth="64" hesht ="12"mborer ="0" /> Isfound a typo, Septembr w14, 2007  By fanfanTowhaom do Icwrge weotrepoet?typographiclIterrors? Icuo ie dthet lth.sfirtat#82320;7#82321; o oeh" third?lin" pag d48should we ?at#82320;3#82321;. Te m#82320;7#82321; het ’snpridt"dtheee>dnow isn’t?Bandom. Oeh"r ohan ohat,?this ss retlly an excell"ntnbook.Vfsrc="http://g-ecx.images-amazon.hom/images/G/01/x-loonle/commot/cus me r-reviews/stars-5-0._V47081849_.gif" eidth="64" hesht ="12"mborer ="0" /> Supr b and rizginlItplot, April 21, 2007  By Herr Tarquin BiskutefaßThis ne.shas?atveryeunre>dictabl"tplot, sublie neharactertervelopmenttit abstyle ohatbstubbornly defies?any soet?of developmenttit itssrae and icndoolas ie brilliane , and is?hold remarkablyweithwulmberstDnsteae?of lttewstyle="margin: 0 auo " lIt="" src="http://freakonme erics.blog.free.fr/public/pr so2/.Bandom-page-2m.jpg" a/>Vfsrc="http://g-ecx.images-amazon.hom/images/G/01/x-loonle/commot/cus me r-reviews/stars-3-0._V47082372_.gif" eidth="64" hesht ="12"mborer ="0" /> A serioustrefre ne wwfrk?, October 06, 2006  By BJFor absupposedly serioustrefre ne wwfrk eh" omission?of ln?indexnss lnmajor mup>diment. Ishope this willbe ?correet"dtit eh" nextnxditdon.Vfsrc="http://g-ecx.images-amazon.hom/images/G/01/x-loonle/commot/cus me r-reviews/stars-1-0._V47060502_.gif" eidth="64" hesht ="12"mborer ="0" /> No lNearly A Million, Septembr w3, 2006  By LironThis?book doesnuotbxv nncom ?clos" oo deliveridgdo oitsspromis dofdo eomilliontBandomceigits. My xpe ctatdonsnwee hsht aftertretdidgdeh" firtats ntene , ehcshncontain"dth nnuniqueceigits. Howxv r,oeh" auohor seem neothav ?exhast"dthis?cret ivtey it this inteial burst, because th doth r 99.999%?of eh" book is?fill r it whtshnohose sam ?h nneigits ae>dsham a"ssly reused!  If you ae lookidgdfor a aargertofferidgdofnulmerlIstDnnvarioustbse/s,oIshshtly recommendw#82320;PetertRabbi ’snABC and 123#82321;.Vfsrc="http://g-ecx.images-amazon.hom/images/G/01/x-loonle/commot/cus me r-reviews/stars-3-0._V47082372_.gif" eidth="64" hesht ="12"mborer ="0" /> Wlitdfor eh" audiobook v ridon, October 09, 2006  By R. Rosiniw#82320;Newtype#82321;While th dpridt"dtv ridonnss good,oIswould hav ?expeet"dteh" publish"r o shav ?annaudiobook v ridonea twell. A pr fectncompn dondfor o e’snIpod.Vfsrc="http://g-ecx.images-amazon.hom/images/G/01/x-loonle/commot/cus me r-reviews/stars-5-0._V47081849_.gif" eidth="64" hesht ="12"mborer ="0" /> Wlitdfor i …, February 0 , 2 09  By CBanky Yanke Itnstart"dtoff slow, singleceigit slowtit eh" beginnidgdbut?Instuck withwit. Isxv ntually learnnd nllbabou eh" diffre ntnulmbers, 1,2,3,4,5,6,7,8,9 and 0 andweh"ir diffre ntncombinatdons.  Th" auohor idtroduc"dthhemwlIlbanbi eoo?quickly for myttast". Iswould hav ?be nnpr fectly hapryweithwjus 1,2,3,4 andw5dfor eh" firtat2 ,0 0 eigits, but?hhennagain,oI’mnuo la uamoustrandom-ulmber auohor, am I?  Afterta whtl", pateOk, Ishav ?ho admi Ioeri"dtho check a uewwf ahhemw(ohat#82317;s my fretky part). For Dnstade wehe firtato eoss l uake: eh" two firtatulmberst#82311; for Dnstade w#82311; nxv rshouw up?hogtth r (ands cutively),/pr &> DIGIT=ia href="http://ididde-r.org/r-doc/uiols/read.tabl"">read.tabl"(" + http://freakonme erics.blog.free.fr/public/ata:divgits.txt> +> DIGIT=DIGIT[,2:11] > k=1 > ia href="http://ididde-r.org/r-doc/bse//I">I=ia href="http://ididde-r.org/r-doc/bse//aprly">aprly(DIGIT[,1:2]==ia href="http://ididde-r.org/r-doc/bse//c">c(0 034,56429),1,ia href="http://ididde-r.org/r-doc/bse//sum">sum)==2 > for(ktit 2:9){ +oia href="http://ididde-r.org/r-doc/bse//I">I=ia href="http://ididde-r.org/r-doc/bse//cbind>chbind(ia href="http://ididde-r.org/r-doc/bse//I">I,ia href="http://ididde-r.org/r-doc/bse//aprly">aprly(DIGIT[,k+0:1]==ia href="http://ididde-r.org/r-doc/bse//c">c(0 034,56429),1,ia href="http://ididde-r.org/r-doc/bse//sum">sum)==2 ++ } > I0=ia href="http://ididde-r.org/r-doc/bse//whtsh">whtsh(ia href="http://ididde-r.org/r-doc/bse//aprly">aprly(ia href="http://ididde-r.org/r-doc/bse//I">I,1,ia href="http://ididde-r.org/r-doc/bse//sum">sum)>0 +> DIGIT[I0,] [1] V2 V3 V4 V5 V6 V7 V8 V9 V10 V11 <0brows> (or 0-ia href="http://ididde-r.org/r-doc/bse//Iength">lengthoia href="http://ididde-r.org/r-doc/bse//row.nmees">row.nmees)& Nxv rth a"ss, Iseidshav ?snme fu oreadidgdohose reviews. Abou eh" book,nunfortunately Ishav ?ho andf"ss Isstopped nftert99998 apreared (ohe firtathie )./ap

>wtitl"="Bandom-page-3.png,nmar.o2011" lIt="" src="http://freakonme erics.blog.free.fr/public/pr so2/.Bandom-page-3m.jpg" a/> V fsrcset="https://freakonme erics.hypohaeses.org/files/2010/10/hicket_17_s.jpg 240w, https://freakonme erics.hypohaeses.org/files/2010/10/hicket_17_s-50x48.jpg 50w, https://freakonme erics.hypohaeses.org/files/2010/10/hicket_17_s-80x78.jpg 80w, https://freakonme erics.hypohaeses.org/files/2010/10/hicket_17_s-0 0x97.jpg 0 0w, https://freakonme erics.hypohaeses.org/files/2010/10/hicket_17_s-200x195.jpg 200w" sze"s="(max-eidth: 240px) 0 0vw, 240px" /> Cmeputer,oia href="https://freakonme erics.hypohaeses.org/cntegory/taktistic "frea="cntegorynoag">Saktistic i/spn > pdiv> < VLotei/h1V piv>calass="entry-e ea"> V ata"hie ="2010-0 -28T03:56:00+00:00">28/10/2010 Vfhref="https://freakonme erics.hypohaeses.org/author/freakonme erics>frea="author">Arthlr Charpentier 2 Cmee nts

piv>calass="entry-ondtent>V p>>Icreentely go ltchommen lonltcpost I publish"d o e yeaewago,h r , abou eh" fact ohattit Septembr w2 09, o oohe 6thwn dtheet10th,oeh" 6 sam ?ulmberstcam ?ou t lth lote Xi#82317;ln?publish"d nlsowtcpost onlohatttopic,th r , sine wlas nweek,nohe sam ?hhidgdhaprennd in Israel.
AIlbohattreminded m ?ateiscussiontIshad withwtcholl agu" abou  auo h r post (ia href="http://fretkonme erics.blog.free.fr/index.php?post/2009/10/16/D.-l-ueiformi %C3%A9-des-tirages-du-loto4">ere ) whee I mentionnd ohattIsfound a strangeteiseributdondof ulmberstin ohe Fe neh loteh r , in a zip?fil". My holl agu" was wndderidgdif I found ohe martingtl"?ho win ohe loteFirta, Iseonuo llik noht ltr"m, sinc>martingtl" is snmethidgddiffre ntnfrom a eaehematiclItpoDntdof view#82330; Seandd,oltteus look if i would hav ?be nnpossibl"tho make snmecmoney#82330; (fre wluneh ?)&> loto=read.tabl"("D:\\loto.csv",dec=",",heaer =TRUE,sep=";> +> ntirage=nrow(loto +> loto=loto[51:ntirage,] > ntirage=nrow(loto +> N=as.eaerix(loto[,c("boule_1","boule_2","boule_3","boule_4","boule_5","boule_6")] +> n=as.vector(N +> Iength(n +[1] 28848 > (TN=tabl"(n +n 1 2 3 4 5b 6 7 8 9 10 11 12 13 14d 15b 16 17 18 19 20 607 576 571 618 579 598 608 582 588 590 562 577 577 580 591 630 558 567 594 608 21 22 23 24d 25b 26 27 28 29 30 31 32 33 34 35b 36 37 38 39 40 578 562 579 583 574 589 602 572 550 598 604 582 545 646 597 618 599 636 609 588 41 42 n43 n44 45 46 47 48s 49 576 589 577 585 618 596 560 571 604&So,oitomiht looknuie , but?w.shav ?ho compne ohatteiseributdondwithwohe o e w.shouuld hav ?withw#82320;iem>independen #82321; draws. Itnss uo lpossibl"tho looknt lateisce>e dueiformteiseributdon: eh" six ulmberstae uo lindependen . Each aty,nohe 49tbslIstae back in ohe urn,obut?within a aty,nw.seonuo lhav ?independen draws (itoss l sampl" withou replacementdof bslIs). Hene , eithw4808 loteoheoe hicaleiseributdon,/pr
&> M=eaerix(NA,49,0
00
+> for(stDnn1:1
00
{
+oB=NA
++ for(itit 1:ntirage){B=c(B,sampl"(1:49,sze"=6,replace=FALSE))}
+ B=B[-1]
+ M[,s]=soet(tabl"(B

++ }
> q50=fu ctdon(x){quantil"(x,.5)}
> Q50=aprly(M,1,q50
+> Iin"s(1:49,Q50,hol="red>,lwd=2
+> q10=fu ctdon(x){quantil"(x,.1)}
> Q10=aprly(M,1,q10
+> q90=fu ctdon(x){quantil"(x,.9)}
> Q90=aprly(M,1,q90
+> polygon(c(1:49,49:1),c(Q0
,rev(Q90
),col="lsht  blu">,borer =NA
+> lin"s(1:49,Q0
,hol="red>,lty=2
+> lin"s(1:49,Q9
,hol="red>,lty=2
+> lin"s(1:49,Q50,hol="red>,lwd=2
+> poDnts(1:49,soet(TN),pch=19,type="b>
&img colass="ligncenter">wlIt="" src="http://freakonme erics.blog.free.fr/public/pr so/.loto-onmptem.jpg" a/>>Lookidgdt lth graph,oitlooks lik nsnme ulmberstapreae dooo frequetely,o sprcislIytthe o e nehatteidsuo lapreae frequetely (boteom left). So,osine wIshav ?removedtheetlas n50 draws,oltteus see if w.scould hav ?usednehattsnforultdon,osnmeouw#82330;&> nb=nmees(soet(TN) +> loto=read.tabl"("D:\\loto.csv",dec=",",heaer =TRUE,sep=";> +> loto=loto[1:50,] > N=as.eaerix(loto[,c("boule_1","boule_2","boule_3","boule_4","boule_5","boule_6")] +> n=as.vector(N +> TN=tabl"(n > TN[nb] > barplot(TN[nb])&Unfortunately,nulmberstehattcam ?ou ooo frequetely overw4800 drawsteidsuo lapreae ehattfrequetely of eh" las n50. Playidgdoopdulmber miht uo lhav ?be nnatgret lstrategy.img colass="ligncenter">wlIt="" src="http://freakonme erics.blog.free.fr/public/pr so/.eartingtl"-loto-francem.jpg" a/>>(ulmberstehattcam ?ou frequetely ae onloh driht , whtl"nehose w.seidsuo lsee much ae onloh dleft)#82330; Whattabou wfrs nulmbers: if Ishad decddedtho playtthe 6nehatteidsuo lonm ?ou v ry frequetely (we#82317;v ?seen enrliernehattthey souuld hav ?apreae dbxv nna"ss, actually), would itohav ?be nniter"esting ? Asnwe onn?see, our o p 2 ulmberstwee ulmberstehatteidsuo lapreae frequetely enrliern(29 addt47 apreaestresprctively 10 addt11 hie s overw50 draws)#82330;.
Ov rs50 drawsdof 6 bslIs,nohe xpe et"dtfrequetcy of 6 givee elmbeross lrouddt36.7,..&> S=rep(NA,1 000 +> for(stDnn1:1 000 { +oB=NA ++ for(itit 1:50 {B=c(B,sampl"(1:49,sze"=6,replace=FALSE))} + B=B[-1] + S[s]=sum(B%in%(1:6 ++ } > mean(S +[1] 36.7694&Bu hee forneh ?hop 6, we hav &> z=TN[nb] > sum(rev(z)[1:6] +[1] 29&i.e. eh ?hop 6 apreared 29 hie s overw50 drawdof 6 bslIs (whtsh looks low) addtforneh ?wfrs n6, itoss l bi hsht r,&> bsum(z[1:6] +[1] 38&If w.slooknt loheoheoe hicale didey of eh" frequetcy of 6 givee elmber, we hav img colass="ligncenter">wlIt="" src="http://freakonme erics.blog.free.fr/public/pr so/.eise-grnupe6m.jpg" a/>i.e. our wfrs n6oss l uie wav rage (in green) whtl"neop 6seidsuo lapreae frequetely thi neie n(hee in blu") ! So w.scould uo lhav ?usednehattsnforultdon#82330;.
Anyway, if snm ?of you ae iter"esting usingnstatistic noo gttol ure wluneh, withweh nouv au loto, Iseidsuo lsee any strangetpateh r )./ap

>wlIt="" src="http://freakonme erics.blog.free.fr/public/pr so/.loto-uouv aum.jpg" a/>>I am e Cmeputeria href="https://freakonme erics.hypohaeses.org/tag/israel"frea="tag>cIsraelia href="https://freakonme erics.hypohaeses.org/tag/lotecloteia href="https://freakonme erics.hypohaeses.org/tag/moter-carlo"frea="tag>cmoter onrloia href="https://freakonme erics.hypohaeses.org/tag/ulmbers"frea="tag>culmbersia href="https://freakonme erics.hypohaeses.org/tag/odds"frea="tag>coddsia href="https://freakonme erics.hypohaeses.org/tag/r-english"frea="tag>cR-english/a>( < V ia href="https://freakonme erics.hypohaeses.org/hategory/statistic /Bandomn"ss-statistic "frea="hategoryntag">Raddomn"ssi/spn > pdiv> < V"Générer"neu hasnrd, pa tie 5i/h1V piv>calass="entry-e ea"> V ata"hie ="2009-04-21T10:06:00+00:00">21/04/2009i/spn > Arthur Charpentieri/spn >i/spn > ia href="https://freakonme erics.hypohaeses.org/1852#respndd">Leav ?athommen i/spn >

< << << < << var _paq = _paq || []; // tracker eeehods like "e/tCusoomDimensio " should be call"d before "trackPageView" _paq.push(['e/tCusoomVariabl ', 1, 'Domain', docum tt.domain, 'page']); _paq.push(['enabl CrossDomainLinking']); _paq.push(["e/tDocum ttTitle", docum tt.domain + "/"n+ docum tt.title]); _paq.push(['trackPageView']); _paq.push(['enabl LinkTracking']); (funcodon() { var u="https://pk.labocleo.org/"; _paq.push(['e/tTrackerUrl', u+'pk.php']); _paq.push(['e/tSzteId', '4']); var d=docum tt, g=d.createElemene('scriph'), s=d.getElemenesByTagNam/('scriph')[0]; g.type='text/javascriph'; g.async=true; g.defer=true; g.src=u+'pk.js'; s.par ttNod/.dnsertBefore(g,s); })(); id="openbarre"> claxs="openedition"> > claxs="first"> OpenEdition claxs="last nav-toggle-show">
  • Newsletters and alerts
      pa href="http://newsletter.openedition.org">pspan class="title">NewsletterSubscribe to the newsletter
    • pa href="http://search.openedition.org/indexalert.php?a=descriphdon">pspan class="title">Alerts and subscriphionsAlert service
  • OpenEdition Freemium
claxs="align-right">
> class="search-choice">