# Analyse de la pratique du sport  ----

## Chargement des données et des packages ----

library(tidyverse)
## -- Attaching packages --------------------------------------- tidyverse 1.3.0 --
## v ggplot2 3.3.3     v purrr   0.3.4
## v tibble  3.0.4     v dplyr   1.0.2
## v tidyr   1.1.2     v stringr 1.4.0
## v readr   1.4.0     v forcats 0.5.0
## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()
library(labelled)
library(questionr)
library(gtsummary)
## #Uighur
library(GGally)
## Registered S3 method overwritten by 'GGally':
##   method from   
##   +.gg   ggplot2
library(effects)
## Loading required package: carData
## lattice theme set by effectsTheme()
## See ?effectsTheme for details.
library(ggeffects)

theme_gtsummary_language("fr", decimal.mark = ",", big.mark = " ")
## Setting `language: fr` theme
theme_gtsummary_mean_sd()

data(hdv2003)

## Étiquettes de variables ------

hdv2003 <- hdv2003 %>%
  set_variable_labels(
    sport = "Pratique du sport ?",
    age = "Âge",
    sexe = "Sexe",
    nivetud = "Niveau d'études",
    relig = "Rapport à la religion",
    heures.tv = "Heures quotidiennes de TV"
  )

## Recodages -------

hdv2003$groupe_age <- hdv2003$age %>% 
  cut(
    include.lowest = TRUE,
    right = FALSE,
    dig.lab = 4,
    breaks = c(18, 25, 45, 60, 97)
  ) %>% 
  fct_recode(
    "18-24" = "[18,25)",
    "25-44" = "[25,45)",
    "45-59" = "[45,60)",
    "60 et +" = "[60,97]"
  )
var_label(hdv2003$groupe_age) <- "Groupe d'âges"

hdv2003$sexe <- hdv2003$sexe %>%
  fct_relevel("Femme")

hdv2003$etudes <- hdv2003$nivetud %>% 
  fct_recode(
    "Primaire" = "N'a jamais fait d'etudes",
    "Primaire" = "A arrete ses etudes, avant la derniere annee d'etudes primaires",
    "Primaire" = "Derniere annee d'etudes primaires",
    "Secondaire" = "1er cycle",
    "Secondaire" = "2eme cycle",
    "Technique/Professionnel" = "Enseignement technique ou professionnel court",
    "Technique/Professionnel" = "Enseignement technique ou professionnel long",
    "Supérieur" = "Enseignement superieur y compris technique superieur"
  ) %>%
  fct_explicit_na("Manquant")
var_label(hdv2003$etudes) <- "Niveau d'études"

## Analyse univariée --------

hdv2003 %>%
  select(sport, groupe_age, sexe, etudes, relig, heures.tv) %>%
  tbl_summary(
    statistic = all_categorical() ~ "{p}% [{n}]",
    digits = all_categorical() ~ c(1, 0)
  )
Caractéristique N = 2 0001
Pratique du sport ?
Non 63,8% [1 277]
Oui 36,1% [723]
Groupe d'âges
18-24 8,5% [169]
25-44 35,3% [706]
45-59 30,4% [609]
60 et + 25,8% [516]
Sexe
Femme 55,0% [1 101]
Homme 45,0% [899]
Niveau d'études
Primaire 23,3% [466]
Secondaire 19,4% [387]
Technique/Professionnel 29,7% [594]
Supérieur 22,0% [441]
Manquant 5,6% [112]
Rapport à la religion
Pratiquant regulier 13,3% [266]
Pratiquant occasionnel 22,1% [442]
Appartenance sans pratique 38,0% [760]
Ni croyance ni appartenance 20,0% [399]
Rejet 4,7% [93]
NSP ou NVPR 2,0% [40]
Heures quotidiennes de TV 2,25 (1,78)
Manquant 5

1 % [n]; Moyenne (ET) 

## Analyse bivariée --------

hdv2003 %>%
  select(sport, groupe_age, sexe, etudes, relig, heures.tv) %>%
  tbl_summary(
    by = "sport",
    percent = "row",
    statistic = all_categorical() ~ "{p}% [{n}]",
    digits = all_categorical() ~ c(1, 0)
  ) %>%
  add_overall(last = TRUE) %>%
  add_p()
