5.3 Mise en œuvre dans R

Pour calculer le khi-deux entre deux variables qualitatives, nous utilisons la fonction de base :

chisq.test(x = ..., y = ...) qui renvoie le nombre de degrés de liberté, les valeurs du khi-deux et de p.

# Importation du csv
dataHLM <- read.csv("data/bivariee/hlm.csv")
# Calcul du khi-deux avec la fonction de base chisq.test
chisq.test(x = dataHLM$Taille, y = dataHLM$Periode)
## 
##  Pearson's Chi-squared test
## 
## data:  dataHLM$Taille and dataHLM$Periode
## X-squared = 63.543, df = 12, p-value = 5.063e-09

Pour la construction du tableau de contingence, deux options sont possibles dépendamment de la structure de votre tableau de données. Premier cas de figure : votre tableau comprend une ligne par observation avec les différentes modalités dans deux colonnes (ici Periode et Taille). Dans la syntaxe ci-dessous, pour chacune des deux variables qualitatives, nous créons un facteur afin de spécifier un intitulé à chaque modalité (factor(levels =c(....), labels = c(..)). Puis, nous utilisons la fonction CrossTable du package gmodels. Pour obtenir les fréquences théoriques, les contributions locales au khi-deux et les déviations, nous spécifions les options suivantes : expected=TRUE, chisq=TRUE, resid=TRUE.

library("gmodels")
# Premiers enregistrements du tableau
head(dataHLM)
##   Periode Taille
## 1       5      1
## 2       5      1
## 3       5      2
## 4       5      1
## 5       5      1
## 6       5      2
# La variable Periode comprend 5 modalités (de 1 à 5)
table(dataHLM$Periode)
## 
##  1  2  3  4  5 
## 56 48 48 47 80
# La variable Taille comprend 4 modalités (de 1 à 4)
table(dataHLM$Taille)
## 
##  1  2  3  4 
## 80 53 99 47
# Création d'un facteur pour les cinq modalités de la période de construction
dataHLM$Periode <- factor(dataHLM$Periode, 
                          levels = c(1,2,3,4,5), 
                          labels = c("<1975", 
                                     "1975-1979", 
                                     "1980-1984", 
                                     "1985-1989", 
                                     "1990-1994"))
# Création d'un facteur pour les quatre modalités de la taille des habitations
dataHLM$Taille <- factor(dataHLM$Taille, 
                         levels = c(1,2,3,4), 
                         labels = c("<25 log.", 
                                    "25-49", 
                                    "50-99", 
                                    "100 et +"))
# Pour construire un tableau de contingence, nous utilisons 
# la fonction CrossTable du package gmodels. 
# Les deux lignes ci-dessous sont mises en commentaire pour ne pas répéter le tableau.
# CrossTable(x=dataHLM$Taille, y=dataHLM$Periode, digits = 2,
#           expected=TRUE, chisq = TRUE, resid = TRUE, format="SPSS")

Deuxième cas de figure : vous disposez déjà d’un tableau de contingence, soit les fréquences observées (\(f_{ij}\)). Nous n’utilisons donc pas la fonction CrossTable, mais directement la fonction chisq.test.

# Importation des données
df1 <-  read.csv("data/bivariee/data_transport.csv", stringsAsFactors = FALSE)
df1 # Visualisation du tableau
##                                     ModeTransport  Homme  Femme
## 1 Automobile, camion ou fourgonnette - conducteur 689400 561830
## 2   Automobile, camion ou fourgonnette - passager  21315  40010
## 3                             Transport en commun 181435 238330
## 4                                          A pied  43715  54360
## 5                                      Bicyclette  24295  13765
## 6                                     Autre moyen   8395   6970
Matrice <- as.matrix(df1[, c("Homme","Femme")])
dimnames(Matrice) <- list(unique(df1$ModeTransport), Sexe=c("Homme","Femme"))
# Notez que vous pouvez saisir vos données directement si vous avez peu d'observations
Femme <- c(689400, 21315, 181435, 43715, 24295, 8395) # Vecteur de valeurs pour les femmes
Homme <- c(561830, 40010, 238330, 54360, 13765, 6970) # Vecteur de valeurs pour les hommes
Matrice <- as.table(cbind(Femme, Homme)) # Création du tableau
# Nom des deux variables et de leurs modalités respectives
dimnames(Matrice) <- list(Transport=c("Automobile (conducteur)",
                                      "Automobile (passager)",
                                      "Transport en commun",                            
                                      "À pied",
                                      "Bicyclette",
                                      "Autre moyen"),
                          Sexe=c("Homme","Femme"))
# Test du khi-deux
test <- chisq.test(Matrice)
print(test)
## 
##  Pearson's Chi-squared test
## 
## data:  Matrice
## X-squared = 29134, df = 5, p-value < 2.2e-16
# Fréquences observées (Fij)
test$observed
##                          Sexe
## Transport                  Homme  Femme
##   Automobile (conducteur) 689400 561830
##   Automobile (passager)    21315  40010
##   Transport en commun     181435 238330
##   À pied                   43715  54360
##   Bicyclette               24295  13765
##   Autre moyen               8395   6970
# Fréquences théoriques (FTij)
round(test$expected,0)
##                          Sexe
## Transport                  Homme  Femme
##   Automobile (conducteur) 643313 607917
##   Automobile (passager)    31530  29795
##   Transport en commun     215820 203945
##   À pied                   50425  47650
##   Bicyclette               19568  18492
##   Autre moyen               7900   7465
# Déviations (Fij - FTij)
round(test$observed-test$expected,0)
##                          Sexe
## Transport                  Homme  Femme
##   Automobile (conducteur)  46087 -46087
##   Automobile (passager)   -10215  10215
##   Transport en commun     -34385  34385
##   À pied                   -6710   6710
##   Bicyclette                4727  -4727
##   Autre moyen                495   -495
# Contributions au khi-deux
round((test$observed-test$expected)^2/test$expected,2)
##                          Sexe
## Transport                   Homme   Femme
##   Automobile (conducteur) 3301.74 3493.98
##   Automobile (passager)   3309.37 3502.05
##   Transport en commun     5478.22 5797.18
##   À pied                   892.81  944.80
##   Bicyclette              1141.71 1208.19
##   Autre moyen               31.04   32.85
# Marges en ligne et en colonne
colSums(Matrice)
##  Homme  Femme 
## 968555 915265
rowSums(Matrice)
## Automobile (conducteur)   Automobile (passager)     Transport en commun 
##                 1251230                   61325                  419765 
##                  À pied              Bicyclette             Autre moyen 
##                   98075                   38060                   15365
# Grand total
sum(Matrice)
## [1] 1883820
# Pourcentages
round(Matrice/sum(Matrice)*100,2)
##                          Sexe
## Transport                 Homme Femme
##   Automobile (conducteur) 36.60 29.82
##   Automobile (passager)    1.13  2.12
##   Transport en commun      9.63 12.65
##   À pied                   2.32  2.89
##   Bicyclette               1.29  0.73
##   Autre moyen              0.45  0.37
# Pourcentages en ligne
round(Matrice/rowSums(Matrice)*100,2)
##                          Sexe
## Transport                 Homme Femme
##   Automobile (conducteur) 55.10 44.90
##   Automobile (passager)   34.76 65.24
##   Transport en commun     43.22 56.78
##   À pied                  44.57 55.43
##   Bicyclette              63.83 36.17
##   Autre moyen             54.64 45.36
# Pourcentages en colonne
round(Matrice/colSums(Matrice)*100,2)
##                          Sexe
## Transport                 Homme Femme
##   Automobile (conducteur) 71.18 58.01
##   Automobile (passager)    2.33  4.37
##   Transport en commun     18.73 24.61
##   À pied                   4.78  5.94
##   Bicyclette               2.51  1.42
##   Autre moyen              0.92  0.76

Pour obtenir les autres mesures d’association (tableau 5.2), nous pourrons utiliser la syntaxe suivante :

df1 <- read.csv("data/bivariee/hlm.csv")
# Fonction pour calculer les autres mesures d'association
AutresMesuresKhi2 <- function(x, y){
  testChi2 <- chisq.test(x, y) # Calcul du khi-deux
  n  <- sum(testChi2$observed)  # Nombre d'observations
  nc <- ncol(testChi2$observed) # Nombre de colonnes
  l <- nrow(testChi2$observed) # Nombre de lignes
  dl <- (nc-1)*(l-1)            # Nombre de degrés de libertés
  chi2 <- testChi2$statistic   # Valeur du khi-deux
  Pchi2 <- testChi2$p.value    # P pour le khi-deux
  
  #Ratio de vraisemblance du khi-deux
  G <- 2*sum(testChi2$observed*log(testChi2$observed/testChi2$expected)) # G2
  PG <- pchisq(G, df=dl, lower.tail = FALSE) # P pour le G22
  
  # khi-deux de Mantel-Haenszel avec le package DescTools
  MHTest <- DescTools::MHChisqTest(testChi2$observed)
  MH <- MHTest$statistic
  PMH <- MHTest$p.value
  
  # Coefficient de correlation polychorique
  df1 <- data.frame("x" = as.factor(x),
                 "y" = as.factor(y))
  polychoricCorr <- correlation::cor_test(df1,"x","y",method = "polychoric")
  polyR <- polychoricCorr$rho
  polyP <- polychoricCorr$p

  # Coefficient Phi et V de Cramer
  phi <- sqrt(chi2/n)
  vc <- sqrt(chi2/(n*min(nc-1,l-1)))
  
  # Tableau pour les résultats
  dfsortie <- data.frame(
        Statistique = c("Khi-deux", 
                        "Ratio de vraisemblance du  khi-deux", 
                        "Khi-deux de Mantel-Haenszel",
                        "Corrélation Polychorique",
                        "Coefficient de Phi",
                        "V de Cramer"), 
        Valeur = round(c(chi2, G, MH, polyR, phi, vc),3), 
        P = round(c(Pchi2, PG, PMH, polyP , NA, NA),10))
  return(dfsortie)
}

dfkhi2 <- AutresMesuresKhi2(df1$Periode, df1$Taille)

# Impression du tableau avec le package stargazer
library(stargazer)
stargazer(dfkhi2, type="text", summary=FALSE, rownames=FALSE, align = FALSE, digits = 3,
          title="Mesures d'association entre les deux variables qualitatives")
Tableau 5.2: Mesures d’association entre deux variables qualitatives
Statistique Valeur P
Khi-deux 63,543 0
Ratio de vraisemblance du khi-deux 67,286 0
Khi-deux de Mantel-Haenszel 48,486 0
Corrélation Polychorique -0,479 0
Coefficient de Phi 0,477
V de Cramer 0,276