Categorical Data Analysis – UNAL

Introducción

En este trabajo se analiza el conjunto de datos Alligator, el cual describe el tipo de alimento consumido por caimanes en diferentes lagos de Florida. A partir de estas observaciones, se ajustan modelos multinomiales que permiten estudiar cómo la combinación entre el lago, el tamaño del caimán y su sexo influye en la probabilidad de consumir cada tipo de alimento (pescado, aves, invertebrados, entre otros). En resumen mi objetivo es estudiar si estas variables influyen en la probabilidad de consumir cada tipo de alimento. Para ello usaré modelos logit multinomiales, pruebas de hipótesis por razón de verosimilitudes y visualizaciones de probabilidades ajustadas

Se comienza con modelos simples que sólo incluyen una categoría (sexo, tamaño o lago) y luego se comparan con modelos más complejos utilizando pruebas de razón de verosimilitud. Esto permite identificar qué variables aportan información significativa. Posteriormente, se estiman las probabilidades ajustadas de cada tipo de alimento y se visualizan mediante gráficos que muestran cómo cambia la dieta según el lugar y el tamaño del caimán.

Este análisis me ayudará a interpretar patrones ecológicos en la alimentación de los caimanes y muestra el uso del modelo logit multinomial como herramienta estadística para estudiar respuestas categóricas con múltiples niveles.

Data

rm(list=ls())
#install.packages("vcdExtra")
library("vcdExtra")
library(nnet)
library(vcd)
library(dplyr)
library(tidyr)
library(ggplot2)
#Cargo -> data
data("Alligator")
#Visualizó mis primeras filas 
Alligator[1:5, ]
##      lake  sex  size    food count
## 1 Hancock male small    fish     7
## 2 Hancock male small  invert     1
## 3 Hancock male small reptile     0
## 4 Hancock male small    bird     0
## 5 Hancock male small   other     5
# Tablas descriptivas
table(Alligator$food)
## 
##    bird    fish  invert   other reptile 
##      16      16      16      16      16
table(Alligator$size)
## 
## large small 
##    40    40
table(Alligator$lake)
## 
##   George  Hancock Oklawaha Trafford 
##       20       20       20       20

Estas tablas me ayudan a ver la frecuencia de cada categoría de alimento, tamaño y nombre lago

Ajuste - Modelo Multinomial Completo

Voy a ajustar un modelo de la variable de interés, es decir la de respuesta. Tipo de alimento. Por medio de sus variables explicativas: - lake - size - sex

Estas las ponderaré con el conteo de los Alligator que observe.

mod <- multinom(food ~ lake + size + sex , data = Alligator, weights = count)
## # weights:  35 (24 variable)
## initial  value 352.466903 
## iter  10 value 270.397070
## iter  20 value 268.958046
## final  value 268.932740 
## converged

Análisis (Que entiendo de esto?)

##Comparación entre los modelos (Pruebas de Razón de Verosimilitud)

Primero se ajustaré el modelo nulo y luego se agregaré las variables una por una.

moda <- multinom( food ~ 1, data=Alligator, weights = count);
## # weights:  10 (4 variable)
## initial  value 352.466903 
## iter  10 value 302.185213
## final  value 302.181462 
## converged
summary(moda)
## Call:
## multinom(formula = food ~ 1, data = Alligator, weights = count)
## 
## Coefficients:
##         (Intercept)
## fish      1.9783455
## invert    1.5459246
## other     0.9007867
## reptile   0.3794897
## 
## Std. Errors:
##         (Intercept)
## fish      0.2959077
## invert    0.3054775
## other     0.3288968
## reptile   0.3599370
## 
## Residual Deviance: 604.3629 
## AIC: 612.3629
modb <- multinom(food ~ 1 +sex, data=Alligator, weights = count)
## # weights:  15 (8 variable)
## initial  value 352.466903 
## iter  10 value 301.233557
## final  value 301.129428 
## converged
brief <- summary(modb)

#Estadísticas Z
abs(brief$coefficients/brief$standard.errors)
##          (Intercept)    sexmale
## fish    3.991278e+00 0.61774897
## invert  3.424194e+00 0.01666029
## other   1.736392e+00 0.14713114
## reptile 6.634074e-05 0.83244081
# P -values aproximados
1-pnorm(abs(brief$coefficients/brief$standard.errors))
##          (Intercept)   sexmale
## fish    3.285911e-05 0.2683704
## invert  3.083130e-04 0.4933538
## other   4.124730e-02 0.4415143
## reptile 4.999735e-01 0.2025801

