
# Matriz de datos
load(url("http://pagines.uab.cat/plopez/sites/pagines.uab.cat.plopez/files/Voto.rda"))
attach(Voto)
Voto=data.frame(Voto,Diarios,Edad,Hijos,Titestud,Indsocec,Ingresos,Actividad,Politica,Sexo,Raza)
Voto=data.frame(na.omit(Voto)) # Se suprimen los casos con valores perdidos de todas las variables


# Paquetes adicionales para la regresión logística
install.packages("fmsb") # Función NagelkerkeR2
library(fmsb)
install.packages("rcompanion")
library(rcompanion)
install.packages("blorr") # Modelos de ARL binaria: https://blorr.rsquaredacademy.com/
library(blorr)

############## Análisis de regresión logística binaria simple

##### Análisis 1: Voto y Diarios 

# Tabla de contingencia entre Voto y Diarios 
library (Deducer)
TC = contingency.tables(row.vars=d(Voto),col.vars=d(Diarios),data=Voto)
TC = tables<-add.chi.squared(TC)
TC = print(tables,prop.r=F,prop.c=T,prop.t=F)

# Para ver la codificación dummy y ver cómo interpretar el modelo
contrasts(Voto$Voto)
contrasts(Voto$Diarios)

# ARL binaria simple Voto según Diarios (cualitativa)
ARLog = glm(Voto ~ Diarios, family=binomial(link='logit'), data=Voto, na.action=na.omit)
summary(ARLog)
# Intervalos de confianza de los coeficientes
confint(ARLog)
# Exponencial de los coeficientes: modelo multiplicativo
exp(coef(ARLog))
exp(confint(ARLog))
exp(cbind(OR = coef(ARLog), confint(ARLog)))
# Resultados generales del modelo con blorr
blr_regress(ARLog, odd_conf_limit=TRUE)
# Para analizar la tabla de desvianza (significación de las variables del modelo)
anova(ARLog,test="Chisq")
# Bondad de ajuste del modelo
NagelkerkeR2(ARLog)
blr_model_fit_stats(ARLog)
nagelkerke(ARLog)

##### Análisis 2: Voto y Educ 

frequencies(Voto[c("Educ")], r.digits=1)

# ARL binaria simple Voto según Educ (cuantitativa)
ARLog = glm(Voto ~ Educ, family=binomial(link='logit'), data=Voto, na.action=na.omit)
summary(ARLog)
confint(ARLog)
exp(cbind(OR = coef(ARLog), confint(ARLog)))
blr_regress(ARLog, odd_conf_limit=TRUE)
anova(ARLog,test="Chisq")
NagelkerkeR2(ARLog)
blr_model_fit_stats(ARLog)

# Log de las razones predichas (linear.predictors) y razones predichas
# Probabilidad estimada (fitted.values): 1/(1+exp(-(1.772-0.212*Educ)))
# Residuos
Predictors=exp(ARLog$linear.predictors)
ARLogres=data.frame(ARLog$linear.predictors,Predictors,ARLog$fitted.values,ARLog$residuals)
frequencies(ARLogres[c("ARLog.linear.predictors","Predictors","ARLog.fitted.values",
                       "ARLog.residuals")], r.digits=1)
save(ARLogres,file="C:/Datos/Votores.rda")

############## Análisis de regresión logística binaria múltiple

##### Análisis 3: Voto según Diarios y Educ 

ARLog = glm(Voto ~ Diarios + Educ, family=binomial(link='logit'), data=Voto, na.action=na.omit)
summary(ARLog)
confint(ARLog)
exp(cbind(OR = coef(ARLog), confint(ARLog)))
blr_regress(ARLog, odd_conf_limit=TRUE)
anova(ARLog,test="Chisq")
NagelkerkeR2(ARLog)
blr_model_fit_stats(ARLog)

# Selección del modelo por etapas
blr_step_aic_both(ARLog)


##### Análisis 4: Voto según múltiples variables

frequencies(Voto[c("Diarios","Edad","Hijos","Titestud","Indsocec","Ingresos","Actividad",
                   "Politica","Sexo","Raza")], r.digits=1)

# Introducción de todas las variables
ARLog = glm(Voto ~ Diarios+Edad+Hijos+Titestud+Indsocec+Ingresos+Actividad+Politica+Sexo+Raza, 
            family=binomial(link='logit'), data=Voto, na.action=na.omit)
summary(ARLog)
NagelkerkeR2(ARLog)

# Selección del modelo
blr_step_aic_both(ARLog)
plot(blr_step_aic_both(ARLog))
ARLog = glm(Voto ~ Diarios+Indsocec+Politica+Edad+Titestud+Actividad, 
            family=binomial(link='logit'), data=Voto, na.action=na.omit)
summary(ARLog)
blr_regress(ARLog, odd_conf_limit=TRUE)
anova(ARLog,test="Chisq")
ARLog0=glm(Voto ~ 1, family=binomial(link='logit'), data=Voto, na.action=na.omit) # Modelo constante
anova(ARLog0,ARLog,test="Chisq") # Significación global del modelo
NagelkerkeR2(ARLog)



