1 Introducción y planteamiento del problema

Una organización busca entender qué factores explican que un empleado decida rotar (dejar su cargo actual) y anticipar, para cada trabajador, la probabilidad de que esto ocurra en el próximo período. Esta información le permitiría a la gerencia intervenir de manera preventiva sobre los empleados con mayor riesgo, en lugar de reaccionar una vez la renuncia ya ocurrió.

El objetivo de este informe es construir un modelo de regresión logística binomial que estime la probabilidad de rotación de un empleado a partir de un conjunto de variables organizacionales, identificar cuáles de esas variables presentan una asociación estadísticamente significativa con la probabilidad de rotación y estimar la magnitud y dirección de dicha asociación, y evaluar qué tan bien discrimina el modelo entre empleados que rotan y empleados que permanecen. Al tratarse de un diseño observacional y transversal, el alcance de estas conclusiones es asociativo, no causal — vuelvo sobre esta distinción de forma explícita al interpretar los coeficientes.

La variable respuesta es Rotación, codificada como “Sí” (el empleado rotó) y “No” (el empleado permaneció). Al ser una variable binaria, el modelo de regresión lineal múltiple no es apropiado: no garantiza que las predicciones queden dentro del rango \([0,1]\), viola los supuestos de normalidad y homoscedasticidad de los residuales, y su ajuste (medido por \(R^2\)) suele ser artificialmente bajo cuando la variable dependiente solo toma dos valores. Por esa razón se recurre al modelo logit, que transforma el problema de estimar una probabilidad en el de estimar el logaritmo de una razón de probabilidades (log-odds), garantizando que la probabilidad estimada esté siempre acotada entre 0 y 1.

2 Selección de variables e hipótesis

Aunque el tamaño muestral (1.470 empleados) permitiría en principio considerar un conjunto amplio de predictores, opté por un modelo parsimonioso de seis variables por dos razones concretas, no solo por preferencia de simplicidad. Primero, una regla empírica estándar en regresión logística (Peduzzi et al., 1996) recomienda al menos 10 eventos de la clase minoritaria por cada variable candidata para que las estimaciones de máxima verosimilitud sean estables; con 237 casos de rotación (la clase minoritaria), el límite razonable de predictores está alrededor de 23-24 — es decir, usar las 24 variables disponibles simultáneamente dejaría el modelo prácticamente en el borde de esa regla, sin margen para interacciones ni términos adicionales. Segundo, varias variables del conjunto de datos miden esencialmente el mismo constructo (por ejemplo, Antigüedad, Antigüedad_Cargo, Años_acargo_con_mismo_jefe y Años_ultima_promoción son todas proxies del tiempo que un empleado lleva vinculado a la empresa, y están correlacionadas entre sí por encima de 0.55), por lo que incluirlas todas añadiría redundancia sin aportar información nueva proporcional al costo en grados de libertad e interpretabilidad.

Antes de seleccionar, es necesario clasificar correctamente el tipo de cada variable, porque varias aparecen codificadas como enteros sin ser cuantitativas en sentido estricto: Educación, Satisfacción_Ambiental, Satisfación_Laboral, Rendimiento_Laboral y Equilibrio_Trabajo_Vida son escalas ordinales tipo Likert (1 a 4 o 1 a 5). En este informe las trato como variables cuasi-cuantitativas de intervalo aproximado, una convención común en la práctica aplicada, y lo declaro explícitamente para que la decisión sea trazable.

Seleccioné 3 variables categóricas y 3 cuantitativas, con un criterio doble: (a) sustento teórico en la literatura de comportamiento organizacional y rotación de personal, y (b) baja colinealidad entre sí, verificada empíricamente antes de comprometerme con el conjunto final.

2.1 Variables categóricas

Horas_Extra (Sí / No). La teoría del agotamiento de recursos en el trabajo (job demands-resources model; Maslach, describiendo el síndrome de burnout) plantea que exigencias laborales sostenidas por encima de la capacidad de recuperación del trabajador erosionan su compromiso y aumentan la intención de abandonar el cargo. Hipótesis: trabajar horas extra aumenta la probabilidad de rotación.

Viaje de Negocios (No_Viaja / Raramente / Frecuentemente). El modelo de conflicto trabajo-familia de Greenhaus y Beutell (1985) explica que cuando las demandas del rol laboral interfieren con las del rol familiar o personal, aumenta la tensión y, con ella, la probabilidad de abandonar el empleo. Viajar con frecuencia es una de las formas más directas de ese conflicto. Hipótesis: se plantea una tendencia ordinal creciente — a mayor frecuencia de viaje, mayor probabilidad de rotación — a contrastar empíricamente, aunque el modelo la trata como variable categórica (cada nivel comparado contra la referencia) y no impone una tendencia lineal entre niveles.

Estado_Civil (Casado / Divorciado / Soltero). El modelo clásico de rotación de Mobley (1977) y buena parte de la literatura posterior asocian los menores costos de salida (menor cantidad de ataduras económicas y familiares) con una mayor probabilidad de dejar un empleo. Los solteros, en promedio, enfrentan menores costos de oportunidad al cambiar de trabajo. Hipótesis: los solteros rotan más que los casados o divorciados.

2.2 Variables cuantitativas (cuasi-cuantitativas en el caso de la primera)

Satisfación_Laboral (escala 1-4). Es uno de los predictores más replicados en la literatura de rotación; el metaanálisis de Griffeth, Hom y Gaertner (2000) la identifica como uno de los correlatos más consistentes de la intención de rotación. Hipótesis: a mayor satisfacción laboral, menor probabilidad de rotación (signo negativo esperado).

Ingreso_Mensual. Desde la teoría de intercambio económico, el salario es la contraprestación directa más visible de la relación laboral. Hipótesis: a mayor ingreso, menor probabilidad de rotación.

Antigüedad_Cargo (años en el cargo actual). La literatura describe dos mecanismos en tensión: por un lado, mayor tiempo en el cargo puede generar arraigo y mayores costos de salida (menor rotación); por otro, puede reflejar estancamiento profesional sin promoción (mayor rotación). Planteo la hipótesis como abierta a contrastar empíricamente con el signo del coeficiente, en lugar de asumir una dirección a priori.

2.3 Verificación de colinealidad entre las variables elegidas

Antes de comprometerme con el modelo final, verifiqué que las tres variables cuantitativas elegidas no estuvieran midiendo esencialmente lo mismo. La correlación de Pearson entre ellas resultó prácticamente nula entre Satisfación_Laboral y las otras dos (\(r \approx 0\)), y moderada-baja entre Ingreso_Mensual y Antigüedad_Cargo (\(r \approx 0.36\)). No existe un umbral universal único a partir del cual una correlación se vuelve “preocupante” — la verificación definitiva la hago más adelante con el Factor de Inflación de Varianza (VIF) sobre el modelo ya ajustado, que sí tiene un criterio de referencia más establecido; esta revisión bivariada previa solo sirve como primer filtro razonable antes de comprometerme con el conjunto final de variables.

3 Carga y preparación de los datos

# Cargo las librerías necesarias para el análisis completo:
# tidyverse: manipulación de datos y graficación (incluye dplyr, ggplot2, readr)
# janitor: estandarización de nombres de columnas
# broom: extracción ordenada (tidy) de resultados de modelos
# car: cálculo del factor de inflación de varianza (VIF)
# caret: partición estratificada de datos y matriz de confusión
# pROC: curva ROC y área bajo la curva (AUC)
# kableExtra: tablas HTML con formato profesional
# scales: formato de ejes y etiquetas (porcentajes, separadores de miles)
# rstatix: cálculo de Cramér's V (tamaño del efecto para tablas de contingencia)
# effsize: cálculo de la d de Cohen (tamaño del efecto para diferencias de medias)
library(tidyverse)
library(janitor)
library(broom)
library(car)
library(caret)
library(pROC)
library(kableExtra)
library(scales)
library(rstatix)
library(effsize)

# Fijo una semilla para que la partición train/test y cualquier proceso
# aleatorio sean exactamente reproducibles por cualquier persona que
# vuelva a ejecutar este documento.
set.seed(123)
# Instalo (solo si no está ya instalado) y cargo paqueteMODELOS, que
# distribuye la base de datos "rotacion" directamente desde el
# repositorio de GitHub del curso. Esto evita depender de un archivo
# local y garantiza que trabajo con la misma versión de los datos
# que el resto del curso.
if (!requireNamespace("paqueteMODELOS", quietly = TRUE)) {
  if (!requireNamespace("devtools", quietly = TRUE)) install.packages("devtools")
  devtools::install_github("centromagis/paqueteMODELOS", force = TRUE)
}
library(paqueteMODELOS)
data("rotacion")

