On s’intéresse au problème suivant.

Monsieur et madame Martin ont 2 enfants dont une fille. Quelle est la probabilité que l’autre soit un garçon ?

L’objectif de ce premier TD est de vous montrer comment simuler ce type de situation afin d’obtenir une évaluation de cette probabilité

Première version

Génération

Prenons un échantillon de \(N=100\) familles de deux enfants et notons \(X_1[N]\) une variable aléatoire valant 0 si le premier enfant de la \(N\)ème famille est un garçon et 1 si c’est une fille. De même, on note \(X_2[N]\) une variable aléatoire valant 0 si le second enfant de la \(N\)ème famille est un garçon et 1 si c’est une fille.

set.seed(42);  # pour fixer la graine du générateur aléatoire et s'assurer qu'on aura bien les mêmes valeurs à chaque fois qu'on génère ce document 
N = 100;
X1=c();        # création d'un tableau vide
X2=c();        # création d'un tableau vide
for(i in (1:N)) {
  if(runif(1)<1/2) {
    X1[i] <-0
  } else {
    X1[i] <-1
  }
  if(runif(1)<1/2) {
    X2[i] <-0
  } else {
    X2[i] <-1    
  }
}

Regardons \(X_1\) et \(X_2\):

X1;
##   [1] 1 0 1 1 1 0 1 0 1 0 1 1 0 0 0 1 0 0 0 1 0 0 0 1 1 0 0 0 1 0 1 1 1 0 1
##  [36] 0 0 0 0 1 1 0 1 0 0 1 0 1 0 1 1 0 1 1 1 1 1 1 0 1 0 1 1 1 0 1 1 1 1 1
##  [71] 0 0 1 0 0 1 1 1 0 0 1 1 1 0 0 1 0 0 1 1 1 0 1 1 0 0 0 0 1 1
X2;
##   [1] 1 1 1 0 1 1 0 1 0 1 0 1 1 1 1 1 1 1 0 1 0 1 1 1 1 0 1 1 0 1 1 1 0 1 0
##  [36] 0 0 1 0 0 0 1 1 0 0 0 1 1 1 1 0 0 1 1 0 1 0 1 0 1 0 1 0 1 1 1 1 1 0 1
##  [71] 0 1 1 1 0 0 0 1 0 1 1 0 0 0 0 0 0 0 1 1 1 0 1 1 0 1 0 0 0 1

Génération à la R

Bon, cette version me donne vraiment des boutons donc voyons comment on écrirait ça en R…

set.seed(42);
N = 100;
X1 = ifelse(runif(N)<1/2,0,1);
X2 = ifelse(runif(N)<1/2,0,1);

Voire en tirant vraiment parti des fonctions du langage:

set.seed(42);
N = 100;
X1 = rbinom(N,1,1/2);
X2 = rbinom(N,1,1/2)

Comptage

Bien, maintenant, comptons le nombre de familles ayant au moins une fille ainsi que celles parmi celles-ci ayant un fils.

familles_avec_au_moins_une_fille = 0
familles_avec_au_moins_une_fille_et_un_fils = 0
for(i in (1:N)) {
  if(X1[i]==1 || X2[i]==1) {
    familles_avec_au_moins_une_fille <- familles_avec_au_moins_une_fille + 1
    if(X1[i]==0 || X2[i]==0) {
      familles_avec_au_moins_une_fille_et_un_fils <- familles_avec_au_moins_une_fille_et_un_fils + 1
    }
  }
}
familles_avec_au_moins_une_fille
## [1] 77
familles_avec_au_moins_une_fille_et_un_fils
## [1] 43

Une estimation de la probabilité précédente est donc:

familles_avec_au_moins_une_fille_et_un_fils/familles_avec_au_moins_une_fille
## [1] 0.5584

Ah, pas mal, c’est presque 1/2…

Comptage à la R

Bon, le code précédent est vraiment horrible, donc réécrivons le “à la R”:

familles_avec_au_moins_une_fille = sum(X1|X2);
familles_avec_au_moins_une_fille
## [1] 77

Maintenant la même chose pour compter les familles avec un fils et une fille:

familles_avec_au_moins_une_fille_et_un_fils = sum((X1|X2)&(!X1|!X2));
familles_avec_au_moins_une_fille_et_un_fils
## [1] 43

ou mieux:

familles_avec_au_moins_une_fille_et_un_fils = sum(xor(X1,X2));
familles_avec_au_moins_une_fille_et_un_fils
## [1] 43

Un peu de ménage

Bon, donc nettoyons un peu tout ça.

eval_proba_garcon_sachant_fille = function(N=100) {
  X1 = rbinom(N,1,1/2);
  X2 = rbinom(N,1,1/2);
  sum(xor(X1,X2))/sum(X1|X2);
}

