Freakonometrics

Depuis le début de la semaine, après avoir déposé les enfants au camp de jour du Musée des Beaux Arts, je viens à pieds à l'université. Chaque fois, je me dis que je pourrais prendre le bus, mais comme aucun ne vient, je commence en marchant, en me disant que si un bus passe, je le prendrais. Et tous les matins, j'arrive au bureau sans m’être fait dépassé par le moindre bus. Et bien sur, j'en ai croisé un paquet qui passaient en sens inverse...

C'est vraiment trop injuste ? En tant que statisticien, je dirais que non, tout simplement à cause du biais de sélection. Si on veut formaliser un peu, on va oublier l'aléa, et raisonner avec des bus qui passent de manière purement déterministe. On va aussi supposer qu'il y en a autant dans un sens que dans l'autre...

J'ai une distance à parcourir, on suppose que les bus sont espacés d'une distance , et qu'ils avancent à une vitesse . Autrement dit, les bus passent tous les secondes (si ma vitesse est exprimée en secondes).

Moi, j'avance à une vitesse avec (oui, on va supposer que je vais moins vite que le bus... ce qui n'est pas forcément une hypothèse faible aux heures de pointes, mais disons que le problème n'a de sens que si aller en bus me permet d'aller plus vite). Le temps que je vais mettre si je fais tous le trajet à pied est . Maintenant, comptons les bus qui passent en face. Je vais croiser tous ceux qui sont déjà sur ma portion de trajet, et il y en a . En plus, je croiserais tous ceux qui vont arriver à l'université, et qui n'y sont pas encore, soit , i.e. le temps qu'il me reste à marcher divisé par le temps qui s'écoule entre deux bus. On a alors un total de bus à croiser, en face.

De mon coté cette fois. Je peux compter tous ceux qui vont arriver à l'université, sauf que cette fois, je devrais enlever tous ceux qui sont déjà sur le trajet, et que je ne croiserais pas (car je vais moins vite qu'eux, par hypothèse). Le nombre de bus qui vont me doubler sera alors cette fois

Si je regarde le ratio du nombre de bus croisés sur le nombre de bus qui m'aura doublé, j'obtiens

Les plus malins auront noté que ce ratio est aussi le ratio entre la vitesse des bus qui me croisent divisé par la vitesse des bus qui me doublent, quand je suis le référentiel (les vitesses sont exprimées par rapport à moi).

Par exemple, si je vais deux fois mois vite que le bus, il y aura trois fois plus de bus qui passent en sens inverse. Et sept fois plus si je vais juste 25% moins vite que le bus. Ce qui n'est pas absurde, quand on pense que je fais le trajet en 20 minutes, et le bus 15...

There was (once again) a nice puzzle in http://www.futilitycloset.com/. Bayes and a good friend are playing Russian roulette. The revolver has six chambers. He puts two bullets in two adjacent chambers, spin the cylinder, hold the gun to his friend's head, and pull the trigger. It clicks. So it is now Bayes's turn: he can choose either to spin the cylinder again or leave it as it is. Which is better? Hopefully, Bayes knows his theorem: if he does spin it, the probability of getting killed is 2 out of 6 (four empty chambers out of six), but if he does not, since his friend is still alive, then the hammer should be next to one of the four cylinders in red, below
So here, there is 3 chance out of 4 to survive, i.e. the probability of getting killed is 1 out of 4 (while it was 1 out of 3 when spinning). So Bayes should not spin. And as always, it is possible to see it is a more general result: more generally, in a revolver with chambers, it there are bullets in adjacent chambers,  if the first player survives, the probability of getting killed is k over , when spinning, while it would be 1 over if we don't. Not spinning is better if and only if

Once again, there was a nice maths puzzle on http://www.futilitycloset.com/ last week (but without further reference). The question was the following, "Five points are located in an equilateral triangle with 10-inch sides (or on its perimeter). What’s the maximum distance between the two closest points?" Actually, this is simply an application of Dirichlet's pigeonholes theorem, as mentioned in the answer of the puzzle, "Connect the midpoints of the triangle’s sides to make four smaller triangles. Because there are five points, two of them must fall within one of these triangles. The maximum distance between these two is 5 inches."

Thus, with Dirichlet's pigeonholes theorem, we know not only the maximal minimum distance, but also where points must be (on corners of inside triangles). Here, there might be two possibilities (with also the different shapes obtained using rotations),

then the minimal distance is . Which can be visualized below,

Based on that pigeonhole theorem, I have the intuition that this result is valid (here we just count the number of inside-triangles), but can we check if this is correct or not ? One idea might be to draw points randomly, and so see where points might end, or at least the maximal distance obtained over millions of random draws... But standard monte carlo might take a while... so we can use two ideas. One idea is from quasi-monte carlo techniques: since we want points to be to be as separate as possible, we do not need to draw randomly in the triangle, but perhaps we can draw randomly points on some grid. A second idea from the latin hypercube technique (and that pigeonholes theorem): instead to generating points randomly in the triangle, perhaps we can draw them in specific regions. For instance, with five points, we know that 4 points have to be in those sub-triangles, and the additional point in any one of those triangle. And because of symmetry, with five points, we can claim that this additional point has to be specifically in one sub-triangle. A random sample with five points within the same sub-triangle will be useless (and a waste of computational time).

With the following code, we define a grid, for a triangle, either upward, or downward, starting from some point (on the left), with a given length, and a given number of subdivision.

TRIANGLES=function(xinf,yinf,l,n,updown="up"){
X=NA;Y=NA
for(i in n:1){
u=xinf+seq(0+(n-i)/2/(n-1),1-(n-i)/2/(n-1),length=i)*l
if(updown=="up") v=rep(yinf+(n-i)*sqrt(3)/2*l/(n-1),i)
if(updown=="down") v=rep(yinf-(n-i)*sqrt(3)/2*l/(n-1),i)
X=c(X,u);Y=c(Y,v)}
return(cbind(X[-1],Y[-1]))}

Here are grid with respectively 20 and 50 points on the lower side. It is then possible to define 4 grids, corresponding to the four sub-triangles,

k=3;st=6
firstgrid=TRIANGLES(0,0,(k-2)/(k-1),k-1)
secondgrid1=TRIANGLES(
firstgrid[1,1],firstgrid[1,2],1*(k-2)/(k-1),st)
secondgrid2=TRIANGLES(
firstgrid[2,1],firstgrid[1,2],1*(k-2)/(k-1),st)
secondgrid3=TRIANGLES(
firstgrid[3,1],firstgrid[3,2],1*(k-2)/(k-1),st)
secondgrid4=TRIANGLES(
firstgrid[3,1],firstgrid[3,2],1*(k-2)/(k-1),st,updown="down")

Then, we just draw randomly five points on that grid, in the four sub-triangles,

N=5
Dmax=0
setpointmax=matrix(0,4,2)
indice=c(1:4,sample(1:4,size=N-4,replace=FALSE))
tindice=table(indice)
indice1=sample(1:nrow(secondgrid1),size=
tindice[1],replace=FALSE)
indice2=sample(1:nrow(secondgrid2),size=
tindice[2],replace=FALSE)
indice3=sample(1:nrow(secondgrid3),size=
tindice[3],replace=FALSE)
indice4=sample(1:nrow(secondgrid4),size=
tindice[4],replace=FALSE)
setpoint=rbind(secondgrid1[indice1,],secondgrid2[indice2,],
secondgrid3[indice3,],secondgrid4[indice4,])

No, we can run a code, where we keep in mind locations of the five points each time we beak a record,

D=min(dist(setpoint,"euclidean"))
if(D>Dmax){Dmax=D
setpointmax=setpoint}

Here are some locations obtained after running the algorithm a few times (with five points)

On the graph below, we can visualize the time it takes before having a record, and the convergence towards 1/2 (which is the true value of the maximal distance)

The convergence is slow... extremely slow... However, we can run the same code for more than five points, e.g. seven points (actually, here sub-triangles are not used here, and it looks like we have been lucky here, since the convergence was rather fast),

Following the two previous posts (here and there), on the number of people that don't get wet while playing with water pistols, consider now an indoor version, in a non-convex room (i.e. player behind wall are now, somehow, protected). In the previous posts, players where playing on a square field, and I briefly mentioned that if the field was a disk, results would have been (roughly) the same: so far, the shape of the field was not an issue. But what if the field is no longer convex,