# clean_names() estandariza los nombres de columna: los pasa a minúsculas,
# reemplaza espacios por guiones bajos y elimina tildes/diéresis, lo que
# evita tener que escribir nombres entre comillas invertidas (backticks)
# en el resto del documento.
rotacion <- rotacion %>% clean_names()

# Reviso la estructura general del conjunto de datos.
glimpse(rotacion)
## Rows: 1,470
## Columns: 24
## $ rotacion                    <chr> "Si", "No", "Si", "No", "No", "No", "No", …
## $ edad                        <dbl> 41, 49, 37, 33, 27, 32, 59, 30, 38, 36, 35…
## $ viaje_de_negocios           <chr> "Raramente", "Frecuentemente", "Raramente"…
## $ departamento                <chr> "Ventas", "IyD", "IyD", "IyD", "IyD", "IyD…
## $ distancia_casa              <dbl> 1, 8, 2, 3, 2, 2, 3, 24, 23, 27, 16, 15, 2…
## $ educacion                   <dbl> 2, 1, 2, 4, 1, 2, 3, 1, 3, 3, 3, 2, 1, 2, …
## $ campo_educacion             <chr> "Ciencias", "Ciencias", "Otra", "Ciencias"…
## $ satisfaccion_ambiental      <dbl> 2, 3, 4, 4, 1, 4, 3, 4, 4, 3, 1, 4, 1, 2, …
## $ genero                      <chr> "F", "M", "M", "F", "M", "M", "F", "M", "M…
## $ cargo                       <chr> "Ejecutivo_Ventas", "Investigador_Cientifi…
## $ satisfacion_laboral         <dbl> 4, 2, 3, 3, 2, 4, 1, 3, 3, 3, 2, 3, 3, 4, …
## $ estado_civil                <chr> "Soltero", "Casado", "Soltero", "Casado", …
## $ ingreso_mensual             <dbl> 5993, 5130, 2090, 2909, 3468, 3068, 2670, …
## $ trabajos_anteriores         <dbl> 8, 1, 6, 1, 9, 0, 4, 1, 0, 6, 0, 0, 1, 0, …
## $ horas_extra                 <chr> "Si", "No", "Si", "Si", "No", "No", "Si", …
## $ porcentaje_aumento_salarial <dbl> 11, 23, 15, 11, 12, 13, 20, 22, 21, 13, 13…
## $ rendimiento_laboral         <dbl> 3, 4, 3, 3, 3, 3, 4, 4, 4, 3, 3, 3, 3, 3, …
## $ anos_experiencia            <dbl> 8, 10, 7, 8, 6, 8, 12, 1, 10, 17, 6, 10, 5…
## $ capacitaciones              <dbl> 0, 3, 3, 3, 3, 2, 3, 2, 2, 3, 5, 3, 1, 2, …
## $ equilibrio_trabajo_vida     <dbl> 1, 3, 3, 3, 3, 2, 2, 3, 3, 2, 3, 3, 2, 3, …
## $ antiguedad                  <dbl> 6, 10, 0, 8, 2, 7, 1, 1, 9, 7, 5, 9, 5, 2,…
## $ antiguedad_cargo            <dbl> 4, 7, 0, 7, 2, 7, 0, 0, 7, 7, 4, 5, 2, 2, …
## $ anos_ultima_promocion       <dbl> 0, 1, 0, 3, 2, 3, 0, 0, 1, 7, 0, 0, 4, 1, …
## $ anos_acargo_con_mismo_jefe  <dbl> 5, 7, 0, 0, 2, 6, 0, 0, 8, 7, 3, 8, 3, 2, …
# Cuento valores faltantes por columna. Un data frame de calidad para
# modelar debe estar libre de NA no gestionados, o al menos con una
# regla explícita de tratamiento si existen.
na_por_columna <- colSums(is.na(rotacion))
sum(na_por_columna)
## [1] 0

El conjunto de datos tiene 1.470 registros y 24 variables, sin valores faltantes en ninguna columna, por lo que no fue necesario imputar ni eliminar observaciones por esta razón. Esto simplifica la etapa de limpieza, pero no exime de revisar la calidad de las categorías (niveles mal escritos, categorías con muy pocas observaciones) antes de modelar.

# Recodifico las variables categóricas de interés como factores, fijando
# explícitamente el nivel de referencia (baseline) de cada una. Esta
# decisión no es trivial: el nivel de referencia es contra el que se
# interpretan todos los demás coeficientes del modelo.
#
# - horas_extra: "No" como referencia (el escenario sin sobrecarga).
# - viaje_de_negocios: "No_Viaja" como referencia (el escenario de
#   menor conflicto trabajo-familia), en lugar del orden alfabético
#   por defecto que dejaría "Frecuentemente" como base.
# - estado_civil: "Casado" como referencia, coincide con el orden
#   alfabético pero lo dejo explícito para que la decisión sea legible.
# - rotacion: fijo "No" como referencia para que "Si" sea la categoría
#   de éxito (Y = 1) en el modelo logístico.
rotacion <- rotacion %>%
  mutate(
    rotacion          = factor(rotacion, levels = c("No", "Si")),
    horas_extra       = factor(horas_extra, levels = c("No", "Si")),
    viaje_de_negocios = factor(viaje_de_negocios,
                                levels = c("No_Viaja", "Raramente", "Frecuentemente")),
    estado_civil      = factor(estado_civil,
                                levels = c("Casado", "Divorciado", "Soltero")),
    # Reescalo el ingreso antes de modelar. La unidad original
    # (ingreso_mensual tal como viene en el dataset) es demasiado
    # pequeña para que un cambio de una unidad produzca un odds ratio
    # distinguible de 1 — no porque el efecto no exista, sino porque
    # la unidad no es relevante en esa escala. Verifico primero el
    # rango real antes de decidir el factor de reescalamiento, en
    # lugar de asumir una denominación (por ejemplo, pesos en
    # millones) sin comprobarla.
    ingreso_miles     = ingreso_mensual / 1000
  )

range(rotacion$ingreso_mensual)
## [1]  1009 19999
mean(rotacion$ingreso_mensual)
## [1] 6502.931
# Verifico que la recodificación no haya generado NA por desajuste
# de escritura entre los niveles originales y los que definí (por
# ejemplo, diferencias de mayúsculas/minúsculas). clean_names() solo
# estandariza los NOMBRES de columna, nunca los VALORES dentro de
# ellas, así que este chequeo es indispensable cada vez que se
# recodifica manualmente a partir de un factor().
recodificacion_check <- rotacion %>%
  summarise(across(c(rotacion, horas_extra, viaje_de_negocios, estado_civil),
                    ~ sum(is.na(.x))))

print(recodificacion_check)
## # A tibble: 1 × 4
##   rotacion horas_extra viaje_de_negocios estado_civil
##      <int>       <int>             <int>        <int>
## 1        0           0                 0            0
stopifnot(
  "La recodificación generó NA: revisa que los niveles en factor(..., levels = ...) coincidan EXACTAMENTE (mayúsculas incluidas) con los valores originales del dataset" =
    sum(recodificacion_check) == 0
)

Verifico explícitamente que la recodificación no haya introducido valores faltantes, en lugar de asumirlo: una discrepancia de mayúsculas/minúsculas o de escritura entre los niveles declarados y los valores reales del dataset convertiría silenciosamente esas observaciones en NA, lo que sesgaría todo el análisis posterior sin que se note a simple vista.

El range() y la media de ingreso_mensual confirman que la variable se mueve en un orden de magnitud de miles (no de millones): dividir entre 1.000, en lugar de entre 1.000.000, es la escala que efectivamente produce una variación observable dentro del rango real de los datos. De ahí la variable derivada ingreso_miles, que es la que utilizo en el modelo en lugar de ingreso_mensual en su unidad original.

4 Análisis univariado

4.1 Variable respuesta: Rotación

# Tabla de frecuencias absolutas y relativas de la variable respuesta.
tabla_rotacion <- rotacion %>%
  count(rotacion) %>%
  mutate(porcentaje = round(n / sum(n) * 100, 1))

tabla_rotacion %>%
  kable(col.names = c("Rotación", "Frecuencia", "Porcentaje (%)"),
        align = "c", caption = "Distribución de la variable Rotación") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)
