Commit 54d16143 authored by 3e6bf7b47a05a05ae3d6af86121dcb5d's avatar 3e6bf7b47a05a05ae3d6af86121dcb5d
parents 21aae8b9 a2fba00c
#+TITLE: Votre titre #+TITLE: À propos du calcul de $\pi$
#+AUTHOR: Votre nom #+AUTHOR: François Févotte
#+DATE: La date du jour #+DATE: 2020-04-06
#+LANGUAGE: fr #+LANGUAGE: fr
# #+PROPERTY: header-args :eval never-export # #+PROPERTY: header-args :eval never-export
...@@ -11,83 +11,67 @@ ...@@ -11,83 +11,67 @@
#+HTML_HEAD: <script type="text/javascript" src="http://www.pirilampo.org/styles/lib/js/jquery.stickytableheaders.js"></script> #+HTML_HEAD: <script type="text/javascript" src="http://www.pirilampo.org/styles/lib/js/jquery.stickytableheaders.js"></script>
#+HTML_HEAD: <script type="text/javascript" src="http://www.pirilampo.org/styles/readtheorg/js/readtheorg.js"></script> #+HTML_HEAD: <script type="text/javascript" src="http://www.pirilampo.org/styles/readtheorg/js/readtheorg.js"></script>
* Quelques explications * En demandant à la lib maths
Ceci est un document org-mode avec quelques exemples de code Mon ordinateur m'indique que $\pi$ vaut /approximativement/ :
python. Une fois ouvert dans emacs, ce document peut aisément être #+begin_src python :results value :exports both :session
exporté au format HTML, PDF, et Office. Pour plus de détails sur from math import *
org-mode vous pouvez consulter https://orgmode.org/guide/. pi
#+end_src
Lorsque vous utiliserez le raccourci =C-c C-e h o=, ce document sera #+RESULTS:
compilé en html. Tout le code contenu sera ré-exécuté, les résultats : 3.141592653589793
récupérés et inclus dans un document final. Si vous ne souhaitez pas
ré-exécuter tout le code à chaque fois, il vous suffit de supprimer
le # et l'espace qui sont devant le ~#+PROPERTY:~ au début de ce
document.
Comme nous vous l'avons montré dans la vidéo, on inclue du code * En utilisant la méthode des aiguilles de Buffon
python de la façon suivante (et on l'exécute en faisant ~C-c C-c~):
#+begin_src python :results output :exports both Mais calculé avec la *méthode* des [[https://fr.wikipedia.org/wiki/Aiguille_de_Buffon][aiguilles de Buffon]], on obtiendrait comme
print("Hello world!") *approximation* :
#+end_src #+BEGIN_SRC python :results value :exports both :session
import numpy as np
np.random.seed(seed=42)
N = 10000
x = np.random.uniform(size=N, low=0, high=1)
theta = np.random.uniform(size=N, low=0, high=pi/2)
2/(sum((x+np.sin(theta))>1)/N)
#+END_SRC
#+RESULTS: #+RESULTS:
: Hello world! : 3.128911138923655
Voici la même chose, mais avec une session python, donc une * Avec un argument "fréquentiel" de surface
persistance d'un bloc à l'autre (et on l'exécute toujours en faisant
~C-c C-c~). Sinon, une méthode plus simple à comprendre et ne faisant pas intervenir d'appel
#+begin_src python :results output :session :exports both à la fonction sinus se base sur le fait que si $X \sim U(0,1)$ et $Y \sim
import numpy U(0,1)$ alors $P\left[X^2 + Y^2 \leq 1\right] = \pi/4$ (voir [[https://fr.wikipedia.org/wiki/M%C3%A9thode_de_Monte-Carlo#D%C3%A9termination_de_la_valeur_de_%CF%80][méthode de Monte
x=numpy.linspace(-15,15) Carlo sur Wikipedia]]). Le code suivant illustre ce fait :
print(x) #+BEGIN_SRC python :results output file :exports both :session :var matplot_lib_filename="montecarlo.png"
#+end_src import matplotlib.pyplot as plt
np.random.seed(seed=42)
N = 1000
x = np.random.uniform(size=N, low=0, high=1)
y = np.random.uniform(size=N, low=0, high=1)
accept = (x*x+y*y) <= 1
reject = np.logical_not(accept)
fig, ax = plt.subplots(1)
ax.scatter(x[accept], y[accept], c='b', alpha=0.2, edgecolor=None)
ax.scatter(x[reject], y[reject], c='r', alpha=0.2, edgecolor=None)
ax.set_aspect('equal')
plt.savefig(matplot_lib_filename)
print(matplot_lib_filename)
#+END_SRC
#+RESULTS: #+RESULTS:
#+begin_example [[file:montecarlo.png]]
[-15. -14.3877551 -13.7755102 -13.16326531 -12.55102041
-11.93877551 -11.32653061 -10.71428571 -10.10204082 -9.48979592 Il est alors aisé d'obtenir une approximation (pas terrible) de $\pi$ en
-8.87755102 -8.26530612 -7.65306122 -7.04081633 -6.42857143 comptant combien de fois, en moyenne, $X^2 + Y^2$ est inférieur à 1 :
-5.81632653 -5.20408163 -4.59183673 -3.97959184 -3.36734694 #+BEGIN_SRC python :results value :exports both :session
-2.75510204 -2.14285714 -1.53061224 -0.91836735 -0.30612245 4*np.mean(accept)
0.30612245 0.91836735 1.53061224 2.14285714 2.75510204 #+END_SRC
3.36734694 3.97959184 4.59183673 5.20408163 5.81632653
6.42857143 7.04081633 7.65306122 8.26530612 8.87755102
9.48979592 10.10204082 10.71428571 11.32653061 11.93877551
12.55102041 13.16326531 13.7755102 14.3877551 15. ]
#+end_example
Et enfin, voici un exemple de sortie graphique:
#+begin_src python :results output file :session :var matplot_lib_filename="./cosxsx.png" :exports results
import matplotlib.pyplot as plt
plt.figure(figsize=(10,5))
plt.plot(x,numpy.cos(x)/x)
plt.tight_layout()
plt.savefig(matplot_lib_filename)
print(matplot_lib_filename)
#+end_src
#+RESULTS: #+RESULTS:
[[file:./cosxsx.png]] : 3.112
Vous remarquerez le paramètre ~:exports results~ qui indique que le code
ne doit pas apparaître dans la version finale du document. Nous vous
recommandons dans le cadre de ce MOOC de ne pas changer ce paramètre
(indiquer ~both~) car l'objectif est que vos analyses de données soient
parfaitement transparentes pour être reproductibles.
Attention, la figure ainsi générée n'est pas stockée dans le document
org. C'est un fichier ordinaire, ici nommé ~cosxsx.png~. N'oubliez pas
de le committer si vous voulez que votre analyse soit lisible et
compréhensible sur GitLab.
Enfin, n'oubliez pas que nous vous fournissons dans les ressources de
ce MOOC une configuration avec un certain nombre de raccourcis
claviers permettant de créer rapidement les blocs de code python (en
faisant ~<p~, ~<P~ ou ~<PP~ suivi de ~Tab~).
Maintenant, à vous de jouer! Vous pouvez effacer toutes ces
informations et les remplacer par votre document computationnel.
#+TITLE: Votre titre #+TITLE: Exercice 02
#+AUTHOR: Votre nom #+AUTHOR: Votre nom
#+DATE: La date du jour #+DATE: La date du jour
#+LANGUAGE: fr #+LANGUAGE: fr
...@@ -11,83 +11,65 @@ ...@@ -11,83 +11,65 @@
#+HTML_HEAD: <script type="text/javascript" src="http://www.pirilampo.org/styles/lib/js/jquery.stickytableheaders.js"></script> #+HTML_HEAD: <script type="text/javascript" src="http://www.pirilampo.org/styles/lib/js/jquery.stickytableheaders.js"></script>
#+HTML_HEAD: <script type="text/javascript" src="http://www.pirilampo.org/styles/readtheorg/js/readtheorg.js"></script> #+HTML_HEAD: <script type="text/javascript" src="http://www.pirilampo.org/styles/readtheorg/js/readtheorg.js"></script>
* Quelques explications
Ceci est un document org-mode avec quelques exemples de code #+BEGIN_SRC python :exports none :session
python. Une fois ouvert dans emacs, ce document peut aisément être import numpy
exporté au format HTML, PDF, et Office. Pour plus de détails sur #+END_SRC
org-mode vous pouvez consulter https://orgmode.org/guide/.
Lorsque vous utiliserez le raccourci =C-c C-e h o=, ce document sera * Les données
compilé en html. Tout le code contenu sera ré-exécuté, les résultats
récupérés et inclus dans un document final. Si vous ne souhaitez pas
ré-exécuter tout le code à chaque fois, il vous suffit de supprimer
le # et l'espace qui sont devant le ~#+PROPERTY:~ au début de ce
document.
Comme nous vous l'avons montré dans la vidéo, on inclue du code #+begin_src python :results output :exports code :session
python de la façon suivante (et on l'exécute en faisant ~C-c C-c~): x = [14.0, 7.6, 11.2, 12.8, 12.5, 9.9, 14.9, 9.4, 16.9, 10.2, 14.9, 18.1, 7.3,
9.8, 10.9,12.2, 9.9, 2.9, 2.8, 15.4, 15.7, 9.7, 13.1, 13.2, 12.3, 11.7, 16.0,
#+begin_src python :results output :exports both 12.4, 17.9, 12.2, 16.2, 18.7, 8.9, 11.9, 12.1, 14.6, 12.1, 4.7, 3.9, 16.9, 16.8,
print("Hello world!") 11.3, 14.4, 15.7, 14.0, 13.6, 18.0, 13.6, 19.9, 13.7, 17.0, 20.5, 9.9, 12.5,
13.2, 16.1, 13.5, 6.3, 6.4, 17.6, 19.1, 12.8, 15.5, 16.3, 15.2, 14.6, 19.1,
14.4, 21.4, 15.1, 19.6, 21.7, 11.3, 15.0, 14.3, 16.8, 14.0, 6.8, 8.2, 19.9,
20.4, 14.6, 16.4, 18.7, 16.8, 15.8, 20.4, 15.8, 22.4, 16.2, 20.3, 23.4, 12.1,
15.5, 15.4, 18.4, 15.7, 10.2, 8.9, 21.0]
#+end_src #+end_src
#+RESULTS: #+RESULTS:
: Hello world!
Voici la même chose, mais avec une session python, donc une * Quelques statistiques
persistance d'un bloc à l'autre (et on l'exécute toujours en faisant
~C-c C-c~). ** Moyenne
#+begin_src python :results output :session :exports both
import numpy #+BEGIN_SRC python :results value :exports both :session
x=numpy.linspace(-15,15) "%.2f" % numpy.mean(x)
print(x) #+END_SRC
#+end_src
#+RESULTS: #+RESULTS:
#+begin_example : 14.11
[-15. -14.3877551 -13.7755102 -13.16326531 -12.55102041
-11.93877551 -11.32653061 -10.71428571 -10.10204082 -9.48979592
-8.87755102 -8.26530612 -7.65306122 -7.04081633 -6.42857143 ** Ecart-type
-5.81632653 -5.20408163 -4.59183673 -3.97959184 -3.36734694
-2.75510204 -2.14285714 -1.53061224 -0.91836735 -0.30612245 #+BEGIN_SRC python :results value :exports both :session
0.30612245 0.91836735 1.53061224 2.14285714 2.75510204 "%.2f" % numpy.std(x, ddof=1)
3.36734694 3.97959184 4.59183673 5.20408163 5.81632653 #+END_SRC
6.42857143 7.04081633 7.65306122 8.26530612 8.87755102
9.48979592 10.10204082 10.71428571 11.32653061 11.93877551
12.55102041 13.16326531 13.7755102 14.3877551 15. ]
#+end_example
Et enfin, voici un exemple de sortie graphique:
#+begin_src python :results output file :session :var matplot_lib_filename="./cosxsx.png" :exports results
import matplotlib.pyplot as plt
plt.figure(figsize=(10,5))
plt.plot(x,numpy.cos(x)/x)
plt.tight_layout()
plt.savefig(matplot_lib_filename)
print(matplot_lib_filename)
#+end_src
#+RESULTS: #+RESULTS:
[[file:./cosxsx.png]] : 4.334094455301447
Vous remarquerez le paramètre ~:exports results~ qui indique que le code ** Médiane
ne doit pas apparaître dans la version finale du document. Nous vous
recommandons dans le cadre de ce MOOC de ne pas changer ce paramètre #+BEGIN_SRC python :results value :exports both :session
(indiquer ~both~) car l'objectif est que vos analyses de données soient "%.2f" % numpy.median(x)
parfaitement transparentes pour être reproductibles. #+END_SRC
Attention, la figure ainsi générée n'est pas stockée dans le document #+RESULTS:
org. C'est un fichier ordinaire, ici nommé ~cosxsx.png~. N'oubliez pas : 14.5
de le committer si vous voulez que votre analyse soit lisible et
compréhensible sur GitLab.
Enfin, n'oubliez pas que nous vous fournissons dans les ressources de ** Valeurs extrêmes
ce MOOC une configuration avec un certain nombre de raccourcis
claviers permettant de créer rapidement les blocs de code python (en #+BEGIN_SRC python :results value :exports both :session
faisant ~<p~, ~<P~ ou ~<PP~ suivi de ~Tab~). numpy.min(x), numpy.max(x)
#+END_SRC
Maintenant, à vous de jouer! Vous pouvez effacer toutes ces
informations et les remplacer par votre document computationnel. #+RESULTS:
| 2.8 | 23.4 |
#+TITLE: Votre titre #+TITLE: Exercice 02-3 : Réaliser un affichage graphique
#+AUTHOR: Votre nom #+AUTHOR: Votre nom
#+DATE: La date du jour #+DATE: La date du jour
#+LANGUAGE: fr #+LANGUAGE: fr
...@@ -11,59 +11,33 @@ ...@@ -11,59 +11,33 @@
#+HTML_HEAD: <script type="text/javascript" src="http://www.pirilampo.org/styles/lib/js/jquery.stickytableheaders.js"></script> #+HTML_HEAD: <script type="text/javascript" src="http://www.pirilampo.org/styles/lib/js/jquery.stickytableheaders.js"></script>
#+HTML_HEAD: <script type="text/javascript" src="http://www.pirilampo.org/styles/readtheorg/js/readtheorg.js"></script> #+HTML_HEAD: <script type="text/javascript" src="http://www.pirilampo.org/styles/readtheorg/js/readtheorg.js"></script>
* Quelques explications #+BEGIN_SRC python :exports none :session
import numpy
Ceci est un document org-mode avec quelques exemples de code import matplotlib.pyplot as plt
python. Une fois ouvert dans emacs, ce document peut aisément être #+END_SRC
exporté au format HTML, PDF, et Office. Pour plus de détails sur
org-mode vous pouvez consulter https://orgmode.org/guide/.
Lorsque vous utiliserez le raccourci =C-c C-e h o=, ce document sera
compilé en html. Tout le code contenu sera ré-exécuté, les résultats
récupérés et inclus dans un document final. Si vous ne souhaitez pas
ré-exécuter tout le code à chaque fois, il vous suffit de supprimer
le # et l'espace qui sont devant le ~#+PROPERTY:~ au début de ce
document.
Comme nous vous l'avons montré dans la vidéo, on inclue du code
python de la façon suivante (et on l'exécute en faisant ~C-c C-c~):
#+begin_src python :results output :exports both
print("Hello world!")
#+end_src
#+RESULTS: #+RESULTS:
: Hello world!
Voici la même chose, mais avec une session python, donc une * Les données
persistance d'un bloc à l'autre (et on l'exécute toujours en faisant
~C-c C-c~). #+begin_src python :results output :exports code :session
#+begin_src python :results output :session :exports both x = [14.0, 7.6, 11.2, 12.8, 12.5, 9.9, 14.9, 9.4, 16.9, 10.2, 14.9, 18.1, 7.3,
import numpy 9.8, 10.9,12.2, 9.9, 2.9, 2.8, 15.4, 15.7, 9.7, 13.1, 13.2, 12.3, 11.7, 16.0,
x=numpy.linspace(-15,15) 12.4, 17.9, 12.2, 16.2, 18.7, 8.9, 11.9, 12.1, 14.6, 12.1, 4.7, 3.9, 16.9, 16.8,
print(x) 11.3, 14.4, 15.7, 14.0, 13.6, 18.0, 13.6, 19.9, 13.7, 17.0, 20.5, 9.9, 12.5,
13.2, 16.1, 13.5, 6.3, 6.4, 17.6, 19.1, 12.8, 15.5, 16.3, 15.2, 14.6, 19.1,
14.4, 21.4, 15.1, 19.6, 21.7, 11.3, 15.0, 14.3, 16.8, 14.0, 6.8, 8.2, 19.9,
20.4, 14.6, 16.4, 18.7, 16.8, 15.8, 20.4, 15.8, 22.4, 16.2, 20.3, 23.4, 12.1,
15.5, 15.4, 18.4, 15.7, 10.2, 8.9, 21.0]
#+end_src #+end_src
#+RESULTS: #+RESULTS:
#+begin_example
[-15. -14.3877551 -13.7755102 -13.16326531 -12.55102041
-11.93877551 -11.32653061 -10.71428571 -10.10204082 -9.48979592
-8.87755102 -8.26530612 -7.65306122 -7.04081633 -6.42857143
-5.81632653 -5.20408163 -4.59183673 -3.97959184 -3.36734694
-2.75510204 -2.14285714 -1.53061224 -0.91836735 -0.30612245
0.30612245 0.91836735 1.53061224 2.14285714 2.75510204
3.36734694 3.97959184 4.59183673 5.20408163 5.81632653
6.42857143 7.04081633 7.65306122 8.26530612 8.87755102
9.48979592 10.10204082 10.71428571 11.32653061 11.93877551
12.55102041 13.16326531 13.7755102 14.3877551 15. ]
#+end_example
Et enfin, voici un exemple de sortie graphique: * Plot
#+begin_src python :results output file :session :var matplot_lib_filename="./cosxsx.png" :exports results
import matplotlib.pyplot as plt
#+begin_src python :results output file :session :var matplot_lib_filename="./plot.png" :exports results
plt.figure(figsize=(10,5)) plt.figure(figsize=(10,5))
plt.plot(x,numpy.cos(x)/x) plt.plot(x)
plt.tight_layout() plt.tight_layout()
plt.savefig(matplot_lib_filename) plt.savefig(matplot_lib_filename)
...@@ -71,23 +45,19 @@ print(matplot_lib_filename) ...@@ -71,23 +45,19 @@ print(matplot_lib_filename)
#+end_src #+end_src
#+RESULTS: #+RESULTS:
[[file:./cosxsx.png]] [[file:./plot.png]]
Vous remarquerez le paramètre ~:exports results~ qui indique que le code * Histogramme
ne doit pas apparaître dans la version finale du document. Nous vous
recommandons dans le cadre de ce MOOC de ne pas changer ce paramètre
(indiquer ~both~) car l'objectif est que vos analyses de données soient
parfaitement transparentes pour être reproductibles.
Attention, la figure ainsi générée n'est pas stockée dans le document #+begin_src python :results output file :session :var matplot_lib_filename="./hist.png" :exports results
org. C'est un fichier ordinaire, ici nommé ~cosxsx.png~. N'oubliez pas plt.figure(figsize=(10,5))
de le committer si vous voulez que votre analyse soit lisible et plt.hist(x)
compréhensible sur GitLab. plt.tight_layout()
Enfin, n'oubliez pas que nous vous fournissons dans les ressources de plt.savefig(matplot_lib_filename)
ce MOOC une configuration avec un certain nombre de raccourcis print(matplot_lib_filename)
claviers permettant de créer rapidement les blocs de code python (en #+end_src
faisant ~<p~, ~<P~ ou ~<PP~ suivi de ~Tab~).
Maintenant, à vous de jouer! Vous pouvez effacer toutes ces #+RESULTS:
informations et les remplacer par votre document computationnel. [[file:./hist.png]]
#+TITLE: Votre titre #+TITLE: Quelques chiffres autour de l'épidémie de Covid-19
#+AUTHOR: Votre nom #+AUTHOR: François Févotte
#+DATE: La date du jour # #+DATE: La date du jour
#+LANGUAGE: fr #+LANGUAGE: fr
# #+PROPERTY: header-args :eval never-export
#+HTML_HEAD: <link rel="stylesheet" type="text/css" href="http://www.pirilampo.org/styles/readtheorg/css/htmlize.css"/> #+HTML_HEAD: <link rel="stylesheet" type="text/css" href="http://www.pirilampo.org/styles/readtheorg/css/htmlize.css"/>
#+HTML_HEAD: <link rel="stylesheet" type="text/css" href="http://www.pirilampo.org/styles/readtheorg/css/readtheorg.css"/> #+HTML_HEAD: <link rel="stylesheet" type="text/css" href="http://www.pirilampo.org/styles/readtheorg/css/readtheorg.css"/>
...@@ -11,83 +10,223 @@ ...@@ -11,83 +10,223 @@
#+HTML_HEAD: <script type="text/javascript" src="http://www.pirilampo.org/styles/lib/js/jquery.stickytableheaders.js"></script> #+HTML_HEAD: <script type="text/javascript" src="http://www.pirilampo.org/styles/lib/js/jquery.stickytableheaders.js"></script>
#+HTML_HEAD: <script type="text/javascript" src="http://www.pirilampo.org/styles/readtheorg/js/readtheorg.js"></script> #+HTML_HEAD: <script type="text/javascript" src="http://www.pirilampo.org/styles/readtheorg/js/readtheorg.js"></script>
* Quelques explications
Ceci est un document org-mode avec quelques exemples de code * Identification du système
python. Une fois ouvert dans emacs, ce document peut aisément être
exporté au format HTML, PDF, et Office. Pour plus de détails sur
org-mode vous pouvez consulter https://orgmode.org/guide/.
Lorsque vous utiliserez le raccourci =C-c C-e h o=, ce document sera #+BEGIN_SRC shell :exports both
compilé en html. Tout le code contenu sera ré-exécuté, les résultats uname -a
récupérés et inclus dans un document final. Si vous ne souhaitez pas #+END_SRC
ré-exécuter tout le code à chaque fois, il vous suffit de supprimer
le # et l'espace qui sont devant le ~#+PROPERTY:~ au début de ce
document.
Comme nous vous l'avons montré dans la vidéo, on inclue du code #+RESULTS:
python de la façon suivante (et on l'exécute en faisant ~C-c C-c~): : Linux aluminium 4.14.0-3-amd64 #1 SMP Debian 4.14.17-1 (2018-02-14) x86_64 GNU/Linux
#+begin_src python :results output :exports both #+BEGIN_SRC python :results value :session :exports both
print("Hello world!") import sys
#+end_src sys.version
#+END_SRC
#+RESULTS: #+RESULTS:
: Hello world! : 3.7.3rc1 (default, Mar 13 2019, 11:01:15)
: [GCC 8.3.0]
* Récupération des données
Voici la même chose, mais avec une session python, donc une On prend les données distribuées par [[https://ourworldindata.org][Our World In Data]], elles-mêmes issues de
persistance d'un bloc à l'autre (et on l'exécute toujours en faisant l'[[https://www.ecdc.europa.eu/][ECDC]] :
~C-c C-c~).
#+begin_src python :results output :session :exports both #+BEGIN_SRC python :results output :session :exports none
import numpy import numpy as np
x=numpy.linspace(-15,15) import matplotlib.pyplot as plt
print(x) import pandas as pd
#+end_src plt.style.use('seaborn')
#+END_SRC
#+RESULTS:
#+BEGIN_SRC python :results value :session :exports both
url = "https://covid.ourworldindata.org/data/ecdc/full_data.csv"
data_raw = pd.read_csv(url, index_col='date')
data_raw.tail(3)
#+END_SRC
#+RESULTS:
: location new_cases new_deaths total_cases total_deaths
: date
: 2020-04-06 Zimbabwe 0 0 9 1
: 2020-04-07 Zimbabwe 0 0 9 1
: 2020-04-08 Zimbabwe 1 0 10 1
* Evolution comparée de l'épidémie dans quelques pays
Je préfère fonder l'analyse sur le nombre de décès, qui me semble plus fiable et
comparable d'un pays à l'autre que, par exemple, le nombre de cas dépistés. Afin
de recaler l'axe temporel, on positionne $t=0$ au moment où le nombre total de
morts dépasse 100.
La logique est la suivante : on considère le nombre $n^c$ (de décès ou de cas)
dans un pays $c$. On suppose que l'évolution de $n^c$ en fonction du temps $t$
suit une loi de croissance exponentielle (au moins dans une phase initiale). On
s'attend ainsi à avoir une expression du type :
\[
n^c(t) = n^c_0 \, e^{\lambda^c \left(t-t^c_0\right)}.
\]
Dans cette expression, la valeur initiale $n^c_0 = n^c(t^c_0)$ est liée à la
taille de la population considérée. La constante de temps $\lambda^c$ en revanche,
n'a pas vraiment de raison de l'être : elle reflète simplement le rythme auquel
l'épidémie se propage.
Si l'on se donne un seuil $\tau$ unique, on peut fixer pour chaque pays
l'instant initial $t^c_0$ de sorte que
\[
n^c(t^c_0) = \tau.
\]
En effectuant un recalage de l'instant initial, c'est à dire en effectuant le
changement de variable temporelle $\tilde{t}^c=t-t^c_0$, on obtient
\[
n^c(\tilde{t}^c) = \tau \, e^{\lambda^c\,\tilde{t}^c},
\]
ce qui signifie que tous les pays deviennent comparables, indépendamment de leur
taille : on s'est ramené à comparer les constantes de temps $\lambda^c$ entre
elles.
En pratique, on prend
\[
t^c_0 = \min \left\{t ; n^c(t) > \tau \right\}.
\]
#+BEGIN_SRC python :results value :exports both :session
# Parameters
column = "total_deaths"
threshold = 100
countries = ["France", "Italy", "Spain", "Belgium", "Germany", "Switzerland"]
total = {}
for country in countries:
total[country] = data_raw.query(f"location == '{country}' & {column} > {threshold}")
total[countries[0]]
#+END_SRC
#+RESULTS: #+RESULTS:
#+begin_example #+begin_example
[-15. -14.3877551 -13.7755102 -13.16326531 -12.55102041 location new_cases new_deaths total_cases total_deaths
-11.93877551 -11.32653061 -10.71428571 -10.10204082 -9.48979592 date
-8.87755102 -8.26530612 -7.65306122 -7.04081633 -6.42857143 2020-03-16 France 924 36 5423 127
-5.81632653 -5.20408163 -4.59183673 -3.97959184 -3.36734694 2020-03-17 France 1210 21 6633 148
-2.75510204 -2.14285714 -1.53061224 -0.91836735 -0.30612245 2020-03-18 France 1097 27 7730 175
0.30612245 0.91836735 1.53061224 2.14285714 2.75510204 2020-03-19 France 1404 69 9134 244
3.36734694 3.97959184 4.59183673 5.20408163 5.81632653 2020-03-20 France 1861 128 10995 372
6.42857143 7.04081633 7.65306122 8.26530612 8.87755102 2020-03-21 France 1617 78 12612 450
9.48979592 10.10204082 10.71428571 11.32653061 11.93877551 2020-03-22 France 1847 112 14459 562
12.55102041 13.16326531 13.7755102 14.3877551 15. ] 2020-03-23 France 1559 112 16018 674
2020-03-24 France 3838 186 19856 860
2020-03-25 France 2446 240 22302 1100
2020-03-26 France 2931 231 25233 1331
2020-03-27 France 3922 365 29155 1696
2020-03-28 France 3809 299 32964 1995
2020-03-29 France 4611 319 37575 2314
2020-03-30 France 2599 292 40174 2606
2020-03-31 France 4376 418 44550 3024
2020-04-01 France 7578 499 52128 3523
2020-04-02 France 4861 509 56989 4032
2020-04-03 France 2116 471 59105 4503
2020-04-04 France 5233 2004 64338 6507
2020-04-05 France 4267 1053 68605 7560
2020-04-06 France 1873 518 70478 8078
2020-04-07 France 3912 833 74390 8911
2020-04-08 France 3777 1417 78167 10328
#+end_example #+end_example
Et enfin, voici un exemple de sortie graphique: #+BEGIN_SRC python :session :results output file :exports results :var pltfile="total.png"
#+begin_src python :results output file :session :var matplot_lib_filename="./cosxsx.png" :exports results plt.figure(figsize=(10,5))
import matplotlib.pyplot as plt for country in total.keys():
x = total[country]['total_deaths']
plt.plot(np.arange(len(x)), x, label=country, linewidth=3)
plt.xlabel(f"Number of days since {column} > {threshold}")
plt.yscale("log")
plt.legend(fontsize='x-large')
plt.tight_layout()
plt.savefig(pltfile)
plt.close()
print(pltfile)
#+END_SRC
#+RESULTS:
[[file:total.png]]
plt.figure(figsize=(10,5))
plt.plot(x,numpy.cos(x)/x)
plt.tight_layout()
plt.savefig(matplot_lib_filename) * Vers un passage du pic ?
print(matplot_lib_filename)
#+end_src On s'intéresse ici à l'atteinte d'un pic du nombre de décès quotidiens liés à
Covid-19. Là encore, l'axe temporel est recalé afin de permettre une
comparaison entre pays : $t=0$ est positionné au moment où le nombre quotidien
de morts dépasse 15.
#+BEGIN_SRC python :session :results value :exports both
threshold = 15
daily = {}
for country in countries:
daily[country] = data_raw.query(f"location == '{country}' & new_deaths > {threshold}")
daily[countries[0]]
#+END_SRC
#+RESULTS:
#+begin_example
location new_cases new_deaths total_cases total_deaths
date
2020-03-14 France 785 18 3661 79
2020-03-16 France 924 36 5423 127
2020-03-17 France 1210 21 6633 148
2020-03-18 France 1097 27 7730 175
2020-03-19 France 1404 69 9134 244
2020-03-20 France 1861 128 10995 372
2020-03-21 France 1617 78 12612 450
2020-03-22 France 1847 112 14459 562
2020-03-23 France 1559 112 16018 674
2020-03-24 France 3838 186 19856 860
2020-03-25 France 2446 240 22302 1100
2020-03-26 France 2931 231 25233 1331
2020-03-27 France 3922 365 29155 1696
2020-03-28 France 3809 299 32964 1995
2020-03-29 France 4611 319 37575 2314
2020-03-30 France 2599 292 40174 2606
2020-03-31 France 4376 418 44550 3024
2020-04-01 France 7578 499 52128 3523
2020-04-02 France 4861 509 56989 4032
2020-04-03 France 2116 471 59105 4503
2020-04-04 France 5233 2004 64338 6507
2020-04-05 France 4267 1053 68605 7560
2020-04-06 France 1873 518 70478 8078
2020-04-07 France 3912 833 74390 8911
2020-04-08 France 3777 1417 78167 10328
#+end_example
Le pic observé pour la France au 04/04 (jour 20 sur le graphique) correspond à la prise en compte
instantanée de tous les décès en EHPAD, dont le décompte est disponible à partir
du 01/02 (/cf./ [[https://www.gouvernement.fr/info-coronavirus/carte-et-donnees][données du gouvernement]])
#+BEGIN_SRC python :session :results output file :exports results :var pltfile="daily.png"
plt.figure(figsize=(10,5))
for country in daily.keys():
x = daily[country]['new_deaths']
plt.plot(np.arange(len(x)), x, label=country, linewidth=3)
plt.xlabel(f"Number of days since daily deaths > {threshold}")
plt.yscale("log")
plt.legend(fontsize='x-large')
plt.tight_layout()
plt.savefig(pltfile)
plt.close()
print(pltfile)
#+END_SRC
#+RESULTS: #+RESULTS:
[[file:./cosxsx.png]] [[file:daily.png]]
Vous remarquerez le paramètre ~:exports results~ qui indique que le code # Local Variables:
ne doit pas apparaître dans la version finale du document. Nous vous # org-babel-python-command: "python3"
recommandons dans le cadre de ce MOOC de ne pas changer ce paramètre # org-confirm-babel-evaluate: nil
(indiquer ~both~) car l'objectif est que vos analyses de données soient # End:
parfaitement transparentes pour être reproductibles.
Attention, la figure ainsi générée n'est pas stockée dans le document
org. C'est un fichier ordinaire, ici nommé ~cosxsx.png~. N'oubliez pas
de le committer si vous voulez que votre analyse soit lisible et
compréhensible sur GitLab.
Enfin, n'oubliez pas que nous vous fournissons dans les ressources de
ce MOOC une configuration avec un certain nombre de raccourcis
claviers permettant de créer rapidement les blocs de code python (en
faisant ~<p~, ~<P~ ou ~<PP~ suivi de ~Tab~).
Maintenant, à vous de jouer! Vous pouvez effacer toutes ces
informations et les remplacer par votre document computationnel.
...@@ -74,25 +74,44 @@ température (en Fahrenheit) et la pression (en psi), et enfin le ...@@ -74,25 +74,44 @@ température (en Fahrenheit) et la pression (en psi), et enfin le
nombre de dysfonctionnements relevés. nombre de dysfonctionnements relevés.
* Inspection graphique des données * Inspection graphique des données
Les vols où aucun incident n'est relevé n'apportant aucune information
sur l'influence de la température ou de la pression sur les
dysfonctionnements, nous nous concentrons sur les expériences où au
moins un joint a été défectueux.
#+begin_src python :results value :session *python* :exports both # Les vols où aucun incident n'est relevé n'apportant aucune information
data = data[data.Malfunction>0] # sur l'influence de la température ou de la pression sur les
data # dysfonctionnements, nous nous concentrons sur les expériences où au
#+end_src # moins un joint a été défectueux.
#+RESULTS: # #+begin_src python :results value :session *python* :exports both
: Date Count Temperature Pressure Malfunction # data = data[data.Malfunction>0]
: 1 11/12/81 6 70 50 1 # data
: 8 2/03/84 6 57 200 1 # #+end_src
: 9 4/06/84 6 63 200 1
: 10 8/30/84 6 70 200 1 # #+RESULTS:
: 13 1/24/85 6 53 200 2 # #+begin_example
: 20 10/30/85 6 75 200 2 # Date Count Temperature Pressure Malfunction
: 22 1/12/86 6 58 200 1 # 0 4/12/81 6 66 50 0
# 1 11/12/81 6 70 50 1
# 2 3/22/82 6 69 50 0
# 3 11/11/82 6 68 50 0
# 4 4/04/83 6 67 50 0
# 5 6/18/82 6 72 50 0
# 6 8/30/83 6 73 100 0
# 7 11/28/83 6 70 100 0
# 8 2/03/84 6 57 200 1
# 9 4/06/84 6 63 200 1
# 10 8/30/84 6 70 200 1
# 11 10/05/84 6 78 200 0
# 12 11/08/84 6 67 200 0
# 13 1/24/85 6 53 200 2
# 14 4/12/85 6 67 200 0
# 15 4/29/85 6 75 200 0
# 16 6/17/85 6 70 200 0
# 17 7/29/85 6 81 200 0
# 18 8/27/85 6 76 200 0
# 19 10/03/85 6 79 200 0
# 20 10/30/85 6 75 200 2
# 21 11/26/85 6 76 200 0
# 22 1/12/86 6 58 200 1
# #+end_example
Très bien, nous avons une variabilité de température importante mais Très bien, nous avons une variabilité de température importante mais
la pression est quasiment toujours égale à 200, ce qui devrait la pression est quasiment toujours égale à 200, ce qui devrait
...@@ -145,26 +164,27 @@ logmodel.summary() ...@@ -145,26 +164,27 @@ logmodel.summary()
#+begin_example #+begin_example
Generalized Linear Model Regression Results Generalized Linear Model Regression Results
============================================================================== ==============================================================================
Dep. Variable: Frequency No. Observations: 7 Dep. Variable: Frequency No. Observations: 23
Model: GLM Df Residuals: 5 Model: GLM Df Residuals: 21
Model Family: Binomial Df Model: 1 Model Family: Binomial Df Model: 1
Link Function: logit Scale: 1.0 Link Function: logit Scale: 1.0000
Method: IRLS Log-Likelihood: -3.6370 Method: IRLS Log-Likelihood: -3.9210
Date: Fri, 20 Jul 2018 Deviance: 3.3763 Date: Thu, 09 Apr 2020 Deviance: 3.0144
Time: 16:56:08 Pearson chi2: 0.236 Time: 12:50:13 Pearson chi2: 5.00
No. Iterations: 5 No. Iterations: 6
Covariance Type: nonrobust
=============================================================================== ===============================================================================
coef std err z P>|z| [0.025 0.975] coef std err z P>|z| [0.025 0.975]
------------------------------------------------------------------------------- -------------------------------------------------------------------------------
Intercept -1.3895 7.828 -0.178 0.859 -16.732 13.953 Intercept 5.0850 7.477 0.680 0.496 -9.570 19.740
Temperature 0.0014 0.122 0.012 0.991 -0.238 0.240 Temperature -0.1156 0.115 -1.004 0.316 -0.341 0.110
=============================================================================== ===============================================================================
#+end_example #+end_example
L'estimateur le plus probable du paramètre de température est 0.0014 # L'estimateur le plus probable du paramètre de température est 0.0014
et l'erreur standard de cet estimateur est de 0.122, autrement dit on # et l'erreur standard de cet estimateur est de 0.122, autrement dit on
ne peut pas distinguer d'impact particulier et il faut prendre nos # ne peut pas distinguer d'impact particulier et il faut prendre nos
estimations avec des pincettes. # estimations avec des pincettes.
* Estimation de la probabilité de dysfonctionnant des joints toriques * Estimation de la probabilité de dysfonctionnant des joints toriques
La température prévue le jour du décollage est de 31°F. Essayons La température prévue le jour du décollage est de 31°F. Essayons
...@@ -187,36 +207,43 @@ print(matplot_lib_filename) ...@@ -187,36 +207,43 @@ print(matplot_lib_filename)
#+RESULTS: #+RESULTS:
[[file:proba_estimate_python.png]] [[file:proba_estimate_python.png]]
Comme on pouvait s'attendre au vu des données initiales, la D'après cette analyse, la probabilité d'échec des joints toriques serait
température n'a pas d'impact notable sur la probabilité d'échec des supérieure à 80% à pour une température de 31°F.
joints toriques. Elle sera d'environ 0.2, comme dans les essais
précédents où nous il y a eu défaillance d'au moins un joint. Revenons # Comme on pouvait s'attendre au vu des données initiales, la
à l'ensemble des données initiales pour estimer la probabilité de # température n'a pas d'impact notable sur la probabilité d'échec des
défaillance d'un joint: # joints toriques. Elle sera d'environ 0.2, comme dans les essais
# précédents où nous il y a eu défaillance d'au moins un joint. Revenons
#+begin_src python :results output :session *python* :exports both # à l'ensemble des données initiales pour estimer la probabilité de
data = pd.read_csv("shuttle.csv") # défaillance d'un joint:
print(np.sum(data.Malfunction)/np.sum(data.Count))
#+end_src # #+begin_src python :results output :session *python* :exports both
# data = pd.read_csv("shuttle.csv")
#+RESULTS: # print(np.sum(data.Malfunction)/np.sum(data.Count))
: 0.06521739130434782 # #+end_src
Cette probabilité est donc d'environ $p=0.065$, sachant qu'il existe # #+RESULTS:
un joint primaire un joint secondaire sur chacune des trois parties du # : 0.06521739130434782
lançeur, la probabilité de défaillance des deux joints d'un lançeur
est de $p^2 \approx 0.00425$. La probabilité de défaillance d'un des # Cette probabilité est donc d'environ $p=0.065$, sachant qu'il existe
lançeur est donc de $1-(1-p^2)^3 \approx 1.2%$. Ça serait vraiment # un joint primaire un joint secondaire sur chacune des trois parties du
pas de chance... Tout est sous contrôle, le décollage peut donc avoir # lançeur, la probabilité de défaillance des deux joints d'un lançeur
lieu demain comme prévu. # est de $p^2 \approx 0.00425$. La probabilité de défaillance d'un des
# lançeur est donc de $1-(1-p^2)^3 \approx 1.2%$. Ça serait vraiment
Seulement, le lendemain, la navette Challenger explosera et emportera # pas de chance... Tout est sous contrôle, le décollage peut donc avoir
avec elle ses sept membres d'équipages. L'opinion publique est # lieu demain comme prévu.
fortement touchée et lors de l'enquête qui suivra, la fiabilité des
joints toriques sera directement mise en cause. Au delà des problèmes # Seulement, le lendemain, la navette Challenger explosera et emportera
de communication interne à la NASA qui sont pour beaucoup dans ce # avec elle ses sept membres d'équipages. L'opinion publique est
fiasco, l'analyse précédente comporte (au moins) un petit # fortement touchée et lors de l'enquête qui suivra, la fiabilité des
problème... Saurez-vous le trouver ? Vous êtes libre de modifier cette # joints toriques sera directement mise en cause. Au delà des problèmes
analyse et de regarder ce jeu de données sous tous les angles afin # de communication interne à la NASA qui sont pour beaucoup dans ce
d'expliquer ce qui ne va pas. # fiasco, l'analyse précédente comporte (au moins) un petit
# problème... Saurez-vous le trouver ? Vous êtes libre de modifier cette
# analyse et de regarder ce jeu de données sous tous les angles afin
# d'expliquer ce qui ne va pas.
# Local Variables:
# org-babel-python-command: "/home/francois/.julia/conda/3/bin/python"
# org-confirm-babel-evaluate: nil
# End:
module2/exo5/freq_temp_python.png

12.3 KB | W: | H:

module2/exo5/freq_temp_python.png

13 KB | W: | H:

module2/exo5/freq_temp_python.png
module2/exo5/freq_temp_python.png
module2/exo5/freq_temp_python.png
module2/exo5/freq_temp_python.png
  • 2-up
  • Swipe
  • Onion skin
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment