# logistic_ex_corr.r

# Chargement des fonctions et packages

library(Design)

##############
# Exercice 5.1

# Data
load("ARTC.rda")
summary(ARTC)

# 1. Modèle de départ
M1 <- glm(C~relevel(A,2)+relevel(R,2)+relevel(T,2),data=ARTC,family=binomial)
summary(M1)

# Une seule variable significative: T

# Qualité globale du modèle
# Test de la déviance
1-pchisq(457.39,404)
#-> p=0.034: modèle insuffisant pour expliquer les données

# Comparaison avec le modèle constant
st <- 467.45-457.39
dl <- 407-404
1-pchisq(st,dl)
# -> p=0.018: le modèle M1 est meilleur que le modèle constant.

# Calcul de probabilité
logit1 <- 1.4356+0.3581+0-0.5083 
logit1
p1 <- 1/(1+exp(-logit1))
p1
# -> 0.7834

# Autres statistiques
M1b <- lrm(C~relevel(A,2)+relevel(R,2)+relevel(T,2),data=ARTC)
M1b


# 2. Amélioration du modèle
# a) suppression des variables non-significatives
summary(glm(C~relevel(A,2)+relevel(T,2),data=ARTC,family=binomial))
summary(glm(C~relevel(T,2),data=ARTC,family=binomial))

logit <- 1.4214-0.6278
p <- 1/(1+exp(-logit))
p

# b) procédure automatique basée sur BIC
step(M1,direction="backward",k=log(408))
# -> Même résultat: une seule variable explicative restante: T

# c) rajouter des interactions: entre A et T
M2 <- glm(C~relevel(A,2)*relevel(T,2)+relevel(R,2),data=ARTC,family=binomial)
summary(M2)

# -> plus rien de significatif: mauvaise idée


##############
# Exercice 5.2

load("TravailDomicile.rda")
summary(TravailDomicile)

# 1. régression logistique
M1 <- glm(relevel(Travail,2)~Sexe+Etat_civil+Enfants,data=TravailDomicile,family=binomial)
summary(M1)

# Test du modèle lui-même
1-pchisq(1594.8,2202)
# -> le modèle est acceptable pour expliquer la variable dépendante

# Test du modèle contre le modèle nul:
1-pchisq(2487.7-1494.8,2207-2202)
# le modèle est fortement préféré au modèle nul

# calcul des odd-ratios
exp(M1$coefficients)

# ajout d'une interaction
M1b <- glm(relevel(Travail,2)~Sexe+Etat_civil*Enfants,data=TravailDomicile,family=binomial)
summary(M1b)


# 2. régression multinomiale

# librairie mlogit
library(mlogit)

# nouveau jeu de données
Data <- mlogit.data(TravailDomicile,varying=NULL,choice="Etat_civil",shape="wide")

# modèle
M2 <- mlogit(Etat_civil~1|Travail+Sexe+Enfants,data=Data,reflevel="celibataire")
summary(M2)

# 3. trois régressions logistiques séparées
M3_divorce <- mlogit(Etat_civil~1|Travail+Sexe+Enfants,data=Data,reflevel="celibataire",alt.subset=c("celibataire","divorce"))
summary(M3_divorce)
M3_marie <- mlogit(Etat_civil~1|Travail+Sexe+Enfants,data=Data,reflevel="celibataire",alt.subset=c("celibataire","marie"))
summary(M3_marie)
M3_veuf <- mlogit(Etat_civil~1|Travail+Sexe+Enfants,data=Data,reflevel="celibataire",alt.subset=c("celibataire","veuf"))
summary(M3_veuf)

# -> Les résultats du point 3. sont similaires à ceux du point 2.