Distribución de la variable Rotación
Rotación Frecuencia Porcentaje (%)
No 1233 83.9
Si 237 16.1
ggplot(tabla_rotacion, aes(x = rotacion, y = n, fill = rotacion)) +
  geom_col(width = 0.55) +
  geom_text(aes(label = paste0(n, " (", porcentaje, "%)")),
            vjust = -0.4, size = 4.2, fontface = "bold") +
  scale_fill_manual(values = c("No" = "#2C7FB8", "Si" = "#D7191C")) +
  scale_y_continuous(limits = c(0, 1400), expand = expansion(mult = c(0, 0.05))) +
  labs(title = "Distribución de empleados según rotación",
       x = NULL, y = "Número de empleados") +
  theme_minimal(base_size = 13) +
  theme(legend.position = "none",
        plot.title = element_text(face = "bold"))

El 16.1 % de los empleados (237 de 1.470) rotó en el período observado, mientras que el 83.9 % restante permaneció en su cargo. Este es un desbalance de clases moderado, no extremo, pero suficiente para que deba tenerlo en cuenta más adelante: un modelo que simplemente prediga “No rota” para todos los empleados ya acertaría en el 83.9 % de los casos sin haber aprendido nada útil, así que la exactitud (accuracy) por sí sola será una métrica engañosa para evaluar este modelo. Volveré sobre esto en la sección de evaluación.

4.2 Variables categóricas seleccionadas

# Construyo un gráfico de distribución para cada una de las tres
# variables categóricas seleccionadas, en un solo panel con facet_wrap
# para facilitar la comparación visual.
categoricas_long <- rotacion %>%
  select(horas_extra, viaje_de_negocios, estado_civil) %>%
  pivot_longer(everything(), names_to = "variable", values_to = "nivel") %>%
  mutate(variable = dplyr::recode(variable,
                            horas_extra = "Horas extra",
                            viaje_de_negocios = "Viaje de negocios",
                            estado_civil = "Estado civil"))

ggplot(categoricas_long, aes(x = nivel, fill = variable)) +
  geom_bar(width = 0.6) +
  geom_text(stat = "count", aes(label = after_stat(count)),
            vjust = -0.4, size = 3.8) +
  facet_wrap(~ variable, scales = "free_x", ncol = 1) +
  labs(title = "Distribución de las variables categóricas seleccionadas",
       x = NULL, y = "Número de empleados") +
  theme_minimal(base_size = 13) +
  theme(legend.position = "none",
        strip.text = element_text(face = "bold", size = 12),
        plot.title = element_text(face = "bold"))

Las tres variables tienen todas sus categorías bien pobladas (ninguna con menos de 100 observaciones), lo que reduce el riesgo de que el modelo logístico presente separación cuasi-perfecta (coeficientes que divergen a valores extremos por falta de datos en alguna combinación de categorías). La mayoría de los empleados no trabaja horas extra, viaja raramente y está casado; esto es relevante porque el nivel de referencia de cada variable en el modelo corresponde, en los tres casos, a la categoría más frecuente.

4.3 Variables cuantitativas seleccionadas

cuantitativas_long <- rotacion %>%
  select(satisfacion_laboral, ingreso_mensual, antiguedad_cargo) %>%
  pivot_longer(everything(), names_to = "variable", values_to = "valor") %>%
  mutate(variable = dplyr::recode(variable,
                            satisfacion_laboral = "Satisfacción laboral (escala 1-4)",
                            ingreso_mensual = "Ingreso mensual",
                            antiguedad_cargo = "Antigüedad en el cargo (años)"))

ggplot(cuantitativas_long, aes(x = valor, fill = variable)) +
  geom_histogram(bins = 20, color = "white") +
  facet_wrap(~ variable, scales = "free", ncol = 1) +
  labs(title = "Distribución de las variables cuantitativas seleccionadas",
       x = NULL, y = "Frecuencia") +
  theme_minimal(base_size = 13) +
  theme(legend.position = "none",
        strip.text = element_text(face = "bold", size = 12),
        plot.title = element_text(face = "bold"))

rotacion %>%
  select(satisfacion_laboral, ingreso_mensual, antiguedad_cargo) %>%
  pivot_longer(everything(), names_to = "Variable", values_to = "valor") %>%
  group_by(Variable) %>%
  summarise(Media = round(mean(valor), 1),
            Mediana = round(median(valor), 1),
            `Desv. estándar` = round(sd(valor), 1),
            Mínimo = min(valor),
            Máximo = max(valor)) %>%
  mutate(Variable = dplyr::recode(Variable,
                            satisfacion_laboral = "Satisfacción laboral",
                            ingreso_mensual = "Ingreso mensual",
                            antiguedad_cargo = "Antigüedad en el cargo")) %>%
  kable(align = "c", caption = "Estadísticos descriptivos de las variables cuantitativas") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)
Estadísticos descriptivos de las variables cuantitativas
Variable Media Mediana Desv. estándar Mínimo Máximo
Antigüedad en el cargo 4.2 3 3.6 0 18
Ingreso mensual 6502.9 4919 4708.0 1009 19999
Satisfacción laboral 2.7 3 1.1 1 4

Satisfación_Laboral se distribuye de forma prácticamente uniforme entre sus cuatro niveles, sin un nivel dominante. Ingreso_Mensual presenta un sesgo marcado hacia la derecha (cola larga de salarios altos), típico de las distribuciones salariales en organizaciones con jerarquías de cargos amplias. Antigüedad_Cargo también está sesgada a la derecha, concentrada en los primeros años del cargo, con una cola de empleados que llevan mucho más tiempo en la misma posición.

5 Análisis bivariado

Analizo ahora la relación de cada una de las seis variables seleccionadas con Rotación, contrastando con la hipótesis planteada en la sección 2.

5.1 Variables categóricas frente a Rotación

# Función auxiliar que construye, para una variable categórica dada,
# la tabla de porcentaje de rotación por nivel.
tabla_pct_rotacion <- function(datos, variable) {
  datos %>%
    count({{ variable }}, rotacion) %>%
    group_by({{ variable }}) %>%
    mutate(porcentaje = round(n / sum(n) * 100, 1)) %>%
    ungroup() %>%
    filter(rotacion == "Si") %>%
    select({{ variable }}, porcentaje)
}

pct_horas_extra   <- tabla_pct_rotacion(rotacion, horas_extra)
pct_viaje         <- tabla_pct_rotacion(rotacion, viaje_de_negocios)
pct_estado_civil  <- tabla_pct_rotacion(rotacion, estado_civil)

pct_horas_extra %>%
  kable(col.names = c("Horas extra", "% que rotó"), align = "c",
        caption = "Porcentaje de rotación según horas extra") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)
Porcentaje de rotación según horas extra
Horas extra % que rotó
No 10.4
Si 30.5
pct_viaje %>%
  kable(col.names = c("Viaje de negocios", "% que rotó"), align = "c",
        caption = "Porcentaje de rotación según frecuencia de viaje") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)
Porcentaje de rotación según frecuencia de viaje
Viaje de negocios % que rotó
No_Viaja 8.0
Raramente 15.0
Frecuentemente 24.9
pct_estado_civil %>%
  kable(col.names = c("Estado civil", "% que rotó"), align = "c",
        caption = "Porcentaje de rotación según estado civil") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)
Porcentaje de rotación según estado civil
Estado civil % que rotó
Casado 12.5
Divorciado 10.1
Soltero 25.5
bivariado_cat_long <- bind_rows(
  pct_horas_extra %>% rename(nivel = horas_extra) %>% mutate(variable = "Horas extra"),
  pct_viaje %>% rename(nivel = viaje_de_negocios) %>% mutate(variable = "Viaje de negocios"),
  pct_estado_civil %>% rename(nivel = estado_civil) %>% mutate(variable = "Estado civil")
)

ggplot(bivariado_cat_long, aes(x = nivel, y = porcentaje, fill = variable)) +
  geom_col(width = 0.55) +
  geom_text(aes(label = paste0(porcentaje, "%")), vjust = -0.4, size = 4) +
  facet_wrap(~ variable, scales = "free_x", ncol = 1) +
  scale_y_continuous(limits = c(0, 40), expand = expansion(mult = c(0, 0.05))) +
  labs(title = "Porcentaje de rotación por categoría",
       x = NULL, y = "% de empleados que rotó") +
  theme_minimal(base_size = 13) +
  theme(legend.position = "none",
        strip.text = element_text(face = "bold", size = 12),
        plot.title = element_text(face = "bold"))

