course-v1:UPSUD+42001: Intro stat avec R
David Sissoko
5 octobre 2018
- Introduction
- Chapitre1:Semaine1:Définitons
- Chapitre2:Semaine1:Représentations graphiques
- Chapitre3:Semaine1:Mesure de position et de dispersion: Principes de base
- Chapitre4: Mesure de position et de dispersion: Pratique
- Lab1:Semaine1:Gestion des données, dataframe, var numériques (quant) et cat (qual)
- Lab2:Semaine1:Indexation critériée d’observations, sélection de variables, graphiques univariés
- Semaine1:Devoir évalué par les pairs
- Chapitre5:Semaine2:Intervalles de confiance
- Chapitre6:Semaine2:Coefficient de corrélation
- Chapitre7:Semaine2:Risque relatif et odds-ratio
- Lab3:Semaine2:Langage R Markdown
- Quizz Exercice 2:Semaine2
- Quiz exercice 2bis (facultatif):Semaine2
- Chapitre8:Semaine3:Tests statistiques le p (petit p) - Test de fisher
- Chapitre9:Semaine3:Test de Neyman et Pearson (tests d’hypothèses)
- Chapitre10:Semaine3:Comparaison de 2 pourcentages
- Chapitre11:Semaine3:Comparaison de 2 moyennes
- Chapitre12:Semaine3:test de nullité d’un coef de corrélation
- Lab4:Semaine3:Tests d’association graphique bivariés
- Quiz exercice3:Semaine3
- Quiz exercice3bis (facultatif):Semaine3
- Semaine3:Devoir évalué par les pairs (2/3)
- Chapitre13:Semaine4:Régression linéaire simple
- Chapitre14:Semaine4:Régression, corrélation et test t
- Chapitre15:Semaine4:Régression linéaire multiple, analyse de variance
- variable explicative catégorielle à plus de 2 classes
- Choix de la modalité de référence
- Estimation de l’effet globale d’une variable catégorielle à plus de 2 classes
- Interaction entre 2 variables explicatives quand le modèle additif n’est plus pertinent
- Analyse de variance
- Condition de validité du modèle de régression
- Chapitre16::Semaine4:Introduction à la régression logistique
- Quiz Exercice 4:Semaine4
- Quiz Exercice 4bis (facultatif):Semaine4
- Chapitre17:Semaine4:Régression logistique multiple
- Quizz
- Lab5:Semaine4:ANOVA Régression linéaire ; régression logistique
- Quiz Exercice4:Semaine4
- Quiz Exercice4 (Facultatif):Semaine4
- Notes et Bookmarks
Introduction
Ce mooc est à l’adresse: https://www.fun-mooc.fr/courses/course-v1:UPSUD+42001+session10/info A la fin de ce document, un chapitre de note et d’adresse intéressante notamment pour l’utilisation de librairie R externes a été créé
Chapitre1:Semaine1:Définitons
| Terme | Définition |
|---|---|
| probabilité d’1 evt | Fréquence d’apparition de l’evt |
| variable aléatoire discrète | c’est 1 var qui ne prend que quelques modlaités (valeurs): par exemple: 2, 3.2, 0.05 |
| variable aléatoire discrète | c’est 1 var dont le nombre de modalité est théroriquement infini (mais peux être compris dans un intervalle) |
| Probilité d’avoir 1 nombre pair lors du tirage d’un dé équilibré à 6 faces | parmi 1,2,3,4,5,6 les modalités paires sont 2,4,6: il y a donc 3 modalités / 6 possibles = 50% de chance |
| Distribution normale | c’est la fréquence des modalités d’une variable aléatoire suite à un nombre important de tirage aléatoire de cette variable: elle a une forme de cloche |
| Variable qualitative | variable discrète dont les modalités permettent de qualifier: exemple: couleur d’un pull |
| Variable quantitative | variable discrète ou quanti qui permet de valoriser: exemple : poids ou age d’un individu |
Chapitre2:Semaine1:Représentations graphiques
Le but est de représenter sur un graphique la distribution de variables qualitatives et quantitatives. Quelques expemples de graphiques:
- les diagrammes en bâtons,
- les camemberts,
- les histogrammes,
- les boîtes à moustaches,
- les diagrammes cartésiens,
- les diagrammes en fagot.
Importer des données et positionner le répertoire de travail
Avant tout, il faut importer les données:
setwd("~/Documents/mooc-fun") #Répertoire par défaut
smp<-read.csv2("smp2.csv")
str(smp) # vérificaiton ducontenu du fichier## 'data.frame': 799 obs. of 26 variables:
## $ age : int 31 49 50 47 23 34 24 52 42 45 ...
## $ prof : Factor w/ 8 levels "agriculteur",..: 3 NA 7 6 8 6 3 2 6 6 ...
## $ duree : int 4 NA 5 NA 4 NA NA 5 4 NA ...
## $ discip : int 0 0 0 0 1 0 0 0 1 0 ...
## $ n.enfant : int 2 7 2 0 1 3 5 2 1 2 ...
## $ n.fratrie : int 4 3 2 6 6 2 3 9 12 5 ...
## $ ecole : int 1 2 2 1 1 2 1 2 1 2 ...
## $ separation : int 0 1 0 1 1 0 1 0 1 0 ...
## $ juge.enfant : int 0 0 0 0 NA 0 1 0 1 0 ...
## $ place : int 0 0 0 1 1 0 1 0 0 0 ...
## $ abus : int 0 0 0 0 0 0 0 0 1 1 ...
## $ grav.cons : int 1 2 2 1 2 1 5 1 5 5 ...
## $ dep.cons : int 0 0 0 0 1 0 1 0 1 0 ...
## $ ago.cons : int 1 0 0 0 0 0 0 0 0 0 ...
## $ ptsd.cons : int 0 0 0 0 0 0 0 0 0 0 ...
## $ alc.cons : int 0 0 0 0 0 0 0 0 1 1 ...
## $ subst.cons : int 0 0 0 0 0 0 1 0 1 0 ...
## $ scz.cons : int 0 0 0 0 0 0 0 0 0 0 ...
## $ char : int 1 1 1 1 1 1 1 1 4 1 ...
## $ rs : int 2 2 2 2 2 1 3 2 3 2 ...
## $ ed : int 1 2 3 2 2 2 3 2 3 2 ...
## $ dr : int 1 1 2 2 2 1 2 2 1 2 ...
## $ suicide.s : int 0 0 0 1 0 0 3 0 4 0 ...
## $ suicide.hr : int 0 0 0 0 0 0 1 0 1 0 ...
## $ suicide.past: int 0 0 0 0 1 0 1 0 1 0 ...
## $ dur.interv : int NA 70 NA 105 NA NA 105 84 78 60 ...A ce stade, on observe dans R, dans l’onglet Environnement le dataframe smp que l’on peux manipuler.
Diagramme baton (bar)
Le but est de représenter la distribution d’une variable: Par exemple: le nombre de détenus par métier
##Tout d'abord, on va générer une table afin d'isoler la variable métier
tablemetier<-table(smp$prof)
##Ensuite, on peux générer le barplot
barplot(tablemetier)
Diagramme camembert (pie)
pie(tablemetier)
histogramme
C’est le grand classique pour représenter la distribution d’une variable aléatoire quantitative continue
hist(smp$age, col = "red", main = "Distribution de l'age des détenus", xlab = "Age", ylab = "Fréquence")
boite à moustache (boxplot)
Utile pour les variables quantitatives
boxplot(smp$age, xlab="Age")
- à l’intérieur de la boîte, vous avez 50 % des données,
- vous avez ensuite une moustache supérieure,
- et entre le bord supérieur de la boîte et la moustache supérieure, vous avez 25 % des données,
- et entre le bord inférieur de la boîte et la moustache inférieure, vous avez aussi 25% des données.
Ce n’est pas exactement cela car on peux voir des modalités au dessus de la moustache supérieure. La moustache supérieur n’est pas égale au max des modalités mais : min(max des données, 1,5 écart-types au dessus du bord supérieur de la boîte).
L’intéret de boxplot est de représenter la distrib d’une var quanti en fonction de sous groupes. Par exemple quelle est la distrib de l’age par rapport à la consommation de substance:
boxplot(smp$age~smp$subst.cons, ylab="Age", xlab="Consommation de substance", main="Distribution de l'age par rapport à la consommation de substance") # subst.cons prend comme valeur 0 si pas de consommation et 1 si consommation
On remarque que les jeunes consomment plus que les vieux
Graphique x/y (plot)
plot(smp$age, smp$n.enfant)
On remarque que en moyenne, plus un détenu est agé, plus son nombre d’enfant est elevé
Les modalités identiques se superposent et ne sont pas visibles. Si on veux voir toutes le smodalités même les doublons, on peux décaler visuellement les valeurs identiques avec jitter
plot(jitter(smp$age, factor = 5), jitter(smp$n.enfant, factor = 5))
Données temporelles
On charge un fichier contenant des données temporelles.
## 'data.frame': 1053 obs. of 3 variables:
## $ NUMERO: int 96 96 96 96 96 96 96 96 157 157 ...
## $ VISIT : int 0 4 7 14 21 28 42 56 0 4 ...
## $ HDRS : int 34 26 12 7 5 1 1 1 27 19 ...On utilise la fonction plotmeans pour représenter l’évolution d’une variable dans le temps
## Installing package into '/home/david/R/x86_64-pc-linux-gnu-library/3.4'
## (as 'lib' is unspecified)##
## Attaching package: 'gplots'## The following object is masked from 'package:stats':
##
## lowessEnsuite on utilise plotmeans
plotmeans(repdat$HDRS~repdat$VISIT, barcol="orange")## Warning in arrows(x, li, x, pmax(y - gap, li), col = barcol, lwd = lwd, :
## zero-length arrow is of indeterminate angle and so skipped
## Warning in arrows(x, li, x, pmax(y - gap, li), col = barcol, lwd = lwd, :
## zero-length arrow is of indeterminate angle and so skipped## Warning in arrows(x, ui, x, pmin(y + gap, ui), col = barcol, lwd = lwd, :
## zero-length arrow is of indeterminate angle and so skipped
## Warning in arrows(x, ui, x, pmin(y + gap, ui), col = barcol, lwd = lwd, :
## zero-length arrow is of indeterminate angle and so skipped
On constate que l’état de santé des patients évolue eu fur et à mesure du temps
Si on veux représenter l’évolution de l’état pour chaque patient sur une seule courbe, on utilise interaction.plot. Cette fonction prend 3 paramétres: le temps (VISIT), la variable qui indique chaque sujet (NUMERO), la variable à représenter (HDRS)
interaction.plot(repdat$VISIT, repdat$NUMERO, repdat$HDRS,legend=FALSE)
legend=FALSE permet d’éviter de faire déborder la légende sur le graphique
Chapitre3:Semaine1:Mesure de position et de dispersion: Principes de base
Comment aggréger les mesures pour synthétiser ? on utilise les mesures de position et de dispersion: moyenne, médiane, écart-type Par exemple: autour de quel age se situe la population des fumeurs de canabis: Pour répondre à cette question, on peux faire une mesure de position. On a également besoin d’évaluer la dispersion: cad la répartition en fonction de sous groupe de position
Position
| Mesures de position | Mesures de dispersion |
|---|---|
| Pour une var catégorielle: il suffit de lister le pourcentage pour chaque catégorie. Pour une var quanti, on peux utiliser la moyenne et la médiane | pour var quanti et cat, on peux utiliser les quartiles ou l’écart-type |
médiane = 50% des obs sont au dessus et 50% au dessous ==> intuitif et robuste (les valeurs extremes ne modifient pas la médianep pusiqu’on coupe l’échantillon en 2) moyenne = centre de gravité ==> simple à calculer contrairement à la médiane, propriété méthématique du barycentre, intérêt comptable pour le calcul de KI
Distribution
les quartiles découpent le jeu d’observation en 4 parties égales. Par exemple pour 20 observations, on aura:
- Q1: 25% des sujets dont les observations sont inférieures au 1er quartile
- Q3: 25% des sujets dont les observations sont supérieures au 3éme quartiles
- 50% des sujets dont les obs se regroupe autour de la médiane
obs<-c(1,2,3,4,5,8,12,34)
boxplot(obs)
L’écart type est la racine carré de la variance: \[et= \sqrt{variance}\] La variance est la moyenne des carrés des écarts à la moyenne: \[variance = (\sum_{i=1}^n (x_i-m)^2)/n\] L’écart type correspond à une inertie et a de bonnes propriétés mathématiques. L’écart-type au carré = variance = moyenne des carrées moins la moyenne au carré: \[et^2=var=(\sum_{i=1}^n x^2_i)/n-(\sum_{i=1}^n x_i)/n)^2\]
ET n’a pas beaucoup de sens , sauf quand la variable a une distribution normale. Dans ce cas, [m - et ; et + m] représnete 2/3 des données (le gros de la troupe)
Chapitre4: Mesure de position et de dispersion: Pratique
Commande R pour les meusres de position et de dispertion d’un dataframe
summary(smp)## age prof duree discip
## Min. :19.0 ouvrier :227 Min. :1.000 Min. :0.000
## 1st Qu.:28.0 sans emploi :222 1st Qu.:4.000 1st Qu.:0.000
## Median :37.0 employe :135 Median :5.000 Median :0.000
## Mean :38.9 artisan : 90 Mean :4.302 Mean :0.232
## 3rd Qu.:48.0 prof.intermediaire: 58 3rd Qu.:5.000 3rd Qu.:0.000
## Max. :83.0 (Other) : 61 Max. :5.000 Max. :1.000
## NA's :2 NA's : 6 NA's :223 NA's :6
## n.enfant n.fratrie ecole separation
## Min. : 0.000 Min. : 0.000 Min. :1.000 Min. :0.0000
## 1st Qu.: 0.000 1st Qu.: 2.000 1st Qu.:1.000 1st Qu.:0.0000
## Median : 1.000 Median : 3.000 Median :2.000 Median :0.0000
## Mean : 1.755 Mean : 4.287 Mean :1.866 Mean :0.4226
## 3rd Qu.: 3.000 3rd Qu.: 6.000 3rd Qu.:2.000 3rd Qu.:1.0000
## Max. :13.000 Max. :21.000 Max. :5.000 Max. :1.0000
## NA's :26 NA's :5 NA's :11
## juge.enfant place abus grav.cons
## Min. :0.0000 Min. :0.0000 Min. :0.0000 Min. :1.000
## 1st Qu.:0.0000 1st Qu.:0.0000 1st Qu.:0.0000 1st Qu.:2.000
## Median :0.0000 Median :0.0000 Median :0.0000 Median :4.000
## Mean :0.2771 Mean :0.2285 Mean :0.2778 Mean :3.643
## 3rd Qu.:1.0000 3rd Qu.:0.0000 3rd Qu.:1.0000 3rd Qu.:5.000
## Max. :1.0000 Max. :1.0000 Max. :1.0000 Max. :7.000
## NA's :5 NA's :7 NA's :7 NA's :4
## dep.cons ago.cons ptsd.cons alc.cons
## Min. :0.0000 Min. :0.0000 Min. :0.0000 Min. :0.0000
## 1st Qu.:0.0000 1st Qu.:0.0000 1st Qu.:0.0000 1st Qu.:0.0000
## Median :0.0000 Median :0.0000 Median :0.0000 Median :0.0000
## Mean :0.3967 Mean :0.1665 Mean :0.2165 Mean :0.1865
## 3rd Qu.:1.0000 3rd Qu.:0.0000 3rd Qu.:0.0000 3rd Qu.:0.0000
## Max. :1.0000 Max. :1.0000 Max. :1.0000 Max. :1.0000
##
## subst.cons scz.cons char rs
## Min. :0.0000 Min. :0.0000 Min. :1.000 Min. :1.000
## 1st Qu.:0.0000 1st Qu.:0.0000 1st Qu.:1.000 1st Qu.:1.000
## Median :0.0000 Median :0.0000 Median :1.000 Median :2.000
## Mean :0.2653 Mean :0.0826 Mean :1.512 Mean :2.057
## 3rd Qu.:1.0000 3rd Qu.:0.0000 3rd Qu.:2.000 3rd Qu.:3.000
## Max. :1.0000 Max. :1.0000 Max. :4.000 Max. :3.000
## NA's :96 NA's :103
## ed dr suicide.s suicide.hr
## Min. :1.000 Min. :1.000 Min. :0.0000 Min. :0.0000
## 1st Qu.:1.000 1st Qu.:1.000 1st Qu.:0.0000 1st Qu.:0.0000
## Median :2.000 Median :2.000 Median :0.0000 Median :0.0000
## Mean :1.866 Mean :2.153 Mean :0.7942 Mean :0.2013
## 3rd Qu.:3.000 3rd Qu.:3.000 3rd Qu.:1.0000 3rd Qu.:0.0000
## Max. :3.000 Max. :3.000 Max. :5.0000 Max. :1.0000
## NA's :107 NA's :111 NA's :41 NA's :39
## suicide.past dur.interv
## Min. :0.0000 Min. : 0.00
## 1st Qu.:0.0000 1st Qu.: 48.00
## Median :0.0000 Median : 60.00
## Mean :0.2841 Mean : 61.89
## 3rd Qu.:1.0000 3rd Qu.: 75.00
## Max. :1.0000 Max. :120.00
## NA's :14 NA's :50On peux observer: * pour les var quanti: le min, le max, les quartiles, la médiane, la moyenne, le nb de données manquantes * pour les var cat: le nb de sujet pour chaque modalité et le nb de données manquantes
Summary n’est pas très pratique quand le nb de variables est important. Pour éviter cela on utilise summary (package de prettyR)
install.packages("prettyR")## Installing package into '/home/david/R/x86_64-pc-linux-gnu-library/3.4'
## (as 'lib' is unspecified)library(prettyR)
describe(smp)## Description of smp##
## Numeric
## mean median var sd valid.n
## age 38.90 37 176.38 13.28 797
## duree 4.30 5 0.75 0.87 576
## discip 0.23 0 0.18 0.42 793
## n.enfant 1.76 1 3.36 1.83 773
## n.fratrie 4.29 3 11.84 3.44 799
## ecole 1.87 2 0.96 0.98 794
## separation 0.42 0 0.24 0.49 788
## juge.enfant 0.28 0 0.20 0.45 794
## place 0.23 0 0.18 0.42 792
## abus 0.28 0 0.20 0.45 792
## grav.cons 3.64 4 2.73 1.65 795
## dep.cons 0.40 0 0.24 0.49 799
## ago.cons 0.17 0 0.14 0.37 799
## ptsd.cons 0.22 0 0.17 0.41 799
## alc.cons 0.19 0 0.15 0.39 799
## subst.cons 0.27 0 0.20 0.44 799
## scz.cons 0.08 0 0.08 0.28 799
## char 1.51 1 0.73 0.85 703
## rs 2.06 2 0.77 0.88 696
## ed 1.87 2 0.76 0.87 692
## dr 2.15 2 0.69 0.83 688
## suicide.s 0.79 0 2.06 1.44 758
## suicide.hr 0.20 0 0.16 0.40 760
## suicide.past 0.28 0 0.20 0.45 785
## dur.interv 61.89 60 386.89 19.67 749
##
## Factor
##
## prof ouvrier sans emploi employe artisan prof.intermediaire autre
## Count 227.00 222.00 135.0 90.00 58.00 31.00
## Percent 28.41 27.78 16.9 11.26 7.26 3.88
##
## prof cadre agriculteur <NA>
## Count 24 6.00 6.00
## Percent 3 0.75 0.75
## Mode ouvrierCommande R pour des calculs de position et de dispersion sur 1 variable cat ou quant
On peux calculer isolémment les valeurs pour les var quant:
- moyenne avec mean
- écart type avec sd (standard deviation)
Attention, pour aboutir à un résultat, il faut rejeter les données manquantes avec l’option na.rm=TRUE
mean(smp$age, na.rm = TRUE)## [1] 38.89962sd(smp$age, na.rm = TRUE)## [1] 13.28098min(smp$age, na.rm = TRUE)## [1] 19max(smp$age, na.rm = TRUE)## [1] 83pour les var cat, on utilise table **Attention pour les données manquantes, on peux utiliser l’option useNA=“always” afin de savoir combien de données manquantes existent dans l’échantillon
table(smp$prof,useNA="always")##
## agriculteur artisan autre
## 6 90 31
## cadre employe ouvrier
## 24 135 227
## prof.intermediaire sans emploi <NA>
## 58 222 6table(smp$prof, deparse.level=2, useNA="always") # deparse.level=2 permet d'afficher le nom de la variable calculée## smp$prof
## agriculteur artisan autre
## 6 90 31
## cadre employe ouvrier
## 24 135 227
## prof.intermediaire sans emploi <NA>
## 58 222 6En complément, on peux utiliser la librairie Knitr qui fournit la fonction kable pour faire des affichages plus sympa des tableaux
library(knitr)
tab<-table(smp$prof,useNA="always")
kable(tab)| Var1 | Freq |
|---|---|
| agriculteur | 6 |
| artisan | 90 |
| autre | 31 |
| cadre | 24 |
| employe | 135 |
| ouvrier | 227 |
| prof.intermediaire | 58 |
| sans emploi | 222 |
| NA | 6 |
On peux également utilisement une extension rmarkdown qui permet de manipuler des dataframe ou des tableaux avec beaucoup de lignes
rmarkdown::paged_table(smp)La fonction datatable du package DT propose également son propre formatage
install.packages("DT")## Installing package into '/home/david/R/x86_64-pc-linux-gnu-library/3.4'
## (as 'lib' is unspecified)DT::datatable(smp)A noter, qu’il existe de nombreuses librairies qui facilitent le formatage des tableaux:
- printr: étend le fonctionnement de knitr
- kableExtra: étend la fonction kable de knitr
- formattable à utiliser en conjonction avec kableExtra
- panter
- xtable
- tables
- sjPlot
Lab1:Semaine1:Gestion des données, dataframe, var numériques (quant) et cat (qual)
Décrire les données
View(smp) #visulaise un dataframe
names(smp) # affiche le nom des variables## [1] "age" "prof" "duree" "discip"
## [5] "n.enfant" "n.fratrie" "ecole" "separation"
## [9] "juge.enfant" "place" "abus" "grav.cons"
## [13] "dep.cons" "ago.cons" "ptsd.cons" "alc.cons"
## [17] "subst.cons" "scz.cons" "char" "rs"
## [21] "ed" "dr" "suicide.s" "suicide.hr"
## [25] "suicide.past" "dur.interv"str(smp) # affiche pour chaque variable son mode représentation : par exemple int ou Factor qui affiche les différents niveaux d'une var cat## 'data.frame': 799 obs. of 26 variables:
## $ age : int 31 49 50 47 23 34 24 52 42 45 ...
## $ prof : Factor w/ 8 levels "agriculteur",..: 3 NA 7 6 8 6 3 2 6 6 ...
## $ duree : int 4 NA 5 NA 4 NA NA 5 4 NA ...
## $ discip : int 0 0 0 0 1 0 0 0 1 0 ...
## $ n.enfant : int 2 7 2 0 1 3 5 2 1 2 ...
## $ n.fratrie : int 4 3 2 6 6 2 3 9 12 5 ...
## $ ecole : int 1 2 2 1 1 2 1 2 1 2 ...
## $ separation : int 0 1 0 1 1 0 1 0 1 0 ...
## $ juge.enfant : int 0 0 0 0 NA 0 1 0 1 0 ...
## $ place : int 0 0 0 1 1 0 1 0 0 0 ...
## $ abus : int 0 0 0 0 0 0 0 0 1 1 ...
## $ grav.cons : int 1 2 2 1 2 1 5 1 5 5 ...
## $ dep.cons : int 0 0 0 0 1 0 1 0 1 0 ...
## $ ago.cons : int 1 0 0 0 0 0 0 0 0 0 ...
## $ ptsd.cons : int 0 0 0 0 0 0 0 0 0 0 ...
## $ alc.cons : int 0 0 0 0 0 0 0 0 1 1 ...
## $ subst.cons : int 0 0 0 0 0 0 1 0 1 0 ...
## $ scz.cons : int 0 0 0 0 0 0 0 0 0 0 ...
## $ char : int 1 1 1 1 1 1 1 1 4 1 ...
## $ rs : int 2 2 2 2 2 1 3 2 3 2 ...
## $ ed : int 1 2 3 2 2 2 3 2 3 2 ...
## $ dr : int 1 1 2 2 2 1 2 2 1 2 ...
## $ suicide.s : int 0 0 0 1 0 0 3 0 4 0 ...
## $ suicide.hr : int 0 0 0 0 0 0 1 0 1 0 ...
## $ suicide.past: int 0 0 0 0 1 0 1 0 1 0 ...
## $ dur.interv : int NA 70 NA 105 NA NA 105 84 78 60 ...summary(smp)## age prof duree discip
## Min. :19.0 ouvrier :227 Min. :1.000 Min. :0.000
## 1st Qu.:28.0 sans emploi :222 1st Qu.:4.000 1st Qu.:0.000
## Median :37.0 employe :135 Median :5.000 Median :0.000
## Mean :38.9 artisan : 90 Mean :4.302 Mean :0.232
## 3rd Qu.:48.0 prof.intermediaire: 58 3rd Qu.:5.000 3rd Qu.:0.000
## Max. :83.0 (Other) : 61 Max. :5.000 Max. :1.000
## NA's :2 NA's : 6 NA's :223 NA's :6
## n.enfant n.fratrie ecole separation
## Min. : 0.000 Min. : 0.000 Min. :1.000 Min. :0.0000
## 1st Qu.: 0.000 1st Qu.: 2.000 1st Qu.:1.000 1st Qu.:0.0000
## Median : 1.000 Median : 3.000 Median :2.000 Median :0.0000
## Mean : 1.755 Mean : 4.287 Mean :1.866 Mean :0.4226
## 3rd Qu.: 3.000 3rd Qu.: 6.000 3rd Qu.:2.000 3rd Qu.:1.0000
## Max. :13.000 Max. :21.000 Max. :5.000 Max. :1.0000
## NA's :26 NA's :5 NA's :11
## juge.enfant place abus grav.cons
## Min. :0.0000 Min. :0.0000 Min. :0.0000 Min. :1.000
## 1st Qu.:0.0000 1st Qu.:0.0000 1st Qu.:0.0000 1st Qu.:2.000
## Median :0.0000 Median :0.0000 Median :0.0000 Median :4.000
## Mean :0.2771 Mean :0.2285 Mean :0.2778 Mean :3.643
## 3rd Qu.:1.0000 3rd Qu.:0.0000 3rd Qu.:1.0000 3rd Qu.:5.000
## Max. :1.0000 Max. :1.0000 Max. :1.0000 Max. :7.000
## NA's :5 NA's :7 NA's :7 NA's :4
## dep.cons ago.cons ptsd.cons alc.cons
## Min. :0.0000 Min. :0.0000 Min. :0.0000 Min. :0.0000
## 1st Qu.:0.0000 1st Qu.:0.0000 1st Qu.:0.0000 1st Qu.:0.0000
## Median :0.0000 Median :0.0000 Median :0.0000 Median :0.0000
## Mean :0.3967 Mean :0.1665 Mean :0.2165 Mean :0.1865
## 3rd Qu.:1.0000 3rd Qu.:0.0000 3rd Qu.:0.0000 3rd Qu.:0.0000
## Max. :1.0000 Max. :1.0000 Max. :1.0000 Max. :1.0000
##
## subst.cons scz.cons char rs
## Min. :0.0000 Min. :0.0000 Min. :1.000 Min. :1.000
## 1st Qu.:0.0000 1st Qu.:0.0000 1st Qu.:1.000 1st Qu.:1.000
## Median :0.0000 Median :0.0000 Median :1.000 Median :2.000
## Mean :0.2653 Mean :0.0826 Mean :1.512 Mean :2.057
## 3rd Qu.:1.0000 3rd Qu.:0.0000 3rd Qu.:2.000 3rd Qu.:3.000
## Max. :1.0000 Max. :1.0000 Max. :4.000 Max. :3.000
## NA's :96 NA's :103
## ed dr suicide.s suicide.hr
## Min. :1.000 Min. :1.000 Min. :0.0000 Min. :0.0000
## 1st Qu.:1.000 1st Qu.:1.000 1st Qu.:0.0000 1st Qu.:0.0000
## Median :2.000 Median :2.000 Median :0.0000 Median :0.0000
## Mean :1.866 Mean :2.153 Mean :0.7942 Mean :0.2013
## 3rd Qu.:3.000 3rd Qu.:3.000 3rd Qu.:1.0000 3rd Qu.:0.0000
## Max. :3.000 Max. :3.000 Max. :5.0000 Max. :1.0000
## NA's :107 NA's :111 NA's :41 NA's :39
## suicide.past dur.interv
## Min. :0.0000 Min. : 0.00
## 1st Qu.:0.0000 1st Qu.: 48.00
## Median :0.0000 Median : 60.00
## Mean :0.2841 Mean : 61.89
## 3rd Qu.:1.0000 3rd Qu.: 75.00
## Max. :1.0000 Max. :120.00
## NA's :14 NA's :50summary(smp$age)## Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
## 19.0 28.0 37.0 38.9 48.0 83.0 2smp$age[1] # affiche la 1ere observation## [1] 31smp$age[1:10] # affiche les obs entre 1 et 10## [1] 31 49 50 47 23 34 24 52 42 45min(smp$age, na.rm = TRUE) # calcul de l'age min en enlevant les valeurs manquantes## [1] 19smp$abus[1:10] # Première valeurs de abus## [1] 0 0 0 0 0 0 0 0 1 1unique(smp$abus[1:10]) # liste de smodalités uniques de la var abus## [1] 0 1head(smp$abus, n=10) # affiche les 10 premières valeurs## [1] 0 0 0 0 0 0 0 0 1 1length(smp$abus) # nb total d'observation (nb ligne du tableau smp)## [1] 799table(smp$abus) # renvoie un tableau d'effecitf associé à chaque modalité de la var##
## 0 1
## 572 220table(smp$abus, useNA = "always") # pou afficher les valeurs manquantes##
## 0 1 <NA>
## 572 220 7Passer d’une var quanti à 1 var quali avec factor
Dans le fichier suite à l’import, abus est uen var cat qui a été importée comme une variable numérique. Cela est confirmé avec la commande summary:
summary(smp$abus)## Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
## 0.0000 0.0000 0.0000 0.2778 1.0000 1.0000 7Si on souhaite la traiter comme une var quali (cat), on peux utiliser pour cela la commande factor. R va alors asoscier des niveaux pour chaque modalité différente de la variable
head(smp$abus)## [1] 0 0 0 0 0 0head(factor(smp$abus)) # association de niveau## [1] 0 0 0 0 0 0
## Levels: 0 1abus.factor <-factor(smp$abus, levels=c(0,1), labels=c("Non", "Oui"))
table(smp$abus, useNA = "always") # sans niveau on a 1 affichage peu intuitif##
## 0 1 <NA>
## 572 220 7table(abus.factor, useNA = "always") # avec les niveaux c'est beaucoup plus clair## abus.factor
## Non Oui <NA>
## 572 220 7A noter que abus.factor a été ajouté à l’espace de travail
Création de niveaux aggrégés avec factor et levels
On peux créer des aggrégations et diminuer le nb de modaliés en créant des catégories. Par exemple, on va créer 5 catégories de détenus à partir du nombre d’enfants qu’ils ont: on peux créer 6 niveaux niveau 1: 1 enfant, niveau 2: 2 enfants, etc. A partir de 6 enfants et plus on pourra aggréger en 1 seul niveau. Tout d’abord on observe la variable:
names(smp) # liste les variables du dataframe## [1] "age" "prof" "duree" "discip"
## [5] "n.enfant" "n.fratrie" "ecole" "separation"
## [9] "juge.enfant" "place" "abus" "grav.cons"
## [13] "dep.cons" "ago.cons" "ptsd.cons" "alc.cons"
## [17] "subst.cons" "scz.cons" "char" "rs"
## [21] "ed" "dr" "suicide.s" "suicide.hr"
## [25] "suicide.past" "dur.interv"head(smp$n.enfant) # 1eres observations## [1] 2 7 2 0 1 3summary(smp$n.enfant) # R considère que cette variable est numérique## Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
## 0.000 0.000 1.000 1.755 3.000 13.000 26table(smp$n.enfant) # on constate que l'on pourrait construire des catégories car la plupart du temps, il y a au max 4 enfants##
## 0 1 2 3 4 5 6 7 8 9 10 11 13
## 214 220 125 101 55 31 7 7 7 2 2 1 1table(smp$n.enfant>4) #715 observations avec 4 enfants au max et seulement 78 obs au delà##
## FALSE TRUE
## 715 58Ensuite, on va créer une nouvelle variable n.enfant.cat traitée comme un facteur R:
factor(smp$n.enfant)## [1] 2 7 2 0 1 3 5 2 1 2 0 <NA> 0 3
## [15] 3 0 2 1 <NA> 1 1 0 0 3 1 3 1 2
## [29] 3 0 1 3 2 1 2 1 1 3 1 0 0 1
## [43] 0 1 2 2 1 1 0 2 0 2 0 1 1 0
## [57] 2 2 0 0 0 0 1 0 1 3 1 0 1 1
## [71] 1 0 0 1 0 1 1 0 0 0 1 0 0 0
## [85] 0 2 1 4 4 0 1 7 1 0 0 1 0 4
## [99] 1 3 2 4 0 1 0 0 1 1 3 1 0 2
## [113] 2 1 0 0 1 3 5 2 3 2 1 1 0 3
## [127] 10 4 2 2 0 1 3 0 1 0 1 0 0 3
## [141] 0 1 1 0 0 11 2 0 0 1 0 1 2 2
## [155] 2 0 1 4 2 0 0 1 0 1 1 0 2 1
## [169] 2 1 0 0 1 5 0 1 0 0 0 3 3 3
## [183] 1 0 5 1 3 1 1 1 0 6 0 0 0 3
## [197] 1 0 2 5 2 3 3 1 1 1 1 7 2 5
## [211] 4 0 2 0 0 1 0 2 3 0 2 3 4 1
## [225] 0 0 1 1 3 1 0 1 3 3 3 2 8 2
## [239] 4 3 0 1 0 1 0 2 0 0 2 9 3 0
## [253] 1 1 2 4 0 1 7 4 0 8 1 1 2 0
## [267] 8 5 2 1 8 1 2 3 3 8 2 1 1 4
## [281] 0 3 0 2 3 1 <NA> <NA> 2 <NA> 3 2 0 0
## [295] 1 5 1 3 3 <NA> 5 4 <NA> <NA> 2 <NA> <NA> <NA>
## [309] <NA> <NA> 2 3 1 2 1 0 2 <NA> 1 3 0 1
## [323] 0 5 7 3 2 3 <NA> 4 4 1 0 0 <NA> 0
## [337] 5 1 <NA> 3 2 1 3 2 1 <NA> 0 0 1 0
## [351] 3 1 <NA> <NA> 0 1 <NA> 3 0 0 1 3 0 5
## [365] 2 1 5 4 0 1 2 0 0 1 0 4 2 0
## [379] 0 5 3 2 2 2 6 0 1 1 0 5 3 1
## [393] 0 0 2 1 1 1 1 0 2 1 1 0 2 4
## [407] 2 2 2 3 3 4 3 0 0 0 1 0 0 4
## [421] 0 1 0 0 0 4 0 1 0 2 3 0 0 3
## [435] 1 2 2 4 0 2 0 5 8 0 2 1 1 3
## [449] 2 4 0 4 0 0 5 1 3 1 1 3 3 1
## [463] 4 4 1 1 9 0 1 1 1 1 3 1 2 5
## [477] 0 2 4 4 0 0 3 10 3 2 2 0 0 2
## [491] 5 2 0 1 0 1 2 4 0 1 2 13 4 3
## [505] 1 1 0 6 2 0 1 4 3 0 1 2 0 0
## [519] 2 3 1 0 3 1 1 5 1 1 1 0 4 3
## [533] 2 0 4 2 1 1 2 1 0 3 1 1 2 3
## [547] 1 3 3 5 0 0 0 0 3 4 0 0 1 5
## [561] 0 0 1 1 3 0 1 1 1 0 4 1 1 0
## [575] 2 4 1 2 1 5 1 1 4 1 1 1 2 0
## [589] 2 0 4 5 0 0 1 0 2 4 1 4 2 0
## [603] 1 0 0 1 0 2 1 1 2 0 1 0 2 1
## [617] 1 0 1 1 1 0 1 1 1 0 0 0 3 1
## [631] 0 0 1 0 1 4 1 1 2 1 2 4 1 2
## [645] 6 1 1 2 1 3 1 0 1 1 2 0 1 0
## [659] 0 2 0 1 0 0 1 0 0 3 3 1 6 <NA>
## [673] 2 1 3 1 1 0 <NA> 2 1 1 1 3 1 2
## [687] 3 4 1 2 <NA> 2 2 0 4 0 0 3 1 3
## [701] 3 0 2 1 1 4 0 3 4 1 1 0 1 2
## [715] 1 2 4 3 2 1 4 2 4 1 <NA> 0 7 2
## [729] 1 1 3 4 3 2 3 2 1 0 2 3 3 4
## [743] 3 0 6 6 3 5 4 4 <NA> 1 3 2 4 3
## [757] 2 3 0 0 8 0 3 0 0 0 5 5 5 1
## [771] 2 5 1 1 5 4 3 0 0 3 0 0 0 2
## [785] 2 3 3 7 2 3 0 1 5 2 1 1 3 0
## [799] 2
## Levels: 0 1 2 3 4 5 6 7 8 9 10 11 13smp$n.enfant.cat<-factor(smp$n.enfant)
table(smp$n.enfant.cat) # on voit que 14 factors ont été crées: on va regrouper en 6 factors: 1, 2, 3, 4, 5, 5+##
## 0 1 2 3 4 5 6 7 8 9 10 11 13
## 214 220 125 101 55 31 7 7 7 2 2 1 1levels(smp$n.enfant.cat)## [1] "0" "1" "2" "3" "4" "5" "6" "7" "8" "9" "10" "11" "13"nlevels(smp$n.enfant.cat)## [1] 13levels(smp$n.enfant.cat)[6:13] # Sélection des derniers niveaux de 6 à 13 ## [1] "5" "6" "7" "8" "9" "10" "11" "13"levels(smp$n.enfant.cat)[6:13]<-"5+" # Création d'une seule modalité de la variable levels pour les niveaux aggrégé appellé "5+"
table(smp$n.enfant.cat) # on vérifie que pour cette variable créée, on a bien 6 niveaux au lieu de 13 : les niveaux 6 à 13 ont été aggrégés en 1 seul niveau##
## 0 1 2 3 4 5+
## 214 220 125 101 55 58Sauvegardes
save(smp, file="smp.sav") #sauvegarde d'un dataframe
## savehistory("commande_history.R") #sauvegarde de l'historique des commandesLab2:Semaine1:Indexation critériée d’observations, sélection de variables, graphiques univariés
On recharge le dataframe générées avec les catégories n.enfant.cat ##Mécanisme de base de l’indexation
load("smp.sav") #Rechargement d'une sauvegarde effectuée avec save
names(smp) # liste des variables du dataframe## [1] "age" "prof" "duree" "discip"
## [5] "n.enfant" "n.fratrie" "ecole" "separation"
## [9] "juge.enfant" "place" "abus" "grav.cons"
## [13] "dep.cons" "ago.cons" "ptsd.cons" "alc.cons"
## [17] "subst.cons" "scz.cons" "char" "rs"
## [21] "ed" "dr" "suicide.s" "suicide.hr"
## [25] "suicide.past" "dur.interv" "n.enfant.cat"smp$age[1] # accès à la 1ere observation de la variable age## [1] 31smp[1,1] # accès à la 1ere observation de la 1ere variable## [1] 31smp[1,"age"] # accès à la 1ere observation de la variable age## [1] 31smp[1,"prof"] # accès à la 1ere observation de la variable prof qui est catégoriel## [1] autre
## 8 Levels: agriculteur artisan autre cadre employe ... sans emploihead(smp$prof)## [1] autre <NA> prof.intermediaire
## [4] ouvrier sans emploi ouvrier
## 8 Levels: agriculteur artisan autre cadre employe ... sans emploihead(smp$prof=="agriculteur") # restriction sur la variable prof pour la modalité agriculteur## [1] FALSE NA FALSE FALSE FALSE FALSEtable(smp$prof=="agriculteur") # tableau des valeurs de la variable prof qui vaut agriculteur##
## FALSE TRUE
## 787 6which(smp$prof=="agriculteur") # numéro des observations pour lesquelles la valeur de prof==agriculteur## [1] 15 312 384 391 439 442smp$age[which(smp$prof=="agriculteur")] # valeurs de la variable age pour lesquels la prof vaut agriculteur## [1] 64 42 37 36 35 79Filtres avec la commande subset
Cette commande permet de filtrer les observations du dataframe directement avec une expression. 3 parametres donc:
- le dataframe
- le filtre à appliquer sur les lignes du dataframe
- la variable que l’on souhaite afficher
subset(smp,prof=="agriculteur",age)## age
## 15 64
## 312 42
## 384 37
## 391 36
## 439 35
## 442 79subset(smp,prof=="agriculteur",1:5) #extension de la sélection à plus de 1 variable. Ici les variables 1 à 5: age, prof, durée, discip, n.enfant## age prof duree discip n.enfant
## 15 64 agriculteur NA 0 3
## 312 42 agriculteur 4 0 3
## 384 37 agriculteur 5 1 2
## 391 36 agriculteur 4 1 3
## 439 35 agriculteur 3 0 0
## 442 79 agriculteur 5 0 5names(smp)[1:5]## [1] "age" "prof" "duree" "discip" "n.enfant"subset(smp,prof=="agriculteur",c(1,3,4,5)) # liste des variables à afficher quand elles ne se suivent pas## age duree discip n.enfant
## 15 64 NA 0 3
## 312 42 4 0 3
## 384 37 5 1 2
## 391 36 4 1 3
## 439 35 3 0 0
## 442 79 5 0 5subset(smp,prof=="agriculteur" & n.enfant>2,c(1,3,4,5)) #combinaison de 2 filtres: le filtre sur la profession ET sur le nombre d'enfant## age duree discip n.enfant
## 15 64 NA 0 3
## 312 42 4 0 3
## 391 36 4 1 3
## 442 79 5 0 5Filtre encore plus complexe: On filtre les valeurs manquantes pour la variable duree
subset(smp,prof=="agriculteur" & n.enfant>2,c(1,3,4,5) & complete.cases(duree))## age prof duree discip n.enfant n.fratrie ecole separation
## 15 64 agriculteur NA 0 3 2 1 0
## 312 42 agriculteur 4 0 3 6 1 0
## 391 36 agriculteur 4 1 3 4 3 1
## 442 79 agriculteur 5 0 5 6 2 0
## juge.enfant place abus grav.cons dep.cons ago.cons ptsd.cons alc.cons
## 15 0 0 0 1 0 0 0 0
## 312 0 0 0 4 1 1 0 1
## 391 1 1 0 2 0 0 1 0
## 442 0 0 0 NA 0 0 0 0
## subst.cons scz.cons char rs ed dr suicide.s suicide.hr suicide.past
## 15 0 0 1 1 1 3 0 0 0
## 312 0 0 2 1 3 2 3 1 0
## 391 1 0 1 NA 3 1 0 0 0
## 442 0 0 1 2 1 1 0 0 0
## dur.interv n.enfant.cat
## 15 80 3
## 312 NA 3
## 391 NA 3
## 442 85 5+Filtre avec table et prop.table
Création d’un tableau d’effectif à partir de la commande table
table.n.enfant.cat<-table(smp$n.enfant.cat) # création d'un tableau d'effectif
sum(table.n.enfant.cat) #calcul de l'effectif## [1] 773100*table.n.enfant.cat/sum(table.n.enfant.cat) # effectifs en proportion##
## 0 1 2 3 4 5+
## 27.684347 28.460543 16.170763 13.065977 7.115136 7.503234100*prop.table(table.n.enfant.cat) # commande équivalente qui calcul directement le tableau de fréquence relative rapporté à l'effecitf total##
## 0 1 2 3 4 5+
## 27.684347 28.460543 16.170763 13.065977 7.115136 7.503234round(prop.table(table.n.enfant.cat), 3)##
## 0 1 2 3 4 5+
## 0.277 0.285 0.162 0.131 0.071 0.075Représentation graphique des fréquences relatives (graphiques univariés)
Variables catégorielles (par exemple n.enfant.cat)
enfant.cat.freq<-100*prop.table(table.n.enfant.cat)
barplot(enfant.cat.freq)#représentation graphique du tableau de fréquences relatives
barplot(enfant.cat.freq, ylim=c(0,30)) # l'échelle des Y est comprise dans un intervalle spécifié
barplot(enfant.cat.freq, ylim=c(0,30), las=1)# Agrément du graphique: les labels en x et y sont orientés horizontalement
Variables numériques (par exemple age)
head(smp$age)## [1] 31 49 50 47 23 34hist(smp$age) #R choisit lui même les classes à regroupées
hist(smp$age, nclass=8) # personnalisation en figeant le nombre de classe à utiliser
hist(smp$age, nclass=8, prob=TRUE, las=1) #covnersion de l'historgramme en historgramme de densité
lines(density(smp$age, na.rm=TRUE)) #utilisation de la densité non paramétrique qui permet de ne pas trop oberver les classes mais uniquement une forme de densité générale
Semaine1:Devoir évalué par les pairs
Chapitre5:Semaine2:Intervalles de confiance
Un intervalle de confiance permet d’estimer la borne min et la borne max de la prévalence d’une modalité sur une population à partir d’un échantillon et selon un niveau de probabilité donné. Par exemple: il y a 95% de chances que la proportion de “déprimés” dans le milieu carcéral soit comprise entre 8% et 64%. Cette estimation étant faite à partir d’un échantillon de 10 détenus tirés au sort la formule de calcul de l’intervalle de confiance 95% est:\[[M - 1.96 * ET ; M + 1.96 * ET]\] M est la modalité à estimer ET est l’écart type de m pour l’échantillon.
smp<-read.csv2("smp2.csv")
library(prettyR)
describe(smp$age)## Description of structure(list(x = c(31L, 49L, 50L, 47L, 23L, 34L, 24L, 52L, 42L, 45L, 31L, NA, 21L, 40L, 64L, 67L, 60L, 63L, NA, 28L, 20L, 30L, 32L, 31L, 26L, 42L, 32L, 40L, 41L, 27L, 24L, 38L, 39L, 36L, 29L, 41L, 36L, 41L, 21L, 21L, 46L, 22L, 21L, 35L, 45L, 38L, 19L, 21L, 27L, 40L, 39L, 47L, 24L, 36L, 39L, 22L, 38L, 37L, 29L, 23L, 36L, 42L, 56L, 28L, 36L, 38L, 43L, 29L, 64L, 25L, 51L, 35L, 30L, 37L, 26L, 36L, 58L, 32L, 30L, 26L, 27L, 23L, 24L, 39L, 43L, 39L, 26L, 44L, 37L, 40L, 24L, 46L, 26L, 38L, 37L, 30L, 39L, 36L, 39L, 28L, 27L, 51L, 48L, 47L, 41L, 35L, 25L, 31L, 44L, 40L, 29L, 34L, 49L, 57L, 33L, 35L, 32L, 34L, 46L, 45L, 31L, 42L, 48L, 34L, 34L, 64L, 50L, 53L, 49L, 53L, 37L, 42L, 55L, 32L, 33L, 40L, 29L, 32L, 23L, 61L, 39L, 30L, 37L, 30L, 39L, 49L, 44L, 40L, 56L, 43L, 27L, 21L, 44L, 50L, 50L, 20L, 37L, 42L, 27L, 22L, 25L, 20L, 21L, 19L, 25L, 24L, 49L, 24L, 26L, 35L, 22L, 24L, 23L, 46L, 26L, 41L, 51L, 20L, 30L, 37L, 49L, 28L, 28L, 51L, 40L, 33L, 25L, 29L, 40L, 43L, 35L, 50L, 44L, 35L, 24L, 43L, 26L, 45L, 42L, 45L, 48L, 45L, 34L, 31L, 40L, 22L, 42L, 38L, 38L, 40L, 46L, 26L, 29L, 25L, 40L, 43L, 28L, 29L, 32L, 28L, 57L, 31L, 71L, 33L, 24L, 22L, 25L, 26L, 52L, 33L, 38L, 39L, 41L, 52L, 33L, 39L, 59L, 33L, 50L, 58L, 23L, 41L, 43L, 42L, 22L, 57L, 41L, 30L, 66L, 49L, 46L, 28L, 59L, 35L, 44L, 83L, 34L, 49L, 60L, 56L, 46L, 62L, 41L, 27L, 53L, 48L, 66L, 66L, 55L, 61L, 43L, 54L, 38L, 51L, 51L, 50L, 56L, 53L, 49L, 41L, 44L, 64L, 42L, 52L, 72L, 43L, 30L, 32L, 43L, 25L, 27L, 25L, 52L, 39L, 42L, 59L, 46L, 62L, 50L, 24L, 43L, 32L, 67L, 28L, 44L, 19L, 20L, 23L, 26L, 28L, 31L, 42L, 57L, 30L, 36L, 53L, 33L, 25L, 22L, 42L, 25L, 32L, 23L, 45L, 48L, 35L, 37L, 38L, 24L, 47L, 61L, 38L, 27L, 27L, 26L, 30L, 47L, 37L, 30L, 41L, 29L, 37L, 28L, 47L, 26L, 50L, 23L, 60L, 37L, 48L, 41L, 28L, 54L, 61L, 33L, 31L, 25L, 66L, 26L, 29L, 29L, 53L, 24L, 48L, 40L, 47L, 40L, 41L, 54L, 25L, 36L, 44L, 32L, 27L, 31L, 34L, 34L, 71L, 20L, 54L, 39L, 50L, 36L, 37L, 43L, 28L, 21L, 35L, 36L, 53L, 36L, 38L, 66L, 62L, 38L, 24L, 49L, 21L, 34L, 29L, 36L, 29L, 33L, 34L, 57L, 65L, 25L, 36L, 31L, 54L, 49L, 42L, 30L, 20L, 23L, 21L, 23L, 39L, 45L, 29L, 21L, 54L, 77L, 23L, 32L, 58L, 49L, 26L, 40L, 51L, 62L, 45L, 41L, 30L, 52L, 20L, 36L, 34L, 35L, 30L, 46L, 79L, 66L, 19L, 41L, 51L, 26L, 56L, 33L, 39L, 72L, 45L, 59L, 21L, 41L, 43L, 55L, 26L, 49L, 29L, 26L, 28L, 77L, 61L, 63L, 30L, 49L, 48L, 45L, 32L, 56L, 48L, 64L, 73L, 33L, 74L, 54L, 27L, 49L, 45L, 27L, 53L, 62L, 54L, 37L, 56L, 60L, 33L, 34L, 32L, 44L, 49L, 46L, 67L, 39L, 59L, 63L, 81L, 38L, 58L, 42L, 73L, 48L, 41L, 28L, 44L, 45L, 46L, 50L, 27L, 56L, 46L, 42L, 25L, 23L, 26L, 19L, 24L, 24L, 32L, 23L, 24L, 33L, 21L, 33L, 41L, 24L, 31L, 19L, 25L, 51L, 39L, 22L, 20L, 30L, 34L, 28L, 20L, 20L, 33L, 24L, 32L, 37L, 25L, 24L, 29L, 19L, 37L, 56L, 49L, 60L, 29L, 22L, 20L, 49L, 33L, 30L, 29L, 25L, 62L, 41L, 33L, 44L, 60L, 24L, 24L, 33L, 27L, 45L, 33L, 44L, 23L, 23L, 35L, 36L, 28L, 24L, 27L, 27L, 28L, 27L, 40L, 52L, 19L, 31L, 21L, 33L, 23L, 30L, 23L, 31L, 48L, 24L, 24L, 26L, 32L, 29L, 38L, 23L, 50L, 26L, 47L, 38L, 24L, 24L, 19L, 25L, 31L, 33L, 26L, 38L, 23L, 37L, 19L, 49L, 33L, 30L, 38L, 30L, 26L, 27L, 21L, 31L, 19L, 26L, 28L, 49L, 35L, 25L, 32L, 27L, 20L, 30L, 25L, 21L, 54L, 27L, 22L, 39L, 21L, 54L, 49L, 23L, 36L, 59L, 50L, 24L, 47L, 42L, 41L, 33L, 46L, 23L, 19L, 39L, 38L, 40L, 39L, 40L, 44L, 26L, 48L, 47L, 23L, 25L, 20L, 45L, 44L, 57L, 39L, 55L, 19L, 34L, 28L, 33L, 19L, 33L, 27L, 46L, 47L, 22L, 27L, 26L, 52L, 56L, 44L, 63L, 34L, 41L, 38L, 37L, 58L, 37L, 24L, 60L, 26L, 21L, 52L, 20L, 37L, 32L, 32L, 58L, 49L, 32L, 37L, 46L, 50L, 44L, 47L, 37L, 38L, 50L, 56L, 30L, 34L, 43L, 55L, 43L, 31L, 55L, 41L, 68L, 45L, 48L, 42L, 71L, 38L, 46L, 65L, 51L, 57L, 57L, 71L, 40L, 43L, 71L, 48L, 34L, 69L, 43L, 35L, 62L, 34L, 51L, 48L, 36L, 44L, 49L, 74L, 19L, 56L, 57L, 65L, 52L, 77L, 29L, 37L, 45L, 40L, 72L, 27L, 56L, 35L, 30L, 37L, 30L, 40L, 54L, 26L, 48L, 83L, 32L, 22L, 48L, 67L, 58L, 37L, 24L, 34L, 39L, 38L, 39L, 56L, 35L, 26L, 70L, 68L, 42L, 41L, 40L, 26L, 50L, 27L, 28L, 44L, 31L, 38L, 71L)), .Names = "x", row.names = c(NA, -799L), class = "data.frame")##
## Numeric
## mean median var sd valid.n
## x 38.9 37 176.38 13.28 797la moyenne (mean)de l’age est: 38.9 ans
l’écart type (sd) de l’age des détenus est: 13.28
l’écart type de la moyenne de l’age est: \[13.28/\sqrt{797}\]
l’intervalle de confiance 95% de l’age des détenus est donc [37.98 : 39.82]
38.9-1.96*13.28/sqrt(797)## [1] 37.9780138.9+1.96*13.28/sqrt(797)## [1] 39.82199L’intervalle de confiance d’un pourcentage peut être estimé à partir de méthodes qui convergent lorsque la taille de l’échantillon est importante. pour 3 individus présentant une modalité commune (par exemple dépression) parmi un échantillon de 10 personnes tirées au sort. Quel est l’intervalle de confiance de la prévalence de la dépression: 30% ?
install.packages("binom")## Installing package into '/home/david/R/x86_64-pc-linux-gnu-library/3.4'
## (as 'lib' is unspecified)library(binom)
binom.confint(3,10,method = "all")## method x n mean lower upper
## 1 agresti-coull 3 10 0.3000000 0.10333842 0.6076747
## 2 asymptotic 3 10 0.3000000 0.01597423 0.5840258
## 3 bayes 3 10 0.3181818 0.07454423 0.5794516
## 4 cloglog 3 10 0.3000000 0.07113449 0.5778673
## 5 exact 3 10 0.3000000 0.06673951 0.6524529
## 6 logit 3 10 0.3000000 0.09976832 0.6236819
## 7 probit 3 10 0.3000000 0.08991347 0.6150429
## 8 profile 3 10 0.3000000 0.08470272 0.6065091
## 9 lrt 3 10 0.3000000 0.08458545 0.6065389
## 10 prop.test 3 10 0.3000000 0.08094782 0.6463293
## 11 wilson 3 10 0.3000000 0.10779127 0.6032219binom.confint(3,10,method = "prop.test")## method x n mean lower upper
## 1 prop.test 3 10 0.3 0.08094782 0.6463293binom.confint(3,10,method = "exact")## method x n mean lower upper
## 1 exact 3 10 0.3 0.06673951 0.6524529binom.confint(300,1000,method = "all")## method x n mean lower upper
## 1 agresti-coull 300 1000 0.3000000 0.2723966 0.3291341
## 2 asymptotic 300 1000 0.3000000 0.2715974 0.3284026
## 3 bayes 300 1000 0.3001998 0.2719448 0.3286787
## 4 cloglog 300 1000 0.3000000 0.2718595 0.3285966
## 5 exact 300 1000 0.3000000 0.2717211 0.3294617
## 6 logit 300 1000 0.3000000 0.2723865 0.3291466
## 7 probit 300 1000 0.3000000 0.2722277 0.3289871
## 8 profile 300 1000 0.3000000 0.2721340 0.3288893
## 9 lrt 300 1000 0.3000000 0.2721419 0.3289000
## 10 prop.test 300 1000 0.3000000 0.2719222 0.3296354
## 11 wilson 300 1000 0.3000000 0.2724068 0.3291239Chapitre6:Semaine2:Coefficient de corrélation
L’objectif est de pouvoir estimer la force de corrélation entre 2 variables. Cette force peut être estimé par le coefficient de corrélation de Pearson: r * Si r = 0 la corrélation est nulle: les variables sous jacentes sont indépendantes (ssi elle suivent une loi normale) * Si r = +/- 1 x et y sont mutuellement déterminées: \[y=ax+b\]
L’interprétation de ce coeff n’est pas facile. Le plus simple est de représenter graphiquement les 2 variables X et y afin de visualiser leur corrélation
plot(jitter(smp$age),jitter(smp$n.enfant))
str(smp)## 'data.frame': 799 obs. of 26 variables:
## $ age : int 31 49 50 47 23 34 24 52 42 45 ...
## $ prof : Factor w/ 8 levels "agriculteur",..: 3 NA 7 6 8 6 3 2 6 6 ...
## $ duree : int 4 NA 5 NA 4 NA NA 5 4 NA ...
## $ discip : int 0 0 0 0 1 0 0 0 1 0 ...
## $ n.enfant : int 2 7 2 0 1 3 5 2 1 2 ...
## $ n.fratrie : int 4 3 2 6 6 2 3 9 12 5 ...
## $ ecole : int 1 2 2 1 1 2 1 2 1 2 ...
## $ separation : int 0 1 0 1 1 0 1 0 1 0 ...
## $ juge.enfant : int 0 0 0 0 NA 0 1 0 1 0 ...
## $ place : int 0 0 0 1 1 0 1 0 0 0 ...
## $ abus : int 0 0 0 0 0 0 0 0 1 1 ...
## $ grav.cons : int 1 2 2 1 2 1 5 1 5 5 ...
## $ dep.cons : int 0 0 0 0 1 0 1 0 1 0 ...
## $ ago.cons : int 1 0 0 0 0 0 0 0 0 0 ...
## $ ptsd.cons : int 0 0 0 0 0 0 0 0 0 0 ...
## $ alc.cons : int 0 0 0 0 0 0 0 0 1 1 ...
## $ subst.cons : int 0 0 0 0 0 0 1 0 1 0 ...
## $ scz.cons : int 0 0 0 0 0 0 0 0 0 0 ...
## $ char : int 1 1 1 1 1 1 1 1 4 1 ...
## $ rs : int 2 2 2 2 2 1 3 2 3 2 ...
## $ ed : int 1 2 3 2 2 2 3 2 3 2 ...
## $ dr : int 1 1 2 2 2 1 2 2 1 2 ...
## $ suicide.s : int 0 0 0 1 0 0 3 0 4 0 ...
## $ suicide.hr : int 0 0 0 0 0 0 1 0 1 0 ...
## $ suicide.past: int 0 0 0 0 1 0 1 0 1 0 ...
## $ dur.interv : int NA 70 NA 105 NA NA 105 84 78 60 ...cor(smp$age, smp$n.enfant,use="complete.obs")## [1] 0.4326039Chapitre7:Semaine2:Risque relatif et odds-ratio
Définition du RR et de l’OD
soit les données d’enquêtes suivantes: ||Enrhumés| |:-:|:-:|:-:| ||Oui|Non| |Tabac: Oui|30|600| |Tabac: Non|30|300|
La question est d’évaluer la prévalence du rhume selon que l’on est fumeur ou pas
RR = (30 / 330) / (30 / 630) #Calcul du Risque Relatif
OD = (30 / 300) / (30 / 600) #Calcul de l'ODDSRATIO
RR## [1] 1.909091OD## [1] 2Pour RR, on peut dire que les fumeurs on 1,9 plus de chance d’être enrhumés que les non fumeurs.
Pour ODDSRATIO, c’est à peu près la même chose surtout si l’échantillon est important
- Si la modalité observée est peu fréquente, alors RR = OD
- Si l’échantillon est faible, alors une étude cas-témoin (50% d’individu avec la modalité et 50% sans la modalité) est possible et c’est OD qui est utilisé car RR n’a plus de sens
Corrélation
La corrélation entre le niveau d’évitement du danger et la dépression dans la population Il faut tout d’abord transformer le niveau d’évitement du danger en variable binaire
smp$ed.bin<-ifelse(smp$ed>2,1,0)
str(smp)## 'data.frame': 799 obs. of 27 variables:
## $ age : int 31 49 50 47 23 34 24 52 42 45 ...
## $ prof : Factor w/ 8 levels "agriculteur",..: 3 NA 7 6 8 6 3 2 6 6 ...
## $ duree : int 4 NA 5 NA 4 NA NA 5 4 NA ...
## $ discip : int 0 0 0 0 1 0 0 0 1 0 ...
## $ n.enfant : int 2 7 2 0 1 3 5 2 1 2 ...
## $ n.fratrie : int 4 3 2 6 6 2 3 9 12 5 ...
## $ ecole : int 1 2 2 1 1 2 1 2 1 2 ...
## $ separation : int 0 1 0 1 1 0 1 0 1 0 ...
## $ juge.enfant : int 0 0 0 0 NA 0 1 0 1 0 ...
## $ place : int 0 0 0 1 1 0 1 0 0 0 ...
## $ abus : int 0 0 0 0 0 0 0 0 1 1 ...
## $ grav.cons : int 1 2 2 1 2 1 5 1 5 5 ...
## $ dep.cons : int 0 0 0 0 1 0 1 0 1 0 ...
## $ ago.cons : int 1 0 0 0 0 0 0 0 0 0 ...
## $ ptsd.cons : int 0 0 0 0 0 0 0 0 0 0 ...
## $ alc.cons : int 0 0 0 0 0 0 0 0 1 1 ...
## $ subst.cons : int 0 0 0 0 0 0 1 0 1 0 ...
## $ scz.cons : int 0 0 0 0 0 0 0 0 0 0 ...
## $ char : int 1 1 1 1 1 1 1 1 4 1 ...
## $ rs : int 2 2 2 2 2 1 3 2 3 2 ...
## $ ed : int 1 2 3 2 2 2 3 2 3 2 ...
## $ dr : int 1 1 2 2 2 1 2 2 1 2 ...
## $ suicide.s : int 0 0 0 1 0 0 3 0 4 0 ...
## $ suicide.hr : int 0 0 0 0 0 0 1 0 1 0 ...
## $ suicide.past: int 0 0 0 0 1 0 1 0 1 0 ...
## $ dur.interv : int NA 70 NA 105 NA NA 105 84 78 60 ...
## $ ed.bin : num 0 0 1 0 0 0 1 0 1 0 ...table(smp$ed.bin,smp$ed, deparse.level = 2, useNA = "always")## smp$ed
## smp$ed.bin 1 2 3 <NA>
## 0 315 155 0 0
## 1 0 0 222 0
## <NA> 0 0 0 107On peux ensuite chercher à calculer la corrélation entre smp\(ed.bin et smp\)dep.cons
## Installing package into '/home/david/R/x86_64-pc-linux-gnu-library/3.4'
## (as 'lib' is unspecified)La fonction twoby2 permet de calculer la corrélation:
twoby2(1-smp$ed.bin, 1-smp$dep.cons)## 2 by 2 table analysis:
## ------------------------------------------------------
## Outcome : 0
## Comparing : 0 vs. 1
##
## 0 1 P(0) 95% conf. interval
## 0 126 96 0.5676 0.5016 0.6312
## 1 135 335 0.2872 0.2481 0.3298
##
## 95% conf. interval
## Relative Risk: 1.9760 1.6456 2.3726
## Sample Odds Ratio: 3.2569 2.3361 4.5408
## Conditional MLE Odds Ratio: 3.2508 2.3037 4.6035
## Probability difference: 0.2803 0.2020 0.3549
##
## Exact P-value: 0
## Asymptotic P-value: 0
## ------------------------------------------------------Lab3:Semaine2:Langage R Markdown
Ce lab présente R Markdown. Toutes les notes de cours sont déjà en R Markdown.
Quizz Exercice 2:Semaine2
Quel est le nombre moyen d’enfants (variable n.enfant) chez les individus diagnostiqués comme dépressifs (dep.cons = 1) (2 chiffres après la virgule) ?
describe(subset(smp,smp$dep.cons==1,select=n.enfant))## Description of subset(smp, smp$dep.cons == 1, select = n.enfant)##
## Numeric
## mean median var sd valid.n
## n.enfant 1.76 1 2.81 1.68 310Autre solution
mean(smp[smp$dep.cons==1,'n.enfant'],na.rm = TRUE)## [1] 1.764516Donner la borne supérieure de l’intervalle interquartile pour la variable duree chez les individus dont l’âge est strictement inférieur à 35 ans (2 chiffres après la virgule). Borne supérieure :
summary(subset(smp,smp$age<35,select=duree))## duree
## Min. :1.000
## 1st Qu.:3.250
## Median :4.000
## Mean :4.005
## 3rd Qu.:5.000
## Max. :5.000
## NA's :129À partir du data frame smp, nous souhaitons recoder l’âge (variable age) en variable catégorielle en considérant 4 intervalles de classe dont les bornes intermédiaires sont définies à partir des 1er, 2ème et 3ème quartiles. Les bornes inférieures et supérieures de la première et dernière classe seront naturellement les valeurs minimale et maximale observées pour la variable age. À l’exception de la première classe dont les deux bornes d’intervalle seront fermées (c’est-à-dire que les bornes seront inclues dans l’intervalle), les bornes inférieures des classes suivantes (2 à 4) seront ouvertes, et les bornes supérieures fermées. Indiquer l’effectif associé à la 3ème classe ainsi constituée.
on cherche tout d’abord à connaître les valeurs min, max de chaque quartile pour la variable age. Les 4 intervalles de classe sont: [Min : Q1] ]Q1 : Median] ]Median : Q3] ]Q3 : Max]
summary(smp$age)## Min. 1st Qu. Median Mean 3rd Qu. Max. NA's
## 19.0 28.0 37.0 38.9 48.0 83.0 2Pou répondre à la question, on peux:
- construire un subset en filtrant
- effectuer un filtre directement sur smp
- recoder une nouvelle variable catégorielle tel que suggéré par l’énoncé puis appliquer un describe sur cette variable
describe(smp[smp$age >37 & smp$age<=48, 'age'])## Description of structure(list(x = c(47L, 42L, 45L, NA, 40L, NA, 42L, 40L, 41L, 38L, 39L, 41L, 41L, 46L, 45L, 38L, 40L, 39L, 47L, 39L, 38L, 42L, 38L, 43L, 39L, 43L, 39L, 44L, 40L, 46L, 38L, 39L, 39L, 48L, 47L, 41L, 44L, 40L, 46L, 45L, 42L, 48L, 42L, 40L, 39L, 39L, 44L, 40L, 43L, 44L, 42L, 46L, 41L, 40L, 40L, 43L, 44L, 43L, 45L, 42L, 45L, 48L, 45L, 40L, 42L, 38L, 38L, 40L, 46L, 40L, 43L, 38L, 39L, 41L, 39L, 41L, 43L, 42L, 41L, 46L, 44L, 46L, 41L, 48L, 43L, 38L, 41L, 44L, 42L, 43L, 43L, 39L, 42L, 46L, 43L, 44L, 42L, 42L, 45L, 48L, 38L, 47L, 38L, 47L, 41L, 47L, 48L, 41L, 48L, 40L, 47L, 40L, 41L, 44L, 39L, 43L, 38L, 38L, 42L, 39L, 45L, 40L, 45L, 41L, 46L, 41L, 39L, 45L, 41L, 43L, 48L, 45L, 48L, 45L, 44L, 46L, 39L, 38L, 42L, 48L, 41L, 44L, 45L, 46L, 46L, 42L, 41L, 39L, 41L, 44L, 45L, 44L, 40L, 48L, 38L, 47L, 38L, 38L, 38L, 39L, 47L, 42L, 41L, 46L, 39L, 38L, 40L, 39L, 40L, 44L, 48L, 47L, 45L, 44L, 39L, 46L, 47L, 44L, 41L, 38L, 46L, 44L, 47L, 38L, 43L, 43L, 41L, 45L, 48L, 42L, 38L, 46L, 40L, 43L, 48L, 43L, 48L, 44L, 45L, 40L, 40L, 48L, 48L, 39L, 38L, 39L, 42L, 41L, 40L, 44L, 38L)), .Names = "x", row.names = c(NA, -211L), class = "data.frame")##
## Numeric
## mean median var sd valid.n
## x 42.48 42 9.94 3.15 209La réponse est 209
Quiz exercice 2bis (facultatif):Semaine2
na.fail part en échec si il recoit des informations contenant NA
na.fail(smp[c(20,221,342,446,531),])## age prof duree discip n.enfant n.fratrie ecole separation
## 20 28 ouvrier 5 0 1 1 2 0
## 221 57 sans emploi 5 0 2 2 2 0
## 342 37 sans emploi 4 0 1 2 1 1
## 446 51 ouvrier 3 0 1 3 3 0
## 531 51 artisan 5 0 4 5 1 0
## juge.enfant place abus grav.cons dep.cons ago.cons ptsd.cons alc.cons
## 20 0 0 0 2 0 0 0 0
## 221 0 0 1 5 1 0 0 0
## 342 0 0 0 1 0 0 0 0
## 446 0 0 0 2 0 0 0 0
## 531 0 0 1 6 1 0 0 0
## subst.cons scz.cons char rs ed dr suicide.s suicide.hr suicide.past
## 20 1 0 1 2 2 2 0 0 1
## 221 0 0 2 1 3 3 4 1 1
## 342 0 0 1 1 2 2 0 0 0
## 446 0 0 1 2 2 2 0 0 0
## 531 0 1 2 3 3 3 3 1 0
## dur.interv ed.bin
## 20 95 0
## 221 50 1
## 342 30 0
## 446 85 0
## 531 68 1La réponse est NON
Pour obtenir un résumé
summary(table(smp[smp$age>=20 & smp$age<=30,'age']))## Number of cases in table: 242
## Number of factors: 1Chapitre8:Semaine3:Tests statistiques le p (petit p) - Test de fisher
Rappel sur le chargement d’un fichier de donnée
Chapitre9:Semaine3:Test de Neyman et Pearson (tests d’hypothèses)
Chapitre10:Semaine3:Comparaison de 2 pourcentages
La prévalence de la dépression est-elle plus elevée chez les détenus qui ont un haut niveau d’évitement que chez les détenus qui ont un bas niveau d’évitement du danger Création d’une variable binaire dans le dataframe smp
Description de la dataframe smp
## 'data.frame': 799 obs. of 27 variables:
## $ age : int 31 49 50 47 23 34 24 52 42 45 ...
## $ prof : Factor w/ 8 levels "agriculteur",..: 3 NA 7 6 8 6 3 2 6 6 ...
## $ duree : int 4 NA 5 NA 4 NA NA 5 4 NA ...
## $ discip : int 0 0 0 0 1 0 0 0 1 0 ...
## $ n.enfant : int 2 7 2 0 1 3 5 2 1 2 ...
## $ n.fratrie : int 4 3 2 6 6 2 3 9 12 5 ...
## $ ecole : int 1 2 2 1 1 2 1 2 1 2 ...
## $ separation : int 0 1 0 1 1 0 1 0 1 0 ...
## $ juge.enfant : int 0 0 0 0 NA 0 1 0 1 0 ...
## $ place : int 0 0 0 1 1 0 1 0 0 0 ...
## $ abus : int 0 0 0 0 0 0 0 0 1 1 ...
## $ grav.cons : int 1 2 2 1 2 1 5 1 5 5 ...
## $ dep.cons : int 0 0 0 0 1 0 1 0 1 0 ...
## $ ago.cons : int 1 0 0 0 0 0 0 0 0 0 ...
## $ ptsd.cons : int 0 0 0 0 0 0 0 0 0 0 ...
## $ alc.cons : int 0 0 0 0 0 0 0 0 1 1 ...
## $ subst.cons : int 0 0 0 0 0 0 1 0 1 0 ...
## $ scz.cons : int 0 0 0 0 0 0 0 0 0 0 ...
## $ char : int 1 1 1 1 1 1 1 1 4 1 ...
## $ rs : int 2 2 2 2 2 1 3 2 3 2 ...
## $ ed : int 1 2 3 2 2 2 3 2 3 2 ...
## $ dr : int 1 1 2 2 2 1 2 2 1 2 ...
## $ suicide.s : int 0 0 0 1 0 0 3 0 4 0 ...
## $ suicide.hr : int 0 0 0 0 0 0 1 0 1 0 ...
## $ suicide.past: int 0 0 0 0 1 0 1 0 1 0 ...
## $ dur.interv : int NA 70 NA 105 NA NA 105 84 78 60 ...
## $ ed.bin : num 0 0 1 0 0 0 1 0 1 0 ...Création d’une table de travail pour croiser les variables ed.b (évitement du danger) et dep.cons (état dépressif)
## smp$dep.cons
## smp$ed.bin 0 1 <NA>
## 0 335 135 0
## 1 96 126 0
## <NA> 51 56 0Stockage de la table dans une variable temporaire et calcul des pourcentages plutôt que des valeurs absolues
Pourcentage de dépression selon que les détenus ont ou pas un haut niveau d’évitement du danger: On constate que les détenus à faible niveau d’évitement du danger (ed.bin=0) sont 2 fois moins déprimés (28,7%) que ceux ayant 1 haut niveau d’évitement du danger (56,7%)
prop.table(tab,1)*100## smp$dep.cons
## smp$ed.bin 0 1
## 0 71.27660 28.72340
## 1 43.24324 56.75676Pourcentage des détenus ayant 1 haut niveau de ed selon qu’ils sont ou pas déprimés: On constate la même chose
prop.table(tab,2)*100## smp$dep.cons
## smp$ed.bin 0 1
## 0 77.72622 51.72414
## 1 22.27378 48.27586Pour s’assurer que cette proportion n’est pas due au hasard, on fait un test chi-2. On constate que le p est très faible.
chisq.test(smp$ed.bin,smp$dep.cons,correct = FALSE)##
## Pearson's Chi-squared test
##
## data: smp$ed.bin and smp$dep.cons
## X-squared = 50.442, df = 1, p-value = 1.228e-12On peux également faire une test exact de Fisher: On constate que p est également très faible.
fisher.test(smp$ed.bin,smp$dep.cons)##
## Fisher's Exact Test for Count Data
##
## data: smp$ed.bin and smp$dep.cons
## p-value = 2.033e-12
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
## 2.303664 4.603460
## sample estimates:
## odds ratio
## 3.250819Chapitre11:Semaine3:Comparaison de 2 moyennes
Existe-t-il en moyenne, une différence significative d’age entre les détenus ayant un haut niveau d’évitement du danger et les autres ? Pour cela, on peux comparer la moyenne d’age des 2 groupes (haut niveau d’évitement et bas niveau d’évitement) Pour savoir si cette comparaison est pertinente, on effectue un test t de Student. Les conditions de validité de ce test sont : (au moins 30 sujets par groupe OU la variable suit une loi normale) ET (la variance de la variable est identique dans chaque groupe)
Comment se répartit l’age chez les détenus: Il y a une vague forme de loi normale
hist(smp$age)
Pour mieux vérifier la normalité, on peux faire une diagramme de normalité
qqnorm(smp$age);qqline(smp$age)
Pour vérifier l’égalité de la variance () ou l’écart type (sd), on utilise by
by(smp$age, smp$ed.bin, sd, na.rm = TRUE)## smp$ed.bin: 0
## [1] 13.38593
## --------------------------------------------------------
## smp$ed.bin: 1
## [1] 13.29636Après avoir vérifié les conditions de valité du test t, on effectue enfin ce test (on vérifie que les variances de chaque groupe sont proches avec l’option var.equal = TRUE)
t.test(smp$age, smp$ed.bin, var.equal = TRUE)##
## Two Sample t-test
##
## data: smp$age and smp$ed.bin
## t = 76.369, df = 1487, p-value < 2.2e-16
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## 37.58791 39.56972
## sample estimates:
## mean of x mean of y
## 38.8996236 0.3208092SI on ne peux pas utiliser de test t (n< 30 ET la variable ne suit pas de loi normale), alors on peux vérifier avec 1 test de Wilcox ou de Man Witney
wilcox.test(smp$age~smp$ed.bin)##
## Wilcoxon rank sum test with continuity correction
##
## data: smp$age by smp$ed.bin
## W = 56770, p-value = 0.06091
## alternative hypothesis: true location shift is not equal to 0Chapitre12:Semaine3:test de nullité d’un coef de corrélation
la distrib de x OU de y doit suivre une loi normale
cor.test(smp$age,smp$rs)##
## Pearson's product-moment correlation
##
## data: smp$age and smp$rs
## t = -6.02, df = 694, p-value = 2.825e-09
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## -0.2922516 -0.1509579
## sample estimates:
## cor
## -0.2227744Si ni x, ni y ne suivent une loi normale, alors on peux faire un test de nullité avec la méthode de spearman
cor.test(smp$age,smp$rs,method ="spearman")## Warning in cor.test.default(smp$age, smp$rs, method = "spearman"): Cannot
## compute exact p-value with ties##
## Spearman's rank correlation rho
##
## data: smp$age and smp$rs
## S = 68743000, p-value = 2.567e-09
## alternative hypothesis: true rho is not equal to 0
## sample estimates:
## rho
## -0.2233474Comparaison d’une moyenne à une référence: comparer la moyenne des ages à un age de référence
t.test(smp$age, mu=24)##
## One Sample t-test
##
## data: smp$age
## t = 31.672, df = 796, p-value < 2.2e-16
## alternative hypothesis: true mean is not equal to 24
## 95 percent confidence interval:
## 37.97618 39.82307
## sample estimates:
## mean of x
## 38.89962Comparaison d’une variable associée à une population mesurée en 2 temps (avant et après): ce sont les tests appariés
Mc Nemar: mcnmar.test(b.debut, b.fin) t de Student pour population appariée: t.test(b.debut, b.fin, paired=TRUE)
Lab4:Semaine3:Tests d’association graphique bivariés
Croisement de 2 variables 2 variables qual (tableau de contigence) 1 var qual et 1 var numerique ==> permet de faire par exemple des moyennes
tab <-table(smp$subst.cons, smp$abus)
tab##
## 0 1
## 0 441 140
## 1 131 80prop.table(tab)##
## 0 1
## 0 0.5568182 0.1767677
## 1 0.1654040 0.1010101prop.table(tab, margin=1) #propostion en ligne: les effectifs sont rapportés aux totaux ligne##
## 0 1
## 0 0.7590361 0.2409639
## 1 0.6208531 0.3791469prop.table(tab, margin=2) #proportion en colonne: les effectifs sont rapportés aux totaux colonne##
## 0 1
## 0 0.7709790 0.6363636
## 1 0.2290210 0.3636364xtab<-xtabs(~ subst.cons + abus, smp)
xtab## abus
## subst.cons 0 1
## 0 441 140
## 1 131 80barplot(xtab)
barplot(xtab, beside=TRUE, main="Proportion des détenus consommant des substances en fonction qu'ils ont subi des abus (0) ou pas (1)")
Test CHI2
chi2<-chisq.test(tab)
chi2##
## Pearson's Chi-squared test with Yates' continuity correction
##
## data: tab
## X-squared = 14.052, df = 1, p-value = 0.0001779Affichage des effectifs observés
chi2$observed##
## 0 1
## 0 441 140
## 1 131 80chi2$expected##
## 0 1
## 0 419.6111 161.38889
## 1 152.3889 58.61111Test de Fisher
fisher.test(tab)##
## Fisher's Exact Test for Count Data
##
## data: tab
## p-value = 0.0002193
## alternative hypothesis: true odds ratio is not equal to 1
## 95 percent confidence interval:
## 1.351339 2.728231
## sample estimates:
## odds ratio
## 1.921985Croisement age et subst.cons
head(smp$age)## [1] 31 49 50 47 23 34table(smp$subst.cons,useNA = "always")##
## 0 1 <NA>
## 587 212 0Décire la variable age en fonction des 2 modalités de la variable subst.cons (0 ou 1)
tapply(smp$age, smp$subst.cons, mean, na.rm=TRUE)## 0 1
## 41.97099 30.36967Pour réaliser un test t de student
t.test(smp$age[smp$subst.cons== 0],smp$age[smp$subst.cons== 1], var.equal = TRUE)#age pour lesquels la conso de substance == 0 et les ages pour lesquels la conso est à 1##
## Two Sample t-test
##
## data: smp$age[smp$subst.cons == 0] and smp$age[smp$subst.cons == 1]
## t = 11.785, df = 795, p-value < 2.2e-16
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## 9.668959 13.533684
## sample estimates:
## mean of x mean of y
## 41.97099 30.36967Ecriture plus rapide (comme avec xtabs)
t.test(age~subst.cons,smp) # age expliqué par la variable subst.cons##
## Welch Two Sample t-test
##
## data: age by subst.cons
## t = 15.24, df = 666.83, p-value < 2.2e-16
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## 10.10664 13.09600
## sample estimates:
## mean in group 0 mean in group 1
## 41.97099 30.36967Notation par formule pour faire une moyenne (mean)
aggregate(age~subst.cons,smp, mean)## subst.cons age
## 1 0 41.97099
## 2 1 30.36967boxplot(age~subst.cons,smp, mean)
boxplot(smp$age, smp$subst.cons)
## Installing package into '/home/david/R/x86_64-pc-linux-gnu-library/3.4'
## (as 'lib' is unspecified)gplots::lmplot2(lm(age~subst.cons, smp))
#gplots::plot.lm2(age~subst.cons, smp) #Cette fonction est dépréciéeEstimation de la durée d’interview selon que la personne est dépressive ou pas la variable durée d’interview : smp\(dur.interv état dépressif: smp\)dep.cons
aggregate(dur.interv~dep.cons,smp, median)## dep.cons dur.interv
## 1 0 60
## 2 1 65Donner la borne inférieure d’un intervalle de confiance à 95 % pour la corrélation linéaire (Pearson) entre les variables durée d’interview (dur.interv) et âge (age) (3 chiffres après la virgule).
cor.test(smp$age, smp$dur.interv)##
## Pearson's product-moment correlation
##
## data: smp$age and smp$dur.interv
## t = 2.3487, df = 745, p-value = 0.0191
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.01408787 0.15650345
## sample estimates:
## cor
## 0.08573358Nous souhaitons vérifier si la durée d’interview (dur.interv) diffère sensiblement selon que les individus ont déjà effectué une tentative de suicide dans le passé ou non (suicide.past) à l’aide d’un test de Wilcoxon. Quel degré de significativité ?
wilcox.test(smp$dur.interv~smp$suicide.past)##
## Wilcoxon rank sum test with continuity correction
##
## data: smp$dur.interv by smp$suicide.past
## W = 41892, p-value = 1.355e-07
## alternative hypothesis: true location shift is not equal to 0Quiz exercice3:Semaine3
Quiz exercice3bis (facultatif):Semaine3
Quelle est la valeur de la statistique de test (c’est-à-dire la valeur du « t », que l’on donnera pour simplifier en valeur absolue) pour un test de Student comparant les durées moyennes d’intervention (dur.interv) selon le diagnostic de dépression (dep.cons). Nous supposerons l’égalité des variances dans les deux groupes (3 chiffres après la virgule).
t.test(smp$dur.interv~smp$dep.cons,var.equal = TRUE)##
## Two Sample t-test
##
## data: smp$dur.interv by smp$dep.cons
## t = -5.2583, df = 747, p-value = 1.9e-07
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -10.457001 -4.771515
## sample estimates:
## mean in group 0 mean in group 1
## 58.92341 66.53767Le fisher test permet d’obtenir l’ODDS RATIO
Semaine3:Devoir évalué par les pairs (2/3)
Chapitre13:Semaine4:Régression linéaire simple
Chargement du fichier
Affichage d’un nuage de points
plot(smp.long$age,smp.long$dur.inter)
Avec jitter on fait apparaître les points qui sont masqués par d’autres points On peux également augmenter le jitter
plot(jitter(smp.long$age,factor = 10),jitter(smp.long$dur.inter,factor = 10))
Dessin de la droite de régression linéraire: durée = a + b * age
plot(jitter(smp.long$age,factor = 10),jitter(smp.long$dur.inter,factor = 10))
abline(lm(smp.long$dur.interv~smp.long$age))
Création du modéle de régression linéaie simple: il introduit du bruit à la droite: durée = a + b * age + bruit Pour savoir si b est statistiquement <> 0, il faut estimer le petit p du modèle
lm2<-lm(dur.interv~age,data = smp.long)
summary(lm2)##
## Call:
## lm(formula = dur.interv ~ age, data = smp.long)
##
## Residuals:
## Min 1Q Median 3Q Max
## -62.470 -14.402 -1.712 12.341 60.055
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 57.04091 2.22028 25.691 <2e-16 ***
## age 0.12625 0.05375 2.349 0.0191 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 19.57 on 745 degrees of freedom
## (52 observations deleted due to missingness)
## Multiple R-squared: 0.00735, Adjusted R-squared: 0.006018
## F-statistic: 5.516 on 1 and 745 DF, p-value: 0.0191Le modèle indique que a=57.05091 et b=0.12625 le p de a est <2e-16: c’est le résultat du test de nullité de a le p de b est 0.0191: c’est le résultat du test de nullité de b qui est p<5% ==> Donc a et b sont statistiquement <> 0
Chapitre14:Semaine4:Régression, corrélation et test t
Quand l’age du détenu augmente alors la durée de l’entretien augmente. Les 2 variables sont effectivement corrélées.
mod1<-lm(smp.long$dur.interv~smp.long$age)
summary(mod1)##
## Call:
## lm(formula = smp.long$dur.interv ~ smp.long$age)
##
## Residuals:
## Min 1Q Median 3Q Max
## -62.470 -14.402 -1.712 12.341 60.055
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 57.04091 2.22028 25.691 <2e-16 ***
## smp.long$age 0.12625 0.05375 2.349 0.0191 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 19.57 on 745 degrees of freedom
## (52 observations deleted due to missingness)
## Multiple R-squared: 0.00735, Adjusted R-squared: 0.006018
## F-statistic: 5.516 on 1 and 745 DF, p-value: 0.0191Le test de nullité présente p tres petit Est-ce que le test de la nullité de la corrélation entre la durée d’interview et l’age des détenus présente la meme résutlat ?
cor.test(smp.long$dur.interv,smp.long$age)##
## Pearson's product-moment correlation
##
## data: smp.long$dur.interv and smp.long$age
## t = 2.3487, df = 745, p-value = 0.0191
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
## 0.01408787 0.15650345
## sample estimates:
## cor
## 0.08573358On constate que le p est très faible et on peux rejeter la nullité de la corrélation. Donc les 2 vaiables sont corrélées et le p du test de nullité de la corrélation = strictement au p du test de non nullité
A noter que cor = b * ET(age)/ET(durée)
0.12625 * sd(smp.long$age, na.rm = TRUE)/sd(smp$dur.interv,na.rm = TRUE)Régrssion linéarie avec un variable binaire: la durée de l’interviwe est-elle liée à l’état dépressif ?
plot(smp.long$dep.cons,jitter(smp.long$dur.inter))
abline(lm(smp.long$dur.interv~smp.long$dep.cons))
La variation de durée des entretiens entre le groupe des déprimés et ceux des non déprimés
mod2<-lm(dur.interv~dep.cons, data =smp.long)
summary(mod2)##
## Call:
## lm(formula = dur.interv ~ dep.cons, data = smp.long)
##
## Residuals:
## Min 1Q Median 3Q Max
## -62.538 -13.923 1.077 12.077 61.077
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 58.9234 0.9041 65.171 < 2e-16 ***
## dep.cons 7.6143 1.4481 5.258 1.9e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 19.33 on 747 degrees of freedom
## (50 observations deleted due to missingness)
## Multiple R-squared: 0.03569, Adjusted R-squared: 0.0344
## F-statistic: 27.65 on 1 and 747 DF, p-value: 1.9e-07chapitre 14: la reg linéaire est une généralisation du test t de student
Régression linéaire entre 2 var quanti mais aussi entre 1 var quanti et une var binaire Par exemple :
Chapitre15:Semaine4:Régression linéaire multiple, analyse de variance
variable à expliquer: durée entretien ==> quantitative variable explicative: age (quant), dépression (binaire), abus de substance (binaire), schizophrnie (binaire) ##régression linéaire multiple
mod3<-lm(dur.interv~age+dep.cons+subst.cons+scz.cons, data = smp)
summary(mod3)##
## Call:
## lm(formula = dur.interv ~ age + dep.cons + subst.cons + scz.cons,
## data = smp)
##
## Residuals:
## Min 1Q Median 3Q Max
## -63.654 -14.522 -1.193 11.482 62.482
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 48.90105 2.62213 18.649 < 2e-16 ***
## age 0.22096 0.05708 3.871 0.000118 ***
## dep.cons 7.38932 1.44783 5.104 4.24e-07 ***
## subst.cons 5.25157 1.74318 3.013 0.002678 **
## scz.cons 2.27256 2.52323 0.901 0.368062
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 19.1 on 742 degrees of freedom
## (52 observations deleted due to missingness)
## Multiple R-squared: 0.05833, Adjusted R-squared: 0.05325
## F-statistic: 11.49 on 4 and 742 DF, p-value: 4.692e-09variable explicative catégorielle à plus de 2 classes
Si on veut utiliser une variable profession comme variable explicative, il faut la recoder en plusieurs variables binaires avec une des modalités définie commme modalité de référence (par exemple agriculteur) ==> R le fait automatiquement
mod4<-lm(dur.interv~age+dep.cons+subst.cons+scz.cons+prof, data = smp)
summary(mod4)##
## Call:
## lm(formula = dur.interv ~ age + dep.cons + subst.cons + scz.cons +
## prof, data = smp)
##
## Residuals:
## Min 1Q Median 3Q Max
## -63.280 -14.164 -1.337 10.959 63.184
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 62.79202 10.20779 6.151 1.26e-09 ***
## age 0.21289 0.05884 3.618 0.000317 ***
## dep.cons 7.36792 1.45840 5.052 5.53e-07 ***
## subst.cons 5.34589 1.76902 3.022 0.002599 **
## scz.cons 2.50439 2.54734 0.983 0.325863
## profartisan -11.48515 9.82936 -1.168 0.243005
## profautre -10.28748 10.33482 -0.995 0.319862
## profcadre -19.29636 10.38568 -1.858 0.063574 .
## profemploye -13.55809 9.76340 -1.389 0.165358
## profouvrier -14.01270 9.72111 -1.441 0.149880
## profprof.intermediaire -13.01926 9.96911 -1.306 0.191977
## profsans emploi -14.27866 9.71782 -1.469 0.142174
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 19.11 on 731 degrees of freedom
## (56 observations deleted due to missingness)
## Multiple R-squared: 0.06595, Adjusted R-squared: 0.05189
## F-statistic: 4.692 on 11 and 731 DF, p-value: 5.825e-07Choix de la modalité de référence
Pour changer la modalité de référence: on choisit celle qui est la plus fréquente
smp$prof<-relevel(smp$prof, ref="ouvrier")
mod5<-lm(dur.interv~age+dep.cons+subst.cons+scz.cons+prof, data = smp)
summary(mod5)##
## Call:
## lm(formula = dur.interv ~ age + dep.cons + subst.cons + scz.cons +
## prof, data = smp)
##
## Residuals:
## Min 1Q Median 3Q Max
## -63.280 -14.164 -1.337 10.959 63.184
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 48.77932 2.83938 17.180 < 2e-16 ***
## age 0.21289 0.05884 3.618 0.000317 ***
## dep.cons 7.36792 1.45840 5.052 5.53e-07 ***
## subst.cons 5.34589 1.76902 3.022 0.002599 **
## scz.cons 2.50439 2.54734 0.983 0.325863
## profagriculteur 14.01270 9.72111 1.441 0.149880
## profartisan 2.52755 2.48989 1.015 0.310381
## profautre 3.72522 3.99637 0.932 0.351567
## profcadre -5.28366 4.25567 -1.242 0.214798
## profemploye 0.45460 2.12659 0.214 0.830785
## profprof.intermediaire 0.99344 2.95809 0.336 0.737089
## profsans emploi -0.26596 1.87727 -0.142 0.887375
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 19.11 on 731 degrees of freedom
## (56 observations deleted due to missingness)
## Multiple R-squared: 0.06595, Adjusted R-squared: 0.05189
## F-statistic: 4.692 on 11 and 731 DF, p-value: 5.825e-07Estimation de l’effet globale d’une variable catégorielle à plus de 2 classes
drop1(mod5,.~.,test="F") # effet globale de la variable profession## Single term deletions
##
## Model:
## dur.interv ~ age + dep.cons + subst.cons + scz.cons + prof
## Df Sum of Sq RSS AIC F value Pr(>F)
## <none> 266846 4395.6
## age 1 4778.4 271624 4406.8 13.0899 0.0003173 ***
## dep.cons 1 9317.1 276163 4419.1 25.5233 5.527e-07 ***
## subst.cons 1 3333.6 270180 4402.8 9.1322 0.0025992 **
## scz.cons 1 352.8 267199 4394.6 0.9666 0.3258633
## prof 7 2295.5 269142 4388.0 0.8983 0.5071556
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Pas d’effet de la profession sur la durée de l’interview (p>5%)
Interaction entre 2 variables explicatives quand le modèle additif n’est plus pertinent
Si il y a un lien entre dépression et abus de substance de façon combinée sur la durée de l’entretien ? On peux le vérifier et le tester avec l’instruction *
mod6<-lm(dur.interv~age+dep.cons*subst.cons+scz.cons, data = smp)
summary(mod6)##
## Call:
## lm(formula = dur.interv ~ age + dep.cons * subst.cons + scz.cons,
## data = smp)
##
## Residuals:
## Min 1Q Median 3Q Max
## -62.032 -14.251 -1.163 11.472 62.313
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 49.51693 2.65788 18.630 < 2e-16 ***
## age 0.21728 0.05711 3.805 0.000154 ***
## dep.cons 6.15780 1.69775 3.627 0.000306 ***
## subst.cons 3.17244 2.29849 1.380 0.167931
## scz.cons 1.97233 2.53094 0.779 0.436059
## dep.cons:subst.cons 4.49688 3.24296 1.387 0.165963
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 19.08 on 741 degrees of freedom
## (52 observations deleted due to missingness)
## Multiple R-squared: 0.06077, Adjusted R-squared: 0.05443
## F-statistic: 9.588 on 5 and 741 DF, p-value: 7.024e-09Il y a une nouvelle ligne dep.cons:subst.cons. Au bout de cette ligne le p vaut 16% ! Donc statistiquement il n’y a pas d’interaction entre dépression et substance sur la durée de l’interview ==> donc quand un détenu est déprimé et consomme des substances, alors la durée est augmentée de 7 + 5 minutes Attention, l’ajout du terme d’interaction empêche d’interpréter les lignes dep et subst
Analyse de variance
Analyse de variance = régression linéaire multiple ou toutes les variables explicatives sont catégorielles (ou qualitatives) C’est un cas particulier qui n’a qu’un intérêt historique L’intérêt réside uniquement quand on veux comparer la moyenne d’une variable quanti à plus de 2 groupes: Exemple: comparer la durée moyenne d’interview en fonction des professions (pas de test t possible car il y a plus de 2 professions)
mod7<-lm(dur.interv~prof, data=smp) #Analyse de variance à 1 facteur
summary(mod7)##
## Call:
## lm(formula = dur.interv ~ prof, data = smp)
##
## Residuals:
## Min 1Q Median 3Q Max
## -61.731 -13.826 -1.731 12.947 58.912
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 61.7315 1.3359 46.211 <2e-16 ***
## profagriculteur 17.0185 9.9071 1.718 0.0863 .
## profartisan 2.0941 2.5033 0.837 0.4031
## profautre 2.4993 4.0755 0.613 0.5399
## profcadre -4.7750 4.3063 -1.109 0.2679
## profemploye 0.3220 2.1742 0.148 0.8823
## profprof.intermediaire 1.3440 3.0096 0.447 0.6553
## profsans emploi -0.6432 1.9168 -0.336 0.7373
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 19.63 on 735 degrees of freedom
## (56 observations deleted due to missingness)
## Multiple R-squared: 0.008295, Adjusted R-squared: -0.001149
## F-statistic: 0.8783 on 7 and 735 DF, p-value: 0.5231drop1(mod7,.~.,test='F') # effet globale de la variable profession## Single term deletions
##
## Model:
## dur.interv ~ prof
## Df Sum of Sq RSS AIC F value Pr(>F)
## <none> 283316 4432.1
## prof 7 2369.9 285686 4424.3 0.8783 0.5231Comment faire pour comparer un pourcentage quand on a plusieurs sous groupes. Par exemple comparer la prévalence de la dépression en fonction des professions:
chisq.test(smp$dep.cons, smp$prof)## Warning in chisq.test(smp$dep.cons, smp$prof): Chi-squared approximation
## may be incorrect##
## Pearson's Chi-squared test
##
## data: smp$dep.cons and smp$prof
## X-squared = 6.7205, df = 7, p-value = 0.4586Condition de validité du modèle de régression
3 conditions principales normalité du terme de bruit variance du bruit ne doit dépendre ni des valeurs de la variable à expliquer ni des valeurs des variables explicatives le bruit ne doit pas avoir de structure de corrélation interne
mod3<-lm(dur.interv~age+dep.cons+subst.cons+scz.cons, data = smp)
hist(resid(mod3), col="grey",main="normalité du bruit")
Chapitre16::Semaine4:Introduction à la régression logistique
Comment expliquer un risque suicidaire élevé par la durée de la peine, l’existence de mesures disciplinaires; l’existence d’abus. Les variables explicatives sont “emmelées”. Un modèle de régression linéaire multiple ne permet pas d’estimer le risque. haut risque de suicide = a + bdurée + cdiscipl + dabus + bruit cette équation ne fonctionne pas car la variable à expliquer est binaire Il faudrait une variable quanti. pour éviter ce problème on remplace la variable binaire à expliquer par: log (prob(HRsuicide=oui)/(1-prob(HRsuicide=oui))= a + bdurée + cdiscipl + dabus
Pour un premier modèle simple ou la seule variable explicative est abus log (prob(HRsuicide=oui)/(1-prob(HRsuicide=oui))= a + b*durée
mod1<-glm(suicide.hr~abus,data=smp, family = "binomial")
summary(mod1)##
## Call:
## glm(formula = suicide.hr ~ abus, family = "binomial", data = smp)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -0.8446 -0.6020 -0.6020 -0.6020 1.8959
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.6161 0.1154 -14.003 < 2e-16 ***
## abus 0.7688 0.1897 4.052 5.07e-05 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 760.21 on 752 degrees of freedom
## Residual deviance: 744.26 on 751 degrees of freedom
## (46 observations deleted due to missingness)
## AIC: 748.26
##
## Number of Fisher Scoring iterations: 4le p est très faible: au risque de 5%, il existe une association statistiquement significative entre des antécédents d’abus dans l’enfance et un haut risque suicidaire pour un détenu incarcéré le coéf b = 0.7688: est facile à interpréter uniquement avec une variable explicative binaire: exp(b) = exp(0.7688) = 2.16 = ODDSRATIO qui associe la variable explicative à la variable à expliquée. Cela ne fonctionne que si la variable à expliquer est suffisament rare (détenus à haut risque suicidaire). Dans ce cas, les antécédents d’abus multiplient par 2 le risque suicidaire
Calcul de l’ODDRATIO directement avec twoby2
exp(0.7688)## [1] 2.157176library(Epi)
twoby2(1-smp$suicide.hr, 1-smp$abus) # calcul de l'ODD-Ratio## 2 by 2 table analysis:
## ------------------------------------------------------
## Outcome : 0
## Comparing : 0 vs. 1
##
## 0 1 P(0) 95% conf. interval
## 0 63 90 0.4118 0.3366 0.4913
## 1 147 453 0.2450 0.2122 0.2810
##
## 95% conf. interval
## Relative Risk: 1.6807 1.3276 2.1276
## Sample Odds Ratio: 2.1571 1.4873 3.1287
## Conditional MLE Odds Ratio: 2.1547 1.4577 3.1764
## Probability difference: 0.1668 0.0837 0.2525
##
## Exact P-value: 1e-04
## Asymptotic P-value: 1e-04
## ------------------------------------------------------Quiz Exercice 4:Semaine4
exp(1.95)## [1] 7.028688Quiz Exercice 4bis (facultatif):Semaine4
Chapitre17:Semaine4:Régression logistique multiple
mod2<-glm(suicide.hr~abus+discip+duree, data =smp, family = "binomial")
summary(mod2)##
## Call:
## glm(formula = suicide.hr ~ abus + discip + duree, family = "binomial",
## data = smp)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.3200 -0.6655 -0.6012 -0.4997 2.0700
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.02462 0.49635 -0.050 0.960439
## abus 0.62289 0.22764 2.736 0.006213 **
## discip 0.52809 0.23767 2.222 0.026287 *
## duree -0.39862 0.11723 -3.400 0.000673 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 555.94 on 549 degrees of freedom
## Residual deviance: 533.26 on 546 degrees of freedom
## (249 observations deleted due to missingness)
## AIC: 541.26
##
## Number of Fisher Scoring iterations: 4Tous les p sont < 5% ==> les 3 variables explicatives sont statistiquement associées à un haut risque suicidaire et ce toute chose égale par ailleurs: antécédents d’abus est statistiquement associé à un haut risque suicidaire
signe des coefficients: 1 abus / 0 pas d’abus : coef positif ==> abus augmente le 1 mesure disciplinaire / 0 absence de mesures: coeff positif ==> mesure disciplinaire durée : de 1 à 5 : coeff négatif ==> 1 durée elevée diminue la probabilité d’être à haut risque suicidaire
amplitude des coefficients Il faut faire exp(coeff) Pour les variables binaires avec 1 comme facteur de risque on peux interpreter directement le exp: on interprete comme les ODDSRATIO: l’existence d’antécédents d’abus dans l’enfance multiplie par 2 le risque d’être à haut risque suicidaire. pour les variables non binaires: la durée : on a 0.7 ==> quand on passe d’un cran dans la durée d’incarcération (1 à 5), le niveau de haut risque suicidaire diminue de 30% (1-0.7)
exp(coefficients(mod2))## (Intercept) abus discip duree
## 0.9756803 1.8643147 1.6956873 0.6712485Ce qui change ou pas entre la régression linéaire multiple et la régression logistique multiple
les variables catégorielles à plus de 2 classes peuvent être utilisées pour connaître l’effet global d’une var catégorielle, on fait drop1(mod1,.~.,test=“Chisq”) pour trouver une synergie entre variable explicative: interaction(smp\(duree*smp\)discipl)
Conditions de validité de la régressison logistique
il faut au moins 5 à 10 evenements par variable explicative. Par exemple : explication de la variable “forme de schizophrenie” : 54 détenus de cette sorte utilisons les variables explicatives: age , trauma, prof. cela correspond à 1+1+ 8 classes recodée en 7 varaibles binaires (8 métiers) = 9 variables explicatives 10 événements par variable explicative: 910 = 90 mais on a 54 détenus 5 évenements par variable explicative: 95 = 45 et on a 54 détenus: donc on est bon
Quizz
Un modèle de régression logistique multiple est réalisé pour étudier les liens entre l’obésité (Oui/Non) en fonction des variables explicatives suivantes : âge (continue) ; sexe (F/H) ; origine géographique (5 classes). La condition de validité d’un tel modèle est elle vérifiée sachant qu’il y a 20 personnes obèse dans l’échantillon : On a 1 + 1 + 4 = 6 variables explicatives et il faut 5 à 10 evts par variable explicatives soit au minimum 6*5 = 30 evts et on a que 20 obèses
Lab5:Semaine4:ANOVA Régression linéaire ; régression logistique
Chargement du fichier smp
Extraction d’un subset
head(subset(smp,prof=="sans emploi" | prof=="cadre" | prof=="prof.intermediaire", c(age,n.enfant, prof)))## age n.enfant prof
## 3 50 2 prof.intermediaire
## 5 23 1 sans emploi
## 11 31 0 prof.intermediaire
## 17 60 2 prof.intermediaire
## 23 32 0 sans emploi
## 27 32 1 sans emploismpb<-subset(smp,prof=="sans emploi" | prof=="cadre" | prof=="prof.intermediaire", c(age,n.enfant, prof))
summary(smpb)## age n.enfant prof
## Min. :19.00 Min. : 0.000 sans emploi :222
## 1st Qu.:27.00 1st Qu.: 0.000 prof.intermediaire: 58
## Median :36.00 Median : 1.000 cadre : 24
## Mean :38.42 Mean : 1.648 ouvrier : 0
## 3rd Qu.:47.25 3rd Qu.: 3.000 agriculteur : 0
## Max. :83.00 Max. :13.000 artisan : 0
## NA's :11 (Other) : 0#Recalcul des niveaux de la variable prof
smpb$prof<-factor(smpb$prof)
table(smp$prof)##
## ouvrier agriculteur artisan
## 227 6 90
## autre cadre employe
## 31 24 135
## prof.intermediaire sans emploi
## 58 222table(smpb$prof)##
## cadre prof.intermediaire sans emploi
## 24 58 222Résumer le nb enfant moyen en fonction de la profession:
aggregate(n.enfant~prof, smpb,mean)## prof n.enfant
## 1 cadre 2.166667
## 2 prof.intermediaire 2.107143
## 3 sans emploi 1.469484boxplot(n.enfant~prof, smpb)
Création d’une ANOVA avec la commande lm
help(lm)
m<-lm(n.enfant~prof, data = smpb)
m #correspond ##
## Call:
## lm(formula = n.enfant ~ prof, data = smpb)
##
## Coefficients:
## (Intercept) profprof.intermediaire profsans emploi
## 2.16667 -0.05952 -0.69718drop1(m, test="F") #tableau d'analyse de variance## Single term deletions
##
## Model:
## n.enfant ~ prof
## Df Sum of Sq RSS AIC F value Pr(>F)
## <none> 947.74 349.96
## prof 2 25.05 972.79 353.60 3.8325 0.02276 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Explication de n.enfant par une vairable numérique
m<-lm(n.enfant~age, data = smpb)
m #affiche intercept cad le terme d'ordonnée à l'origine et age qui représente la pente##
## Call:
## lm(formula = n.enfant ~ age, data = smpb)
##
## Coefficients:
## (Intercept) age
## -1.00902 0.06849summary(m) #affichage des tests associés##
## Call:
## lm(formula = n.enfant ~ age, data = smpb)
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.2646 -0.9087 -0.2511 0.5708 9.0094
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.009022 0.262696 -3.841 0.00015 ***
## age 0.068488 0.006357 10.773 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.546 on 291 degrees of freedom
## (11 observations deleted due to missingness)
## Multiple R-squared: 0.2851, Adjusted R-squared: 0.2827
## F-statistic: 116.1 on 1 and 291 DF, p-value: < 2.2e-16Création du modele à partir d’un subset
m<-lm(n.enfant~age, smp, subset = prof=="sans emploi" | prof=="cadre" | prof=="prof.intermediaire")
summary(m)##
## Call:
## lm(formula = n.enfant ~ age, data = smp, subset = prof == "sans emploi" |
## prof == "cadre" | prof == "prof.intermediaire")
##
## Residuals:
## Min 1Q Median 3Q Max
## -4.2646 -0.9087 -0.2511 0.5708 9.0094
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.009022 0.262696 -3.841 0.00015 ***
## age 0.068488 0.006357 10.773 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.546 on 291 degrees of freedom
## (17 observations deleted due to missingness)
## Multiple R-squared: 0.2851, Adjusted R-squared: 0.2827
## F-statistic: 116.1 on 1 and 291 DF, p-value: < 2.2e-16coef(m)## (Intercept) age
## -1.00902159 0.06848829coef(m)[2] #pente## age
## 0.06848829coef(m)["age"] #pente = coef de régression assoicé à la variable age## age
## 0.06848829confint(m) # intervalles de confiance à 95% par défaut## 2.5 % 97.5 %
## (Intercept) -1.52604619 -0.49199700
## age 0.05597581 0.08100076anova(m) # tableau d'analyse de variance associé à la régression## Analysis of Variance Table
##
## Response: n.enfant
## Df Sum Sq Mean Sq F value Pr(>F)
## age 1 277.35 277.35 116.05 < 2.2e-16 ***
## Residuals 291 695.44 2.39
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1predict(m, data.frame(age=c(20,30,40)), interval = "confidence") #prédiction à partir du modèle à partir de valeurs non observés avec ajout des intervalles de confiance via l'option intervall## fit lwr upr
## 1 0.3607441 0.06588459 0.6556037
## 2 1.0456270 0.83652266 1.2547313
## 3 1.7305099 1.55212966 1.9088901fit valeurs prédites lwr et upr correspondent aux intervalles de confiance de la prédiction
Régression logistique avec variable binaire
smp$n.enfant.bin<- ifelse(smp$n.enfant>2,1,0)
table(smp$n.enfant)##
## 0 1 2 3 4 5 6 7 8 9 10 11 13
## 214 220 125 101 55 31 7 7 7 2 2 1 1table(smp$n.enfant.bin)##
## 0 1
## 559 214commange glm
help(glm)
m<-glm(n.enfant.bin~age,smp,family=binomial("logit")) #création d'un modèle de régresion logistique avec comme échelle de lien le logit
summary(m) #affichage d'un résumé des coeffs de régression ##
## Call:
## glm(formula = n.enfant.bin ~ age, family = binomial("logit"),
## data = smp)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.8551 -0.7525 -0.5326 0.8763 2.1301
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -3.827089 0.312803 -12.235 <2e-16 ***
## age 0.069487 0.007016 9.904 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 912.06 on 772 degrees of freedom
## Residual deviance: 794.16 on 771 degrees of freedom
## (26 observations deleted due to missingness)
## AIC: 798.16
##
## Number of Fisher Scoring iterations: 4Quiz Exercice4:Semaine4
Pour réaliser une ANOVA à un facteur entre les variables y (réponse) et z (facteur) dans le data frame d, nous utilisons la commande
Nous considérons l’âge (age) des individus ayant 4 enfants ou plus (n.enfant) et dont la catégorie socio-professionnelle (prof) figure parmi les modalités suivantes : « sans emploi », « ouvrier », « cadre » et « employé ». Pour ce sous-ensemble de l’échantillon du data frame smp, le rapport entre les deux variances les plus extrêmes dans ces 4 groupes est :
mysubset<-subset(smp,prof=="sans emploi" | prof=="cadre" | prof=="ouvrier" | prof=="employe" & n.enfant>=4, c(age,prof)) #extraction du subset
table(mysubset$prof)##
## ouvrier agriculteur artisan
## 227 0 0
## autre cadre employe
## 0 24 16
## prof.intermediaire sans emploi
## 0 222head(mysubset)## age prof
## 4 47 ouvrier
## 5 23 sans emploi
## 6 34 ouvrier
## 9 42 ouvrier
## 10 45 ouvrier
## 16 67 ouvrieraggregate(age~prof,mysubset,var)## prof age
## 1 ouvrier 136.6297
## 2 cadre 159.0362
## 3 employe 135.7167
## 4 sans emploi 199.8942199.8942/135.7167## [1] 1.472878>1
Nous souhaitons réaliser une ANOVA à un facteur en considérant l’âge (age) comme variable réponse, et la taille de la fratrie (n.fratrie) recodée en 3 classes (0-2, 3-4, 5+) comme variable explicative. Les bornes des intervalles sont inclues pour chacune des trois classes. Indiquer le résultat du test F de Fisher-Snedecor d’égalité des moyennes :
smp$n.fratrie.cat<-factor(smp$n.fratrie)
table(smp$n.fratrie.cat) # on voit que 21 factors ont été crées: on va regrouper en 3 factors: (0-2, 3-4, 5+)##
## 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 16 17 18
## 76 86 111 129 90 75 63 41 33 31 22 14 10 4 3 3 1 2
## 19 20 21
## 3 1 1levels(smp$n.fratrie.cat)## [1] "0" "1" "2" "3" "4" "5" "6" "7" "8" "9" "10" "11" "12" "13"
## [15] "14" "16" "17" "18" "19" "20" "21"nlevels(smp$n.fratrie.cat)## [1] 21levels(smp$n.fratrie.cat)[1:3]<-"0-2"
levels(smp$n.fratrie.cat)[2:3]<-"3-4"
levels(smp$n.fratrie.cat)[3:18]<-"5+"
table(smp$n.fratrie.cat)##
## 0-2 3-4 5+
## 273 219 307m<-glm(n.fratrie.cat~age, smp,family=binomial("logit")) #création d'un modèle de régresion logistique avec comme échelle de lien le logit
summary(m) #affichage d'un résumé des coeffs de régression ##
## Call:
## glm(formula = n.fratrie.cat ~ age, family = binomial("logit"),
## data = smp)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.5735 -1.4184 0.8602 0.9194 1.1104
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 1.113956 0.233148 4.778 1.77e-06 ***
## age -0.011496 0.005588 -2.057 0.0397 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1021.8 on 796 degrees of freedom
## Residual deviance: 1017.6 on 795 degrees of freedom
## (2 observations deleted due to missingness)
## AIC: 1021.6
##
## Number of Fisher Scoring iterations: 4anova(m)## Analysis of Deviance Table
##
## Model: binomial, link: logit
##
## Response: n.fratrie.cat
##
## Terms added sequentially (first to last)
##
##
## Df Deviance Resid. Df Resid. Dev
## NULL 796 1021.8
## age 1 4.2214 795 1017.6Nous nous intéressons à la relation entre la variable séparation (separation) et l’âge (age) des individus, que l’on modélise à l’aide d’une régression logistique. Donner la borne inférieure de l’intervalle de confiance à 95 % pour l’odds-ratio (3 chiffres après la virgule).
m<-glm(separation~age, data = smp, family = "binomial")
summary(m)##
## Call:
## glm(formula = separation ~ age, family = "binomial", data = smp)
##
## Deviance Residuals:
## Min 1Q Median 3Q Max
## -1.1529 -1.0567 -0.9597 1.2748 1.5691
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) 0.187954 0.224813 0.836 0.4031
## age -0.012936 0.005536 -2.337 0.0195 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1073.4 on 787 degrees of freedom
## Residual deviance: 1067.9 on 786 degrees of freedom
## (11 observations deleted due to missingness)
## AIC: 1071.9
##
## Number of Fisher Scoring iterations: 4fitted(m)## 1 2 3 4 5 6 7
## 0.4469338 0.3903316 0.3872576 0.3965057 0.4726329 0.4373618 0.4694097
## 8 9 10 11 13 14 15
## 0.3811365 0.4120820 0.4027130 0.4469338 0.4790858 0.4183640 0.3452553
## 16 17 18 20 21 22 23
## 0.3365360 0.3570442 0.3481854 0.4565452 0.4823150 0.4501335 0.4437384
## 24 25 26 27 28 29 30
## 0.4469338 0.4629713 0.4120820 0.4437384 0.4183640 0.4152195 0.4597566
## 31 32 33 34 35 36 37
## 0.4694097 0.4246727 0.4215151 0.4310059 0.4533374 0.4152195 0.4310059
## 38 39 40 41 42 43 44
## 0.4152195 0.4790858 0.4790858 0.3996053 0.4758584 0.4790858 0.4341812
## 45 46 47 48 49 50 51
## 0.4027130 0.4246727 0.4855457 0.4790858 0.4597566 0.4183640 0.4215151
## 52 53 54 55 56 57 58
## 0.3965057 0.4694097 0.4310059 0.4215151 0.4758584 0.4246727 0.4278363
## 59 60 61 62 63 64 65
## 0.4533374 0.4726329 0.4310059 0.4120820 0.3690087 0.4565452 0.4310059
## 66 67 68 69 70 71 72
## 0.4246727 0.4089515 0.4533374 0.3452553 0.4661891 0.3841925 0.4341812
## 73 74 75 76 77 78 79
## 0.4501335 0.4278363 0.4629713 0.4310059 0.3630053 0.4437384 0.4501335
## 80 81 82 83 84 85 86
## 0.4629713 0.4597566 0.4726329 0.4694097 0.4215151 0.4089515 0.4215151
## 87 88 89 90 91 92 93
## 0.4629713 0.4058284 0.4278363 0.4183640 0.4694097 0.3996053 0.4629713
## 94 95 96 97 98 99 100
## 0.4246727 0.4278363 0.4501335 0.4215151 0.4310059 0.4215151 0.4565452
## 101 102 103 104 105 106 107
## 0.4597566 0.3841925 0.3934144 0.3965057 0.4152195 0.4341812 0.4661891
## 108 109 110 111 112 113 114
## 0.4469338 0.4058284 0.4183640 0.4533374 0.4373618 0.3903316 0.3660018
## 115 116 117 118 119 120 121
## 0.4405477 0.4341812 0.4437384 0.4373618 0.3996053 0.4027130 0.4469338
## 122 123 124 126 127 128 129
## 0.4120820 0.3934144 0.4373618 0.3452553 0.3872576 0.3780900 0.3903316
## 130 131 132 133 134 135 136
## 0.3780900 0.4278363 0.4120820 0.3720259 0.4437384 0.4405477 0.4183640
## 137 138 139 140 141 142 143
## 0.4533374 0.4437384 0.4726329 0.3540800 0.4215151 0.4501335 0.4278363
## 144 145 146 147 148 149 150
## 0.4501335 0.4215151 0.3903316 0.4058284 0.4183640 0.3690087 0.4089515
## 151 152 153 154 155 156 157
## 0.4597566 0.4790858 0.4058284 0.3872576 0.3872576 0.4823150 0.4278363
## 158 159 160 161 162 163 164
## 0.4120820 0.4597566 0.4758584 0.4661891 0.4823150 0.4790858 0.4855457
## 165 166 167 168 169 170 171
## 0.4661891 0.4694097 0.3903316 0.4694097 0.4629713 0.4341812 0.4758584
## 172 173 174 175 176 177 178
## 0.4694097 0.4726329 0.3996053 0.4629713 0.4152195 0.3841925 0.4823150
## 179 180 181 182 183 184 185
## 0.4501335 0.4278363 0.3903316 0.4565452 0.4565452 0.3841925 0.4183640
## 186 187 188 189 190 191 192
## 0.4405477 0.4661891 0.4533374 0.4183640 0.4089515 0.4341812 0.3872576
## 193 194 195 196 197 198 199
## 0.4058284 0.4341812 0.4694097 0.4089515 0.4629713 0.4027130 0.4120820
## 200 201 202 203 204 205 206
## 0.4027130 0.3934144 0.4027130 0.4373618 0.4469338 0.4183640 0.4758584
## 207 208 209 210 211 212 213
## 0.4120820 0.4246727 0.4246727 0.4183640 0.3996053 0.4629713 0.4533374
## 214 215 216 217 218 219 220
## 0.4661891 0.4183640 0.4089515 0.4565452 0.4533374 0.4437384 0.4565452
## 221 222 223 224 225 226 227
## 0.3660018 0.4469338 0.3250819 0.4405477 0.4694097 0.4758584 0.4661891
## 228 229 230 231 232 233 234
## 0.4629713 0.3811365 0.4405477 0.4246727 0.4215151 0.4152195 0.3811365
## 235 236 237 238 239 240 241
## 0.4405477 0.4215151 0.3600193 0.4405477 0.3872576 0.3630053 0.4726329
## 242 243 244 245 246 247 248
## 0.4152195 0.4089515 0.4120820 0.4758584 0.3660018 0.4152195 0.4501335
## 249 250 251 252 253 254 255
## 0.3394304 0.3903316 0.3996053 0.4565452 0.3600193 0.4341812 0.4058284
## 256 257 258 259 260 261 262
## 0.2919880 0.4373618 0.3903316 0.3570442 0.3690087 0.3996053 0.3511271
## 263 264 265 266 267 268 269
## 0.4152195 0.4597566 0.3780900 0.3934144 0.3394304 0.3394304 0.3720259
## 270 271 272 273 274 275 276
## 0.3540800 0.4089515 0.3750530 0.4246727 0.3841925 0.3841925 0.3872576
## 277 278 279 280 281 282 283
## 0.3690087 0.3780900 0.3903316 0.4152195 0.4058284 0.3452553 0.4120820
## 284 285 286 287 288 289 290
## 0.3811365 0.3222501 0.4089515 0.4501335 0.4437384 0.4089515 0.4661891
## 291 292 293 294 295 296 297
## 0.4597566 0.4661891 0.3811365 0.4215151 0.4120820 0.3600193 0.3996053
## 298 299 300 301 302 303 304
## 0.3511271 0.3872576 0.4694097 0.4089515 0.4437384 0.3365360 0.4565452
## 305 306 307 308 309 310 311
## 0.4058284 0.4855457 0.4823150 0.4726329 0.4629713 0.4565452 0.4469338
## 312 313 314 315 316 317 318
## 0.4120820 0.3660018 0.4501335 0.4310059 0.3780900 0.4405477 0.4661891
## 319 320 321 322 323 324 325
## 0.4758584 0.4120820 0.4661891 0.4437384 0.4726329 0.4027130 0.3934144
## 326 327 328 329 330 331 332
## 0.4341812 0.4278363 0.4246727 0.4694097 0.3965057 0.3540800 0.4246727
## 333 334 335 336 337 338 339
## 0.4597566 0.4597566 0.4629713 0.4501335 0.3965057 0.4278363 0.4501335
## 340 341 342 343 344 345 346
## 0.4152195 0.4533374 0.4278363 0.4565452 0.3965057 0.4629713 0.3872576
## 347 348 349 350 351 352 353
## 0.4726329 0.3570442 0.4278363 0.3934144 0.4152195 0.4565452 0.3750530
## 355 356 357 358 359 360 361
## 0.4405477 0.4469338 0.4661891 0.3394304 0.4629713 0.4533374 0.4533374
## 362 363 364 365 366 367 368
## 0.3780900 0.4694097 0.3934144 0.4183640 0.3965057 0.4183640 0.4152195
## 369 370 371 372 373 374 375
## 0.3750530 0.4661891 0.4310059 0.4058284 0.4437384 0.4597566 0.4469338
## 376 377 378 379 380 381 382
## 0.4373618 0.4373618 0.3250819 0.4823150 0.3750530 0.4215151 0.3872576
## 383 384 385 386 387 388 389
## 0.4310059 0.4278363 0.4089515 0.4565452 0.4790858 0.4341812 0.4310059
## 390 391 392 393 394 395 396
## 0.3780900 0.4310059 0.4246727 0.3394304 0.3511271 0.4246727 0.4694097
## 397 398 399 400 401 402 403
## 0.3903316 0.4790858 0.4373618 0.4533374 0.4310059 0.4533374 0.4405477
## 404 405 406 407 408 409 410
## 0.4373618 0.3660018 0.3423369 0.4661891 0.4310059 0.4469338 0.3750530
## 411 412 413 414 415 416 417
## 0.3903316 0.4120820 0.4501335 0.4823150 0.4726329 0.4790858 0.4726329
## 418 419 420 421 422 423 424
## 0.4215151 0.4027130 0.4533374 0.4790858 0.3750530 0.3082888 0.4726329
## 425 426 427 428 429 430 431
## 0.4437384 0.3630053 0.3903316 0.4629713 0.4183640 0.3841925 0.3511271
## 432 433 434 435 436 437 438
## 0.4027130 0.4152195 0.4501335 0.3811365 0.4823150 0.4310059 0.4373618
## 439 440 441 442 443 444 445
## 0.4341812 0.4501335 0.3996053 0.3027992 0.3394304 0.4855457 0.4152195
## 446 447 448 449 450 451 452
## 0.3841925 0.4629713 0.3690087 0.4405477 0.4215151 0.3222501 0.4027130
## 453 454 455 456 457 458 459
## 0.3600193 0.4790858 0.4152195 0.4089515 0.3720259 0.4629713 0.3903316
## 460 461 462 463 464 465 466
## 0.4533374 0.4629713 0.4565452 0.3082888 0.3540800 0.3481854 0.4501335
## 467 468 469 470 471 472 474
## 0.3903316 0.3934144 0.4027130 0.4437384 0.3690087 0.3934144 0.3194313
## 475 476 477 478 479 480 481
## 0.4405477 0.3166256 0.3750530 0.4597566 0.3903316 0.4027130 0.4597566
## 482 483 484 486 487 488 489
## 0.3780900 0.3511271 0.3750530 0.3690087 0.3570442 0.4405477 0.4373618
## 490 491 492 493 494 495 496
## 0.4437384 0.4058284 0.3903316 0.3996053 0.3365360 0.4215151 0.3600193
## 497 498 499 500 501 502 503
## 0.3481854 0.2973652 0.4246727 0.3630053 0.4120820 0.3194313 0.3934144
## 504 505 506 507 508 509 510
## 0.4152195 0.4565452 0.4058284 0.4027130 0.3996053 0.3872576 0.4597566
## 511 512 513 514 515 516 517
## 0.3690087 0.3996053 0.4120820 0.4661891 0.4726329 0.4629713 0.4855457
## 518 519 520 521 522 523 524
## 0.4694097 0.4694097 0.4437384 0.4726329 0.4694097 0.4405477 0.4790858
## 525 526 527 528 529 530 531
## 0.4405477 0.4152195 0.4694097 0.4469338 0.4855457 0.4661891 0.3841925
## 532 533 534 535 536 537 538
## 0.4215151 0.4758584 0.4823150 0.4501335 0.4373618 0.4565452 0.4823150
## 539 540 541 542 543 544 545
## 0.4823150 0.4405477 0.4694097 0.4437384 0.4278363 0.4661891 0.4694097
## 546 547 548 549 550 551 552
## 0.4533374 0.4855457 0.4278363 0.3690087 0.3903316 0.3570442 0.4533374
## 553 554 555 556 557 558 559
## 0.4758584 0.4823150 0.3903316 0.4405477 0.4501335 0.4533374 0.4661891
## 560 561 562 563 564 565 566
## 0.3511271 0.4152195 0.4405477 0.4058284 0.3570442 0.4694097 0.4694097
## 567 568 569 570 571 572 573
## 0.4405477 0.4597566 0.4027130 0.4405477 0.4058284 0.4726329 0.4726329
## 574 575 576 577 578 579 580
## 0.4341812 0.4310059 0.4565452 0.4694097 0.4597566 0.4597566 0.4565452
## 581 582 583 584 585 586 587
## 0.4597566 0.4183640 0.3811365 0.4855457 0.4469338 0.4790858 0.4405477
## 588 589 590 591 592 593 594
## 0.4726329 0.4501335 0.4726329 0.4469338 0.3934144 0.4694097 0.4694097
## 595 596 597 598 599 600 601
## 0.4629713 0.4437384 0.4533374 0.4246727 0.4726329 0.3872576 0.4629713
## 602 603 604 605 606 607 608
## 0.3965057 0.4246727 0.4694097 0.4694097 0.4855457 0.4661891 0.4469338
## 609 610 611 612 613 614 615
## 0.4405477 0.4629713 0.4246727 0.4726329 0.4278363 0.4855457 0.3903316
## 616 617 618 619 620 621 622
## 0.4405477 0.4501335 0.4246727 0.4501335 0.4629713 0.4597566 0.4790858
## 623 624 625 626 627 628 629
## 0.4469338 0.4855457 0.4629713 0.4565452 0.3903316 0.4341812 0.4661891
## 630 631 632 633 634 635 636
## 0.4437384 0.4597566 0.4823150 0.4501335 0.4661891 0.4790858 0.3750530
## 637 638 639 640 641 642 643
## 0.4597566 0.4758584 0.4215151 0.4790858 0.3750530 0.3903316 0.4726329
## 644 645 646 647 648 649 650
## 0.4310059 0.3600193 0.3872576 0.4694097 0.3965057 0.4120820 0.4152195
## 651 652 653 654 655 656 657
## 0.4405477 0.3996053 0.4726329 0.4855457 0.4215151 0.4246727 0.4183640
## 658 659 660 661 662 663 664
## 0.4215151 0.4183640 0.4058284 0.4629713 0.3934144 0.3965057 0.4726329
## 665 666 667 668 669 670 671
## 0.4661891 0.4823150 0.4027130 0.4058284 0.3660018 0.4215151 0.3720259
## 672 673 674 675 676 677 678
## 0.4855457 0.4373618 0.4565452 0.4405477 0.4855457 0.4405477 0.4597566
## 679 680 681 682 683 684 685
## 0.3996053 0.3965057 0.4758584 0.4597566 0.4629713 0.3811365 0.3690087
## 686 687 688 689 690 691 692
## 0.4058284 0.3481854 0.4373618 0.4152195 0.4246727 0.4278363 0.3630053
## 693 694 695 696 697 698 699
## 0.4278363 0.4694097 0.3570442 0.4629713 0.4790858 0.3811365 0.4823150
## 700 701 703 704 705 706 707
## 0.4278363 0.4437384 0.3630053 0.3903316 0.4437384 0.4278363 0.3996053
## 708 709 710 711 712 713 714
## 0.3872576 0.4058284 0.3965057 0.4278363 0.4246727 0.3872576 0.3690087
## 715 716 717 718 719 720 721
## 0.4501335 0.4373618 0.4089515 0.3720259 0.4089515 0.4469338 0.3720259
## 722 723 724 726 727 728 729
## 0.4152195 0.3336537 0.4027130 0.4120820 0.3250819 0.4246727 0.3996053
## 730 731 732 733 734 735 736
## 0.3423369 0.3841925 0.3660018 0.3660018 0.3250819 0.4183640 0.4089515
## 737 738 739 740 741 742 743
## 0.3250819 0.3934144 0.4373618 0.3307838 0.4089515 0.4341812 0.3511271
## 744 745 746 747 748 749 750
## 0.4373618 0.3841925 0.3934144 0.4310059 0.4058284 0.3903316 0.3166256
## 751 752 753 754 755 756 757
## 0.4855457 0.3690087 0.3660018 0.3423369 0.3811365 0.3082888 0.4533374
## 758 759 760 761 762 763 764
## 0.4278363 0.4027130 0.4183640 0.3222501 0.4597566 0.3690087 0.4341812
## 765 766 767 768 769 770 771
## 0.4501335 0.4278363 0.4501335 0.4183640 0.3750530 0.4629713 0.3934144
## 772 773 774 775 776 777 778
## 0.2919880 0.4437384 0.4758584 0.3934144 0.3365360 0.3630053 0.4278363
## 779 780 781 782 783 784 785
## 0.4694097 0.4373618 0.4215151 0.4246727 0.4215151 0.3690087 0.4341812
## 786 787 788 789 792 793 794
## 0.4629713 0.3279265 0.3336537 0.4120820 0.4629713 0.3872576 0.4597566
## 796 797 798 799
## 0.4058284 0.4469338 0.4246727 0.3250819Quiz Exercice4 (Facultatif):Semaine4
À partir d’un modèle de régression linéaire dont les résultats ont été enregistrés dans une variable appelée mod, l’instruction fitted(mod) permet de fournir :
Les valeurs prédites de la variable réponse pour les valeurs observées de la variable explicative
Pour obtenir un intervalle de confiance à 90 % pour les coefficients d’un modèle de régression logistique stocké dans une variable appelée mod, on utilise la commande :
confint(mod, level=0.90)