diff --git a/module3/exo3/exercice_fr.Rmd b/module3/exo3/exercice_fr.Rmd index 6a18651cee02876d2d304c209a7f8be2885ca5d5..a387fe084536a2195b0d6b97c47ce0ef359f0b13 100644 --- a/module3/exo3/exercice_fr.Rmd +++ b/module3/exo3/exercice_fr.Rmd @@ -1,8 +1,10 @@ --- -title: "Sujet 6 : Autour du Paradoxe de Simpson" +title: 'Sujet 6 : Autour du Paradoxe de Simpson' author: "Clair Ch" date: "20 avril 2020" -output: html_document +output: + pdf_document: default + html_document: default --- @@ -10,7 +12,7 @@ output: html_document knitr::opts_chunk$set(echo = TRUE) ``` -# Instructions +## Instructions Voici les instructions du Sujet 6 : Autour du Paradoxe de Simpson @@ -28,16 +30,19 @@ Voici les instructions du Sujet 6 : Autour du Paradoxe de Simpson 4. Déposez votre étude dans FUN -# Préparation des données +## Préparation des données -## Téléchargement +#### Téléchargement Les données sont disponibles sur le [Gitlab](https://gitlab.inria.fr/learninglab/mooc-rr/mooc-rr-ressources/blob/master/module3/Practical_session/Subject6_smoking.csv) du MOOC ["Recherche Reproductible : principes méthodologiques pour une science transparente"](https://www.fun-mooc.fr/courses/course-v1:inria+41016+session01bis/about) de l'Inria. On peut les récupérer au format .csv à cette adresse : ```{r} -data_url<-"https://gitlab.inria.fr/learninglab/mooc-rr/mooc-rr-ressources/-/raw/master/module3/Practical_session/Subject6_smoking.csv?inline=false" +data_url<-"https://gitlab.inria.fr/learninglab/mooc-rr/mooc-rr-ressources/-/raw/master/ +module3/Practical_session/Subject6_smoking.csv?inline=false" ``` +Elles regroupent, par femme, les habitudes de tabagisme (fumeuse/non fumeuse), le statut lors de la deuxième étude de 1995 (vivante ou morte) et l'âge lors de la première étude (effectuée entre 1972 et 1974). + Pour nous protéger contre une éventuelle disparition ou modification du serveur du Gitlab du MOOC, nous faisons une copie locale de ce jeux de données que nous préservons avec notre analyse. Il est inutile et même risqué 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. @@ -48,7 +53,8 @@ if (!file.exists(data_file)) { } ``` -## Lecture +#### Lecture + ```{r} data<-read.csv(data_file) @@ -63,6 +69,7 @@ tail(data) Nous avons les données pour 1314 femmes, avec dans l'ordre : si la personne fume ou non (colonne `Smoker` : `Yes` ou `No`), si elle est morte ou vivante au moment de la deuxième étude en 1995 (colonne `Status` : `Dead` ou `Alive`) et son âge lors de la première étude faite entre 1972 et 1974 (colonne `Age`). + Y a-t'il des points manquants dans nos données ? ```{r} @@ -72,6 +79,7 @@ data[na_records,] Aucune donnée manquante, parfait ! + Vérifions la classe de chaque colonne ```{r} @@ -79,15 +87,15 @@ class(data$Smoker) class(data$Status) class(data$Age) ``` -La données de la colonne `Age` sont bien de classe numeric -Les données des colonnes `Smoker`et `Status`sont de classe factor, ça posera peut-être problème par la suite ? + +Les données de la colonne `Age` sont bien de classe numeric, les autres colonnes de type factor. -# Décès et tabagisme +## Décès et tabagisme Nous voulons étudier le taux de mortalité pour chaque groupe (fumeuses/non fumeuses) -## Inspection des données +#### Inspection des données Faisons un tableau représentant le nombre de femmes par habitude de tabagisme (fumeuses/non fumeuses) et statut (vivantes/mortes) @@ -96,12 +104,14 @@ t<-table(data$Smoker,data$Status) t ``` + Vérifions que la somme de chaque catégorie correspond bien au nombre total de femmes de l'étude ```{r} sum(t)==nrow(data) ``` + Nous pouvons voir qu'il y a plus de femmes vivantes que de femmes mortes et plus de femmes qui ne fument pas : ```{r} @@ -110,6 +120,7 @@ alive dead<-sum(t[,2]) #femmes mortes dead ``` + Il y a également plus de femmes non fumeuses que fumeuses : ```{r} @@ -119,56 +130,388 @@ smoke<-sum(t[2,]) #femmes fumeuses smoke ``` -## Calcul du taux de mortalité +#### Calcul du taux de mortalité Calculons dans chaque groupe (fumeuses/non fumeuses) le taux de mortalité (rapport entre le nombre de femmes décédées dans un groupe et le nombre total de femmes dans ce groupe) taux de mortalité pour les femmes non fumeuses : + ```{r} nsmoke_morta<-t[1,2]/nsmoke nsmoke_morta ``` + taux de mortalité pour les femmes fumeuses : + ```{r} smoke_morta<-t[2,2]/smoke smoke_morta ``` + Les femmes non fumeuses présentent un taux de mortalité plus élevé ! + De combien est-il plus élevé ? ```{r} nsmoke_morta/smoke_morta ``` -## Calcul de l'intervalle de confiance à 95% du taux de mortalité +Arrondi au centième, les femmes non fumeuses présent donc un taux de mortalité `r round(nsmoke_morta/smoke_morta,digit=2)` fois plus élevé que les femmes fumeuses dans cette étude. C'est surprenant car contre-intruitif ! + + +On peut faire un test statistique pour avoir une idée si cette différence n'est pas due au hasard (test Chi-square de Pearson sur les proportions) : + +```{r} +chisq.test(t) +``` + +La p-value est en dessous de 0.05, il y a donc moins de 5% de chance que cette différence du taux de mortalité soit due au hasard. Mais est-elle vraiment due aux habitudes de tabagisme ? + +#### Calcul de l'intervalle de confiance à 95% du taux de mortalité Pour calculer les intervalles de confiance à 95% du taux de mortalité par catégorie (fumeuses/non fumeuses), nous allons suivre les instructions de l'article [How to Determine the Confidence Interval for a Population Proportion](https://www.dummies.com/education/math/statistics/how-to-determine-the-confidence-interval-for-a-population-proportion/) du site [dummies.com](https://www.dummmies.com) -mettre formule mathématique +La formule donnée pour calculé l'intervalle de confiance à 95% est la suivante : +$$1,96\sqrt{\frac{\beta(1-\beta)}{n}}$$ +$\beta$ représente ici le taux de mortalité des femmes fumeuses ou non fumeuses et $n$ le nombre total de femmes fumeuses ou non fumeuses + +Créons la fonction correspondante : + +```{r} +CI95<-function(morta,n){ + 1.96*sqrt((morta*(1-morta))/n) +} +``` + + +Calculons les intervalles de confiance pour les femmes non fumeuses et les femmes fumeuses : + +```{r} +nsmoke_CI<-CI95(nsmoke_morta,nsmoke) +smoke_CI<-CI95(smoke_morta,smoke) +``` + +Arrondi au centième, les femmes fumeuses ont un taux de mortalité de `r round(smoke_morta,digit=2)`$\pm$ `r round(smoke_CI,digit=2)` + +Arrondi au centième, les femmes non fumeuses ont un taux de mortalité de `r round(nsmoke_morta,digit=2)`$\pm$ `r round(nsmoke_CI,digit=2)` + +#### Représentation graphique + +Pour faire les graphiques, nous utiliserons la librairie `ggplot2`. +Ce code permet d'installer la libraire si nécessaire et de la loader. + +```{r} +if(!require(ggplot2)){ + install.packages("ggplot2") + library(ggplot2) +} +``` + +Pour faire le graphe, il faut d'abord mettre les données de taux de mortalité et d'intervalle de confiance dans une data frame, j'en profite pour les passer en pourcentage : + +```{r} +df_morta<-data.frame(Smoker=c("No","Yes"), Mortality=c(nsmoke_morta*100,smoke_morta*100), CI95=c(nsmoke_CI*100,smoke_CI*100)) +df_morta +``` + +Avec l'aide d'un [tutoriel pour faire +des barplots](http://www.sthda.com/french/wiki/ggplot2-barplots-guide-de-demarrage-rapide-logiciel-r-et-visualisation-de-donnees#barplot-du-comptage-des-observations), faisons un graphe représentant le taux de mortalité en fonction des habitudes de tabagisme : + +```{r} +plot_morta<-ggplot(df_morta,aes(x=Smoker,y=Mortality))+ + geom_bar(stat="identity",fill="steelblue",width=0.5)+ + theme_minimal()+ + ylim(c(0,50))+ #limite de l'axe y + geom_errorbar(aes(ymin=Mortality-CI95,ymax=Mortality+CI95),width=.2,position=position_dodge(.9))+ #barres d'erreurs + labs(title="Taux de Mortalité en fonction des habitudes de tabagisme", subtitle="Taux de Mortalité en pourcentage +/- intervalle de confiance à 95%") +plot_morta +``` + +Les femmes non fumeuses semblent avoir un taux de mortalité plus élevé que les femmes fumeuses. Si on considère ce résultat comme tel, on serait tenté de conclure que la cigarette est bonne pour la santé ! + +## Décés et Tabagisme en fonction de la classe d'âge + +#### Ajout de la tranche d'âge + +Rajoutons une colonne `Tranche` dans le fichier `data`donnant la tranche d'âge (`18-34`: de 18 ans inclus à 35 ans non inclus ; `35-54`: de 35 ans inclus à 55 ans non inclus ; `55-64`: de 55 ans inclus à 65 non inclus ; `>65`: plus de 65 ans inclus), mettons le message `ERROR`si une tranche d'âge n'est pas trouvée + +```{r} +for(i in 1:nrow(data)){ + if(data$Age[i]>=65){data$Tranche[i]<-">65" + }else if(data$Age[i]>=55){data$Tranche[i]<-"55-64" + }else if(data$Age[i]>=35){data$Tranche[i]<-"35-54" + }else if(data$Age[i]>=18){data$Tranche[i]<-"18-34" + }else(data$Tranche[i]<-"ERROR") +} +``` + + +Vérifions que ça ait bien marché : + +```{r} +head(data) +tail(data) +``` + + +Est-ce qu'il y a eu des erreurs ? + +```{r} +sum(data$Tranche=="ERROR") +``` + +Aucune erreur, tout va bien ! + +#### Calcul du taux de mortalité par tranche d'âge et tabagisme + +Faisons une fonction permettant de calculer le nombre de femmes par tranche d'âge selon leur habitude de tabagisme (fumeuse/non fumeuse) et éventuellement leur statut (dead/alive) + +```{r} +comptage<-function(smoker,statut,tranche){ + if(statut=="NULL") + {sum(data$Smoker==smoker&data$Tranche==tranche) + }else{sum(data$Smoker==smoker&data$Status==statut&data$Tranche==tranche)} +} +``` + + +Calculons le nombre de femmes totales et mortes par tranche d'âge et catégories + +```{r} +tranches<-c("18-34","35-54","55-64",">65") + +nsmoke_age<-mapply(comptage,smoker="No",statut="NULL",tranche=tranches) #femmes non fumeuses +smoke_age<-mapply(comptage,smoker="Yes",statut="NULL",tranche=tranches) #femmes fumeuses +nsmoke_age_dead<-mapply(comptage,smoker="No",statut="Dead",tranche=tranches) #femmes non fumeuses mortes +smoke_age_dead<-mapply(comptage,smoker="Yes",statut="Dead",tranche=tranches) #femmes fumeuses mortes + +#nommer les tranches d'âge dans les vecteurs : +names(nsmoke_age)<-tranches +names(smoke_age)<-tranches +names(nsmoke_age_dead)<-tranches +names(smoke_age_dead)<-tranches +``` + + +Vérifions que les nombres de femmes fumeuses et non fumeuses correspondent bien au nombre de femmes totales + +```{r} +sum(smoke_age,nsmoke_age)==nrow(data) +``` + +Vérifions également que le nombre de femmes mortes fumeuses et non fumeuses correspondent bien au nombre de femmes totales mortes : + +```{r} +sum(smoke_age_dead,nsmoke_age_dead)==nrow(data$Status=="Dead") +``` + + +Nous pouvons maintenant calculer le taux de mortalité pour les femmes non fumeuses par tranche d'âge + +```{r} +nsmoke_age_morta<-nsmoke_age_dead/nsmoke_age +nsmoke_age_morta +``` + +Et pour les femmes fumeuses par tranche d'âge : + +```{r} +smoke_age_morta<-smoke_age_dead/smoke_age +smoke_age_morta +``` + + +Mettons les taux de mortalité en pourcentage au centième près selon les tranches d'âges et habitude de tabagisme dans une data frame pour pouvoir les comparer : + +```{r} +df_morta_age<-data.frame(Age=tranches,"Mortality non Smoker"=round(nsmoke_age_morta*100,digit=2),"Mortality Smoker"=round(smoke_age_morta*100,digit=2)) +df_morta_age +``` + +On peut voir tout d'abord que de façon générale, le taux de mortalité est plus important pour les catégories d'âges les plus "avancées", ce qui est un résultat attendu. Ici la mortalité des femmes qui fument semble plus importante que celle des femmes qui ne fument pas dans toutes les tranches d'âge (à part celle des femmes de plus de 65 ans). On a donc une différence quand on regarde le taux de mortalité total et le taux de mortalité par tranche d'âge. Le taux de mortalité plus élevé chez les femmes non fumeuses n'était peut-être pas du au tabagisme... + +#### Calcul de l'intervalle de confiance à 95% + +De la même manière que pour le taux de mortalité global par tranche d'âge, nous calculons l'intervalle de confiance à 95% en prenant en compte la tranche d'âge : + +```{r} +nsmoke_age_IC<-mapply(CI95, morta=nsmoke_age_morta,n=nsmoke_age) +smoke_age_IC<-mapply(CI95,morta=smoke_age_morta,n=nsmoke_age) +``` + + +#### Graphe + +Mettons au propres les données pour faire le graphe: la 1ère colonne repésente l'habitude de tabagisme, la deuxième la tranche d'âge, la troisième le taux de mortalité (en pourcentage au centième près), la quatrième l'intervalle de confiance à 95% (en pourcentage au centième près) + +```{r} +df_morta_age2<-data.frame(Smoker=c(rep("No",4),rep("Yes",4)),Age=rep(c("18-34","35-54","55-64",">65"),2), Mortality=c(nsmoke_age_morta,smoke_age_morta),CI95=c(nsmoke_age_IC,smoke_age_IC)) +df_morta_age2$CI95<-round(df_morta_age2$CI95*100,digits=2) +df_morta_age2$Mortality<-round(df_morta_age2$Mortality*100,digits=2) +df_morta_age2 +``` + +Pour que les tranches d'âges soient dans le bon ordre dans le graphe, nous ordonnons les niveaux de la colonne `Age` + +```{r} +df_morta_age2$Age<-ordered(df_morta_age2$Age, levels=c("18-34","35-54","55-64",">65")) +``` + + +Faisons maintenant le graphe, toujours grâce à la libraire `ggplot2` : + +```{r} +if(!require(ggplot2)){ + install.packages("ggplot2") + library(ggplot2) +} +``` + + +```{r} +plot_morta_age<-ggplot(df_morta_age2,aes(x=Age,y=Mortality,fill=Smoker))+ + geom_bar(stat="identity",position="dodge")+ + scale_fill_brewer(palette="Paired")+ + theme_minimal()+ + geom_errorbar(aes(ymin=Mortality-CI95, ymax=Mortality+CI95), width=.2, position=position_dodge(.9))+ #barres d'erreur + labs(title="Taux de Mortalité par catégorie d'âge en fonction des habitudes de tabagisme", subtitle="Taux de Mortalité en pourcentage +/- intervalle de confiance à 95%") +plot_morta_age +``` + +Ici il apparait que la mortalité est plus importante pour les femmes fumeuses que non fumeuses, la cigarette ne semble au final pas bon pour la santé ! + +#### Paradoxe + +Le paradoxe entre le premier résultat (si on regarde le taux de mortalité de façon globale, il est plus élevé pour les femmes non fumeuses) et ce dernier (si on classe les femmes par catégories d'âges, on a un taux de mortalité plus élevé pour les femmes fumeuses) pourrait provenir du fait que les catégories d'âges ne sont pas représentées de manière égales entre les femmes fumeuses et non fumeuses : les femmes "jeunes" auraient plus tendance à fumer, et au contraire les femmes "âgées" auraient tendance à moins fumer et seraient plus représentées dans le groupe des non fumeuses. Ainsi le taux de mortalité "naturelle" (de vieillesse) est globalement plus important pour le groupe des non fumeuses. + +Nous pouvons vérifier cette hypothèse en analysant la distribution des âges par habitude de tabagisme : + +```{r} +summary(data$Age[data$Smoker=="No"]) +summary(data$Age[data$Smoker=="Yes"]) +``` + +On voit que l'étendu des âges est sensiblement la même entre les femmes non fumeuses et fumeuses (de 18 à environ 89-90 ans). Par contre l'âge médian est plus élevé chez les femmes non fumeuses que non fumeuses (48.40 ans contre 43.10 ans). C'est en concordance avec notre théorie. + + +Comparons maintenant les histogrammes (distribution des âges) pour avoir une meilleure idée de la distribution des âges : + +```{r} +par(mfrow=c(1,2)) #pour avoir les histogrammes côte à côte +hist(data$Age[data$Smoker=="No"],main="Femmes non fumeuses",xlim=c(18,90),ylim=c(0,80),xlab="Age") +hist(data$Age[data$Smoker=="Yes"], main="Femmes fumeuses",xlim=c(18,90),ylim=c(0,80),xlab="Age") +``` + +En comparant ces histogrammes, on voit bien que la distribution des âges n'est pas du tout la même entre les femmes fumeuses et non fumeuses : il y a très peu de femmes "âgées" dans les femmes fumeuses. Ainsi beaucoup moins de femmes meurent de mort "naturelle" (vieillesse) dans le groupe des femmes fumeuses, ce qui abaisse le taux de mortalité global. Pour savoir si la cigarette a en effet un effet sur le taux de mortalité, il faut prendre en compte cette variable âge pour s'affranchir du taux de mortalité du à la vieillesse. + + +## Régression logistique + +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. + +La régression logistique permet d'évaluer et caractériser les relations entre deux variables, ici l'âge et la mortalité. On peut trouver une explication dans ce [document](https://gitlab.inria.fr/learninglab/mooc-rr/mooc-rr-ressources/blob/master/module2/exo5/challenger_fr.pdf) du MOOC "Recherche reproductible" et une mise en pratique dans [celui-ci](https://app-learninglab.inria.fr/moocrr/gitlab/moocrr-session3/moocrr-reproducibility-study/blob/master/challenger.pdf) + +#### Création de la variable Death + +Ajoutons une variable `Death`prenant la valeur de 0 si la personne est décédée et 1 si la personne est vivante + +```{r} +for(i in 1:nrow(data)){ + if(data$Status[i]=="Dead") + {data$Death[i]<-0} + else(data$Death[i]<-1) +} +``` + +On vérifie que ça a bien fonctionné + +```{r} +head(data) +``` + +#### Calcul de la régression logistique + +Faisons une régression logistique de la mortalité en fonction de l'âge selon que la personne fume ou non : + +```{r} +logistic_reg<-glm(Death~Age+Smoker,data=data,family=binomial) +summary(logistic_reg) +``` + +On voit qu'il y a un effet négatif de l'âge sur la mortalité (estimate de -0.099+/-0.005) qui est signicatif (p<0.05). +Il y a aussi un effet négatif de l'habitude de tabagisme (estimate de -0.27+/-0.14) mais qui n'est pas significatif (p>0.05). + +Cette régression logistique ne montre aucune évidence d'un effet du tabagisme sur la mortalité dans ce jeu de données. + + +#### Graphes + +Pour faire le graphe, nous utiliserons la librairie `ggplot2` ```{r} -nsmoke_CI<-1.96*sqrt((nsmoke_morta*(1-nsmoke_morta))/nsmoke) -smoke_CI<-1.96*sqrt((smoke_morta*(1-smoke_morta))/smoke) +if(!require(ggplot2)){ + install.packages("ggplot2") + library(ggplot2) +} ``` +Faisons le graphe de la régression logistique de la mortalité (Death vaut 0 is la personne est décédée et 1 si la personne est morte) en fonction des habitudes de tabagisme et de l'âge, avec en grisé l'intervalle de confiance à 95% +```{r} +ggplot(data,aes(x=Age,y=Death,color=Smoker)) + geom_point(alpha=.3,size=1) + + theme_bw() + + geom_smooth(method = "glm", + method.args = list(family = "binomial"))+ + labs(title="Logistic regression") +``` + +Sur ce graphe, il ne semble pas globalement que le fait de fumer augmente ou non la probabilité d'être décedée 20 ans plus tard. Cependant la courbe des femmes non fumeuses semble faire un coude et être légèrement au dessus de celles des fumeuses pour un âge entre 25 et 65 ans. Re-faire une analyse sur cette tranche d'âge serait intéressant. Le tabagisme n'a peut-être aucun effet dans nos données sur le femmes les plus jeunes et les plus âgées (les plus âgées vont mourir de vieillesse, et peut-être que les femmes les plus jeunes ne fument pas depuis assez longtemps pour que ça ait un effet évident ou qu'elles sont plus protégées d'un effet néfaste du tabac). -## Représentation graphique +#### Régression logistique sur les femmes entre 25 et 65 ans -Pour faire les graphiques, nous utiliserons la librairie `ggplot2` +Pour cette partie nous utiliserons la librarie `dplyr ```{r} -library(ggplot2) +if(!require(dplyr)){ + install.packages("dplyr") + library(dplyr) +} ``` +Créons un deuxième jeu de données ne comprenant que les femmes entre 25 et 65 ans. +```{r} +data2<-filter(data,Age>25&Age<65) +``` -# Truc + +Faisons une régression logistique sur ces données : ```{r} -hist(data$Age) -summary(data$Age) +logistic_reg2<-glm(Death~Age+Smoker,data=data2,family=binomial) +summary(logistic_reg2) ``` +Ici à la fois l'âge et l'habitude de tabagisme ont un effet négatif significatif. Il semble donc que le fait de fumer affecte la probabilité de décéder pour les femmes qui avaient entre 25 et 65 ans. + + +Faisons le nouveau graphe : + +```{r} +ggplot(data2,aes(x=Age,y=Death,color=Smoker)) + geom_point(alpha=.3,size=1) + + theme_bw() + + geom_smooth(method = "glm", + method.args = list(family = "binomial"))+ + labs(title="Logistic regression 2") +``` + +Sur ce nouveau graphe il apparait que le fait de fumer augmente légèrement la probabilité d'être décédée 20 ans plus tard pour les personnes qui avaient entre 25 et 65 ans. + + +## Conclusion + +Les femmes fumeuses qui avaient entre 25 et 65 ans au moment de la première étude semblent avoir une probabilité plus forte de mourir que les femmes non fumeuses. Le fait de fumer augmenterait donc la probabilité de mourir dans les 20 ans si on est d'âge "jeune-moyen". Avec ce jeu de données, il semble donc que le fait de fumer pour les femmes "âgées" ou "très jeunes" n'ait pas d'incidence sur la probabilité de mourir 20 ans plus tard. On peut faire l'interprétation que si la femme est âgée, 20 ans plus tard elle va très probalement mourir qu'elle fume ou non (pour les fumeuses à cause du tabagisme et/ou vieillesse et pour les non fumeuses de vieillesse), et pour les plus jeunes soit elles ne fument pas depuis assez longtemps pour que ça ait une incidence sur leur probabilité de mourir, soit elles sont "protégées"" par leur "jeunesse". Il serait intéressant de faire une étude à intervalles de temps plus rapproché (par exemple tous les 5 ans), pour limiter l'effet de la mort naturelle sur les personnes les plus âgées et en prenant en compte depuis combien de temps la personne fume pour voir si cela a un effet sur la probabilité de mourir. + +Le tabagisme semble donc avoir un effet néfaste sur la mortalité qui est dépendant de l'âge. \ No newline at end of file