# Aplico la prueba chi-cuadrado de independencia a cada variable
# categórica frente a Rotación. La hipótesis nula es que la variable
# categórica y la rotación son independientes; un valor p pequeño
# lleva a rechazarla, es decir, a concluir que sí existe asociación.
chi_horas   <- chisq.test(table(rotacion$horas_extra, rotacion$rotacion))
chi_viaje   <- chisq.test(table(rotacion$viaje_de_negocios, rotacion$rotacion))
chi_civil   <- chisq.test(table(rotacion$estado_civil, rotacion$rotacion))

# Un valor p pequeño solo dice que la asociación es poco probable bajo
# independencia, no dice qué tan grande es esa asociación — con 1.470
# observaciones casi cualquier asociación real, por pequeña que sea,
# termina siendo "significativa". Por eso acompaño cada prueba con
# Cramér's V, que cuantifica la magnitud de la asociación en una
# escala de 0 (independencia total) a 1 (asociación perfecta),
# independiente del tamaño de la muestra.
v_horas <- rstatix::cramer_v(table(rotacion$horas_extra, rotacion$rotacion))
v_viaje <- rstatix::cramer_v(table(rotacion$viaje_de_negocios, rotacion$rotacion))
v_civil <- rstatix::cramer_v(table(rotacion$estado_civil, rotacion$rotacion))

tibble(
  Variable = c("Horas extra", "Viaje de negocios", "Estado civil"),
  `Estadístico chi-cuadrado` = round(c(chi_horas$statistic, chi_viaje$statistic, chi_civil$statistic), 2),
  `Grados de libertad` = c(chi_horas$parameter, chi_viaje$parameter, chi_civil$parameter),
  `Valor p` = signif(c(chi_horas$p.value, chi_viaje$p.value, chi_civil$p.value), 3),
  `Cramér's V` = round(c(v_horas, v_viaje, v_civil), 3)
) %>%
  kable(align = "c", caption = "Pruebas chi-cuadrado de independencia y tamaño del efecto (Cramér's V)") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)
Pruebas chi-cuadrado de independencia y tamaño del efecto (Cramér’s V)
Variable Estadístico chi-cuadrado Grados de libertad Valor p Cramér’s V
Horas extra 87.56 1 0.0e+00 0.246
Viaje de negocios 24.18 2 5.6e-06 0.128
Estado civil 46.16 2 0.0e+00 0.177

Interpreto Cramér’s V con las convenciones habituales, que no son un estándar estricto sino una referencia orientativa: valores cercanos a 0.1 se consideran una asociación pequeña, alrededor de 0.3 moderada, y por encima de 0.5 grande.

Las tres variables categóricas muestran una asociación estadísticamente significativa con la rotación de forma bivariada (valor p muy por debajo de 0.05 en los tres casos, con tamaños del efecto que reporto en la tabla de Cramér’s V más abajo). Horas_Extra exhibe la diferencia más marcada: cerca del triple de rotación entre quienes trabajan horas extra frente a quienes no. Viaje_de_Negocios sigue una tendencia creciente consistente con la hipótesis planteada: la proporción de rotación aumenta a medida que aumenta la frecuencia de viaje. Estado_Civil es consistente con la hipótesis de Mobley: los solteros rotan considerablemente más que casados y divorciados, mientras que estos dos últimos grupos son similares entre sí.

Es importante ser preciso sobre el alcance de esta evidencia: estas son asociaciones bivariadas, sin controlar por las demás variables. Una variable puede resultar significativa de forma aislada y dejar de serlo — o cambiar de magnitud — cuando se controla simultáneamente por las demás en el modelo conjunto de la sección siguiente. Retomo esta comparación explícitamente al interpretar los coeficientes del modelo multivariado.

5.2 Variables cuantitativas frente a Rotación

cuant_bivariado_long <- rotacion %>%
  select(rotacion, satisfacion_laboral, ingreso_mensual, antiguedad_cargo) %>%
  pivot_longer(-rotacion, names_to = "variable", values_to = "valor") %>%
  mutate(variable = dplyr::recode(variable,
                            satisfacion_laboral = "Satisfacción laboral",
                            ingreso_mensual = "Ingreso mensual",
                            antiguedad_cargo = "Antigüedad en el cargo"))

ggplot(cuant_bivariado_long, aes(x = rotacion, y = valor, fill = rotacion)) +
  geom_boxplot(width = 0.5, outlier.alpha = 0.4) +
  facet_wrap(~ variable, scales = "free_y", ncol = 1) +
  scale_fill_manual(values = c("No" = "#2C7FB8", "Si" = "#D7191C")) +
  labs(title = "Distribución de las variables cuantitativas según rotación",
       x = "Rotación", y = NULL) +
  theme_minimal(base_size = 13) +
  theme(legend.position = "none",
        strip.text = element_text(face = "bold", size = 12),
        plot.title = element_text(face = "bold"))

# Aplico la prueba t de Welch (no asume varianzas iguales entre los dos
# grupos, una precaución razonable cuando los tamaños de grupo son muy
# distintos, como ocurre aquí: 1233 vs 237 empleados).
t_satisfaccion <- t.test(satisfacion_laboral ~ rotacion, data = rotacion)
t_ingreso      <- t.test(ingreso_mensual ~ rotacion, data = rotacion)
t_antiguedad   <- t.test(antiguedad_cargo ~ rotacion, data = rotacion)

# La d de Cohen cuantifica la magnitud de la diferencia de medias en
# desviaciones estándar conjuntas, independientemente del tamaño de la
# muestra. Con n = 1.470 casi cualquier diferencia real resulta
# "significativa" en el sentido del valor p; la d de Cohen permite
# distinguir una diferencia significativa pero trivial de una que
# además es sustantivamente relevante.
d_satisfaccion <- effsize::cohen.d(satisfacion_laboral ~ rotacion, data = rotacion)$estimate
d_ingreso      <- effsize::cohen.d(ingreso_mensual ~ rotacion, data = rotacion)$estimate
d_antiguedad   <- effsize::cohen.d(antiguedad_cargo ~ rotacion, data = rotacion)$estimate

tibble(
  Variable = c("Satisfacción laboral", "Ingreso mensual", "Antigüedad en el cargo"),
  `Media (No rota)` = round(c(t_satisfaccion$estimate[1], t_ingreso$estimate[1], t_antiguedad$estimate[1]), 1),
  `Media (Sí rota)` = round(c(t_satisfaccion$estimate[2], t_ingreso$estimate[2], t_antiguedad$estimate[2]), 1),
  `Valor p (t de Welch)` = signif(c(t_satisfaccion$p.value, t_ingreso$p.value, t_antiguedad$p.value), 3),
  `d de Cohen` = round(abs(c(d_satisfaccion, d_ingreso, d_antiguedad)), 3)
) %>%
  kable(align = "c", caption = "Prueba t de Welch y tamaño del efecto (d de Cohen) por grupo de rotación") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)
Prueba t de Welch y tamaño del efecto (d de Cohen) por grupo de rotación
Variable Media (No rota) Media (Sí rota) Valor p (t de Welch) d de Cohen
Satisfacción laboral 2.8 2.5 0.000105 0.283
Ingreso mensual 6832.7 4787.1 0.000000 0.440
Antigüedad en el cargo 4.5 2.9 0.000000 0.442

Interpreto la d de Cohen con las convenciones usuales (también orientativas, no un corte estricto): alrededor de 0.2 es un efecto pequeño, 0.5 moderado, 0.8 o más grande.

Satisfación_Laboral es una escala ordinal de 1 a 4 que en este informe trato como cuasi-cuantitativa; la prueba t asume una escala de intervalo, supuesto defendible pero no automático. Como verificación de sensibilidad, repito la comparación con una prueba no paramétrica que no depende de ese supuesto:

# Wilcoxon-Mann-Whitney no asume que la variable sea de intervalo ni
# que su distribución sea normal; solo requiere que sea al menos
# ordinal, lo que la hace más apropiada en sentido estricto para
# Satisfación_Laboral que la prueba t.
wilcox_satisfaccion <- wilcox.test(satisfacion_laboral ~ rotacion, data = rotacion)

tibble(
  Prueba = c("t de Welch", "Wilcoxon-Mann-Whitney"),
  `Valor p` = signif(c(t_satisfaccion$p.value, wilcox_satisfaccion$p.value), 3)
) %>%
  kable(align = "c", caption = "Satisfacción laboral: comparación entre prueba paramétrica y no paramétrica") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)
Satisfacción laboral: comparación entre prueba paramétrica y no paramétrica
Prueba Valor p
t de Welch 1.05e-04
Wilcoxon-Mann-Whitney 7.96e-05

