#+TITLE:  Autour du Paradoxe de Simpson
#+AUTHOR: Adam Taheraly
#+DATE:   \daily
#+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/readtheorg.css"/>
#+HTML_HEAD: <script src="https://ajax.googleapis.com/ajax/libs/jquery/2.1.3/jquery.min.js"></script>
#+HTML_HEAD: <script src="https://maxcdn.bootstrapcdn.com/bootstrap/3.3.4/js/bootstrap.min.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>

* Chargement des librairies R    :noexport:
#+begin_src R :results output :session *R* :exports both
library(devtools) # Obtention des informations sur la session R
library(dplyr) # Manipulation des données
library(tidyverse) # Collection de packages pour la data science
library(knitr) # Génération de rapports dynamique
library(forcats)
library(broman) # Function intéressante pour R
library(RColorBrewer) # Couleur des figures
library(ggplot2) # Création des figures
library(orgutils) # Utilitaire orgmode
library(xtable) # Jolie tableau
library(kableExtra) # Création de tableau
#+end_src


* Mise en place des données :noexport:
** Chargement
#+begin_src R :results output :session *R* :exports both
simpson_data <- read.csv("Subject6_smoking.csv")
#+end_src

Les données sont sous la forme :

| Libellé | Description                   |
|---------+-------------------------------|
| Smoker  | Fume (Yes)  ou non (No)       |
| Status  | En vie (Alive) ou mort (Dead) |
| Age     | Age du sujet                  |


** Vérification
Vérification du bon chargement des données en visualisant le début et
la fin des données.
#+begin_src R :results output :session *R* :exports both
head(simpson_data)
tail(simpson_data)
#+end_src

Vérification de la présence de données manquante.
#+begin_src R :results output :session *R* :exports both
na_records <- apply(simpson_data, 1, function(x) any(is.na(x)))
simpson_data[na_records,]
#+end_src


Vérification de la classe des données et du nombre de data.
#+begin_src R :results output :session *R* :exports both
nrow(simpson_data)
class(simpson_data$Smoker)
class(simpson_data$Status)
class(simpson_data$Age)
#+end_src


* Analyse des données
** Question 1
#+begin_src R :results output :session *R* :exports both
dead_or_alive <- simpson_data%>%
  group_by(Smoker)%>%
summarise(Number = n(), Alive = sum(Status == "Alive"), Dead = sum(Status == "Dead"), mortality = (sum(Status == "Dead")/n())*100)
dead_or_alive$mortality <- round(dead_or_alive$mortality, 2)
#+end_src

#+begin_src R :results table :colnames yes :session *R* :exports results
xtable(dead_or_alive)
#+end_src


#+begin_src R :results output graphics :file (org-babel-temp-file "figure" ".png") :exports both :width 600 :height 400 :session *R* 
mortality <- ggplot(dead_or_alive, aes(x = Smoker, y = mortality, fill = Smoker)) + geom_bar(stat = "identity") + 
    labs(x = "Smoker", y = "Mortality rate", title = "Mortality rate by behaviour")
mortality
#+end_src

Le résultat est surprenant car on observe un plus fort taux de
mortalité chez les non-fumeur que chez les fumeur.

** Question 2
On recrée un data-frame pour ne pas interagir avec les résultats
précédent.
#+begin_src R :results output :session *R* :exports both
simpson_data_q2 <- simpson_data
simpson_data_q2$Age <- as.numeric(simpson_data_q2$Age)
#+end_src



On regroupe les ages en catégories
#+begin_src R :results output :session *R* :exports both
simpson_data_q2$Age <- cut(simpson_data_q2$Age, c(18, 34, 54, 64, 200), include.lowest = TRUE, labels = c("18-34", "35-54", "55-64", "plus de 65"))
table(simpson_data_q2$Age)
#+end_src

#+begin_src R :results output :session *R* :exports both
dead_or_alive2 <- simpson_data_q2%>%
  group_by(Smoker, Age)%>%