library(sp)
plot(0:2,0:2,col="white",xlab="",ylab="")
MAP=Polygon(cbind(c(0,0,1,1,2,2,0),
c(0,2,2,1,1,0,0)))
polygon(MAP@coords,col="light blue")

touch=function(x1,y1,x2,y2,n=251){
X=seq(x1,x2,length=n)
Y=seq(y1,y2,length=n)
sum(point.in.polygon(X,Y,MAP@coords[,1],
MAP@coords[,2], mode.checked=FALSE)==0)==0
}
 
NOTWETnc=function(n,p){
sx=runif(50)*2;sy=runif(50)*2
IN=which(point.in.polygon(sx,sy,MAP@coords[,1],
MAP@coords[,2], mode.checked=FALSE)==1)
Sx=sx[IN];Sy=sy[IN]
Sx=Sx[1:n];Sy=Sy[1:n]
IN=IN[1:n]
MI=matrix(NA,n,n)
for(i in 1:(n-1)){
for(j in (i+1):(n)){
MI[j,i]=MI[i,j]=touch(Sx[i],Sy[i],Sx[j],Sy[j])
}}
(d=as.matrix(dist(cbind(Sx,Sy),
method = "euclidean",upper=TRUE)))
diag(d)=999999
dpossible=d
dpossible[MI==FALSE]=999999
dmin=apply(dpossible,2,which.min)
#whonotwet=( (1:n) %notin% names(table(dmin)) )
notwet=n-length(table(dmin))
return(notwet)}
 
NOTWET=function(n){
x=runif(n)
y=runif(n)
(d=as.matrix(dist(cbind(x,y),
method = "euclidean",upper=TRUE)))
diag(d)=999999
dmin=apply(d,2,which.min)
notwet=n-length(table(dmin))
return(notwet)}
 
NSim=10000
Nnc=Vectorize(NOTWETnc)(n=rep(11,NSim))
Nc=Vectorize(NOTWET)(n=rep(11,NSim))
T=table(Nc)
Tn=table(Nnc)
plot(as.numeric(names(Tn)),
Tn/NSim,type="b",col="blue")
lines(as.numeric(names(T)),
T/NSim,type="b",col="red",pch=4)

Je suis tombé un peu par hasard sur un vieux papier de Peter Enis, Enis (1973), qui revient sur une relation que l'on prend toujours pour acquis (moi le premier) à savoir . Mais il faut des hypothèses complémentaires pour que cette relation soit vraie. Le fait que existe par exemple ne suffit pas ! (je laisse les amateurs de résultats contre-intuitifs le temps de réfléchir à cette affirmation)

Ce qui est pénible avec les paradoxes, en tant que scientifique, c'est qu'ils nous obligent à remettre en cause notre intuition, sans cesse. Car mine de rien, on fonctionne tous à l'intuition, c'est un peu comme ça qu'on monte un programme de recherche: on a l'intuition que l'on pourra obtenir des résultats, mais tout reste à prouver.

Prenons la question suivante "les garçons on-t-ils plus de sœurs que les filles ?" en moyenne. N'importe qui avec un peu de bon sens répondrait "oui", l'idée naturelle est que la fille ne peut pas se compter elle même comme sœur, donc les filles devraient en compter un peu moins... Et comme toujours la formulation est importante. On ne demande pas, par exemple "les garçons ont-t-ils plus chance d'avoir au moins une sœur que les filles ?" (auquel cas la réponse serait différente).

Maintenant réfléchissons un peu, et faisons un petit modèle où on suppose que la taille de la fratrie ne dépend pas des sexes de la fratrie. Et comptons, pour mieux voir ce qui se passe. Et commençons par un cas simple, par exemple avec trois enfants,

• (G,G,G) - personne n'a de sœur
• (G,G,F), (G,F,G) ou (F,G,G) - la fille n'a pas de sœur, par contre, il y a deux garçons avec une sœur
• (G,F,F), (F,G,F) ou (F,F,G) - dans chaque cas, il y a deux filles avec une sœur, et un garçon avec deux sœurs,
• (F,F,F) - ce qui fait trois filles avec chacune deux sœurs.

• filles et garçons - chacun des garçons a sœurs, alors que chacune des filles a sœurs,