Para que es la Estadistica Z?

Cómo se interpreta el p- value, para que se hace esto?

otra vez? se resta standars errors? why?

##Modelos con v.size and v.lake

modc <- multinom(food ~ 1+size, data=Alligator, weights = count)
## # weights:  15 (8 variable)
## initial  value 352.466903 
## iter  10 value 294.669601
## final  value 294.606678 
## converged
modd <- multinom(food ~ 1+lake, data=Alligator, weights = count)
## # weights:  25 (16 variable)
## initial  value 352.466903 
## iter  10 value 283.282977
## iter  20 value 280.584605
## final  value 280.583843 
## converged
mode <- multinom(food ~ 1+size+lake, data=Alligator, weights = count)
## # weights:  30 (20 variable)
## initial  value 352.466903 
## iter  10 value 271.161880
## iter  20 value 270.060188
## final  value 270.040139 
## converged
modf <- multinom(food ~ 1+size+lake+sex, data=Alligator, weights = count)
## # weights:  35 (24 variable)
## initial  value 352.466903 
## iter  10 value 270.397070
## iter  20 value 268.958046
## final  value 268.932740 
## converged

que significa cada una de estas cositas?

##Comparaciones

anova(moda,modb, test="Chisq")
## Likelihood ratio tests of Multinomial Models
## 
## Response: food
##     Model Resid. df Resid. Dev   Test    Df LR stat.   Pr(Chi)
## 1       1       316   604.3629                                
## 2 1 + sex       312   602.2589 1 vs 2     4 2.104069 0.7166248
anova(moda,modc, test="Chisq")  # Size relevante
## Likelihood ratio tests of Multinomial Models
## 
## Response: food
##      Model Resid. df Resid. Dev   Test    Df LR stat.     Pr(Chi)
## 1        1       316   604.3629                                  
## 2 1 + size       312   589.2134 1 vs 2     4 15.14957 0.004400846
anova(moda,modd, test="Chisq")  # Lake relevante
## Likelihood ratio tests of Multinomial Models
## 
## Response: food
##      Model Resid. df Resid. Dev   Test    Df LR stat.      Pr(Chi)
## 1        1       316   604.3629                                   
## 2 1 + lake       304   561.1677 1 vs 2    12 43.19524 2.092387e-05
anova(modc,mode, test="Chisq")  # Lake relevante con Size
## Likelihood ratio tests of Multinomial Models
## 
## Response: food
##             Model Resid. df Resid. Dev   Test    Df LR stat.      Pr(Chi)
## 1        1 + size       312   589.2134                                   
## 2 1 + size + lake       300   540.0803 1 vs 2    12 49.13308 1.982524e-06
anova(modd,mode, test="Chisq")  # Size relevante con Lake
## Likelihood ratio tests of Multinomial Models
## 
## Response: food
##             Model Resid. df Resid. Dev   Test    Df LR stat.      Pr(Chi)
## 1        1 + lake       304   561.1677                                   
## 2 1 + size + lake       300   540.0803 1 vs 2     4 21.08741 0.0003042796
anova(mode,modf, test="Chisq")  # Sexo NO agrega información
## Likelihood ratio tests of Multinomial Models
## 
## Response: food
##                   Model Resid. df Resid. Dev   Test    Df LR stat.   Pr(Chi)
## 1       1 + size + lake       300   540.0803                                
## 2 1 + size + lake + sex       296   537.8655 1 vs 2     4 2.214799 0.6963208

Esto me indica que las variables size and lake explican el tipo de alimento; el sexo no es significativo

Modelo final seleccionado