Ambas pruebas coinciden en la conclusión (asociación estadísticamente significativa), lo que da tranquilidad de que tratar Satisfación_Laboral como cuasi-cuantitativa no está distorsionando el resultado.

Las diferencias de medias son estadísticamente significativas en las tres variables cuantitativas, con tamaños del efecto que van de pequeños a moderados según la tabla anterior. Los empleados que rotaron tienen, en promedio, menor satisfacción laboral, menor ingreso mensual y menor antigüedad en el cargo que quienes permanecieron. Esto es consistente, de forma bivariada, con la hipótesis de Satisfación_Laboral e Ingreso_Mensual, y con una de las dos direcciones teóricas que había dejado abiertas para Antigüedad_Cargo (arraigo/costo de salida, en vez de estancamiento profesional) — aunque, igual que con las categóricas, esta es evidencia bivariada; confirmo si se sostiene al controlar por las demás variables en la siguiente sección.

6 Estimación del modelo de regresión logística

6.1 Partición de los datos en entrenamiento y prueba

Antes de estimar el modelo separo los datos en un conjunto de entrenamiento (para ajustar los coeficientes) y uno de prueba (para evaluar qué tan bien generaliza el modelo a datos que no vio durante el ajuste). Evaluar el modelo sobre los mismos datos con los que se estimó tiende a sobreestimar su capacidad predictiva real.

Uso caret::createDataPartition() en lugar de un muestreo aleatorio simple porque estratifica por la variable respuesta: garantiza que la proporción de 16.1 % de rotación se mantenga aproximadamente igual tanto en el conjunto de entrenamiento como en el de prueba. Con un desbalance de clases como el de este conjunto de datos, un muestreo aleatorio simple podría, por azar, dejar muy pocos casos de rotación en el conjunto de prueba y volver inestable la evaluación.

indice_entrenamiento <- createDataPartition(rotacion$rotacion, p = 0.70, list = FALSE)
train <- rotacion[indice_entrenamiento, ]
test  <- rotacion[-indice_entrenamiento, ]

# Verifico que la proporción de rotación se haya mantenido similar en
# ambos conjuntos, confirmando que la estratificación funcionó.
tibble(
  Conjunto = c("Entrenamiento", "Prueba"),
  n = c(nrow(train), nrow(test)),
  `% Rotación` = round(c(mean(train$rotacion == "Si"), mean(test$rotacion == "Si")) * 100, 1)
) %>%
  kable(align = "c", caption = "Tamaño y balance de clases en cada conjunto") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)
Tamaño y balance de clases en cada conjunto
Conjunto n % Rotación
Entrenamiento 1030 16.1
Prueba 440 16.1

6.2 Ajuste del modelo

# Estimo el modelo logit por máxima verosimilitud con glm(), family =
# binomial(link = "logit"). El modelo incluye las tres variables
# categóricas y las tres cuantitativas seleccionadas y justificadas
# en la sección 2, estimado únicamente con los datos de entrenamiento.
modelo_logit <- glm(
  rotacion ~ horas_extra + viaje_de_negocios + estado_civil +
             satisfacion_laboral + ingreso_miles + antiguedad_cargo,
  family = binomial(link = "logit"),
  data = train
)

summary(modelo_logit)
## 
## Call:
## glm(formula = rotacion ~ horas_extra + viaje_de_negocios + estado_civil + 
##     satisfacion_laboral + ingreso_miles + antiguedad_cargo, family = binomial(link = "logit"), 
##     data = train)
## 
## Coefficients:
##                                 Estimate Std. Error z value Pr(>|z|)    
## (Intercept)                     -2.06355    0.54487  -3.787 0.000152 ***
## horas_extraSi                    1.43869    0.19307   7.452 9.23e-14 ***
## viaje_de_negociosRaramente       1.13912    0.46136   2.469 0.013547 *  
## viaje_de_negociosFrecuentemente  1.94738    0.48424   4.021 5.78e-05 ***
## estado_civilDivorciado          -0.02691    0.27171  -0.099 0.921106    
## estado_civilSoltero              1.09627    0.21194   5.172 2.31e-07 ***
## satisfacion_laboral             -0.31939    0.08455  -3.778 0.000158 ***
## ingreso_miles                   -0.09290    0.02741  -3.389 0.000701 ***
## antiguedad_cargo                -0.12126    0.03296  -3.679 0.000234 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 909.69  on 1029  degrees of freedom
## Residual deviance: 738.39  on 1021  degrees of freedom
## AIC: 756.39
## 
## Number of Fisher Scoring iterations: 5

6.3 Interpretación de los coeficientes mediante razones de probabilidad (odds ratios)

Los coeficientes del modelo logit están en la escala del logaritmo de odds, que no es directamente interpretable en términos prácticos. Los exponencio para obtener razones de probabilidad (odds ratios), acompañadas de su intervalo de confianza al 95 %.

# broom::tidy() con exponentiate = TRUE devuelve directamente los odds
# ratios y sus intervalos de confianza, en lugar de tener que aplicar
# exp() manualmente a cada coeficiente.
tabla_or <- tidy(modelo_logit, conf.int = TRUE, exponentiate = TRUE) %>%
  mutate(across(where(is.numeric), ~ round(.x, 3)))

tabla_or %>%
  kable(col.names = c("Término", "Odds ratio", "Error estándar", "Estadístico z",
                       "Valor p", "IC 95% inferior", "IC 95% superior"),
        align = "c", caption = "Odds ratios e intervalos de confianza del modelo logit") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)
Odds ratios e intervalos de confianza del modelo logit
Término Odds ratio Error estándar Estadístico z Valor p IC 95% inferior IC 95% superior
(Intercept) 0.127 0.545 -3.787 0.000 0.040 0.350
horas_extraSi 4.215 0.193 7.452 0.000 2.894 6.175
viaje_de_negociosRaramente 3.124 0.461 2.469 0.014 1.364 8.554
viaje_de_negociosFrecuentemente 7.010 0.484 4.021 0.000 2.906 19.925
estado_civilDivorciado 0.973 0.272 -0.099 0.921 0.565 1.645
estado_civilSoltero 2.993 0.212 5.172 0.000 1.984 4.560
satisfacion_laboral 0.727 0.085 -3.778 0.000 0.615 0.857
ingreso_miles 0.911 0.027 -3.389 0.001 0.861 0.959
antiguedad_cargo 0.886 0.033 -3.679 0.000 0.829 0.944

Interpreto cada odds ratio como el factor multiplicativo que aplica a los odds de rotar, manteniendo las demás variables constantes:

  • Horas extra (Sí vs. No): los odds de rotar de un empleado que trabaja horas extra son varias veces más altos que los de uno que no las trabaja. Es, con claridad, el efecto individual más fuerte del modelo.
  • Viaje de negocios: viajar frecuentemente multiplica los odds de rotar frente a no viajar; viajar raramente también los aumenta, aunque en menor magnitud y con un valor p más cercano al límite de significancia convencional (0.05), por lo que interpreto ese coeficiente en particular con algo más de cautela.
  • Estado civil: ser soltero (frente a casado) multiplica los odds de rotar; ser divorciado no muestra una diferencia estadísticamente distinguible de ser casado, controlando por las demás variables del modelo.
  • Satisfacción laboral: cada punto adicional en la escala de satisfacción (1 a 4) multiplica los odds de rotar por un factor menor a 1, es decir, los reduce. Es un efecto protector.
  • Ingreso mensual (en miles): el modelo usa ingreso_miles en lugar del ingreso en su unidad original, precisamente para que el odds ratio sea interpretable en una escala económicamente relevante. El odds ratio reportado corresponde al cambio en los odds de rotar por cada 1.000 unidades monetarias adicionales de ingreso mensual, y es menor a 1: confirma el efecto protector esperado.
  • Antigüedad en el cargo: cada año adicional en el cargo multiplica los odds de rotar por un factor menor a 1, confirmando la dirección negativa observada en el análisis bivariado.

6.4 Verificación del supuesto de ausencia de multicolinealidad

# car::vif() calcula el Factor de Inflación de Varianza. Para
# predictores categóricos con más de dos niveles se reporta el
# GVIF junto con GVIF^(1/(2*Df)), que es la cantidad comparable
# directamente con el umbral usual de referencia (valores por
# encima de 5, o de forma más conservadora de 2.5-3, sugieren
# preocupación).
vif(modelo_logit)
##                         GVIF Df GVIF^(1/(2*Df))
## horas_extra         1.038505  1        1.019071
## viaje_de_negocios   1.027648  2        1.006841
## estado_civil        1.045188  2        1.011110
## satisfacion_laboral 1.029280  1        1.014534
## ingreso_miles       1.082913  1        1.040631
## antiguedad_cargo    1.079651  1        1.039062

