Prérequis : calcul de moyennes et de ratios, techniques de présentations graphiques simples, éventuellement régression logistique

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.

Les données sont disponibles dans ce fichier CSV. Vous trouverez sur chaque ligne si la personne fume ou non, si elle est vivante ou décédée au moment de la seconde étude, et son âge lors du premier sondage.

Cet exercice peut être réalisé indifféremment en R ou en Python. Votre mission si vous l’acceptez :

Partie 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 ?

Partie 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.

Partie 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).

Partie 4 :

Déposez votre étude dans FUN.

Partie 1

Dans un premier temps, il s’agit d’importer les données a analyser, de representer dans un tableau le nombre total de femmes vivantes et décédée sur la période en fonction de leur habitude de tabagisme. Ensuite, nous calculerons le taux de mortalite dans chaque groupe (fumeuses / non fumeuses), nous representerons ces donnees en utilisant un graphique adapte pour finir par commenter ces resultats.

Pakage pouvant etre utile

library(epiDisplay)

## Loading required package: foreign

## Loading required package: survival

## Loading required package: MASS

## Loading required package: nnet

library(epiR)

## Package epiR 1.0-14 is loaded

## Type help(epi.about) for summary information

## Type browseVignettes(package = 'epiR') to learn how to use epiR for applied epidemiological analyses

##

library(prettyR)
library(knitr)
library(kableExtra)
library(ggplot2)

## 
## Attaching package: 'ggplot2'

## The following object is masked from 'package:epiDisplay':
## 
##     alpha

library(mfp)
# Ajouter des packages ICI

Importer les donnees

url <- "https://gitlab.inria.fr/learninglab/mooc-rr/mooc-rr-ressources/-/raw/master/module3/Practical_session/Subject6_smoking.csv?inline=false"
data <- read.csv(url)

Premier apercu du jeu de donne :

head(data)

##   Smoker Status  Age
## 1    Yes  Alive 21.0
## 2    Yes  Alive 19.3
## 3     No   Dead 57.5
## 4     No  Alive 47.1
## 5    Yes  Alive 81.4
## 6     No  Alive 36.8

Colonne 1 represente le statut : fumeur / non fumer Colonne 2 represente le statut : mort / vivant Colonne 3 reoresente l ’age.

Representation graphique

Nous allons calculer les pourcentages de femmes vivantes / mortes en fonction de leur status fumeuse / non fumeuse. Ceci nous permettra d avoir une idee de la distributipon des differentes populations.

mort_fum <- sum(data$Status[data$Smoker=="Yes"]=="Dead")
mort_nonfum <- sum(data$Status[data$Smoker=="No"]=="Dead")
vivant_fum <-sum(data$Status[data$Smoker=="Yes"]=="Alive")
vivant_nonfum <- sum(data$Status[data$Smoker=="No"]=="Alive")
nb_fum <- sum(data$Smoker=="Yes")
nb_nonfum <- sum(data$Smoker=="No")
nb_viv <- sum(data$Status=="Alive")
nb_deces <-sum(data$Status=="Dead")
tot <- sum(data$Status=="Alive"|data$Status=="Dead")

tab_2 <- data.frame(cbind(rbind(mort_fum, vivant_fum, nb_fum)), 
                    rbind(mort_nonfum,vivant_nonfum, nb_nonfum),
                    rbind(nb_deces,nb_viv,tot))

row.names(tab_2) <- c("Decedees", "Vivantes","Total")
names(tab_2) <- c("Fumeuses", "Non Fumeuses", "Total")

kable(tab_2,align = "c")

	Fumeuses	Non Fumeuses	Total
Decedees	139	230	369
Vivantes	443	502	945
Total	582	732	1314

tab_3 <- data.frame(cbind(rbind((sum(data$Status[data$Smoker=="Yes"]=="Dead")/
                                   sum(data$Smoker=="Yes"))*100,
                                (sum(data$Status[data$Smoker=="Yes"]=="Alive")/
                                   sum(data$Smoker=="Yes"))*100)), 
                    rbind((sum(data$Status[data$Smoker=="No"]=="Dead")/
                             sum(data$Smoker=="No"))*100,
                          (sum(data$Status[data$Smoker=="No"]=="Alive")/
                             sum(data$Smoker=="No"))*100))