Caractéristique Non, N = 1 2771 Oui, N = 7231 Total, N = 2 0001 p-value2
Groupe d'âges <0,001
18-24 34,3% [58] 65,7% [111] 100,0% [169]
25-44 50,8% [359] 49,2% [347] 100,0% [706]
45-59 72,2% [440] 27,8% [169] 100,0% [609]
60 et + 81,4% [420] 18,6% [96] 100,0% [516]
Sexe <0,001
Femme 67,8% [747] 32,2% [354] 100,0% [1 101]
Homme 59,0% [530] 41,0% [369] 100,0% [899]
Niveau d'études <0,001
Primaire 89,3% [416] 10,7% [50] 100,0% [466]
Secondaire 69,8% [270] 30,2% [117] 100,0% [387]
Technique/Professionnel 63,6% [378] 36,4% [216] 100,0% [594]
Supérieur 42,2% [186] 57,8% [255] 100,0% [441]
Manquant 24,1% [27] 75,9% [85] 100,0% [112]
Rapport à la religion 0,14
Pratiquant regulier 68,4% [182] 31,6% [84] 100,0% [266]
Pratiquant occasionnel 66,7% [295] 33,3% [147] 100,0% [442]
Appartenance sans pratique 62,2% [473] 37,8% [287] 100,0% [760]
Ni croyance ni appartenance 59,9% [239] 40,1% [160] 100,0% [399]
Rejet 64,5% [60] 35,5% [33] 100,0% [93]
NSP ou NVPR 70,0% [28] 30,0% [12] 100,0% [40]
Heures quotidiennes de TV 2,48 (1,90) 1,83 (1,43) 2,25 (1,78) <0,001
Manquant 2 3 5

1 % [n]; Moyenne (ET)

2 Pearson's Chi-squared test; Welch Two Sample t-test 

hdv2003 %>%
  ggbivariate(
    outcome = "sport",
    explanatory = c("groupe_age", "sexe", "etudes", "relig", "heures.tv")
  )
## Warning: Removed 5 rows containing non-finite values (stat_boxplot).

## Calcul de la régression logistique -----

mod <- glm(
  sport ~ groupe_age + sexe + etudes + relig + heures.tv,
  family = binomial(),
  data = hdv2003
)

mod %>% 
  tbl_regression(
    intercept = TRUE,
    exponentiate = TRUE,
    add_estimate_to_reference_row = TRUE
  )
Caractéristique OR1 95% CI1 p-value
(Intercept) 0,44 0,23 – 0,82 0,011
Groupe d'âges
18-24 1,00 —
25-44 0,66 0,42 – 1,03 0,069
45-59 0,35 0,21 – 0,55 <0,001
60 et + 0,27 0,16 – 0,45 <0,001
Sexe
Femme 1,00 —
Homme 1,56 1,27 – 1,92 <0,001
Niveau d'études
Primaire 1,00 —
Secondaire 2,62 1,79 – 3,87 <0,001
Technique/Professionnel 2,90 2,01 – 4,23 <0,001
Supérieur 6,76 4,65 – 9,96 <0,001
Manquant 8,84 4,67 – 17,0 <0,001
Rapport à la religion
Pratiquant regulier 1,00 —
Pratiquant occasionnel 0,98 0,68 – 1,42 >0,9
Appartenance sans pratique 0,99 0,71 – 1,40 >0,9
Ni croyance ni appartenance 0,80 0,55 – 1,18 0,3
Rejet 0,68 0,39 – 1,19 0,2
NSP ou NVPR 0,92 0,40 – 2,03 0,8
Heures quotidiennes de TV 0,89 0,83 – 0,95 <0,001

1 OR = rapport de cotes, CI = intervalle de confiance 

ggcoef_model(mod, exponentiate = TRUE)

## Effets marginaux ------

mod %>% allEffects() %>% plot()

# plot(allEffects(mod))

ggeffect(mod) %>% 
  plot() %>%
  cowplot::plot_grid(plotlist = .)