modfi <- multinom(food ~ lake + size , data = Alligator, weights = count)
## # weights:  30 (20 variable)
## initial  value 352.466903 
## iter  10 value 271.161880
## iter  20 value 270.060188
## final  value 270.040139 
## converged
summary(modfi)
## Call:
## multinom(formula = food ~ lake + size, data = Alligator, weights = count)
## 
## Coefficients:
##         (Intercept) lakeHancock lakeOklawaha lakeTrafford sizesmall
## fish      2.0929906  -0.6951466    0.6529211  -1.08764377 0.6308607
## invert    0.5439613  -2.3535280    1.5901541   0.03435528 2.0890800
## other     0.1887161   0.1310280    0.6585231   0.42873394 0.9624274
## reptile  -1.2215892   0.5476069    3.1118287   1.84765700 0.2796424
## 
## Std. Errors:
##         (Intercept) lakeHancock lakeOklawaha lakeTrafford sizesmall
## fish      0.6622265   0.7812640     1.201951    0.8416809 0.6424896
## invert    0.7311028   0.9337166     1.227642    0.8631658 0.6939177
## other     0.7902850   0.8919688     1.368459    0.9383104 0.7127188
## reptile   1.2076995   1.3715178     1.585926    1.3173360 0.8062363
## 
## Residual Deviance: 540.0803 
## AIC: 580.0803
# P-values aproximados
brief<-summary(modfi)
1-pnorm(abs(brief$coefficients/brief$standard.errors))
##          (Intercept) lakeHancock lakeOklawaha lakeTrafford   sizesmall
## fish    0.0007873958 0.186794240   0.29349001   0.09813906 0.163074492
## invert  0.2284296349 0.005857712   0.09760991   0.48412570 0.001303838
## other   0.4056322928 0.441606440   0.31518183   0.32386383 0.088450171
## reptile 0.1558883781 0.344846872   0.02487233   0.08037243 0.364352336

que significa esto?

Pobabilidades Predichas: Modelo final

newdat <- expand.grid(
  lake = levels(Alligator$lake),
  size = levels(Alligator$size)
)

probs <- predict(modfi, newdata = newdat, type = "probs")

probs_df <- cbind(newdat, as.data.frame(probs)) %>%
  pivot_longer(cols = levels(Alligator$food),
               names_to = "food",
               values_to = "prob")

ggplot(probs_df,
       aes(x = lake, y = prob, color = food, group = food)) +
  geom_line() +
  geom_point() +
  facet_grid(. ~ size) +
  ylab("Predicted probability") +
  theme_bw()

Este gráfico muestra cómo cambia la dieta según el lago y el tamaño

Tabla tipo Agresti y Gráfico de Mosaic

allitable <- xtabs(count ~ lake + sex + size + food, data = Alligator)

# Tabla estructurada
structable(food ~ lake + sex + size, allitable)
##                       food bird fish invert other reptile
## lake     sex    size                                     
## George   female large         0    8      1     1       0
##                 small         0    3      9     1       1
##          male   large         1    9      0     2       0
##                 small         2   13     10     2       0
## Hancock  female large         2    3      0     3       1
##                 small         2   16      3     3       2
##          male   large         1    4      0     2       0
##                 small         0    7      1     5       0
## Oklawaha female large         1    0      1     0       0
##                 small         0    3      9     2       1
##          male   large         0   13      7     0       6
##                 small         0    2      2     1       0
## Trafford female large         0    0      1     0       0
##                 small         1    2      4     4       1
##          male   large         3    8      6     5       6
##                 small         0    3      7     1       1
# Gráfico de mosaic
mosaic(~ food + lake + size, allitable,
       shade = TRUE,
       labeling = labeling_values)

que es mosaic?

El mosaic permite ver asociación entre dieta y variables explicativas visualmente.

Modelo incluyendo el sexo (solo diagnóstico)

modh <- multinom(food ~ lake + size +sex, data = Alligator, weights = count)
## # weights:  35 (24 variable)
## initial  value 352.466903 
## iter  10 value 270.397070
## iter  20 value 268.958046
## final  value 268.932740 
## converged
newdat <- expand.grid(
  lake = levels(Alligator$lake),
  size = levels(Alligator$size),
  sex  = levels(Alligator$sex)
)

probs <- predict(modh, newdata = newdat, type = "probs")
probs_df <- cbind(newdat, as.data.frame(probs)) |> pivot_longer(
  cols = levels(Alligator$food),
  names_to = "food",
  values_to = "prob")

ggplot(probs_df,
       aes(x = lake, y = prob, color = food, group = food)) +
  geom_line() +
  geom_point() +
  facet_grid(size ~ sex) +
  ylab("Predicted probability") +
  theme_bw()

Se observa que el sexo no cambia apreciablemente las probabilidades

Conclusión Final

  • El tipo de alimento consumido depende significativamente de lago y tamaño del caimán.
  • El sexo no tiene influencia estadísticamente significativa.
  • Los caimanes grandes y de ciertos lagos tienden a consumir diferentes tipos de presas, lo cual evidencia posibles factores ecológicos como disponibilidad de alimento o competencia dentro del hábitat.