Finalement, elle est toute simple cette fonction. Réévaluons notre probabilité (notez que je n’ai pas besoin de préciser N car je lui ai donnée une valeur par défaut lors de la définition de ma fonction:

eval_proba_garcon_sachant_fille()
## [1] 0.7273

Ooups, c’est vraiment plus grand que tout à l’heure… :) Essayons plusieurs fois pour voir:

for(i in (1:10)) {
  print(eval_proba_garcon_sachant_fille());
}
## [1] 0.7465
## [1] 0.7067
## [1] 0.7042
## [1] 0.8148
## [1] 0.6027
## [1] 0.6986
## [1] 0.6806
## [1] 0.6818
## [1] 0.6351
## [1] 0.6707

Ah, oui, c’est quand même carrément plus grand que 1/2 comme l’intuition le suggère. Ce n’est pas très pratique ces séries de nombres. Essayons quelque chose de plus “graphique”.

X = c();
for(i in (1:30)) {
  X[i] = eval_proba_garcon_sachant_fille();
}
plot(X,ylim = c(0,1))

plot of chunk unnamed-chunk-13

On voit bien que pour les différentes populations tirées au hasard, le ratio se ballade autour de 0.66.

Et si on essayait avec une population plus grande directement…

eval_proba_garcon_sachant_fille(100000);
## [1] 0.6679

Bon, et bien si ça devait être un nombre qui tombe rond, je dirais 2/3 du coup… :)

Calcul de la probabilité

Analyse des différents cas

Ici, on a indépendance et équiprobabilité donc il suffit de regarder les différents cas possibles

composition F G
F * X
G X .

Il y a \(3\) situations où la famille comporte une fille et sur ces \(3\) possibilités, il y en a \(2\) où la famille comporte un garçon. La probabilité qu’une famillle comporte un garçon sachant que la famille comporte une fille est donc de \(2/3\).

Un calcul formel

Formellement, si on note \(Y_1\) la variable aléatoire associée au sexe du premier enfant et \(Y_2\) la variable aléatoire associée au sexe du second enfant, on cherche à calculer \(P((Y_1=G \text{ ou } Y_2=G)|(Y_1=F \text{ ou } Y_2=F))\). On cherche donc à calculer une probabilité conditionnelle.

Rappelons qu’on a \(P(A|B)=P(A\text{ et }B)/P(B)\). Dans notre cas, pour utiliser cette formule, on définira \(A\) par la proposition \(\text{Il y a au moins un garçon}\) et \(B\) par la proposition \(\text{Il y a au moins une fille}\).

On a \(P(B) = 1-P(!B) = 1- 1/4 = 3/4\) et \(P(A\text{ et }B)=P(\text{il y a exactement une fille et un garçon})=1/2\). On retrouve bien \(P(A|B)=2/3\). Mais le petit tableau précédent est finalement bien plus édifiant, non ?

Il est intéressant de remarquer que dans le programme précédent, on a défini \(2n\) variables aléatoires (tous les \(X_1[i]\) et \(X_2[i]\)) et pas juste \(2\) (\(Y_1\) et \(Y_2\)) comme dans la preuve précédente. Notre programme simule une population de 100 familles en les considérent comme indépendantes les unes des autres et suivant la même loi. De même quand on tire \(X_1\) et \(X_2\), on suppose que le sexe du second enfant est indépendant de celui du premier et qu’il y a équiprobabilité des deux sexes.

Sauriez vous refaire le code et l’analyse précédente en supposant qu’on n’a plus équiprobabilité mais par exemple 1 chance sur 3 d’avoir une fille et 2 chanes sur 3 d’avoir un garçon ?

Ce qu’il faut retenir

  1. Quand on cherche à évaluer une probabilité, il faut se méfier de son intuition… :)
  2. Écrire une simulation est une bonne façon d’y voir plus clair et de se débarasser d’intuitions trompeuses.
  3. Si ici, il était facile de calculer la probabilité cherchée, sur des cas plus complexes il peut vite devenir très difficile voire impossible de faire un tel calcul. Une simulation est une façon très efficace d’en obtenir un estimateur.
  4. La valeur par simulation reste cependant une variable aléatoire. On verra plus tard comment évaluer l’erreur faite par un tel estimateur.
  5. Vous pouvez programmer avec R comme avec n’importe quel langage mais vous avez tout intérêt à apprendre à utiliser ce qui fait sa force:
    • Des fonctions/opérateurs de haut niveau (il est possible de travailler directement avec des vecteurs, voire avec des tableaux tout entier) qui permettent d’éviter l’usage de boucles.
    • Toutes les fonctions de générations, de tests et d’analyse statistiques, de représentation graphiques adaptées font partie intégrante du langage.
    • De la programmation littérale qui permet de faciliter la reproduction de votre travail.