Modelos de Elección Discreta

Author

Orlando Joaqui-Barandica

Published

September 8, 2025

1 Resumen de la sesión

En esta sesión introducimos modelos de elección discreta para variables dependientes cualitativas. Empezamos con casos binarios y el modelo Logit (formulación, interpretación y validación) y avanzamos a elecciones con m>2 alternativas mediante Logit Multinomial. Cerraremos con criterios para decidir entre binario, multinomial u ordinal y simulaciones reproducibles (código oculto).

2 Objetivos de aprendizaje

Distinguir entre variables binarias, ordinales y nominales y seleccionar el modelo adecuado.
Formular e interpretar un Logit (odds, OR, efectos marginales y ajuste).
Entender la teoría de utilidad aleatoria y la especificación del Logit Multinomial.
Reconocer cuándo usar modelos ordinales (líneas paralelas) vs multinomiales.

3 Marco conceptual

Los modelos de elección discreta estudian la probabilidad de que un individuo elija una opción dentro de un conjunto finito. Ejemplos:

1. Adopción de un crédito (sí/no)
1. Compra online (sí/no)
1. Modo de transporte (bus/auto/bici)
1. Elección de marca (m>2, nominal)
1. Satisfacción en escalas Likert (ordinal).

3.1 Elección binaria y modelo Logit

Sea \((Y_i \in {0,1})\). El Logit modela \((p_i=Pr(Y_i=1\mid X_i))\) con

