You are currently viewing Une énigme pour moi…

Une énigme pour moi…

Bonjour les gens ! En effectuant quelques recherches pour mon prochain livre à paraître uniquement sur mathweb.fr, je suis tombé sur le problème suivant:

On choisit au hasard un nombre entre 0 et 1, et on répète cela jusqu’à ce que la somme des nombres choisis dépasse 1.

Quelle est le nombre moyen de nombres choisis ?

Le problème peut s’implémenter en Python à l’aide de la fonction suivante:

from random import random

def moyenne(n):
    T = 0
    for k in range(n):
        s, i = 0, 0

        while s < 1:
            s += random()
            i += 1
            
        T += i
    
    return T / n

print( moyenne( 10**9 ) )

On obtient alors:

>>> moyenne(10**7)
2.7186881
>>> moyenne(10**8)
2.71829368
>>> moyenne(10**9)
2.718286834

La moyenne semble être égale à “e” (ce n’est qu’une conjecture).

Le problème est que je n’arrivais pas à établir la moindre formule…

La solution de Philippe Rackette

Fort heureusement, il y a parmi les visiteurs et abonné·e·s de ce site de bien meilleurs probabilistes que moi ! C’est ainsi que l’un d’eux m’a proposé la solution suivante (le lien vers le PDF a été mis en commentaire, mais il va sans doute arriver un jour où le fichier correspondant sera effacé, donc j’ai préféré le retranscrire ici.

Posons \(S_{k+1} = P(X>k+1)\) pour tout entier naturel \(k\). Alors,$$S_{k+1}=\int_{x_1+\cdots+x_{k+1}<1,\ x_i \geqslant0} 1\text{d}x_1\cdots \text{d}x_{k+1}$$En effet, la densité de probabilité est homogène et la probabilité qu’il faille plus de k+1 nombres nécessite que la somme des k+1 premiers nombres soit strictement inférieure à 1.

On peur alors écrire pour k non nul:$$S_{k+1}=\int_{x_{k+1}=0}^{x_{k+1}=1^-}\left( \int_{x_1+\ \cdots\ +x_k<1-x_{k+1},\ x_i\geqslant0}1\text{d}x_1\cdots\text{d}x_k \right)\text{d}x_{k+1}$$

En effectuant les changements de variables:$$u_i = \frac{1}{1-x_{k+1}}x_i$$pour i variant de 1 à k, on obtient:$$S_{k+1}=\int_{x_{k+1}=0}^{x_{k+1}=1^-}\left( \int_{u_1+\ \cdots\ +u_k<1,\ u_i\geqslant0}(1-x_{k+1})^k\text{d}u_1\cdots\text{d}u_k \right)\text{d}x_{k+1}$$que l’on peut aussi écrire: $$S_{k+1}=\int_{x_{k+1}=0}^{x_{k+1}=1^-}\left( \int_{u_1+\ \cdots\ +u_k<1,\ u_i\geqslant0}1\text{d}u_1\cdots\text{d}u_k \right)(1-x_{k+1})^k\text{d}x_{k+1}$$c’est-à-dire:$$S_{k+1}=\int_{x_{k+1}=0}^{x_{k+1}=1^-}S_k(1-x_{k+1})^k\text{d}x_{k+1}$$

Comme \(S_k\) est indépendant de \(x_{k+1}\), on a alors:$$S_{k+1}=S_k\int_{x_{k+1}=0}^{x_{k+1}=1^-}(1-x_{k+1})^k\text{d}x_{k+1}$$

Or,$$\int_0^1(1-x)^k\text{d}x=\left[-\frac{(1-x)^{k+1}}{k+1}\right]_0^1=\frac{1}{k+1}$$ d’où:$$S_{k+1}=\frac{1}{k+1}S_k.$$

De plus, \(S_1=1\) car l’événement “plus d’un nombre est nécessaire pour obtenir une somme strictement plus grande que 1” est certain. On a alors:$$\begin{align}S_k&=\frac{1}{k}S_{k-1}\\& = \frac{1}{k}\times\frac{1}{k-1}S_{k-2}\\&=\cdots\\&=\frac{1}{k}\times\frac{1}{k-1}\times\cdots\times\frac{1}{k-i}S_{k-i-1}\\&=\cdots\\&=\frac{1}{k}\times\frac{1}{k-1}\times\cdots\times\frac{1}{k-(k-2)}S_{k-(k-2)-1}\\&=\frac{1}{k}\times\frac{1}{k-1}\times\cdots\times\frac{1}{2}S_{1}\\& = \frac{1}{k!}\end {align}$$

Ainsi,$$P(X\leqslant k)=1-P(X>k)=1-S_k=1-\frac{1}{k!}$$

Par suite, on a:$$\begin{align}P(X=k) & =P(X\leqslant k)-P(X\leqslant k-1)\\&=\left(1-\frac{1}{k!}\right)-\left(1-\frac{1}{(k-1)!}\right)\\&=\frac{1}{(k-1)!}-\frac{1}{k!}\\&=\frac{k-1}{k!} \end{align}$$

On peut alors calculer le nombre moyen que l’on cherchait en déterminant l’espérance de X:$$\begin{align} E(X) & =\sum_{k\geqslant1} k \times P(X=k)\\& = \sum_{k\geqslant1}k\times\frac{k-1}{k!}\\&=\sum_{k\geqslant1}\frac{k-1}{(k-1)!}\\&=\sum_{p\geqslant0}\frac{1}{p!}\\&=\text{e}.\end{align}$$

L’énigme est ainsi résolue! Me voilà rassuré… La conjecture de départ était donc bonne, mais c’est bien mieux de voir pourquoi!

En mathématiques, tout est une question de point de vue: je n’avais pas le bon dès le départ et je ne pouvais donc pas avoir la bonne idée !

Cet article a 16 commentaires

  1. Est-ce que la valeur 1/6 pour P(X=2) est “confirmée” par la simulation ? Cela me semble peu. Je la verrais plutôt à 1/2 selon le “raisonnement” suivant : après un premier tirage aléatoire, qu’on sait d’avance ne pas être le dernier, la somme cumulée est “probablement autour de 1/2” donc il y a une chance sur deux pour qu’au tirage suivant la somme cumulée dépasse 1 et une chance sur deux qu’elle ne dépasse pas 1. En continuant… Si après 2 tirages, 1 n’a toujours pas était atteint, le somme cumulée est “probablement autour de 2/3” et il y a 1 chance sur 3 de ne pas dépasser la somme cumulée 1 au 3e tirage. Si après k tirages, la somme cumulée 1 n’a toujours pas était atteinte, elle est “probablement autour de k/(k+1)” et il y a une chance sur k+1 de ne pas dépasser la somme cumulée 1 au k+1 ème tirage. En admettant (ou expliquant) cela on peut en déduire que l’espérance est bien 1+1/2+1/3!+1/4!+1/5!+… = e

    1. Le raisonnement exposé dans l’article est loin d’être correct je pense… Donc 1/6, en effet, paraît peu. C’est d’ailleurs cette valeur qui indique que mon raisonnement est sans doute faux. Pour ce qui est de ce que vous expliquez, il faudrait un spécialiste des probabilités là… car tout ça est obscure pour moi (si Mickaël Launay venait à passer par ici…:-))

  2. On pourrait commencer par essayer de montrer que pour tout entier k supérieur à 1, \(S_{k+1}=\frac{1}{k+1}S_k\) avec \(S_k=\int_{x_1+…x_k \leq, 1 x_i \geq 0} 1 dx_1 … dx_k\) (la densité de probabilité est homogène). Au rang 1, cela signifie que la moitié inférieure gauche du carré unitaire a une aire de 1/2, au rang 2, que la pyramide (O,I,J,K) du repère canonique de l’espace a un volume de 1/6. \(S_k\) correspond à l’espérance de \(X_1+…X_k\) sachant que \(X_1+…X_k\) est inférieur à 1 où les \(X_k\) suivent une loi uniforme sur [0,1[… et \(S_k=\int_{x_1+…+x_k \leq 1, x_i \geq 0} 1 dx_1 … dx_k\) représente \(P(X > k)=1-P(X \leq k)\)

  3. Christophe Bal

    Le code suivant confirme l’intuition de P. Rackette.

    from random import random
    
    moy = 0
    rep = 10**6
    
    for _ in range(rep):
        N1 = random()
        N2 = random()
    
        if N1 + N2 > 1:
            moy += 1
            
    moy /= rep
    
    print(moy)
    
  4. Ce dernier code donne :

    
    0.499667
    0.500387
    0.499749
    

    qui “montre” que la probabilité que la somme des deux premiers nombres pris au hasard semble être égale à \(\frac{1}{2}\) (précision apportée pour les visiteurs).

  5. Christophe Bal

    Une idée géométrique… À creuser je pense.

    Pour \(X_1\) et \(X_2\), l’ensemble des points \((X_1 ; X_2)\) est dans un carré de côté 1, et les points tels que \(X_1 + X_2 \leqslant 1\) correspondent à la moitié basse du carré d’aire \(\frac{1}{2}\) d’où le 1/2 de proba.

    Pour \(X_1,\ X_2, X_3\), l’ensemble des points \((X_1 ; X_2 ; X_3)\) est dans un cube de côté 1 et les points tels que \(X_1 + X_2 + X_3 \leqslant 1\) correspondent à une pyramide de volume \(\frac{1}{2}\times\frac{1}{3}\) d’où le 1/6 de proba.

    Il faudrait continuer avec de n simplexes droits.

  6. Cette façon de voir est en effet intéressante… mais pour la généraliser à la dimension \(n\), c’est une autre histoire… 🙂

  7. Christophe Bal

    Voici comment généraliser.

    Prenons \(n\) vecteurs unitaires \(e_1,\ , e_2 ,\ \ldots ,\ e_n\) formant un repère orthonormé.

    Il est connu que le volume du polytope formé par ces vecteurs est \( | det(e_1 , e_2 , … , e_n) |\) soit ici 1 .

    Il y a n! façons d’ordonner les vecteurs. Pour un ordre choisi \(e_{k_1} ,\ e_{k_2} ,\ \ldots ,\ e_{k_n}\) , on a un repère affine \( (e_{k_1} ,\ e_{k_2} ,\ \ldots ,\ e_{k_n}) \) d’origine \( e_{k_1} \) avec des abus évidents de notations.

    Dans ce repère, nous avons un volume représentant de façon naturelle les points \( (X_{k_1} ; X_{k_2} ; \ldots ; X_{k_n}) \) tels que \( X_{k_1} + X_{k_2} + \ldots + X_{k_n} \leqslant 1\).

    On a donc qu’un tel volume aura \( \frac{| det(e_1 , e_2 , … , e_n) |}{n!} = \frac{1}{n!} \) pour volume.

  8. CHRISTOPHE BAL

    Ma réponse n’est pas claire dans le choix des repères (ce n’est pas bon en fait). J’ai du mal à exprimer ce choix. L’idée est d’aller chercher les simplexes dans chaque “coin” afin de couvrir le polytope.

    1. Perso, je suis largué 🙂 Même avec ces explications, je ne vois pas comment obtenir une espérance de \(e\) au final.

  9. Christophe Bal

    Je vais essayer de dégager du temps pour rédiger un truc au propre.

    Je vous tiendrais au courant.

    1. Si vous rédigez un document au format PDF (c’est sans doute plus simple que sur un site internet), n’hésitez pas à partager le lien: peut-être que d’autres personnes seraient intéressées 🙂

  10. LAU BOULANGER

    Bonjour Stéphane,
    en ce qui concerne ta première interrogation c’est effectivement un problème d’indépendance.Si on pose X ta variable aléatoire, il faut calculer l’espérance de cette variable aléatoire.
    P(X=1)=0
    P(X=2)=0.5, il faut calculer une intégrale double sur le domaine 0<x<1 et 1-x<y2) c’est à dire P(X1+X2+X3>1|X1+Y<1)
    Donc calculer l'aire du domaine 0<x<1, 0<y<1-x, 1-x-y<z2).
    On obtient 1/6.
    On continue ainsi, j’imagine que l’on peut faire une récurrence, mais j’ai la flemme de rédiger .
    Ensuite l’espérance n’est rien d’autre que la somme des inverse des factorielles, du coup on trouve le bon résultat.
    Bravo pour ton travail sur LaTeX, en tout cas.
    Par contre si tu veux intégrer des résultats statistiques dans LaTeX, tu pourrais regarder Sweave pour R en utilisant l’IDE Rstudio.
    Bon dimanche
    amicalement
    Laurent

  11. LAU BOULANGER

    Désolé je n’avais pas vu la réponse de Philippe !

  12. CHRISTOPHE BAL

    Merci Philippe.

    Il resterait à passer via un argument géométrique au lieu de calculs d’intégrales pour que ce soit plus joli (je n’ai toujours pas le temps de taper ceci proprement mais c’est promis je le ferais d’ici mi juin et je posterais ceci ici).

Laisser un commentaire