#! /usr/bin/env python3
# -*- coding: utf-8; mode: python -*-
r""" Correction d'un exercice de maths avec Python, issu des annales pour les oraux du concours Centrale-Supélec.
- *Sujet :* PC 1, http://www.concours-centrale-supelec.fr/CentraleSupelec/MultiY/C2015/PC-Mat2-2015-27.pdf
- *Dates :* séances de révisions, vendredi 10 et samedi 11 juin (2016).
- *Auteur :* Lilian Besson, pour le cours d'informatique pour tous en prépa MP au Lycée Lakanal (http://perso.crans.org/besson/infoMP/),
- *Licence :* MIT Licence (http://lbesson.mit-license.org).
-----------------------------------------------------------------------------
"""
from __future__ import print_function, division # Python 2 compatibility
import matplotlib.pyplot as plt
# %% Question 0.
[docs]def q0():
r""" Préliminaires :
- Comme pour tous ces exercices, on commence par **importer** les modules requis, ici uniquement `matplotlib <http://matplotlib.org/users/beginner.html>`_ (pour la question `Q5 <#PC_Mat2_2015_27.q5>`_) :
>>> import matplotlib.pyplot as plt
"""
pass
# %% Question 1.
[docs]def q1():
r""" Question 1.
Question purement mathématique.
- Les variables aléatoires :math:`X_n` suivent toutes une loi de Bernoulli de paramètre :math:`p = \frac{1}{2}` (:math:`X_n \thicksim \mathcal{B}(\frac{1}{2})`), et sont indépendantes,
- Et la variable aléatoire :math:`S_n = \sum_{i=1}^{n} X_i` suit donc une loi binomiale de paramètres :math:`(n, p)`, par définition (:math:`S_n \thicksim \mathcal{Bin}(n, p)`.
- Donc on peut utiliser le cours directement :
.. math::
\mathbb{E}(S_n) = n \mathbb{E}(X_1) = n p = n \frac{1}{2} = \frac{n}{2}, \\
\mathbb{V}(S_n) = n \mathbb{V}(X_1) = n p (1-p) = n \frac{1}{4} = \frac{n}{4}.
"""
pass
# %% Question 2.
[docs]def q2():
r""" Question 2.
- On va noter :math:`f(n)` la fonction définie par :math:`\lfloor \frac{n}{2} \rfloor` si :math:`n` est impair (si :math:`n = 2 k + 1`, :math:`f(n) = k`), et :math:`\lfloor \frac{n-1}{2} \rfloor = \frac{n}{2} - 1` si :math:`n` est pair (si :math:`n = 2 k`, :math:`f(n) = k - 1`),
- :math:`f : \mathbb{N}^{*} \to \mathbb{N}` est bien définie, et vérifie :math:`\forall n > 0, f(n) = \min {i \in \mathbb{N} : i < \frac{n}{2}}`.
- En fait, :math:`f(n) = \lfloor \frac{n-1}{2} \rfloor` comme me l'avait fait remarqué un élève.
- Ainsi, comme :math:`S_n` est à support discret dans :math:`[0, n] \cap \mathbb{N}`, on peut écrire :math:`\mathbb{P}(S_n < \frac{n}{2})` comme une somme finie :
.. math::
\mathbb{P}(S_n < \frac{n}{2}) &= \sum_{i = 0}^{i < \frac{n}{2}} \mathbb{P}(S_n = i) \\
&= \sum_{i = 0}^{f(n)} \mathbb{P}(S_n = i) \\
&= \sum_{i = 0}^{f(n)} {n \choose k} p^{k} (1-p)^{n - k} \\
&= \sum_{i = 0}^{f(n)} {n \choose k} (\frac{1}{2})^{k} (\frac{1}{2})^{n - k} \\
&= (\frac{1}{2})^{n} \sum_{i = 0}^{f(n)} {n \choose k}.
- Cette dernière expression semble satisfaisante : c'est une expression sommatoire explicite, et l'énoncé demandait de ne pas chercher à la simplifier.
"""
pass
# %% Question 3.
[docs]def q3():
r""" Question 3.
Ici, pour écrire une fonction pour calculer :math:`{n \choose k}`, `plusieurs approches <https://en.wikipedia.org/wiki/Binomial_coefficient#Binomial_coefficient_in_programming_languages>`_ sont possibles :
1. on écrit une fonction :py:func:`factorial` (de :math:`\mathbb{N} \to \mathbb{N}`) qui calcule :math:`n!`, et on utilise :math:`{n \choose k} = \frac{n!}{k! (n-k)!}`. Ça marche, mais ce sera assez lent, et très redondant en calculs (cf. :py:func:`binom_factorial`).
2. ou on écrit une fonction :py:func:`binom` directement, en se basant sur la propriété élémentaire :math:`{n \choose k} = (\frac{n-k+1}{k}) {n \choose k-1}`.
C'est la deuxième approché, itérative et facile à écrire, que je recommande.
- Note: on pense à l'optimisation gratuite qui consiste à remplacer :math:`k` par :math:`\min(k, n-k)`.
Ça donne quelque chose comme ça : ::
def binom(n, k):
if k < 0 or k > n: # Cas de base (n k) = 0
return 0
if k == 0 or k == n: # Cas de base (n k) = 1
return 1
k = min(k, n - k) # Utilise la symétrie !
produit = 1
for i in range(k):
produit *= (n - i)
for i in range(k):
produit //= (i + 1)
return int(produit) # Force à être entier !
.. seealso::
En fait, on peut aussi utiliser :py:func:`scipy.special.binom` du module :py:mod:`scipy.special`, qui fonctionne très bien (elle renvoit un flottant attention) :
>>> from scipy.special import binom as binom2
>>> binom2(10, 2)
45.0
"""
pass
[docs]def binom(n, k):
r""" Calcule efficacement en :math:`\mathcal{O}(\min(k, n-k))` multiplication et division entières le nombre ``binom(n, k)`` :math:`= {n \choose k}` (et en :math:`\mathcal{O}(1)` en mémoire !).
- Exemples :
>>> binom(6, -3)
0
>>> binom(6, 0)
1
>>> binom(6, 1)
6
>>> binom(6, 2)
15
>>> binom(6, 3)
20
>>> binom(6, 4)
15
>>> binom(6, 5)
6
>>> binom(6, 6)
1
"""
if k < 0 or k > n: # Cas de base (n k) = 0
return 0
if k == 0 or k == n: # Cas de base (n k) = 1
return 1
k = min(k, n - k) # Utilise la symétrie !
produit = 1
for i in range(k):
produit *= (n - i)
for i in range(k):
produit //= (i + 1)
return int(produit) # Force à être entier !
from math import factorial
[docs]def binom_factorial(n, k):
r""" Coefficient binomial calculé en utilisant des factorielles (via :py:func:`math.factorial`) :
.. math: {n \choose k} = \frac{n!}{k! (n-k)!}.
"""
return factorial(n) // (factorial(k) * factorial(n - k))
# %% Question 4.
[docs]def q4():
r""" Question 4.
- La question 2 donne la formule, cf. `Q2 <#PC_Mat2_2015_27.q2>`_,
- Et la question précédente donne une fonction pour calcul :math:`{n \choose k}`, :py:func:`binom`,
- Il s'agit donc juste de combiner les deux : ::
def sn2(n):
n2 = ((n - 1) // 2) # f(n)
somme = sum(binom(n, k) for k in range(1, n2 + 1))
return somme / 2**n
.. hint:: Pour éviter des problèmes de conversions entre entiers (``int``, en précision arbitraire) et flottants (``float``, en précision à 16 chiffres), il vaut mieux calculer ``somme / 2**n`` (division entre deux gros entiers) que ``(1/2)**n * somme`` (petit flottant fois gros entier... = 0 !).
"""
pass
[docs]def sn2(n):
r""" Calcule la formule obtenue en `Q2 <#PC_Mat2_2015_27.q2>`_, en sommant les valeurs de :math:`{n \choose k}` pour :math:`k = 1, \dots, f(n)` (en utilisant :py:func:`binom`).
- Exemples :
>>> sn2(1)
0.0
>>> sn2(2)
0.0
>>> sn2(3)
0.375
>>> sn2(4)
0.25
>>> sn2(5)
0.46875
>>> sn2(50) # doctest: +ELLIPSIS
0.44386...
>>> sn2(500) # doctest: +ELLIPSIS
0.48216...
>>> sn2(1000) # doctest: +ELLIPSIS
0.48738...
>>> sn2(2000) # Ca commence à prendre du temps # doctest: +ELLIPSIS
0.49108...
"""
# n2 = (n // 2) if n % 2 == 1 else (n // 2) - 1 # f(n)
n2 = ((n - 1) // 2) # f(n)
# somme = sum(binom(n, k) for k in range(1, n2 + 1))
somme = sum(binom_factorial(n, k) for k in range(1, n2 + 1))
# return (1 / 2)**n * somme # FIXED Marche pas, cf. explication plus haut
return somme / 2**n
# %% Question 5.
[docs]def q5():
r""" Question 5.
- Avec `Q4 <#PC_Mat2_2015_27.q4>`_, on dispose de la fonction :py:func:`sn2`, et donc on peut calculer :math:`u_{20k} =` ``sn2(20k)``, très facilement.
- On crée une liste de valeurs de :math:`u_{20k}` pour :math:`k = 1,\dots,100`,
- Enfin, on utilise simplement la fonction :py:func:`matplotlib.pyplot.plot` (``plt.plot(X, Y)``) pour afficher ces points : ::
# On créé deux listes
valeurs_k = [ 20 * k for k in range(1, 100 + 1) ]
valeurs_u = [ sn2(k) for k in valeurs_k ]
plt.figure() # Nouvelle figure
plt.plot(valeurs_k, valeurs_u) # On affiche les points (k, sn2(20*k))
plt.title("Valeurs de u(20k) pour k = 1 .. 100")
plt.show() # On affiche la figure
- Cela devrait donner une figure comme ça :
.. image:: PC_Mat2_2015_27__figure1.png
- On peut donc raisonnablement conjecturer que la suite converge vers :math:`\frac{1}{2}`, ce qu'on va monter en `Q7 <#PC_Mat2_2015_27.q7>`_ et `Q8 <#PC_Mat2_2015_27.q8>`_ (de deux façons différentes), mais la convergence semble très lente.
"""
# On créé deux listes
valeurs_k = [20 * k for k in range(1, 100 + 1)]
valeurs_u = [sn2(k) for k in valeurs_k]
plt.figure() # Nouvelle figure
# plt.hold(True)
# On affiche les points (k, sn2(20*k))
plt.plot(valeurs_k, valeurs_u, 'b+-', label=r'$u_{20k}$')
# plt.plot(valeurs_k, [0.50] * len(valeurs_k), 'r:', label=r'Limite $\frac{1}{2}$')
plt.title(r"Valeurs de $u_{20k}$ pour $k = 1,...,100$")
# plt.legend()
plt.savefig("PC_Mat2_2015_27__figure1.png")
plt.show() # On affiche la figure
# %% Question 6.
[docs]def q6():
r""" Question 6.
- Par symétrie, on obtient directement que :math:`\mathbb{P}(S > \frac{n}{2}) = \mathbb{P}(S < \frac{n}{2}) = u_n`.
- Soit :math:`n \geq 1`. De façon immédiate, on peut décomposer l'intervalle entier :math:`[0, n] \cap \mathbb{N}` en trois morceaux disjoints :
.. math::
[0, n] \cap \mathbb{N} &= { i \in \mathbb{N} : 0 \leq i \leq n} \\
&= \{ 0 \leq i < \frac{n}{2} \} \cup \{ i = \frac{n}{2} \} \cup \{ \frac{n}{2} < i \leq n \}.
- Donc en considérant les probabilités pour la v.a. :math:`S_n`, on a directement (par la formule de l'union, cas disjoint) :
.. math::
1 &= \mathbb{P}(S \in [0, n] \cap \mathbb{N}) \\
&= \mathbb{P}(S \in \{ 0 \leq i < \frac{n}{2} \} \cup \{ i = \frac{n}{2} \} \cup \{ \frac{n}{2} < i \leq n \}) \\
&= \mathbb{P}(S \in \{ 0 \leq i < \frac{n}{2} \}) + \mathbb{P}(S \in \{ i = \frac{n}{2} \}) + \mathbb{P}(S \in \{ \frac{n}{2} < i \leq n \}) \\
&= \mathbb{P}(S < \frac{n}{2}) + \mathbb{P}(S_n = \frac{n}{2}) + \mathbb{P}(S > \frac{n}{2}) \\
&= u_n + \mathbb{P}(S_n = \frac{n}{2}) + \mathbb{P}(S > \frac{n}{2}).
- Et de même, on procède pareil pour montrer, trivialement, que :
.. math::
v_n &= \mathbb{P}(S \geq \frac{n}{2}) \\
&= \mathbb{P}(S_n = \frac{n}{2}) + \mathbb{P}(S > \frac{n}{2}) \\
&= \mathbb{P}(S_n = \frac{n}{2}) + \mathbb{P}(S < \frac{n}{2}) \\
&= \mathbb{P}(S_n = \frac{n}{2}) + u_n.
- **Conclusion :** on a bien montré que :math:`\forall n \geq 1, v_n = u_n + \mathbb{P}(S_n = \frac{n}{2})`.
"""
pass
# %% Question 7.
[docs]def q7():
r""" Question 7.
- D'après les calculs faits en `Q6 <#PC_Mat2_2015_27.q6>`_, :math:`u_n \to \frac{1}{2}` *ssi* :math:`\mathbb{P}(S_n = \frac{n}{2}) \to 0` pour :math:`n \to +\infty`.
- Mais, pour :math:`n` impair, on a déjà que :math:`\frac{n}{2} \not\in \mathbb{N}` donc :math:`\mathbb{P}(S_n = \frac{n}{2}) = 0`.
- Et pour :math:`n = 2k` pair, on peut écrire exactement la valeur de :math:`\mathbb{P}(S_n = \frac{n}{2}) = \mathbb{P}(S_n = k) = (\frac{1}{2})^n {n \choose k}`.
- On étudie la limite de :math:`\frac{1}{4^k}{2k \choose k}`, et sans même avoir besoin d'un équivalent (ie. formule de Stirling) on devrait (sa)voir que ça tend vers :math:`0`. Mais soyons rigoureux :
.. note::
*Preuve :* (avec Stirling). On a :
.. math::
{2k \choose k} &= \frac{(2k)!}{(k!)^2} \\
&\thicksim \frac{\sqrt{2 \pi 2k}(\frac{2k}{\mathrm{e}})^{2k}}{(\sqrt{2 \pi k}(\frac{k}{\mathrm{e}})^{k})^2} \\
&\thicksim \frac{2\sqrt{\pi k} 2^{2k} (\frac{k}{\mathrm{e}})^{2k}}{2 \pi k(\frac{k}{\mathrm{e}})^{2k}} \\
&\thicksim \frac{\sqrt{\pi k} 4^k}{\pi k} \\
&\thicksim \frac{4^k}{\sqrt{\pi k}}.
Et donc :math:`\frac{1}{4^k}{2k \choose k} \thicksim \frac{1}{\sqrt{\pi k}}`
pour :math:`k \to +\infty`, ce qui prouve que :math:`\mathbb{P}(S_{2k} = \frac{2k}{2}) \to 0` aussi.
- **Conclusion :** que :math:`n` soit pair ou impair, on a bien montré que :math:`\mathbb{P}(S_n = \frac{n}{2}) \to 0` pour :math:`n \to +\infty`.
- Et donc, cela prouve le résultat observé empiriquement en `Q5 <#PC_Mat2_2015_27.q5>`_ :
.. math::
\lim_{n \to +\infty} u_n = \frac{1}{2}.
.. note::
On a une estimation de la vitesse de convergence de :math:`\mathbb{P}(S_n = \frac{n}{2})` vers :math:`0` (en :math:`\mathcal{O}(\frac{1}{\sqrt{n}}`),
donc cela donne directement un ordre de grandeur de convergence de :math:`u_n` vers sa limite !
Cool non ?
En fait, c'est une conséquence immédiate du `théorème central limite <https://fr.wikipedia.org/wiki/Théorème_central_limite>`_, qui est LE résultat central en statistiques (cf. `ce point là <https://fr.wikipedia.org/wiki/Th%C3%A9or%C3%A8me_central_limite#Convergence_vers_la_limite>`_).
"""
pass
# %% Question 8.
[docs]def q8():
r""" Question 8.
.. hint:: On remarque que l'énoncé nous aide en suggérant ici de se restreindre aux :math:`n` pairs, ce qui est exactement ce qu'on a fait en `Q7 <#PC_Mat2_2015_27.q7>`_.
- De même qu'au dessus, on obtient que :
.. math:: w_n := \mathbb{P}(S_{2n} = n) = \frac{1}{4^n}{2n \choose n}
- Donc on peut calculer explicitement :math:`\ln(\frac{w_{n+1}}{w_n})` :
.. math::
\ln(\frac{w_{n+1}}{w_n}) &= \ln( \frac{1}{4} {2(n+1) \choose (n+1)}/{2n \choose n} ) \\
&= - \ln(4) + \ln( \frac{(2(n+1))!}{((n+1)!)^2} / \frac{(2n)!}{(n!)^2} ) \\
&= - \ln(4) + \ln( \frac{(2n+2)(2n+1)}{(n+1)^2} ) \\
&= - 2\ln(2) + \ln( 2\frac{2n+1}{n+1} ) \\
&= - 2\ln(2) + \ln(2) + \ln(\frac{2n+1}{n+1}) \\
&= - \ln(2) + \ln(\frac{(2n+1)}{n+1}) \\
&= \ln(\frac{(2n+1)}{2(n+1)}) \\
&= -\ln(\frac{2n+2}{2n+1}) = -\ln(\frac{(2n+1) + 1}{2n+1}) \\
&= -\ln(1 + \frac{1}{2n+1}).
- Ensuite, en notant :math:`x = \frac{1}{2n+1}`, :math:`x \to 0` si :math:`n \to \infty`, et comme :math:`\ln(1 + x) \thicksim x`, on obtient :
.. math::
\ln(\frac{w_{n+1}}{w_n}) &= -\ln(1 + \frac{1}{2n+1}) \\
&\thicksim - \frac{1}{2n+1} \to 0.
- Dès lors, on peut sommer la série de terme général :math:`\ln(\frac{w_{n+1}}{w_n})`,
pour obtenir, par comparaison avec la série divergente de terme général de signe constant (négatif) :math:`- \frac{1}{2n+1}`, une série divergente vers :math:`-\infty`.
- Mais on reconnaît une somme téléscopique : :math:`\sum_{n=0}^{N} \ln(\frac{w_{n+1}}{w_n}) = \ln(w_{N+1}) - \ln(w_{0})`, donc ça montre que :math:`\ln(w_{N+1}) \to -\infty` pour :math:`N \to +\infty`.
- Ainsi on obtient directement que :math:`w_n \to 0` pour :math:`n \to +\infty`, c'est ce qui était demandé.
- Enfin on obtient comme voulu que :math:`u_n \to \frac{1}{2}`.
.. note::
Ici aussi, cet équivalent :math:`\ln(\frac{w_{n+1}}{w_n}) \thicksim - \frac{1}{2n+1}` donne une indication
sur la vitesse de convergence de :math:`w_n \to 0` (en :math:`\mathcal{O}(\sqrt{n})`),
qui est la même que celle obtenue en `Q7 <#PC_Mat2_2015_27.q7>`_.
Cool.
"""
pass
# %% Fin de la correction.
if __name__ == '__main__':
from doctest import testmod
print("Test automatique de toutes les doctests ecrites dans la documentation (docstring) de chaque fonction :")
testmod(verbose=True)
# testmod()
print("\nPlus de details sur ces doctests peut etre trouve dans la documentation de Python:\nhttps://docs.python.org/3/library/doctest.html (en anglais)")
# Fin de PC_Mat2_2015_27.py