I. INTRODUCCIÓN
En el presente estudio, se analizará la base de datos Trabajo contiene la información de 11 personas de las cuales se medirán las siguientes variables como PAU (nota en Selectividad), bach (nota media bachiller) y el prestigio del instituto al que asistieron teniendo como variable respuesta la admisión de los estudiantes , las variables se presentan a continuación:
Y : Admisión de los estudiantes ** 1: Admitido ** 0: No admitido X1 : Pau: Nota en Selectividad X2 : bach : Nota media bachiller X3 : Prestigio.
Para alcanzar nuestros objetivos de análisis, aplicaremos tanto la regresión logística clásica como la regresión logística bayesiana. Estos dos enfoques nos permitirán explorar y modelar de manera efectiva la probabilidad de admisión de los estudiantes en función de las variables predictoras mencionadas. La regresión logística clásica se basa en suposiciones y métodos tradicionales, mientras que la regresión logística bayesiana utiliza un enfoque más probabilístico y permite incorporar incertidumbre en nuestras inferencias.
A través de este estudio, buscamos proporcionar una comprensión más profunda de cómo las notas en Selectividad, el desempeño en el bachillerato y el prestigio del instituto pueden influir en la admisión de estudiantes en un contexto específico. Además, la comparación entre los enfoques clásicos y bayesianos en el análisis de datos nos permitirá evaluar la robustez y flexibilidad de estos modelos en este contexto particular.
II. REGRESIÓN LOGÍSTICA CLÁSICA
CONSIDERANDO LOS SIGUIENTES DATOS
# Cargar los datos (reemplaza 'datos.csv' con el nombre de tu archivo de datos)
datos<-read.table("ADMITIDO-EV2.txt",header = T,sep = "\t")
attach(datos)
head(datos,11)## Admitido Pau bach Prestigio
## 1 0 380 3.61 3
## 2 1 660 3.67 3
## 3 1 800 4.00 1
## 4 1 640 3.19 4
## 5 0 520 2.93 4
## 6 1 760 3.00 2
## 7 1 560 2.98 1
## 8 0 400 3.08 2
## 9 1 540 3.39 3
## 10 1 700 3.92 2
## 11 0 380 3.05 1A. Estimar los parámetros con la ecuación de regresión logística binaria clásica.
# Ajustar el modelo de regresión logística binaria
modelo_logistico <- glm(Admitido ~ Pau + bach + Prestigio , data = datos, family = binomial)## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred##
## Call:
## glm(formula = Admitido ~ Pau + bach + Prestigio, family = binomial,
## data = datos)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -3.203e+02 1.287e+06 0.000 1
## Pau 3.678e-01 7.229e+02 0.001 1
## bach 5.520e+01 2.644e+05 0.000 1
## Prestigio -1.420e+01 1.117e+05 0.000 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1.4421e+01 on 10 degrees of freedom
## Residual deviance: 3.9958e-10 on 7 degrees of freedom
## AIC: 8
##
## Number of Fisher Scoring iterations: 25INTERPRETACIÓN - Se observa que, las variables explicativas (Pau , bach y Prestigio) no son significativas (p-value > 5%), indicando que no afecta a la variable respuesta de condición de admitido o no a dicha institucion.
El valor del intercepto es -3.203e+02. Este valor representa la estimación del logaritmo de la razón de condición de admisión cuando todas las demás variables son iguales a cero.
El coeficiente para la variable “Pau” es 3.678e-01. Esto significa que un aumento de una unidad en la nota en Selectividad (“Pau”) se asocia con un aumento en el logaritmo de la razón de estudiantes admitidos de aproximadamente 0.3678 unidades, manteniendo constantes las otras variables.
El coeficiente para la variable “bach” es 5.520e+01. Esto significa que un aumento de una unidad en la nota media del bachillerato (“bach”) se asocia con un aumento en el logaritmo de la razón de o estudiantes admitidos de aproximadamente 55.20 unidades, manteniendo constantes las otras variables.
El coeficiente para la variable “Prestigio” es -1.420e+01. Esto significa que un aumento de una unidad en la medida de prestigio del instituto (“Prestigio”) se asocia con una disminución en el logaritmo de estudiantes admitidos de aproximadamente 14.20 unidades, manteniendo constantes las otras variables.
B. Determinar la bondad de ajuste del modelo.
## ResourceSelection 0.3-6 2023-06-27## Warning in hoslem.test(Admitido, fitted(modelo_logistico), g = 10): The data
## did not allow for the requested number of bins.##
## Hosmer and Lemeshow goodness of fit (GOF) test
##
## data: Admitido, fitted(modelo_logistico)
## X-squared = 1.9979e-10, df = 3, p-value = 1INTERPRETACIÓN
- Los resultados de la prueba de Hosmer-Lemeshow indican que el modelo de regresión logística se ajusta muy bien a los datos, ya que el valor p es 1 y la estadística de prueba (X-squared) es extremadamente pequeña.
C. Hallar el coeficiente OR e interpretar.
## Waiting for profiling to be done...## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred
## Warning: glm.fit: fitted probabilities numerically 0 or 1 occurred## OR 2.5 % 97.5 %
## (Intercept) 8.144029e-140 0.00000e+00 Inf
## Pau 1.444521e+00 1.99955e-17 3.125409e+17
## bach 9.424154e+23 NA Inf
## Prestigio 6.833188e-07 0.00000e+00 NA#Pseudo-R2
library(DescTools)
PseudoR2(modelo_logistico,which = c("McFadden","CoxSnell","Nagelkerke"))## McFadden CoxSnell Nagelkerke
## 1.0000000 0.7304398 1.0000000INTERPRETACIÓN - El intervalo de confianza para el intercepto va desde 0 hasta infinito. Esto sugiere que, dado el valor del coeficiente del intercepto, la probabilidad de ser admitido puede variar ampliamente y no está limitada por un rango
Para la variable “Pau” OR es de aproximadamente 1.444521e+00. Esto indica que, para un aumento unitario en la variable “Pau”, la razón de posibilidades (odds) de ser admitido aumenta en un factor de aproximadamente 1.44.
Para la variable “bach”. El OR es de aproximadamente 9.424154e+23. Esto indica que, para un aumento unitario en la variable “bach”, la razón de posibilidades (odds) de ser admitido aumenta en un factor extremadamente grande.
Para la variable “Prestigio”. El OR es de aproximadamente 6.833188e-07. Esto indica que, para un aumento unitario en la variable “Prestigio”, la razón de posibilidades (odds) de ser admitido disminuye en un factor de aproximadamente 0.0000006833.
D.Realizar predicciones en base a los siguientes valores: X (Pau= 1000, Bach=4.5 y Prestigio =1) e interpretar.
#Tabla de clasificacion
y.pred<-ifelse(modelo_logistico$fitted.values > 0.43, yes = 1, no = 0)
tab_clasif <-table(Admitido, y.pred,dnn = c("observ", "predic"))
tab_clasif## predic
## observ 0 1
## 0 4 0
## 1 0 7## [1] 1#Prediccion
nuevo<-data.frame(Pau=c(1000) , bach=c(4.5) , Prestigio=c(1))
pred<-predict(modelo_logistico,newdata=nuevo,type="response")
pred## 1
## 1INTERPRETACIÓN
La tabla de clasificación muestra la coincidencia entre las observaciones reales y las predicciones del modelo: En la fila 0 (Observación real igual a 0), el modelo hizo 4 predicciones correctas (predijo 0 cuando la observación real era 0) y 0 predicciones incorrectas. En la fila 1 (Observación real igual a 1), el modelo hizo 7 predicciones correctas (predijo 1 cuando la observación real era 1) y 0 predicciones incorrectas. Esto permitió clacular el porcentaje de acierto .
El valor PGA de 1 indica que todas las predicciones del modelo fueron correctas en este conjunto de datos específico.
II.REGRESIÓN LOGÍSTICA BAYESIANA
E.Encontrar los parámetros posteriores para una regresión logística bayesiana asumiendo que siguen una distribución normal: Intercepto -> N (0,9) y para los Bi -> N (0,10)
#REGRESION LOGISTICA BAYESIANA ASUMIENDO NORMALIDAD
# Instala y carga el paquete "rstan" si aún no lo has hecho
library(rstan)## Loading required package: StanHeaders##
## rstan version 2.26.23 (Stan version 2.26.1)## For execution on a local, multicore CPU with excess RAM we recommend calling
## options(mc.cores = parallel::detectCores()).
## To avoid recompilation of unchanged Stan programs, we recommend calling
## rstan_options(auto_write = TRUE)
## For within-chain threading using `reduce_sum()` or `map_rect()` Stan functions,
## change `threads_per_chain` option:
## rstan_options(threads_per_chain = 1)## Do not specify '-march=native' in 'LOCAL_CPPFLAGS' or a Makevars filelibrary(StanHeaders)
# Definir el modelo en formato Stan
modelo_stan <- "
data {
int<lower=0> N; // Número de observaciones
int<lower=0, upper=1> y[N]; // Variable de respuesta binaria
vector[N] x1; // Variable predictora x1
}
parameters {
real alpha; // Intercepto
real beta1; // Coeficiente para x1
}
model {
// Prior para los parámetros
alpha ~ normal(0, 9);
beta1 ~ normal(0, 10);
// Modelo logístico
for (i in 1:N) {
y[i] ~ bernoulli_logit(alpha + beta1 * x1[i]);
}
}
"
# Convertir los datos en formato necesario para Stan
datos_stan <- list(
N = nrow(datos),
y = Admitido,
x1 = Pau , x2 = bach , x3 = Prestigio
)
# Compilar el modelo
modelo_compilado <- stan_model(model_code = modelo_stan)
# Ajustar el modelo utilizando MCMC
ajuste <- sampling(modelo_compilado, data = datos_stan, chains = 4, iter = 2000)##
## SAMPLING FOR MODEL 'anon_model' NOW (CHAIN 1).
## Chain 1:
## Chain 1: Gradient evaluation took 6.5e-05 seconds
## Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0.65 seconds.
## Chain 1: Adjust your expectations accordingly!
## Chain 1:
## Chain 1:
## Chain 1: Iteration: 1 / 2000 [ 0%] (Warmup)
## Chain 1: Iteration: 200 / 2000 [ 10%] (Warmup)
## Chain 1: Iteration: 400 / 2000 [ 20%] (Warmup)
## Chain 1: Iteration: 600 / 2000 [ 30%] (Warmup)
## Chain 1: Iteration: 800 / 2000 [ 40%] (Warmup)
## Chain 1: Iteration: 1000 / 2000 [ 50%] (Warmup)
## Chain 1: Iteration: 1001 / 2000 [ 50%] (Sampling)
## Chain 1: Iteration: 1200 / 2000 [ 60%] (Sampling)
## Chain 1: Iteration: 1400 / 2000 [ 70%] (Sampling)
## Chain 1: Iteration: 1600 / 2000 [ 80%] (Sampling)
## Chain 1: Iteration: 1800 / 2000 [ 90%] (Sampling)
## Chain 1: Iteration: 2000 / 2000 [100%] (Sampling)
## Chain 1:
## Chain 1: Elapsed Time: 0.43 seconds (Warm-up)
## Chain 1: 0.115 seconds (Sampling)
## Chain 1: 0.545 seconds (Total)
## Chain 1:
##
## SAMPLING FOR MODEL 'anon_model' NOW (CHAIN 2).
## Chain 2:
## Chain 2: Gradient evaluation took 1.5e-05 seconds
## Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0.15 seconds.
## Chain 2: Adjust your expectations accordingly!
## Chain 2:
## Chain 2:
## Chain 2: Iteration: 1 / 2000 [ 0%] (Warmup)
## Chain 2: Iteration: 200 / 2000 [ 10%] (Warmup)
## Chain 2: Iteration: 400 / 2000 [ 20%] (Warmup)
## Chain 2: Iteration: 600 / 2000 [ 30%] (Warmup)
## Chain 2: Iteration: 800 / 2000 [ 40%] (Warmup)
## Chain 2: Iteration: 1000 / 2000 [ 50%] (Warmup)
## Chain 2: Iteration: 1001 / 2000 [ 50%] (Sampling)
## Chain 2: Iteration: 1200 / 2000 [ 60%] (Sampling)
## Chain 2: Iteration: 1400 / 2000 [ 70%] (Sampling)
## Chain 2: Iteration: 1600 / 2000 [ 80%] (Sampling)
## Chain 2: Iteration: 1800 / 2000 [ 90%] (Sampling)
## Chain 2: Iteration: 2000 / 2000 [100%] (Sampling)
## Chain 2:
## Chain 2: Elapsed Time: 0.39 seconds (Warm-up)
## Chain 2: 0.108 seconds (Sampling)
## Chain 2: 0.498 seconds (Total)
## Chain 2:
##
## SAMPLING FOR MODEL 'anon_model' NOW (CHAIN 3).
## Chain 3:
## Chain 3: Gradient evaluation took 9e-06 seconds
## Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 0.09 seconds.
## Chain 3: Adjust your expectations accordingly!
## Chain 3:
## Chain 3:
## Chain 3: Iteration: 1 / 2000 [ 0%] (Warmup)
## Chain 3: Iteration: 200 / 2000 [ 10%] (Warmup)
## Chain 3: Iteration: 400 / 2000 [ 20%] (Warmup)
## Chain 3: Iteration: 600 / 2000 [ 30%] (Warmup)
## Chain 3: Iteration: 800 / 2000 [ 40%] (Warmup)
## Chain 3: Iteration: 1000 / 2000 [ 50%] (Warmup)
## Chain 3: Iteration: 1001 / 2000 [ 50%] (Sampling)
## Chain 3: Iteration: 1200 / 2000 [ 60%] (Sampling)
## Chain 3: Iteration: 1400 / 2000 [ 70%] (Sampling)
## Chain 3: Iteration: 1600 / 2000 [ 80%] (Sampling)
## Chain 3: Iteration: 1800 / 2000 [ 90%] (Sampling)
## Chain 3: Iteration: 2000 / 2000 [100%] (Sampling)
## Chain 3:
## Chain 3: Elapsed Time: 0.329 seconds (Warm-up)
## Chain 3: 0.086 seconds (Sampling)
## Chain 3: 0.415 seconds (Total)
## Chain 3:
##
## SAMPLING FOR MODEL 'anon_model' NOW (CHAIN 4).
## Chain 4:
## Chain 4: Gradient evaluation took 8e-06 seconds
## Chain 4: 1000 transitions using 10 leapfrog steps per transition would take 0.08 seconds.
## Chain 4: Adjust your expectations accordingly!
## Chain 4:
## Chain 4:
## Chain 4: Iteration: 1 / 2000 [ 0%] (Warmup)
## Chain 4: Iteration: 200 / 2000 [ 10%] (Warmup)
## Chain 4: Iteration: 400 / 2000 [ 20%] (Warmup)
## Chain 4: Iteration: 600 / 2000 [ 30%] (Warmup)
## Chain 4: Iteration: 800 / 2000 [ 40%] (Warmup)
## Chain 4: Iteration: 1000 / 2000 [ 50%] (Warmup)
## Chain 4: Iteration: 1001 / 2000 [ 50%] (Sampling)
## Chain 4: Iteration: 1200 / 2000 [ 60%] (Sampling)
## Chain 4: Iteration: 1400 / 2000 [ 70%] (Sampling)
## Chain 4: Iteration: 1600 / 2000 [ 80%] (Sampling)
## Chain 4: Iteration: 1800 / 2000 [ 90%] (Sampling)
## Chain 4: Iteration: 2000 / 2000 [100%] (Sampling)
## Chain 4:
## Chain 4: Elapsed Time: 0.345 seconds (Warm-up)
## Chain 4: 0.111 seconds (Sampling)
## Chain 4: 0.456 seconds (Total)
## Chain 4:## $summary
## mean se_mean sd 2.5% 25%
## alpha -12.9587227 0.2282981776 5.31690169 -24.639951332 -16.30488132
## beta1 0.0258001 0.0004381881 0.01022597 0.009041819 0.01814719
## lp__ -3.6309627 0.0311921915 1.01239381 -6.259997745 -4.02609174
## 50% 75% 97.5% n_eff Rhat
## alpha -12.53646473 -9.05291771 -4.20182889 542.3909 1.011068
## beta1 0.02496367 0.03252334 0.04798775 544.6126 1.011580
## lp__ -3.32392429 -2.89839874 -2.60952088 1053.4337 1.003770
##
## $c_summary
## , , chains = chain:1
##
## stats
## parameter mean sd 2.5% 25% 50%
## alpha -12.31260974 4.935022101 -22.500297874 -15.63173691 -12.03762415
## beta1 0.02456887 0.009474092 0.008446578 0.01737288 0.02371102
## lp__ -3.54950503 0.963872446 -6.194878136 -3.91095758 -3.27642355
## stats
## parameter 75% 97.5%
## alpha -8.74675391 -3.64268057
## beta1 0.03107663 0.04422018
## lp__ -2.87503901 -2.60962705
##
## , , chains = chain:2
##
## stats
## parameter mean sd 2.5% 25% 50%
## alpha -13.29053280 5.60170178 -25.418583931 -16.94763901 -12.8029775
## beta1 0.02649818 0.01081698 0.008923131 0.01824144 0.0255255
## lp__ -3.67636386 1.02521539 -6.335130361 -4.08910906 -3.3417660
## stats
## parameter 75% 97.5%
## alpha -9.25712097 -4.35924540
## beta1 0.03393242 0.05006758
## lp__ -2.91440306 -2.62956526
##
## , , chains = chain:3
##
## stats
## parameter mean sd 2.5% 25% 50%
## alpha -12.75503341 5.24532172 -23.723293758 -16.02285538 -12.44441608
## beta1 0.02534364 0.01001717 0.008243277 0.01793467 0.02476782
## lp__ -3.63121542 0.96963430 -5.914004946 -4.08785638 -3.34456822
## stats
## parameter 75% 97.5%
## alpha -8.99931915 -3.64868485
## beta1 0.03186906 0.04737876
## lp__ -2.88922837 -2.61147141
##
## , , chains = chain:4
##
## stats
## parameter mean sd 2.5% 25% 50%
## alpha -13.47671498 5.39236110 -25.291391134 -16.86731941 -12.86403164
## beta1 0.02678973 0.01040636 0.009981899 0.01907589 0.02566218
## lp__ -3.66676638 1.08284205 -6.670972160 -4.02983042 -3.31924826
## stats
## parameter 75% 97.5%
## alpha -9.46465544 -4.73981820
## beta1 0.03363699 0.04934504
## lp__ -2.93056661 -2.59628825## ci_level: 0.8 (80% intervals)## outer_level: 0.95 (95% intervals)
INTERPRETACIÓN
Los resultados muestran estadísticas para los parámetros del modelo, en particular para alpha y beta1. La media (mean) y el intervalo de credibilidad del 95% (2.5% y 97.5%) proporcionan una estimación puntual y un rango de incertidumbre para estos parámetros.
El diagnóstico de convergencia (Rhat) para ambos parámetros es cercano a 1, lo que indica que las cadenas de Markov utilizadas en el muestreo parecen haber convergido adecuadamente.
-La log probabilidad posterior (lp__) es una medida de la bondad de ajuste del modelo a los datos. Cuanto mayor sea este valor, mejor se ajusta el modelo a los datos observados.
III. CONCLUSIÓN
En general, estos modelos proporcionan herramientas valiosas para analizar y predecir la probabilidad de admisión en función de las características de los solicitantes. Los resultados permiten entender cómo estas variables influyen en la toma de decisiones de admisión y ofrecen una visión más completa de la incertidumbre asociada con las estimaciones de los parámetros. La elección entre el modelo clásico y el bayesiano depende de las preferencias y requisitos específicos del análisis y de la cantidad de información disponible.