Ningún valor se acerca al umbral de preocupación, lo que confirma en el modelo ajustado la baja colinealidad que ya había verificado de forma bivariada antes de seleccionar las variables. Los valores de VIF/GVIF no evidencian problemas relevantes de multicolinealidad entre los predictores incluidos, lo que reduce la preocupación por inflación de las varianzas asociada a redundancia lineal entre variables. Esto no implica, por sí solo, que las pruebas de significancia individual sean válidas en un sentido absoluto — su confiabilidad depende también de la especificación correcta del modelo, la independencia de las observaciones y otros supuestos que discuto en la sección de limitaciones — pero sí descarta específicamente la multicolinealidad como explicación de algún coeficiente no significativo.

6.5 Observaciones influyentes (distancia de Cook)

Antes de pasar a la significancia global, reviso si hay observaciones individuales que estén influyendo desproporcionadamente sobre el ajuste del modelo — un empleado atípico que, por sí solo, mueva de forma importante los coeficientes estimados.

# La distancia de Cook, adaptada a GLM, mide cuánto cambiarían los
# coeficientes estimados si se eliminara esa observación del ajuste.
# No hay un punto de corte universal; una referencia común (no un
# estándar estricto) es 4/n, que uso aquí solo como guía visual.
cooksd <- cooks.distance(modelo_logit)
umbral_cook <- 4 / nrow(train)

cooks_df <- tibble(indice = seq_along(cooksd), distancia_cook = cooksd) %>%
  mutate(influyente = distancia_cook > umbral_cook)

ggplot(cooks_df, aes(x = indice, y = distancia_cook, color = influyente)) +
  geom_point(size = 1.4, alpha = 0.7) +
  geom_hline(yintercept = umbral_cook, linetype = "dashed", color = "grey40") +
  scale_color_manual(values = c("FALSE" = "#2C7FB8", "TRUE" = "#D7191C")) +
  labs(title = "Distancia de Cook por observación (conjunto de entrenamiento)",
       x = "Índice de la observación", y = "Distancia de Cook",
       color = paste0("Por encima\ndel umbral (", round(umbral_cook, 4), ")")) +
  theme_minimal(base_size = 13) +
  theme(plot.title = element_text(face = "bold"))

n_influyentes <- sum(cooks_df$influyente)

De las 1030 observaciones del conjunto de entrenamiento, 71 quedan por encima del umbral de referencia (4/n). Dado que ese umbral es solo orientativo — con una muestra de este tamaño es normal que un número pequeño de observaciones lo supere sin que eso implique necesariamente un problema — no las elimino automáticamente del modelo. Si ese número fuera una proporción grande del total, o si uno o dos puntos concentraran una distancia de Cook muy por encima del resto (órdenes de magnitud, no solo por encima del umbral), ahí sí valdría la pena inspeccionar esos registros puntuales antes de confiar en el modelo.

6.6 Significancia global del modelo

# Contraste global: H0: todos los coeficientes (excepto el intercepto)
# son cero, frente a Ha: al menos uno es distinto de cero. El
# estadístico de prueba es la diferencia entre la devianza nula y la
# devianza residual, que se distribuye chi-cuadrado con los grados de
# libertad correspondientes a la diferencia de parámetros estimados.
estadistico_chi2 <- modelo_logit$null.deviance - modelo_logit$deviance
gl <- modelo_logit$df.null - modelo_logit$df.residual
valor_p_global <- pchisq(estadistico_chi2, df = gl, lower.tail = FALSE)

tibble(
  `Devianza nula` = round(modelo_logit$null.deviance, 2),
  `Devianza residual` = round(modelo_logit$deviance, 2),
  `Chi-cuadrado` = round(estadistico_chi2, 2),
  `Grados de libertad` = gl,
  `Valor p` = signif(valor_p_global, 3)
) %>%
  kable(align = "c", caption = "Prueba de significancia global del modelo") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)
Prueba de significancia global del modelo
Devianza nula Devianza residual Chi-cuadrado Grados de libertad Valor p
909.69 738.39 171.3 8 0

El valor p de la prueba global es prácticamente cero, por lo que rechazo \(H_0\): el conjunto de variables seleccionadas, en su totalidad, sí aporta capacidad explicativa significativa sobre la rotación, más allá de lo que aportaría un modelo con solo el intercepto.

7 Evaluación del modelo

7.1 Matriz de confusión y métricas de clasificación

# Genero las probabilidades predichas por el modelo sobre el conjunto
# de prueba (datos que el modelo no vio al momento de estimarse).
test <- test %>%
  mutate(prob_predicha = predict(modelo_logit, newdata = test, type = "response"))

# Clasifico con el punto de corte convencional de 0.5 como primer
# acercamiento.
test <- test %>%
  mutate(clase_predicha_05 = factor(if_else(prob_predicha > 0.5, "Si", "No"),
                                     levels = c("No", "Si")))

matriz_05 <- confusionMatrix(test$clase_predicha_05, test$rotacion, positive = "Si")
matriz_05
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction  No  Si
##         No 359  61
##         Si  10  10
##                                           
##                Accuracy : 0.8386          
##                  95% CI : (0.8009, 0.8718)
##     No Information Rate : 0.8386          
##     P-Value [Acc > NIR] : 0.5316          
##                                           
##                   Kappa : 0.1602          
##                                           
##  Mcnemar's Test P-Value : 2.958e-09       
##                                           
##             Sensitivity : 0.14085         
##             Specificity : 0.97290         
##          Pos Pred Value : 0.50000         
##          Neg Pred Value : 0.85476         
##              Prevalence : 0.16136         
##          Detection Rate : 0.02273         
##    Detection Prevalence : 0.04545         
##       Balanced Accuracy : 0.55687         
##                                           
##        'Positive' Class : Si              
## 
# Visualizo la matriz de confusión como un mapa de calor con las
# frecuencias impresas dentro de cada celda, evitando que el texto se
# empalme con las etiquetas de los ejes.
matriz_df <- as.data.frame(matriz_05$table)

ggplot(matriz_df, aes(x = Reference, y = Prediction, fill = Freq)) +
  geom_tile(color = "white", linewidth = 1) +
  geom_text(aes(label = Freq), color = "black", size = 7, fontface = "bold") +
  scale_fill_gradient(low = "#DEEBF7", high = "#2C7FB8") +
  labs(title = "Matriz de confusión (punto de corte = 0.5)",
       x = "Valor observado", y = "Valor predicho", fill = "Frecuencia") +
  theme_minimal(base_size = 13) +
  theme(plot.title = element_text(face = "bold"))

Con el punto de corte convencional de 0.5, el modelo alcanza una exactitud alta, pero esa exactitud es engañosa por el desbalance de clases que señalé en la sección de análisis univariado. La sensibilidad (la capacidad de detectar correctamente a los empleados que sí van a rotar) resulta baja: el modelo tiende a clasificar como “No rota” a muchos empleados que en realidad sí rotaron, precisamente porque el punto de corte de 0.5 favorece a la clase mayoritaria. La especificidad, en cambio, es muy alta: el modelo es muy bueno identificando a quienes van a permanecer.

Para el objetivo de negocio de este problema (identificar preventivamente a empleados en riesgo de rotar), una sensibilidad baja es una limitación seria: significa dejar pasar sin intervención a una parte importante de los empleados que efectivamente van a renunciar.

7.2 Curva ROC y área bajo la curva (AUC)

# pROC::roc() construye la curva ROC comparando, para cada punto de
# corte posible, la tasa de verdaderos positivos (sensibilidad) contra
# la tasa de falsos positivos (1 - especificidad).
curva_roc <- roc(test$rotacion, test$prob_predicha, levels = c("No", "Si"), direction = "<")
auc_valor <- auc(curva_roc)

# ci.auc() calcula el intervalo de confianza del AUC por bootstrap
# (por defecto, 2000 remuestreos). Reportar el AUC sin su intervalo
# da la falsa impresión de que es un valor exacto y no una estimación
# sujeta a la variabilidad de la muestra de prueba.
ic_auc <- ci.auc(curva_roc)

ggroc(curva_roc, color = "#D7191C", linewidth = 1.1) +
  geom_abline(slope = 1, intercept = 1, linetype = "dashed", color = "grey50") +
  labs(title = paste0("Curva ROC del modelo (AUC = ", round(auc_valor, 4), ")"),
       subtitle = paste0("IC 95% (bootstrap): [", round(ic_auc[1], 4), ", ", round(ic_auc[3], 4), "]"),
       x = "Especificidad", y = "Sensibilidad") +
  theme_minimal(base_size = 13) +
  theme(plot.title = element_text(face = "bold"))

