###### Opción 1: análisis de una tabla de frecuencias partiendo de una matriz de datos

# Data frame Penalty.rda
load(url("http://pagines.uab.cat/plopez/sites/pagines.uab.cat.plopez/files/Penalty.rda"))

# Conversión de la matriz de casos (data frame) a tabla
TC=table(Penalty$PENALTY,Penalty$MURDER,Penalty$VICTIM,dnn=c("PENALTY","MURDER","VICTIM"))
TC

# Análisis log-lineal

# Modelo saturado
MS=glm(Freq ~ PENALTY*MURDER*VICTIM, data=TC, family=poisson(link="log"))
summary(MS)
names(MS)
# Función step: selección de modelo
step(MS, direction="backward", test="Chisq")

# Modelo saturado: cambio de categoría de referencia
MS=glm(Freq ~ relevel(PENALTY,2)*MURDER*VICTIM, poisson, data=TC)
summary(MS)
step(MS, direction="backward", test="Chisq")


###### Opción 2: creación y análisis de una tabla de frecuencias

# Introducción de la tabla como Data Frame
PENALTY=c("Death row","Prison")
MURDER=c("White","White","Black","Black")
VICTIM=c("White","White","White","White","Black","Black","Black","Black")
FRE=c(72,2074,48,239,0,111,11,2209)
TCF=data.frame(PENALTY,MURDER,VICTIM,FRE)
TCF

T1=xtabs(FRE ~ PENALTY + MURDER, data=TCF)
round(prop.table((T1),2)*100,1)
summary(T1)

T2=xtabs(FRE ~ PENALTY + VICTIM, data=TCF)
round(prop.table((T2),2)*100,1)
summary(T2)

T3=xtabs(FRE ~ VICTIM + MURDER, data=TCF)
round(prop.table((T3),2)*100,1)
summary(T3)

T4=xtabs(FRE ~ PENALTY + MURDER + VICTIM, data=TCF)
ftable(T4)
round(prop.table(T4,c(3,2))*100,1)

# Análisis log-lineal
# Modelo saturado
MS=glm(Fre ~ PENALTY*MURDER*VICTIM, data=TCF, family=poisson(link="log"))
summary(MS)
# Modelo saturado: cambio de categoría de referencia
levels(TCF$PENALTY)
MS=glm(Fre ~ relevel(PENALTY,2)*MURDER*VICTIM, poisson, data=TCF)
summary(MS)
names(MS)
# Función step: selección de modelo
step(MS, direction="backward", test="Chisq")
# Modelo homogéneo
MH=glm(Fre ~ relevel(PENALTY, 2) + MURDER + VICTIM + 
       relevel(PENALTY,2):MURDER + relevel(PENALTY, 2):VICTIM + MURDER:VICTIM,
       poisson, data=TCF)
summary(MH)