names(tab_3) <- c("Fumeuses", "Non fumeuses")
row.names(tab_3) <- c("Decedees", "Vivantes")
kable(tab_3, align = "c")

	Fumeuses	Non fumeuses
Decedees	23.88316	31.42076
Vivantes	76.11684	68.57923

barplot(as.matrix(tab_3), col = c("black","white"), space = 1.5,width = 0.5) 
legend("center",xpd=NA, legend = c("Decedees", "Vivantes"), fill = c("black","white"))

On voit que, le proportion decedée chez les non fumeuses est superieur a la proportion de femme décédées chez les fumeuses. Ceci est-il significatif ? Pour verifier cela nous allons calculer les intervalles de confiance (IC).

Calcul des intervalles de confiance des proportions inspire de lien

#  Fonction Intervalle de confiance :

IC = function(x) {
  y <- x/100
  inf = y-(1.96*sqrt((y*(1-y))/n))
  sup = y+(1.96*sqrt((y*(1-y))/n))
  print(c(x,inf*100,sup*100))
}

Chez les fumeuses :

x <- tab_3[1,1]
n <- as.numeric(tab_2[3,1])

IC(x)

## [1] 23.88316 20.41914 27.34719

Chez les non fumeuses :

x <- tab_3[1,2];
n <- tab_2[3,2];

IC(x)

## [1] 31.42077 28.05793 34.78360

L’intervalle de confiance a 95% nous montre que la proportion de femmes fumeuses decedes est compris entre [20.4;27.34] alors que la proportion de femmes non fumeuses decedees est compris entre [28.05;34.74]. A premiere vu, on voit qu il y aurait plus de deces chez les femmes non fumeuses que chez les femmes fumeuses dans la periode a l etude ce qui est a l encontre des idees recues.

Partie 2

Dans la partie 2 de l analyse il s’agit d ajouter une cagorie liee a l age. # colonne 3 du tableau.

Dans un premier temps je vais ajouter une colonne au tableau pour ranger chaque age dans une classe d’age 18-34 ans, 34-54 ans, 55-64 ans, plus de 65 ans.

data$cl_age <- as.factor(ifelse(data$Age>=18&data$Age<35,"[18;34]",
                                ifelse(data$Age>=35&data$Age<54,"[35;54]",
                                       ifelse(data$Age>=55&data$Age<65,"[55;64]","[65;["))))
graphics.off()

tab1(data$cl_age);

## data$cl_age : 
##         Frequency Percent Cum. percent
## [18;34]       416    31.7         31.7
## [35;54]       420    32.0         63.6
## [55;64]       236    18.0         81.6
## [65;[         242    18.4        100.0
##   Total      1314   100.0        100.0

Verification

head(data)

##   Smoker Status  Age  cl_age
## 1    Yes  Alive 21.0 [18;34]
## 2    Yes  Alive 19.3 [18;34]
## 3     No   Dead 57.5 [55;64]
## 4     No  Alive 47.1 [35;54]
## 5    Yes  Alive 81.4   [65;[
## 6     No  Alive 36.8 [35;54]

tail(data)

##      Smoker Status  Age  cl_age
## 1309     No  Alive 42.1 [35;54]
## 1310    Yes  Alive 35.9 [35;54]
## 1311     No  Alive 22.3 [18;34]
## 1312    Yes   Dead 62.1 [55;64]
## 1313     No   Dead 88.6   [65;[
## 1314     No  Alive 39.1 [35;54]

J ai donc ajoute a droite du tableau une categorie cl_age qui range les ages de la colonne 3 dans differentes classe (18-34 ans, 34-54 ans, 55-64 ans, plus de 65 ans.).

tab_5 <- ftable(data[,c(1,2,4)]);
tab_5 <- round(prop.table(tab_5,2)*100,1);
tab_5 <- rbind(tab_5, tab1(data$cl_age)$output.table[,1]);