El modelo alcanza un AUC de 0.7704 (IC 95%: 0.7112 - 0.8296) sobre el conjunto de prueba. Según las guías convencionales de discriminación de modelos de clasificación, esto se interpreta como una capacidad discriminativa aceptable a buena: bien por encima de 0.5 (una clasificación equivalente al azar) pero por debajo del rango de 0.9+ que se consideraría excelente. El intervalo de confianza, calculado por remuestreo bootstrap sobre el conjunto de prueba, recuerda que esta estimación está sujeta a la variabilidad propia de evaluar el modelo sobre una única partición train/test — una validación cruzada repetida daría una estimación más estable del rango plausible del AUC, algo que señalo como limitación más adelante y no desarrollo aquí por alcance de esta actividad.

7.3 Análisis de sensibilidad del punto de corte

El punto de corte de 0.5 y el óptimo de Youden son solo dos puntos entre infinitos posibles. Antes de fijar uno, reviso cómo se comportan sensibilidad, especificidad y precisión a lo largo de una rejilla de umbrales, para entender el trade-off completo en vez de comparar únicamente dos valores.

# coords() con una rejilla explícita de umbrales permite leer cómo
# evolucionan las métricas de clasificación en todo el rango de
# cortes posibles, no solo en un punto.
grilla_cortes <- coords(
  curva_roc,
  x = seq(0.10, 0.60, by = 0.05),
  input = "threshold",
  ret = c("threshold", "sensitivity", "specificity", "ppv", "npv")
)

grilla_cortes %>%
  mutate(across(where(is.numeric), ~ round(.x, 3))) %>%
  kable(col.names = c("Punto de corte", "Sensibilidad", "Especificidad", "Precisión (PPV)", "VPN"),
        align = "c", caption = "Métricas de clasificación a lo largo de una rejilla de puntos de corte") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)
Métricas de clasificación a lo largo de una rejilla de puntos de corte
Punto de corte Sensibilidad Especificidad Precisión (PPV) VPN
0.10 0.859 0.547 0.268 0.953
0.15 0.718 0.672 0.297 0.925
0.20 0.606 0.778 0.344 0.911
0.25 0.549 0.840 0.398 0.906
0.30 0.465 0.873 0.412 0.894
0.35 0.380 0.892 0.403 0.882
0.40 0.282 0.924 0.417 0.870
0.45 0.183 0.954 0.433 0.859
0.50 0.141 0.973 0.500 0.855
0.55 0.113 0.978 0.500 0.851
0.60 0.113 0.989 0.667 0.853
grilla_larga <- grilla_cortes %>%
  select(threshold, sensitivity, specificity) %>%
  pivot_longer(-threshold, names_to = "metrica", values_to = "valor") %>%
  mutate(metrica = dplyr::recode(metrica, sensitivity = "Sensibilidad", specificity = "Especificidad"))

ggplot(grilla_larga, aes(x = threshold, y = valor, color = metrica)) +
  geom_line(linewidth = 1) +
  geom_point(size = 2) +
  geom_vline(xintercept = 0.5, linetype = "dashed", color = "grey50") +
  scale_color_manual(values = c("Sensibilidad" = "#D7191C", "Especificidad" = "#2C7FB8")) +
  scale_y_continuous(labels = scales::percent_format(accuracy = 1)) +
  labs(title = "Sensibilidad y especificidad según el punto de corte",
       subtitle = "La línea vertical marca el corte convencional de 0.5",
       x = "Punto de corte", y = "Valor de la métrica", color = NULL) +
  theme_minimal(base_size = 13) +
  theme(plot.title = element_text(face = "bold"), legend.position = "bottom")

La tabla y el gráfico muestran con claridad el trade-off: a medida que el punto de corte baja, la sensibilidad sube y la especificidad baja, de forma continua a lo largo de toda la rejilla, no solo entre dos valores puntuales. Esto es lo que permite elegir un corte de forma informada en vez de arbitraria.

7.4 Selección de un punto de corte de referencia

Dentro de ese trade-off, el índice de Youden (que maximiza sensibilidad + especificidad − 1) es un criterio estándar, reproducible y defendible como punto de partida cuando no se dispone de una estimación explícita de los costos relativos de cada tipo de error.

corte_optimo <- coords(curva_roc, "best", best.method = "youden",
                        ret = c("threshold", "sensitivity", "specificity"))
corte_optimo %>%
  kable(align = "c", caption = "Punto de corte según el índice de Youden") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)
Punto de corte según el índice de Youden
threshold sensitivity specificity
0.1634485 0.7183099 0.7235772

Es importante ser preciso sobre lo que Youden sí y no optimiza: maximiza la suma de sensibilidad y especificidad por igual, pero no incorpora los costos reales y potencialmente distintos de un falso negativo (un empleado que rota y el modelo no detecta) frente a un falso positivo (un empleado que se interviene preventivamente sin que fuera a rotar). Si la organización pudiera cuantificar esos costos — por ejemplo, el costo de reemplazar a un empleado clave frente al costo de una conversación o revisión salarial preventiva —, el punto de corte óptimo en términos de negocio podría ser distinto (probablemente más bajo que el de Youden, dado que el costo de un falso negativo suele ser mayor en problemas de retención de talento). La tabla de la rejilla anterior es precisamente la herramienta que permitiría a la gerencia elegir ese corte una vez tenga esa información de costos; Youden es el punto de partida razonable en ausencia de ella, no necesariamente el óptimo final.

test <- test %>%
  mutate(clase_predicha_optima = factor(
    if_else(prob_predicha > corte_optimo$threshold, "Si", "No"),
    levels = c("No", "Si")
  ))

matriz_optima <- confusionMatrix(test$clase_predicha_optima, test$rotacion, positive = "Si")
matriz_optima
## Confusion Matrix and Statistics
## 
##           Reference
## Prediction  No  Si
##         No 267  20
##         Si 102  51
##                                           
##                Accuracy : 0.7227          
##                  95% CI : (0.6784, 0.7641)
##     No Information Rate : 0.8386          
##     P-Value [Acc > NIR] : 1               
##                                           
##                   Kappa : 0.3014          
##                                           
##  Mcnemar's Test P-Value : 2.244e-13       
##                                           
##             Sensitivity : 0.7183          
##             Specificity : 0.7236          
##          Pos Pred Value : 0.3333          
##          Neg Pred Value : 0.9303          
##              Prevalence : 0.1614          
##          Detection Rate : 0.1159          
##    Detection Prevalence : 0.3477          
##       Balanced Accuracy : 0.7209          
##                                           
##        'Positive' Class : Si              
## 

Al bajar el punto de corte desde 0.5 hasta el valor de Youden, la sensibilidad aumenta de forma sustancial a costa de una especificidad algo menor, confirmando en un punto concreto el trade-off que ya se veía en la rejilla completa. La elección final de dónde operar sobre ese trade-off es, en el fondo, una decisión de negocio y no puramente estadística: depende de cómo la organización valore el costo de un falso negativo frente a un falso positivo, algo que este análisis no cuantifica pero que la tabla de sensibilidad del punto anterior deja lista para que la gerencia decida con esa información.

8 Predicción de un caso hipotético

Construyo el perfil de un empleado hipotético para ilustrar cómo se usaría el modelo en la práctica: un empleado que trabaja horas extra, viaja frecuentemente, está soltero, reporta baja satisfacción laboral (nivel 1 de 4), tiene un ingreso mensual notoriamente por debajo del promedio de la organización, y lleva apenas un año en su cargo actual.

# En lugar de asumir una cifra de ingreso y afirmar que está "por
# debajo de la media", calculo la media real de ingreso_mensual y
# construyo el perfil hipotético como una fracción explícita de esa
# media, de modo que la afirmación quede demostrada por el propio
# código y no solo enunciada en el texto.
media_ingreso <- mean(train$ingreso_mensual)
ingreso_hipotetico <- round(media_ingreso * 0.5)

empleado_hipotetico <- tibble(
  horas_extra         = factor("Si", levels = levels(train$horas_extra)),
  viaje_de_negocios   = factor("Frecuentemente", levels = levels(train$viaje_de_negocios)),
  estado_civil        = factor("Soltero", levels = levels(train$estado_civil)),
  satisfacion_laboral = 1,
  ingreso_miles       = ingreso_hipotetico / 1000,
  antiguedad_cargo    = 1
)

