---
title: "Sujet 6 : Autour du Paradoxe de Simpson"
author: "Louis Hognon"
date: "19/04/2020"
output: html_document

Consignes : 
1. Représentez dans un tableau le nombre total de femmes vivantes et décédées sur la période en fonction de leur habitude de tabagisme. Calculez dans chaque groupe (fumeuses / non fumeuses) le taux de mortalité (le rapport entre le nombre de femmes décédées dans un groupe et le nombre total de femmes dans ce groupe). Vous pourrez proposer une représentation graphique de ces données et calculer des intervalles de confiance si vous le souhaitez. En quoi ce résultat est-il surprenant ?

2. Reprenez la question 1 (effectifs et taux de mortalité) en rajoutant une nouvelle catégorie liée à la classe d'âge. On considérera par exemple les classes suivantes : 18-34 ans, 34-54 ans, 55-64 ans, plus de 65 ans. En quoi ce résultat est-il surprenant ? Arrivez-vous à expliquer ce paradoxe ? De même, vous pourrez proposer une représentation graphique de ces données pour étayer vos explications.

3. Afin d'éviter un biais induit par des regroupements en tranches d'âges arbitraires et non régulières, il est envisageable d'essayer de réaliser une régression logistique. Si on introduit une variable Death valant 1 ou 0 pour indiquer si l'individu est décédé durant la période de 20 ans, on peut étudier le modèle Death ~ Age pour étudier la probabilité de décès en fonction de l'âge selon que l'on considère le groupe des fumeuses ou des non fumeuses. Ces régressions vous permettent-elles de conclure sur la nocivité du tabagisme ? Vous pourrez proposer une représentation graphique de ces régressions (en n'omettant pas les régions de confiance).

---


#Préparation des données #### 


Les données **autour du Paradoxe de Simpson** se trouvent sont acessibles via le site de [gitlab INRIA](https://gitlab.inria.fr/learninglab/mooc-rr/mooc-rr-ressources/blob/master/module3/Practical_session/Subject6_smoking.csv)

L'URL de ce lien est :

```{r}
data_url= "https://gitlab.inria.fr/learninglab/mooc-rr/mooc-rr-ressources/-/raw/master/module3/Practical_session/Subject6_smoking.csv?inline=false"
```


Pour **eviter une éventuelle disparition des données du serveur gitlab**, nous faisons une **copie locale de ce jeux de données** que nous préservons  Il est inutile et même risquée de télécharger les données à chaque exécution, car dans le cas d'une panne nous pourrions remplacer nos données par un fichier défectueux. Pour cette raison, nous téléchargeons les données seulement si la copie locale n'existe pas.
```{r}
data_file = "simpsons.csv"
if (!file.exists(data_file)) {
    download.file(data_url, data_file, method="auto")
}
```


## Téléchargement des données ###

Nous téléchargons les données via la commande **read.csv**, appliquée à data_url et data_file, dans le cas d'une éventuelle disparition des données du serveur gitlab. 


```{r}
data_exo6 = read.csv(data_url)
data_exo6 = read.csv(data_file)
```



# Description de la base de donnée ###


Avant de commencer les analyses sur la base de données, nous allons d'abord réaliser une description de celle-ci. Puis vérifier si certains paramètres : type de données, données manquantes, données abberantes, 

## Afficher la structure de la base de donnée 

Avec la commande **str** nous allons vérifier le type de données dans chaque  colonne ainsi que les premières valeurs.

```{r}
str(data_exo6)
```


<u>On retrouve bien nos 3 colonnes</u>: 

Smoker : si la personne fume ou non
Status : si elle est vivante ou décédée au moment de la seconde étude
Age : son âge lors du premier sondage

Pour rappel : *En 1972-1974, à Whickham, une ville du nord-est de l'Angleterre, située à environ 6,5 kilomètres au sud-ouest de Newcastle upon Tyne, un sondage d'un sixième des électeurs a été effectué afin d'éclairer des travaux sur les maladies thyroïdiennes et cardiaques (Tunbridge et al. 1977). Une suite de cette étude a été menée vingt ans plus tard (Vanderpump et al. 1995). Certains des résultats avaient trait au tabagisme et cherchaient à savoir si les individus étaient toujours en vie lors de la seconde étude. Par simplicité, nous nous restreindrons aux femmes et parmi celles-ci aux 1314 qui ont été catégorisées comme "fumant actuellement" ou "n'ayant jamais fumé". Il y avait relativement peu de femmes dans le sondage initial ayant fumé et ayant arrêté depuis (162) et très peu pour lesquelles l'information n'était pas disponible (18). La survie à 20 ans a été déterminée pour l'ensemble des femmes du premier sondage.*



## Vérifications : données abbérantes 

On utilise la commande summary afin de voir s'il y a des données abbérantes dans chacune des colones. On précise la colone à inspecter grâce à la commande $.

```{r}
summary(data_exo6$Smoker)
summary(data_exo6$Status)
summary(data_exo6$Age)


```

Il n'y a pas de données abbérantes, cad, des réponses autre que No ou Yes, Alive ou Dead et concernant l'âge pas d'âge < 0 et > 100 ans. 


## Vérificaitons : données manquantes 


```{r}
na_records = apply(data_exo6, 1, function (x) any(is.na(x)))
data_exo6[na_records,]
```

Il n'y a aucune données manquantes 



# Exercices 


## Partie 1 

### **Consgines : ** Représentez dans un tableau le nombre total de femmes vivantes et décédées sur la période en fonction de leur habitude de tabagisme. 

On utilise la fonction table avec les variables demandées : ici Smokers et status


```{r}
tableau_mortalite= table(data_exo6$Smoker,data_1$Status)
tableau_mortalite
```


### **Consignes :**Calculez dans chaque groupe (fumeuses / non fumeuses) le taux de mortalité (le rapport entre le nombre de femmes décédées dans un groupe et le nombre total de femmes dans ce groupe). 

On utilise la fonciton prop.table pour avoir les taux de mortalité dans chaque groupe. 
La fonction round permet d'avoir un arrondi du résultat


```{r}
round(((prop.table((tableau_mortalite)))*100),2)
```


On obtient alors un pourcentage de mortalité de : 

- 17.50 % : Groupe Non Fumeuse

- 10.58 % : groupe Fumeuse 

Afin de ne pas utiliser manuellement les données écrites ci-dessus, nous allons les conserver des vecteurs

```{r}
tableau = round(((prop.table(table(data_exo6$Smoker,data_exo6$Status)))*100),2)

#creation d'un data frame à partir du tableau

export_data_fram = as.data.frame(tableau)

# selection des colonnes et variables qui nous interesse et creation d'un nouveau data fram

new_data_frame = as.data.frame(export_data_fram[3:4,3])

# attribution des valeurs de mortalite pour les groupes à des vecteurs 

mortalite_non_fum = new_data_frame[1,]

mortalite_fumeuse = new_data_frame[2,]

```



###**Consignes : ** Vous pourrez proposer une représentation graphique de ces données et calculer des intervalles de confiance si vous le souhaitez. 



```{r}
plot(data_exo6$Smoker,data_exo6$Status, main = "Statut vital en fonction du statut tabagique", xlab= "Statut Tabagique", ylab="Statut Vital")


```


Pour le calcul des intervals de confiance nous allons avoir besoin d'utiliser le package binom.

```{r}
install.packages("binom")
library(binom)
```

Puis on utilise la fonction binom.confint en precisant la proportion de mortalité dans le groupe fumeuse et ensuite non fumeuse en fonction de la population totale 
```{r}
# Tout d'abord je vais attribuer à la variable N le nombre total de femme dans la population 

N = length(data_exo6$Age)

# On utilise la fonction binom.confint, et ici on choisit la methode exact

IC_mort_fume = binom.confint(mortalite_fumeuse, N, methods = "exact")

IC_mort_non_fum = binom.confint(mortalite_non_fum, N, methods = "exact")


IC_mort_fume
IC_mort_non_fum

```




### **Consignes : **En quoi ce résultat est-il surprenant ?

Ces données parraissent surprennantes puisqu'elles sont contradictoires avec la [littérature](https://www.em-consulte.com/rmr/article/157296). Les femmes qui fumaîent semblent avoir un taux de mortalité amoindrie en comparaison avec celles qui ne fumaient pas. 







## Partie 2 


### **Consignes : ** Reprenez la question 1 (effectifs et taux de mortalité) en rajoutant une nouvelle catégorie liée à la classe d'âge. On considérera par exemple les classes suivantes : 18-34 ans, 34-54 ans, 55-64 ans, plus de 65 ans.



Nous allons créer une nouvelle variable que l'on nomme age_classe correspondant aux différentes catégories d'âges. 


```{r}
# D'abord recuperer l'âge maximal sur la variable Age

max = max(data_exo6$Age)

data_exo6$Age_classe = cut(data_exo6$Age, c(18,34.01,54,64,max),right = TRUE,
                           include.lowest = TRUE,labels = c("18-34 ans", "34-54 ans", "55-64 ans", "> 65 ans"))

table(data_exo6$Age_classe)

```



Puis nous allons créer un nouveau tableau en rentrant la variable Age_classe


```{r}
new_tableau_mortalite= table(data_exo6$Smoker,data_1$Status,data_exo6$Age_classe)
new_tableau_mortalite

```


Puis nous allons regarder le pourcentage de décès dans chacun des groupes

```{r}
round(((prop.table((new_tableau_mortalite)))*100),2)


```


Réalisation d'une figure de la mortalité selon l'âge et le statut tabagique


Pour cela nous utilison le package ggplot2


```{r}
library(ggplot2)


d<-data.frame(Age=c("18-34 ans", "34-54 ans", "55-64 ans", "> 65 ans"),Non_Fum=c(0.46,1.45,3.04,12.56),Fum=c(0.38,3.12,3.88,3.20))
d2<-data.frame(Age=rep(d$Age,2),C=c(d$Non_Fum,d$Fum),type=c(rep("Non_Fum",nrow(d)),rep("Fum",nrow(d))))
ggplot(data=d2,aes(x=Age,y=C,fill=type))+geom_bar(binwidth=0.5,position="dodge")

```



###En quoi ce résultat est-il surprenant ? Arrivez-vous à expliquer ce paradoxe ? De même, vous pourrez proposer une représentation graphique de ces données pour étayer vos explications.


On remarque que si l'on raissone par classe d'âge on obtient des différences de mortalité entre les groupes fumeuses et non fumeuses. Cependant, on observe que les femmes non fumeuses décèdent **3 fois plus** dans la tranche d'âge 64 et plus. 

Cela pourrait s'expliquer par une moyenne d'âge plus élevée

Nous pouvons le vérifier en :

- calculant la moyenne d'âge dans le groupe fumeuse et non fumeuse avec la fonction tapply
- realiser un boxplot avec la fonction plot
- réaliser un test t de student en supposant que la normalité et la variance sont respectée 

```{r}
tapply(data_exo6$Age,data_exo6$Smoker, mean)

plot(data_exo6$Smoker,data_exo6$Age)

# Verfication normalité distribution et normalité de la variance
hist(data_exo6$Age)
by(data_exo6$Age,data_exo6$Smoker, sd, na.rm = TRUE)

t.test(data_exo6$Age~data_exo6$Smoker, var.equal=TRUE)
```

Il y a une différence significative d'âge entre les deux groupes. Le groupe non fumeuse est significativement plus âgée. De plus il y a un pourcentage de femme non fumeuse importante et par conséquent cela renforce les chiffres de mortalité à la fin. Plus de femme âgée et non fumeuse en comparaison avec peu de femme fumeuse et moins jeunes.



## Partie 3 #####


### **Conignes :** Afin d'éviter un biais induit par des regroupements en tranches d'âges arbitraires et non régulières, il est envisageable d'essayer de réaliser une régression logistique. Si on introduit une variable Death valant 1 ou 0 pour indiquer si l'individu est décédé durant la période de 20 ans, on peut étudier le modèle Death ~ Age pour étudier la probabilité de décès en fonction de l'âge selon que l'on considère le groupe des fumeuses ou des non fumeuses. Ces régressions vous permettent-elles de conclure sur la nocivité du tabagisme ? Vous pourrez proposer une représentation graphique de ces régressions (en n'omettant pas les régions de confiance).



On commence par créer une nouvelle variable Death avec la fonction recode 
Alive = 0 
Dead = 1 
et qui devra être numérique

```{r}
data_exo6$Death <- recode(data_exo6$Status,
               "Alive" = "0",
               "Dead" = "1")
data_exo6$Death <- as.numeric(as.character(data_exo6$Death))


class(data_exo6$Death)
```


Pour utiliser la variable smokers dans 


```{r}
data_exo6$Smoke.unit <- recode(data_exo6$Smoker,
               "Yes" = "1",
               "No" = "0")
data_exo6$Smoke.unit <- as.numeric(as.character(data_exo6$Smoke.unit))

class(data_exo6$Smoke.unit)
```



Realisation de la regression logistique avec la fonction glm



```{r}

#on va utiliser utiliser une regression logistique multiple avec fonction glm

mod1 = glm(Death~Age, data= data_exo6, family="binomial")
summary(mod1)

mod2= glm(Death~Age+Smoke.unit, data= data_exo6, family="binomial")
summary(mod2)

```