summarise(Number = n(), Alive = sum(Status == "Alive"), Dead = sum(Status == "Dead"))%>%
ungroup()%>%
group_by(Smoker)%>%
mutate(mortality = (Dead/sum(Number))*100)
dead_or_alive$mortality <- round(dead_or_alive$mortality, 2)
#+end_src

#+begin_src R :results table :colnames yes :session *R* :exports results
xtable(dead_or_alive2)
#+end_src

#+begin_src R :results output graphics :file (org-babel-temp-file "figure" ".png") :exports both :width 600 :height 400 :session *R* 
mortality2 <- ggplot(dead_or_alive2, aes(x = Smoker, y = mortality, fill = Age))+ geom_bar(stat = "identity", position = position_dodge()) + 
    labs(x = "Smoker", y = "Mortality rate", title = "Mortality rate by behaviour")
mortality2
#+end_src

On observe une diminution du taux de mortalité chez les moins de 65
ans parmis les non fumeur par rapport au fumeur. Cependant, chez les
plus de 65 ans, la tendence s'inverse avec un plus grande mortalité
chez les fumeurs par rappor au non fumeur (3 fois plus).

L'inversion de ces proportions peut s'expliquer par un enrichissement
en personnes agées chez les non fumeurs.

** Question 3
On recrée un data-frame pour ne pas interagir avec les résultats
précédent.

#+begin_src R :results output :session *R* :exports both
simpson_data_q3 <- simpson_data
simpson_data_q3$Death <-  factor(simpson_data_q3$Status)
simpson_data_q3$Death <-  fct_recode(simpson_data_q3$Death, "1"="Alive", "0"="Dead")
simpson_data_q3$Death <- as.numeric(as.character(simpson_data_q3$Death))
simpson_data_q3$Age <- as.numeric(as.character(simpson_data_q3$Age))
#+end_src


#+begin_src R :results output :session *R* :exports both
head(simpson_data_q3)
#+end_src


#+begin_src R :results output :session *R* :exports both
reg_log_total <- ggplot(simpson_data_q3, aes(x=Age,y=Death)) + geom_point(alpha=.3,size=3) + 
    theme_bw() + 
    geom_smooth(method = "glm", 
        method.args = list(family = "binomial"),fullrange = TRUE)

reg_log_fumeur <- ggplot(simpson_data_q3[simpson_data_q3$Smoker == "Yes",], aes(x=Age,y=Death)) + geom_point(alpha=.3,size=3) + 
    theme_bw() + 
    geom_smooth(method = "glm", 
        method.args = list(family = "binomial"),fullrange = TRUE)

reg_log_non_fumeur <- ggplot(simpson_data_q3[simpson_data_q3$Smoker == "No",], aes(x=Age,y=Death)) + geom_point(alpha=.3,size=3) + 
    theme_bw() + 
    geom_smooth(method = "glm", 
        method.args = list(family = "binomial"),fullrange = TRUE)
#+end_src


#+begin_src R :results output graphics :file (org-babel-temp-file "figure" ".png") :exports both :width 600 :height 400 :session *R* 
reg_log_total
#+end_src


#+begin_src R :results output graphics :file (org-babel-temp-file "figure" ".png") :exports both :width 600 :height 400 :session *R* 
reg_log_fumeur
#+end_src


#+begin_src R :results output graphics :file (org-babel-temp-file "figure" ".png") :exports both :width 600 :height 400 :session *R* 
reg_log_non_fumeur
#+end_src

Passer 60 ans, bien que les courbes soient semblable, on voit que l'intervalle de confiance chez les fumeurs
est plus grands que chez les non fumeur. Cela s'explique bien par une
plus faible population de personne agée chez les fumeurs.

Etant donné l'augmentation de la mortalités observé dans la question
2, on peut conclure à une mortalité supérieur dans le cas d'une
consommation de tabac chez les moins de 65 ans. Cependant, il est
necessaire de faire une analyse avec des populations plus homogènes
chez les + de 65 ans.