## Warning in rbind(tab_5, tab1(data$cl_age)$output.table[, 1]): number of columns
## of result is not a multiple of vector length (arg 2)

tab_5 <- as.table(tab_5);

row.names(tab_5) <- c("Non fumeuses vivantes ", "Non fumeuses decedees", 
                      "Fumeuses vivantes", "Fumeuses decedees" ,"Effectifs");
colnames(tab_5) <- c("[18;34]", "[35;54]", "[55;64]", "[65;[");

par(mar = c(5,5,5,13))
barplot(tab_5[c(1:4),], col = c("gray","black","gray","black"), density = c(30,30,40,100),
        angle = c(70,70,0,0))
xmin <- par("usr")[1]
xmax <- par("usr")[2]
ymin <- par("usr")[3]
ymax <- par("usr")[4]
par(xpd=TRUE)
lambda <- 0.025
legend(((1 + lambda) * par("usr")[2] - lambda * par("usr")[1]),50, 
       legend = c("Non fumeuses vivantes", "Non fumeuses decedees", 
                  "Fumeuses vivantes", "Fumeuses decedees"), 
       fill = c("gray","black","gray","black"),density = c(30,30,50,100),angle = (c(70,70,0,0)))

La representation ci contre nous permet de representer la repartion des effectifs dans les differentes classes d’age. Il apparait maintenant plus clair qu il y a une biais de distribution des effectifs en fonction de leur tranche d’age. En effet on remarque que les femmes non fumeuse decedee on un effectifs plus grand dans la categorie [65;[ par comparaison au femme fumeuse decedees de la meme categorie. Il va donc de soit de penser que le biais de repartition dans les differentes classe d’age peu avoir un impact sur le taux de femme decedee en fonction de leur statut fumeuse/non fumeuse. C’est le paradoxe des simpsons.

Partie 3

Dans un premier temps on va recoder la variable deces 0 si non decede et 1 si decede.

data$deces <- as.factor(ifelse(data$Status=="Alive",0,1))
head(data)

##   Smoker Status  Age  cl_age deces
## 1    Yes  Alive 21.0 [18;34]     0
## 2    Yes  Alive 19.3 [18;34]     0
## 3     No   Dead 57.5 [55;64]     1
## 4     No  Alive 47.1 [35;54]     0
## 5    Yes  Alive 81.4   [65;[     0
## 6     No  Alive 36.8 [35;54]     0

On voit donc que la variable presente colonne 5 est bien recodee a la norme voulue. On peut faire une premiere regression pour comparer statut fumeur/ non fumeur et statut dece / vivant.

#Age vs statut decede
reg_log_total <- ggplot(data, aes(x=Age,y=deces)) +
  geom_point() +
  geom_smooth()

#Age vs statut fumeur
reg_log_fumeur <- ggplot(data[data$Smoker == "Yes",], aes(x=Age,y=deces)) +
  geom_point() +
  geom_smooth()

#Age vs statut NON fumeur
reg_log_non_fumeur <- ggplot(data[data$Smoker == "No",], aes(x=Age,y=deces)) +
  geom_point() +
  geom_smooth()

reg_log_total

## `geom_smooth()` using method = 'loess' and formula 'y ~ x'

reg_log_fumeur

## `geom_smooth()` using method = 'loess' and formula 'y ~ x'

reg_log_non_fumeur

## `geom_smooth()` using method = 'loess' and formula 'y ~ x'

Les ci-joint montre une repartition semblable entre reg_log_total et reg_log_non_fumeur. Cela peut s’expliquer par une faible population de personne agee chez les fumeurs (cf question 2). De ce fait l’augmentation de la mortalite observe chez les non fumeurs peut en parti etre explique par le biais de repartition des ages chez les deux groupes fumeur et non fumeur… Il serait donc necessaire de refaire une analyse chez des populations plus homogenes …

Sujet 6 : Autour du Paradoxe de Simpson

Olivier Messina

25/04/2020