\[ p_i=\frac{1}{1+e^{-X_i'\beta}} \quad \Rightarrow \quad \log\frac{p_i}{1-p_i}=X_i'\beta. \]

Se estima por máxima verosimilitud; los coeficientes se interpretan como odds ratios \((\exp(\beta_j))\); los efectos en probabilidad mediante efectos marginales.

Important

Cuándo usar Logit (binario): cuando la respuesta tiene 2 categorías excluyentes y nos interesa modelar probabilidades acotadas en \((0,1)\) con elasticidad no lineal. Evita el MLP por predicciones fuera de \((0,1)\), heterocedasticidad y linealidad estricta.

3.2 De binario a multinomial

Cuando \((m>2)\) sin orden, partimos de utilidad aleatoria:

El individuo (i) elige la alternativa (j) con mayor utilidad

\[(U_{ij}=X_i'\beta_j+\varepsilon_{ij})\].

Con errores i.i.d. se deriva el Logit Multinomial (base category logit):

\[Pr(Y_i=j\mid X_i)=\frac{\exp(X_i'\beta_j)}{\sum_{h=0}^{m-1}\exp(X_i'\beta_h)},\]

Relación con múltiples logits binarios y eficiencia de estimar simultáneamente (vs. pares por separado). Considere el supuesto IIA (independencia de alternativas irrelevantes) y sus implicaciones empíricas.

Note

¿Multinomial u ordinal? Si las categorías tienen orden (p. ej., baja, media, alta), prefiera logit/probit ordinal (líneas paralelas). Usar un modelo nominal sobre datos ordinales pierde eficiencia; usar un ordinal sobre datos nominales introduce sesgo.

4 Codigo Modelo Logit - Probit

Un investigador está interesado en cómo las variables, como el GRE (puntaje del examen de Graduados), el GPA (promedio de calificaciones) y el prestigio de la institución de pregrado, tienen efecto en la admisión a la escuela de postgrado.

La variable de respuesta, admitir / no admitir, es una variable binaria.

Code

######################################
#Modelo Logit - Probit
######################################



mydata <- read.csv("Datos/binary.csv")
head(mydata)

Code

summary(mydata)

     admit             gre             gpa             rank      
 Min.   :0.0000   Min.   :220.0   Min.   :2.260   Min.   :1.000  
 1st Qu.:0.0000   1st Qu.:520.0   1st Qu.:3.130   1st Qu.:2.000  
 Median :0.0000   Median :580.0   Median :3.395   Median :2.000  
 Mean   :0.3175   Mean   :587.7   Mean   :3.390   Mean   :2.485  
 3rd Qu.:1.0000   3rd Qu.:660.0   3rd Qu.:3.670   3rd Qu.:3.000  
 Max.   :1.0000   Max.   :800.0   Max.   :4.000   Max.   :4.000

Code

sapply(mydata, sd)

      admit         gre         gpa        rank 
  0.4660867 115.5165364   0.3805668   0.9444602

Code

xtabs(~admit + rank, data = mydata)

     rank
admit  1  2  3  4
    0 28 97 93 55
    1 33 54 28 12

Code

# MODELO LOGIT 

mydata$rank <- factor(mydata$rank)
mylogit <- glm(admit ~ gre + gpa + rank, data = mydata, family = "binomial")

summary(mylogit)


Call:
glm(formula = admit ~ gre + gpa + rank, family = "binomial", 
    data = mydata)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-1.6268  -0.8662  -0.6388   1.1490   2.0790  

Coefficients:
             Estimate Std. Error z value Pr(>|z|)    
(Intercept) -3.989979   1.139951  -3.500 0.000465 ***
gre          0.002264   0.001094   2.070 0.038465 *  
gpa          0.804038   0.331819   2.423 0.015388 *  
rank2       -0.675443   0.316490  -2.134 0.032829 *  
rank3       -1.340204   0.345306  -3.881 0.000104 ***
rank4       -1.551464   0.417832  -3.713 0.000205 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 499.98  on 399  degrees of freedom
Residual deviance: 458.52  on 394  degrees of freedom
AIC: 470.52

Number of Fisher Scoring iterations: 4

Tanto gre como gpa son estadísticamente significativos, al igual que los tres términos para rango. Los coeficientes de regresión logística dan el cambio en las probabilidades del logaritmo del resultado para un aumento de una unidad en la variable predictora.

Por cada cambio de unidades en gre, las probabilidades log de admisión (versus no admisión) aumentan en 0.002.
Para un aumento de una unidad en gpa, las probabilidades de registro de admisión a la escuela de posgrado aumenta en 0,804.
Las variables indicadoras de rango tienen una interpretación ligeramente diferente.
Por ejemplo, haber asistido a una institución de pregrado con un rango de 2, frente a una institución con un rango de 1,cambia las probabilidades de registro de admisión en -0.675.
Debajo de la tabla de coeficientes se encuentran los índices de ajuste, incluidos los residuales nulos y de desviación y el AIC.
Más adelante mostramos un ejemplo de cómo puede usar estos valores para ayudar a evaluar el ajuste del modelo.

Code

confint(mylogit) # Usando Max. Verosimilitud

                    2.5 %       97.5 %
(Intercept) -6.2716202334 -1.792547080
gre          0.0001375921  0.004435874
gpa          0.1602959439  1.464142727
rank2       -1.3008888002 -0.056745722
rank3       -2.0276713127 -0.670372346
rank4       -2.4000265384 -0.753542605

Code

confint.default(mylogit) #Usando errores estándar

                    2.5 %       97.5 %
(Intercept) -6.2242418514 -1.755716295
gre          0.0001202298  0.004408622
gpa          0.1536836760  1.454391423
rank2       -1.2957512650 -0.055134591
rank3       -2.0169920597 -0.663415773
rank4       -2.3703986294 -0.732528724

Code

exp(coef(mylogit)) # Odds ratios

(Intercept)         gre         gpa       rank2       rank3       rank4 
  0.0185001   1.0022670   2.2345448   0.5089310   0.2617923   0.2119375

Code

exp(cbind(OR = coef(mylogit), confint(mylogit))) # Odds ratios - Int. Conf.

                   OR       2.5 %    97.5 %
(Intercept) 0.0185001 0.001889165 0.1665354
gre         1.0022670 1.000137602 1.0044457
gpa         2.2345448 1.173858216 4.3238349
rank2       0.5089310 0.272289674 0.9448343
rank3       0.2617923 0.131641717 0.5115181
rank4       0.2119375 0.090715546 0.4706961

Ahora podemos decir que para un aumento de una unidad en gpa, las probabilidades de ser admitido en la escuela de postgrado (versus no ser admitido) aumentan en un factor de 2.23.

También puede usar las probabilidades pronosticadas para ayudarlo a comprender el modelo. Las probabilidades pronosticadas se pueden calcular para variables predictoras tanto categóricas como continuas. Para crear probabilidades pronosticadas, primero necesitamos crear un nuevo marco de datos con los valores que queremos que las variables independientes adopten para crear nuestras predicciones.

Code

newdata1 <- with(mydata, data.frame(gre = mean(gre), gpa = mean(gpa), rank = factor(1:4)))

newdata1$rankP <- predict(mylogit, newdata = newdata1, type = "response")
newdata1

En el resultado anterior, vemos que la probabilidad prevista de ser aceptado en un programa de postgrado es de 0.52 para estudiantes de las instituciones de pregrado de mayor prestigio (rango = 1) y 0.18 para estudiantes de las instituciones de menor rango (rango = 4), con gre y gpa por sus valores medios. Podemos hacer algo muy similar para crear una tabla de probabilidades predichas que varíe el valor de gre y rank. Vamos a trazarlos, por lo que crearemos 100 valores de gre entre 200 y 800, en cada valor de rango (es decir, 1, 2, 3 y 4).

Code

newdata2 <- with(mydata, data.frame(gre = rep(seq(from = 200, to = 800, length.out = 100),4), 
                                    gpa = mean(gpa), rank = factor(rep(1:4, each = 100))))

El código para generar las probabilidades predichas (la primera línea a continuación) es el mismo que antes, excepto que también vamos a pedir errores estándar para que podamos trazar un intervalo de confianza. Obtenemos las estimaciones en la escala del enlace y transformamos los valores pronosticados y los límites de confianza en probabilidades.

Code

newdata3 <- cbind(newdata2, predict(mylogit, newdata = newdata2, type = "link",
                                    se = TRUE))
newdata3 <- within(newdata3, {
  PredictedProb <- plogis(fit)
  LL <- plogis(fit - (1.96 * se.fit))
  UL <- plogis(fit + (1.96 * se.fit))
})


head(newdata3)

Code

library(ggplot2)

ggplot(newdata3, aes(x = gre, y = PredictedProb)) + 
  geom_ribbon(aes(ymin = LL,ymax = UL, fill = rank), alpha = 0.2) +
  geom_line(aes(colour = rank),size = 1)

4.1 Hosmer-Lemeshow

Code

library(ResourceSelection)
hoslem.test(mylogit$y,fitted(mylogit))


    Hosmer and Lemeshow goodness of fit (GOF) test

data:  mylogit$y, fitted(mylogit)
X-squared = 11.085, df = 8, p-value = 0.1969

No rechazo Ho, indicando que no hay evidencia de mal ajuste. Esto es bueno, ya que aquí sabemos que el modelo está correctamente especificado.

4.2 Modelo Probit

Code

# MODELO PROBIT

myprobit <- glm(admit ~ gre + gpa + rank, data = mydata, family =  binomial(link = "probit"))

summary(myprobit)


Call:
glm(formula = admit ~ gre + gpa + rank, family = binomial(link = "probit"), 
    data = mydata)

Deviance Residuals: 
    Min       1Q   Median       3Q      Max  
-1.6163  -0.8710  -0.6389   1.1560   2.1035  

Coefficients:
             Estimate Std. Error z value Pr(>|z|)    
(Intercept) -2.386836   0.673946  -3.542 0.000398 ***
gre          0.001376   0.000650   2.116 0.034329 *  
gpa          0.477730   0.197197   2.423 0.015410 *  
rank2       -0.415399   0.194977  -2.131 0.033130 *  
rank3       -0.812138   0.208358  -3.898 9.71e-05 ***
rank4       -0.935899   0.245272  -3.816 0.000136 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

(Dispersion parameter for binomial family taken to be 1)

    Null deviance: 499.98  on 399  degrees of freedom
Residual deviance: 458.41  on 394  degrees of freedom
AIC: 470.41

Number of Fisher Scoring iterations: 4

Code

#Hosmer-Lemeshow

library(ResourceSelection)
hoslem.test(myprobit$y,fitted(myprobit))


    Hosmer and Lemeshow goodness of fit (GOF) test

data:  myprobit$y, fitted(myprobit)
X-squared = 12.609, df = 8, p-value = 0.126

Code

# Comparación de Modelos

library(memisc)
mtable(myprobit, mylogit)


Calls:
myprobit: glm(formula = admit ~ gre + gpa + rank, family = binomial(link = "probit"), 
    data = mydata)
mylogit: glm(formula = admit ~ gre + gpa + rank, family = "binomial", 
    data = mydata)

============================================
                   myprobit      mylogit    
--------------------------------------------
  (Intercept)       -2.387***    -3.990***  
                    (0.674)      (1.140)    
  gre                0.001*       0.002*    
                    (0.001)      (0.001)    
  gpa                0.478*       0.804*    
                    (0.197)      (0.332)    
  rank: 2/1         -0.415*      -0.675*    
                    (0.195)      (0.316)    
  rank: 3/1         -0.812***    -1.340***  
                    (0.208)      (0.345)    
  rank: 4/1         -0.936***    -1.551***  
                    (0.245)      (0.418)    
--------------------------------------------
  Log-likelihood  -229.207     -229.259     
  N                400          400         
============================================
  Significance: *** = p < 0.001;   
                ** = p < 0.01;   
                * = p < 0.05

Code

library(marginaleffects)


avg_slopes(mylogit)

Esto me permite cuantificar el cambio promedio en la probabilidad de admisión (admit) asociado con cambios unitarios en las variables explicativas (gre, gpa y rank), proporcionando una interpretación más clara de los efectos directos de cada variable en la probabilidad de admisión.

En este bloque de código, estoy ajustando un modelo de regresión logística (logit) para predecir la probabilidad de admisión (admit) en función del puntaje GRE, el GPA, el rango de la institución de pregrado (rank) y la interacción entre GRE y GPA.

Luego, utilizo la función plot_slopes para visualizar cómo varía la pendiente de la probabilidad de admisión en función de GRE condicionado al GPA y viceversa. Esto me permite entender mejor la interacción entre GRE y GPA y su efecto en la probabilidad de admisión.

Code

mylogit2 <- glm(admit ~ gre + gpa + gre*gpa + rank, data = mydata, family = "binomial")

# Graficar la pendiente del efecto de GRE en la probabilidad de admisión condicionada al GPA
plot_slopes(mylogit2, variables="gre", condition = "gpa")

Code

# Graficar la pendiente del efecto de GPA en la probabilidad de admisión condicionada al GRE
plot_slopes(mylogit2, variables="gpa", condition = "gre")

Code

# Ahora con la variable rank
mylogit2 <- glm(admit ~ gre + gpa + rank*gpa + rank, data = mydata, family = "binomial")

plot_slopes(mylogit2, variables="rank", condition = "gpa")+
  lims(y=c(-1,1))

Code

mylogit2 <- glm(admit ~ gre + gpa + rank*gre + rank, data = mydata, family = "binomial")

plot_slopes(mylogit2, variables="rank", condition = "gre")+
  lims(y=c(-1,1))

4.3 Modelo Multinomial

Ejemplo 1. Las elecciones ocupacionales de las personas pueden verse influenciadas por las ocupaciones de sus padres y su propio nivel educativo. Podemos estudiar la relación entre la elección de ocupación y el nivel educativo y la ocupación del padre. Las elecciones ocupacionales serán la variable de resultado que consta de categorías de ocupaciones.

Ejemplo 2. Un biólogo puede estar interesado en las elecciones de alimentos que hacen los caimanes. Los caimanes adultos pueden tener preferencias diferentes a las de los jóvenes. La variable de resultado aquí serán los tipos de alimentos, y las variables predictoras podrían ser el tamaño de los caimanes y otras variables ambientales.

Ejemplo 3. Los estudiantes que ingresan a la escuela secundaria eligen programas entre programa general, programa vocacional y programa académico. Su elección podría modelarse utilizando su puntaje de escritura y su estatus económico social.

El conjunto de datos contiene variables sobre 200 estudiantes. La variable de resultado es ‘prog’ tipo de programa. Las variables predictoras son el status socioeconómico, ‘ses’ una variable categórica de tres niveles, y la puntuación de escritura, ‘write’ una variable continua.

Code

######################################
#Modelo Multinomial
######################################



library(foreign) 
library(nnet) 
library(ggplot2) 
library(reshape2)


ml  <-  read.dta ( "Datos/hsb.dta" )

head(ml)

A continuación utilizamos la función multinom del nnet paquete para estimar un modelo de regresión logística multinomial. Hay otras funciones en otros paquetes de R capaces de realizar regresión multinomial.

Primero, debemos elegir el nivel de nuestro resultado que deseamos utilizar como línea de base y especificarlo en la función relevel. Luego, ejecutamos nuestro modelo usando multinom. El paquete multinom no incluye el cálculo del valor p para los coeficientes de regresión, por lo que calculamos los valores p utilizando las pruebas de Wald (aquí pruebas z).

Code

ml$prog2 <- relevel(ml$prog, ref = "academic")

test <- multinom(prog2 ~ ses + write, data = ml)

# weights:  15 (8 variable)
initial  value 219.722458 
iter  10 value 179.982880
final  value 179.981726 
converged

Code

summary(test)

Call:
multinom(formula = prog2 ~ ses + write, data = ml)

Coefficients:
         (Intercept)  sesmiddle    seshigh      write
general     2.852198 -0.5332810 -1.1628226 -0.0579287
vocation    5.218260  0.2913859 -0.9826649 -0.1136037

Std. Errors:
         (Intercept) sesmiddle   seshigh      write
general     1.166441 0.4437323 0.5142196 0.02141097
vocation    1.163552 0.4763739 0.5955665 0.02221996

Residual Deviance: 359.9635 
AIC: 375.9635

Code

# Valores z
z <- summary(test)$coefficients/summary(test)$standard.errors
z

         (Intercept)  sesmiddle   seshigh     write
general     2.445214 -1.2018081 -2.261334 -2.705562
vocation    4.484769  0.6116747 -1.649967 -5.112689

Code

# Valores p
p <- (1 - pnorm(abs(z), 0, 1)) * 2
p

          (Intercept) sesmiddle    seshigh        write
general  0.0144766100 0.2294379 0.02373856 6.818902e-03
vocation 0.0000072993 0.5407530 0.09894976 3.176045e-07

Code

# Algunas interpretaciones del modelo "test"

summary(test)

Call:
multinom(formula = prog2 ~ ses + write, data = ml)

Coefficients:
         (Intercept)  sesmiddle    seshigh      write
general     2.852198 -0.5332810 -1.1628226 -0.0579287
vocation    5.218260  0.2913859 -0.9826649 -0.1136037

Std. Errors:
         (Intercept) sesmiddle   seshigh      write
general     1.166441 0.4437323 0.5142196 0.02141097
vocation    1.163552 0.4763739 0.5955665 0.02221996

Residual Deviance: 359.9635 
AIC: 375.9635

Un aumento de una unidad en la variable write se asocia con una disminución en las probabilidades logarítmicas de estar en un programa general versus un programa académico en una cantidad de 0,058.
Un aumento de una unidad en la variable write se asocia con la disminución en las probabilidades logarítmicas de estar en un programa vocacional versus un programa académico. por la cantidad de .1136
Las probabilidades logarítmicas de estar en un programa general frente a un programa académico disminuirán en 1.163 si se pasa de ses=“low”a ses=“high”
Las probabilidades logarítmicas de estar en un programa general frente a un programa académico disminuirán en 0.533 si se pasa de ses=“low”a ses=“middle”, aunque este coeficiente no es significativo.
Las probabilidades logarítmicas de estar en un programa vocacional frente a un programa académico disminuirán en 0.983 si se pasa de ses=“low”a ses=“high”
Las probabilidades logarítmicas de estar en un programa vocacional frente a un programa académico aumentarán en 0.291 si se pasa de ses=“low”a ses=“middle”,aunque este coeficiente no es significativo.

Code

# Odds ratio


exp(coef(test))

         (Intercept) sesmiddle   seshigh     write
general     17.32582 0.5866769 0.3126026 0.9437172
vocation   184.61262 1.3382809 0.3743123 0.8926116

Code

# Probabilidades predichas

pp <- fitted(test)

pp

     academic    general   vocation
1   0.1482764 0.33824545 0.51347812
2   0.1202017 0.18062830 0.69916997
3   0.4186747 0.23680818 0.34451707
4   0.1726885 0.35083838 0.47647308
5   0.1001231 0.16893737 0.73093949
6   0.3533566 0.23779760 0.40884583
7   0.1562562 0.19735040 0.64639339
8   0.1001231 0.16893737 0.73093949
9   0.2331292 0.22039757 0.54647324
10  0.1699402 0.20255307 0.62750669
11  0.2777727 0.37620656 0.34602073
12  0.2917502 0.23360371 0.47464612
13  0.1071687 0.30821953 0.58461179
14  0.2888779 0.22953567 0.48158639
15  0.1482764 0.33824545 0.51347812
16  0.2777727 0.37620656 0.34602073
17  0.3126251 0.37708945 0.31028546
18  0.3293898 0.23309322 0.43751702
19  0.3293898 0.23309322 0.43751702
20  0.6324598 0.20043415 0.16710605
21  0.1998583 0.21215265 0.58798910
22  0.2888779 0.22953567 0.48158639
23  0.3306609 0.37639623 0.29294288
24  0.2777727 0.37620656 0.34602073
25  0.1726885 0.35083838 0.47647308
26  0.3966333 0.23772084 0.36564585
27  0.3676935 0.37276244 0.25954407
28  0.2888779 0.22953567 0.48158639
29  0.2292007 0.36934060 0.40145871
30  0.3865787 0.36985032 0.24357101
31  0.2888779 0.22953567 0.48158639
32  0.1996911 0.36131499 0.43899395
33  0.2141410 0.36565296 0.42020603
34  0.3126251 0.37708945 0.31028546
35  0.2331292 0.22039757 0.54647324
36  0.1998583 0.21215265 0.58798910
37  0.1202017 0.18062830 0.69916997
38  0.7921407 0.14065656 0.06720272
39  0.2292007 0.36934060 0.40145871
40  0.3504532 0.23404073 0.41550609
41  0.1001231 0.16893737 0.73093949
42  0.4852936 0.23070135 0.28400507
43  0.4823541 0.34543275 0.17221315
44  0.4408612 0.23532266 0.32381609
45  0.2777727 0.37620656 0.34602073
46  0.4632141 0.35150976 0.18527617
47  0.3293898 0.23309322 0.43751702
48  0.4632141 0.35150976 0.18527617
49  0.2509706 0.22391068 0.52511877
50  0.3937910 0.23421301 0.37199598
51  0.2331292 0.22039757 0.54647324
52  0.2331292 0.22039757 0.54647324
53  0.2777727 0.37620656 0.34602073
54  0.5202215 0.33179483 0.14798364
55  0.4604652 0.23018005 0.30935477
56  0.3676935 0.37276244 0.25954407
57  0.4630986 0.23327985 0.30362151
58  0.7657822 0.15267894 0.08153887
59  0.6102774 0.29135200 0.09837060
60  0.2917502 0.23360371 0.47464612
61  0.5049270 0.22479322 0.27027983
62  0.2160825 0.21646512 0.56745234
63  0.6876322 0.18315578 0.12921205
64  0.5698425 0.21322354 0.21693394
65  0.3937910 0.23421301 0.37199598
66  0.3293898 0.23309322 0.43751702
67  0.5202215 0.33179483 0.14798364
68  0.3937910 0.23421301 0.37199598
69  0.4604652 0.23018005 0.30935477
70  0.4827571 0.22774111 0.28950178
71  0.7210140 0.17103773 0.10794824
72  0.5698425 0.21322354 0.21693394
73  0.4630986 0.23327985 0.30362151
74  0.4247955 0.36195269 0.21325179
75  0.5698425 0.21322354 0.21693394
76  0.4632141 0.35150976 0.18527617
77  0.4604652 0.23018005 0.30935477
78  0.4604652 0.23018005 0.30935477
79  0.2888779 0.22953567 0.48158639
80  0.3293898 0.23309322 0.43751702
81  0.3676935 0.37276244 0.25954407
82  0.5049270 0.22479322 0.27027983
83  0.6876322 0.18315578 0.12921205
84  0.2777727 0.37620656 0.34602073
85  0.6324598 0.20043415 0.16710605
86  0.3126251 0.37708945 0.31028546
87  0.5049270 0.22479322 0.27027983
88  0.1998583 0.21215265 0.58798910
89  0.6110367 0.20362512 0.18533822
90  0.3088400 0.23158545 0.45957451
91  0.6876322 0.18315578 0.12921205
92  0.7210140 0.17103773 0.10794824
93  0.5291960 0.22405052 0.24675343
94  0.4158872 0.23343321 0.35067961
95  0.6685375 0.18724728 0.14421521
96  0.6110367 0.20362512 0.18533822
97  0.4604652 0.23018005 0.30935477
98  0.4604652 0.23018005 0.30935477
99  0.3293898 0.23309322 0.43751702
100 0.2331292 0.22039757 0.54647324
101 0.4604652 0.23018005 0.30935477
102 0.4604652 0.23018005 0.30935477
103 0.5268859 0.22136713 0.25174694
104 0.2331292 0.22039757 0.54647324
105 0.3937910 0.23421301 0.37199598
106 0.5572058 0.31650507 0.12628915
107 0.4632141 0.35150976 0.18527617
108 0.5049270 0.22479322 0.27027983
109 0.6110367 0.20362512 0.18533822
110 0.5268859 0.22136713 0.25174694
111 0.6102774 0.29135200 0.09837060
112 0.5049270 0.22479322 0.27027983
113 0.5049270 0.22479322 0.27027983
114 0.5049270 0.22479322 0.27027983
115 0.6110367 0.20362512 0.18533822
116 0.5698425 0.21322354 0.21693394
117 0.8158376 0.12901647 0.05514591
118 0.6876322 0.18315578 0.12921205
119 0.6110367 0.20362512 0.18533822
120 0.6685375 0.18724728 0.14421521
121 0.5752548 0.30836653 0.11637863
122 0.6110367 0.20362512 0.18533822
123 0.7200349 0.16949969 0.11046543
124 0.8158376 0.12901647 0.05514591
125 0.6110367 0.20362512 0.18533822
126 0.7921407 0.14065656 0.06720272
127 0.6110367 0.20362512 0.18533822
128 0.6110367 0.20362512 0.18533822
129 0.7200349 0.16949969 0.11046543
130 0.7657822 0.15267894 0.08153887
131 0.6110367 0.20362512 0.18533822
132 0.6324598 0.20043415 0.16710605
133 0.3676935 0.37276244 0.25954407
134 0.7921407 0.14065656 0.06720272
135 0.6685375 0.18724728 0.14421521
136 0.7921407 0.14065656 0.06720272
137 0.5049270 0.22479322 0.27027983
138 0.8370332 0.11788753 0.04507931
139 0.2610484 0.37464158 0.36431000
140 0.7210140 0.17103773 0.10794824
141 0.6876322 0.18315578 0.12921205
142 0.5202215 0.33179483 0.14798364
143 0.5698425 0.21322354 0.21693394
144 0.5698425 0.21322354 0.21693394
145 0.6500011 0.19291318 0.15708576
146 0.6685375 0.18724728 0.14421521
147 0.6685375 0.18724728 0.14421521
148 0.6596792 0.26469707 0.07562378
149 0.5698425 0.21322354 0.21693394
150 0.7921407 0.14065656 0.06720272
151 0.6110367 0.20362512 0.18533822
152 0.8267369 0.12338166 0.04988149
153 0.7200349 0.16949969 0.11046543
154 0.6685375 0.18724728 0.14421521
155 0.6102774 0.29135200 0.09837060
156 0.6685375 0.18724728 0.14421521
157 0.7210140 0.17103773 0.10794824
158 0.8267369 0.12338166 0.04988149
159 0.6596792 0.26469707 0.07562378
160 0.8267369 0.12338166 0.04988149
161 0.8726772 0.09748693 0.02983585
162 0.8559136 0.10735910 0.03672735
163 0.5752548 0.30836653 0.11637863
164 0.7921407 0.14065656 0.06720272
165 0.6110367 0.20362512 0.18533822
166 0.6110367 0.20362512 0.18533822
167 0.8559136 0.10735910 0.03672735
168 0.8267369 0.12338166 0.04988149
169 0.8559136 0.10735910 0.03672735
170 0.7921407 0.14065656 0.06720272
171 0.5572058 0.31650507 0.12628915
172 0.8158376 0.12901647 0.05514591
173 0.6596792 0.26469707 0.07562378
174 0.8726772 0.09748693 0.02983585
175 0.8559136 0.10735910 0.03672735
176 0.6596792 0.26469707 0.07562378
177 0.7508140 0.15740972 0.09177625
178 0.6596792 0.26469707 0.07562378
179 0.6102774 0.29135200 0.09837060
180 0.6876322 0.18315578 0.12921205
181 0.8370332 0.11788753 0.04507931
182 0.6110367 0.20362512 0.18533822
183 0.8559136 0.10735910 0.03672735
184 0.6110367 0.20362512 0.18533822
185 0.8726772 0.09748693 0.02983585
186 0.8726772 0.09748693 0.02983585
187 0.8043127 0.13477968 0.06090766
188 0.8726772 0.09748693 0.02983585
189 0.8267369 0.12338166 0.04988149
190 0.5049270 0.22479322 0.27027983
191 0.8559136 0.10735910 0.03672735
192 0.8267369 0.12338166 0.04988149
193 0.6110367 0.20362512 0.18533822
194 0.8043127 0.13477968 0.06090766
195 0.8370332 0.11788753 0.04507931
196 0.8559136 0.10735910 0.03672735
197 0.6864017 0.18143036 0.13216795
198 0.7508140 0.15740972 0.09177625
199 0.7200349 0.16949969 0.11046543
200 0.6685375 0.18724728 0.14421521

A continuación, si queremos examinar los cambios en la probabilidad predicha asociados con una de nuestras dos variables, podemos crear pequeños conjuntos de datos variando una variable mientras mantenemos la otra constante.

Primero haremos esto manteniendo write su media y examinando las probabilidades predichas para cada nivel de ses.

Code

dses <- data.frame(ses = c("low", "middle", "high"), write = mean(ml$write))
dses

Code

predict(test, newdata = dses, "probs")

   academic   general  vocation
1 0.4396845 0.3581917 0.2021238
2 0.4777488 0.2283353 0.2939159
3 0.7009007 0.1784939 0.1206054

Otra forma de entender el modelo es utilizando las probabilidades predichas para observar las probabilidades promediadas predichas para diferentes valores de la variable predictora continua write dentro de cada nivel de ses.

Code

dwrite <- data.frame(ses = rep(c("low", "middle", "high"),
                               each = 41), write = rep(c(30:70), 3))

## almacene las probabilidades predichas para cada valor de ses
pp.write <- cbind(dwrite, predict(test, newdata = dwrite, type = "probs", se = TRUE))

## calcular las probabilidades medias dentro de cada nivel de ses 
by(pp.write[, 3:5], pp.write$ses, colMeans)

pp.write$ses: high
 academic   general  vocation 
0.6164315 0.1808037 0.2027648 
------------------------------------------------------------ 
pp.write$ses: low
 academic   general  vocation 
0.3972977 0.3278174 0.2748849 
------------------------------------------------------------ 
pp.write$ses: middle
 academic   general  vocation 
0.4256198 0.2010864 0.3732938

Code

# Grafico probabilidades predichas
# Pasar la base a formato long para utilizar ggplot

lpp <- melt(pp.write, id.vars = c("ses", "write"), value.name = "probability")

ggplot(lpp, aes(x = write, y = probability, colour = ses)) +
  geom_line() +
  facet_grid(variable ~ ., scales = "free")

Code

library(marginaleffects)

avg_slopes(test)

En este bloque de código, estoy calculando los efectos promedio marginales en un modelo multinomial para entender cómo las variables independientes ses (nivel socioeconómico) y write (puntuación de escritura) influyen en la probabilidad de estar en diferentes grupos de la variable dependiente (academic, general, vocation).

Estos resultados permiten interpretar el cambio promedio en la probabilidad de estar en cada grupo, condicionado a cambios en las variables independientes, proporcionando una comprensión detallada de los efectos específicos.

Resultados del modelo multinomial

Group: academic

ses: high - low (Estimate: 0.23181, Std. Error: 0.09367, z-value: 2.475, Pr(>|z|): 0.0133)
Un nivel socioeconómico alto comparado con un nivel bajo aumenta la probabilidad de estar en el grupo ‘academic’ en 0.23181 unidades en la escala logit. Este efecto es estadísticamente significativo (p-valor = 0.0133).
ses: middle - low (Estimate: 0.03576, Std. Error: 0.08406, z-value: 0.425, Pr(>|z|): 0.6705)
Un nivel socioeconómico medio comparado con un nivel bajo tiene un efecto positivo muy pequeño y no significativo en la probabilidad de estar en el grupo ‘academic’.
write: dY/dX (Estimate: 0.01714, Std. Error: 0.00288, z-value: 5.959, Pr(>|z|): <0.001)
Un incremento en una unidad en la puntuación de escritura aumenta la probabilidad de estar en el grupo ‘academic’ en 0.01714. Este efecto es altamente significativo (p-valor < 0.001).

Group: general

ses: high - low (Estimate: -0.16006, Std. Error: 0.08686, z-value: -1.843, Pr(>|z|): 0.0654)
Un nivel socioeconómico alto comparado con un nivel bajo disminuye la probabilidad de estar en el grupo ‘general’ en 0.16006 unidades en la escala logit. Este efecto es marginalmente significativo (p-valor = 0.0654).
ses: middle - low (Estimate: -0.12447, Std. Error: 0.07993, z-value: -1.557, Pr(>|z|): 0.1194)
Un nivel socioeconómico medio comparado con un nivel bajo tiene un efecto negativo y no significativo en la probabilidad de estar en el grupo ‘general’.
write: dY/dX (Estimate: -0.00283, Std. Error: 0.00282, z-value: -1.002, Pr(>|z|): 0.3164)
Un incremento en una unidad en la puntuación de escritura disminuye la probabilidad de estar en el grupo ‘general’ en 0.00283. Este efecto no es significativo (p-valor = 0.3164).

Group: vocation

ses: high - low (Estimate: -0.07175, Std. Error: 0.07432, z-value: -0.965, Pr(>|z|): 0.3343)
Un nivel socioeconómico alto comparado con un nivel bajo disminuye la probabilidad de estar en el grupo ‘vocation’ en 0.07175 unidades en la escala logit. Este efecto no es significativo (p-valor = 0.3343).
ses: middle - low (Estimate: 0.08870, Std. Error: 0.06982, z-value: 1.270, Pr(>|z|): 0.2039)
Un nivel socioeconómico medio comparado con un nivel bajo aumenta la probabilidad de estar en el grupo ‘vocation’ en 0.08870 unidades en la escala logit. Este efecto no es significativo (p-valor = 0.2039).
write: dY/dX (Estimate: -0.01431, Std. Error: 0.00251, z-value: -5.698, Pr(>|z|): <0.001)
Un incremento en una unidad en la puntuación de escritura disminuye la probabilidad de estar en el grupo ‘vocation’ en 0.01431. Este efecto es altamente significativo (p-valor < 0.001).

5 Taller

5.1 Objetivo

Modelar una elección con ≥3 categorías nominales y responder preguntas de política/gestión con probabilidades predichas, efectos marginales, etc.

5.2 Fuentes sugeridas

World Bank Data (DataBank), OECD Data, Our World in Data, UN Data, FAOSTAT, IMF Data.

Elijan un dataset con ≥ 300 observaciones (país-año, hogares, empresas, etc) y variables claras. Importante que la frecuencia en las clases sea balanceada.

5.3 Ejemplos de preguntas

Transporte

¿Qué factores explican elegir bus/metro/bici/auto? Y: modo. X: tiempo, costo, clima, ingresos. Contrafactual: −10% en tiempo del metro → ¿nuevos shares?

¿Cómo cambian las elecciones Uber/taxi/bus según hora del día? Y: servicio. X: hora pico, tarifa dinámica, lluvia. Contrafactual: +15% tarifa Uber.

Educación

¿Qué perfila la elección de programa académico/general/vocacional? Y: programa. X: puntaje, SES, género, escuela previa. Contrafactual: +10 puntos en “write”.

¿Qué determina escoger universidad pública/privada/externa? Y: tipo de institución. X: ingreso hogar, becas, distancia. Contrafactual: beca del 50% en privada.

Salud

¿Qué lleva a elegir EPS A/EPS B/privado? Y: proveedor. X: prima, tiempo de espera, edad, crónicos. Contrafactual: reducir espera en EPS B a 20 min.

En urgencias, ¿por qué eligen clínica X/clínica Y/hospital? Y: centro. X: distancia, saturación, seguro. Contrafactual: abrir nueva sede de la clínica X más cerca.

Finanzas

¿Qué define elegir tarjeta débito/crédito/fintech para pagar? Y: método. X: ingresos, comisiones, cashback, edad. Contrafactual: +2% cashback en fintech.

¿Qué influye en tipo de cuenta (básica/premium/pyme)? Y: producto. X: tarifas, uso mensual, score. Contrafactual: eximir cuota de manejo en “premium”.

Consumo/Marketing

¿Qué explica elegir marca A/B/C de café? Y: marca. X: precio, sello orgánico, origen, reseñas. Contrafactual: certificación orgánica para B.

Delivery: Rappi/UberEats/DidiFood Y: plataforma. X: costo envío, tiempo, cupones, rating. Contrafactual: cupón 30% en Rappi.

Energía/Ambiente

¿Qué determina el combustible de cocina gas/electricidad/biomasa? Y: combustible. X: ingreso, subsidio, zona, precio relativo. Contrafactual: bajar precio de electricidad 15%.

5.4 Pasos (mínimos)

Pregunta y alternativas Define la decisión a modelar y las categorías nominales (≥3). Explica por qué no es binario ni ordinal. Fija la categoría base.
Datos y limpieza Fuente, tamaño (≥300 obs. sugerido), definiciones, recodificaciones a factor, outliers/NA, y una tabla de frecuencias por categoría.
Modelo (MNL) Especifica y ~ x1 + x2 + … (incluye 1 interacción simple si tiene sentido). Explica la plausibilidad del IIA (antes de estimar): ¿hay “clones”?
Estimación y lectura Ajusta con nnet::multinom. Reporta coeficientes, odds ratios (OR) y p-valores (Wald). Interpreta por categoría vs. baseline.
Validación/diagnóstico Pseudo-R², log-likelihood, tasa de acierto (matriz de confusión).
Predicciones Construye perfiles y reporta probabilidades predichas P(y=j | x); grafica P vs. un predictor clave con facetas por categoría.