tibble(
  `Ingreso del empleado hipotético` = ingreso_hipotetico,
  `Ingreso medio de la organización` = round(media_ingreso, 1),
  `¿Está por debajo de la media?` = ingreso_hipotetico < media_ingreso
) %>%
  kable(align = "c", caption = "Verificación del ingreso hipotético frente a la media real") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)
Verificación del ingreso hipotético frente a la media real
Ingreso del empleado hipotético Ingreso medio de la organización ¿Está por debajo de la media?
3215 6430.6 TRUE
probabilidad_hipotetica <- predict(modelo_logit, newdata = empleado_hipotetico, type = "response")

tibble(
  `Probabilidad estimada de rotar` = round(probabilidad_hipotetica, 3),
  `Punto de corte de referencia` = round(corte_optimo$threshold, 3),
  `Decisión sugerida` = if_else(probabilidad_hipotetica > corte_optimo$threshold,
                                 "Intervenir preventivamente", "No requiere intervención prioritaria")
) %>%
  kable(align = "c", caption = "Predicción para el empleado hipotético") %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE)
Predicción para el empleado hipotético
Probabilidad estimada de rotar Punto de corte de referencia Decisión sugerida
0.843 0.163 Intervenir preventivamente

El perfil hipotético combina, de manera deliberada, todos los factores de riesgo identificados como significativos en el modelo (horas extra, viaje frecuente, soltería, baja satisfacción, bajo ingreso relativo, poca antigüedad en el cargo), por lo que es razonable esperar una probabilidad de rotación estimada considerablemente por encima del punto de corte de referencia, lo que llevaría a recomendar una intervención preventiva sobre este empleado.

9 Conclusiones y recomendaciones

Factores que aumentan el riesgo de rotación (odds ratio mayor a 1): trabajar horas extra es, con claridad, el factor de mayor magnitud. Le siguen la frecuencia de viajes de negocio y el estado civil de soltería.

Factores que reducen el riesgo de rotación (odds ratio menor a 1): mayor satisfacción laboral, mayor ingreso mensual y mayor antigüedad en el cargo actual.

Con base en estos hallazgos, una posible estrategia de retención debería priorizar:

  1. Gestión de horas extra: redistribuir la carga de trabajo, contratar personal adicional en las áreas con mayor sobrecarga sostenida, o compensar de forma más atractiva el tiempo extra trabajado, en lugar de depender de él de manera estructural.
  2. Racionalización de los viajes de negocio: evaluar si algunos viajes frecuentes pueden sustituirse por reuniones virtuales, o introducir compensaciones específicas (tiempo de descanso adicional, bonificaciones) para los roles con mayor exigencia de viaje.
  3. Atención a empleados solteros y a quienes llevan poco tiempo en su cargo actual, como dos factores de riesgo independientes: cada uno aumenta el riesgo de rotación por separado, según el modelo. No incluí una interacción entre Estado_Civil y Antigüedad_Cargo, así que no puedo afirmar que la combinación de ambos (“soltero y con poca antigüedad”) sea más riesgosa que la suma de sus efectos individuales — eso sería una hipótesis razonable para una iteración futura del modelo, no una conclusión que los datos actuales respalden. Con esa salvedad, planes de desarrollo de carrera más visibles y mentoría temprana son acciones razonables para ambos factores.
  4. Revisión de brechas salariales: dado que el ingreso mensual tiene un efecto protector significativo, conviene identificar si existen cargos o niveles específicos con salarios rezagados frente al mercado.
  5. Seguimiento activo de la satisfacción laboral: dado que es un predictor consistente tanto en este modelo como en la literatura de rotación, una encuesta de clima con periodicidad regular permitiría anticipar caídas en la satisfacción antes de que se traduzcan en renuncias.

Es importante ser cuidadoso con el alcance de estas conclusiones: el modelo identifica asociaciones estadísticas, no relaciones causales. Por ejemplo, no puedo afirmar que reducir las horas extra de un empleado específico reducirá su probabilidad de rotar en la magnitud exacta implicada por el odds ratio; solo puedo afirmar que, en los datos observados, los empleados con más horas extra rotan más, controlando por las demás variables del modelo.

10 Buenas prácticas y limitaciones del análisis

  • Desbalance de clases: con solo 16.1 % de casos positivos, la exactitud del modelo no es una métrica confiable por sí sola; prioricé la curva ROC, el AUC y un análisis explícito de sensibilidad del punto de corte.
  • Ausencia de valores faltantes: verificado explícitamente antes de modelar; no fue necesario aplicar ninguna estrategia de imputación.
  • Multicolinealidad: no se evidenciaron valores de VIF/GVIF preocupantes, ni correlaciones bivariadas altas entre los predictores finalmente incluidos.
  • Separación cuasi-perfecta: no se observan indicios evidentes de separación en el ajuste (todas las categorías tienen un tamaño de celda razonable y los errores estándar del modelo son razonables); una evaluación formal con herramientas como brglm2::detect_separation() podría incorporarse como diagnóstico adicional, pero no se realizó en este informe.
  • Validación: el AUC y las métricas de clasificación se calcularon sobre una única partición entrenamiento/prueba (70/30). Esa estimación depende, en alguna medida, de qué observaciones cayeron en cada conjunto. Una validación cruzada estratificada (por ejemplo, 10 particiones repetidas varias veces) daría una estimación más estable del desempeño esperado del modelo, a costa de mayor tiempo de cómputo; no se implementó aquí por alcance de la actividad.
  • Calibración: evalué la capacidad discriminativa del modelo (su habilidad para ordenar correctamente casos positivos por encima de negativos, vía AUC), pero no evalué formalmente su calibración — es decir, si una probabilidad predicha de 0.3 corresponde efectivamente a una frecuencia observada de rotación cercana al 30 %. Un modelo puede discriminar razonablemente bien y, al mismo tiempo, estar mal calibrado; herramientas como el Brier score o una curva de calibración complementarían este análisis.
  • Linealidad del logit: la regresión logística asume que el logaritmo de los odds varía linealmente con cada predictor cuantitativo. No verifiqué formalmente este supuesto para Ingreso_Mensual ni Antigüedad_Cargo (ambas con distribuciones sesgadas identificadas en el análisis univariado); pruebas como Box-Tidwell o una comparación de AIC contra una especificación con splines (splines::ns()) permitirían confirmar si una forma no lineal ajusta mejor, en particular para Antigüedad_Cargo, donde el efecto podría no ser monótono en todo su rango.
  • Métricas adicionales de clasificación: dado el desbalance de clases, una curva de precisión-recall (PR-AUC) sería un complemento informativo al AUC de la curva ROC, que no se incluyó en este informe.
  • Causalidad: los resultados deben leerse como asociaciones estadísticas controladas por las demás variables del modelo, no como efectos causales; el diseño es observacional y transversal (una sola medición en el tiempo), no experimental.
  • Reproducibilidad: el uso de set.seed(123) garantiza que la partición entrenamiento/prueba sea idéntica en cualquier ejecución posterior de este documento.
  • Alcance de las variables incluidas: trabajar con 6 de las 24 variables disponibles se justificó por una regla de eventos-por-variable y por evitar redundancia entre predictores correlacionados (sección 2), no por una regla general de que “menos variables es mejor”. Esto implica que parte de la variabilidad de la rotación queda sin explicar por este modelo en particular (evidenciado en un AUC que, sin ser bajo, tampoco es excelente); no se comparó formalmente contra un modelo con el conjunto completo de variables.

11 Referencias

  • Peduzzi, P., Concato, J., Kemper, E., Holford, T. R., & Feinstein, A. R. (1996). A simulation study of the number of events per variable in logistic regression analysis. Journal of Clinical Epidemiology.
  • Mobley, W. H. (1977). Intermediate linkages in the relationship between job satisfaction and employee turnover. Journal of Applied Psychology.
  • Greenhaus, J. H., & Beutell, N. J. (1985). Sources of conflict between work and family roles. Academy of Management Review.
  • Griffeth, R. W., Hom, P. W., & Gaertner, S. (2000). A meta-analysis of antecedents and correlates of employee turnover: Update, moderator tests, and research implications for the next millennium. Journal of Management.
  • Hosmer, D. W., Lemeshow, S., & Sturdivant, R. X. (2013). Applied Logistic Regression (3.ª ed.). Wiley.
  • James, G., Witten, D., Hastie, T., & Tibshirani, R. (2021). An Introduction to Statistical Learning (2.ª ed.). Springer.
  • Gujarati, D. N. (2009). Econometría (5.ª ed.). McGraw-Hill.