Rotación de Cargo
Análisis Mediante Regresión Logística Binomial
Leidy Torres Ramírez y Carolina Marín Ortíz
Preparación de los datos
La empresa cuenta con datos históricos de sus trabajadores (edad, antigüedad, ingreso, satisfacción, horas extra, entre otros) y quiere comprender y prever los factores que influyen en la rotación de empleados. La gerencia planea estimar un modelo de regresión logística binomial que permita estimar la probabilidad de que un empleado rote e identificar qué factores inciden en mayor medida en esa probabilidad, para diseñar acciones de retención.
Librerías y funciones auxiliares
# Instalación de paquetes
# install.packages("devtools")
# devtools::install_github("centromagis/paqueteMODELOS", force = TRUE)
library(paqueteMODELOS)
library(dplyr)
library(tidyr)
library(purrr)
library(tibble)
library(ggplot2)
library(scales)
library(patchwork)
library(knitr)
library(kableExtra)
library(broom)
col_azul <- "#018ABE"
col_naranja <- "#FF7F00"
theme_set(theme_minimal(base_size = 11) +
theme(strip.text = element_text(face = "bold"),
plot.title = element_text(face = "bold")))
set.seed(123) # reproducibilidad
# ---- Funciones auxiliares ----
estilo_tabla <- function(x, caption, digits = 3, col.names = NA, ...) {
kable(x, caption = caption, digits = digits, col.names = col.names, align = "c", ...) %>%
kable_styling(bootstrap_options = c("striped", "hover", "condensed"),
full_width = FALSE, position = "center", font_size = 12)
}
# Formato de p-valores y estrellas de significancia
fmt_p <- function(p) {
ifelse(is.na(p), NA_character_,
ifelse(p < 0.001, "<0.001", formatC(p, format = "f", digits = 3)))
}
# Un "*" solo en una celda de la tabla se convierte en viñeta y la celda queda vacía.
signif_stars <- function(p) {
as.character(cut(p, breaks = c(-Inf, 0.001, 0.01, 0.05, 0.1, Inf),
labels = c("\\*\\*\\*", "\\*\\*", "\\*", ".", " ")))
}
f_num <- function(x, d = 2) formatC(x, format = "f", digits = d)
# Etiquetas legibles para los términos del modelo
etq_termino <- function(term) {
mapa <- c(
"(Intercept)" = "Intercepto",
"horas_extraSi" = "Horas extra: Sí (vs. No)",
"estado_civilDivorciado" = "Estado civil: Divorciado (vs. Casado)",
"estado_civilSoltero" = "Estado civil: Soltero (vs. Casado)",
"viaje_negociosRaramente" = "Viaje: Raramente (vs. no viaja)",
"viaje_negociosFrecuentemente" = "Viaje: Frecuentemente (vs. no viaja)",
"ingreso_miles" = "Ingreso mensual (por cada 1.000)",
"antiguedad_cargo" = "Antigüedad en el cargo (por año)",
"distancia_casa" = "Distancia a casa (por unidad)"
)
ifelse(term %in% names(mapa), unname(mapa[term]), term)
}
# Asimetría (coeficiente de Fisher-Pearson) y conteo de atípicos (regla 1.5*RIC)
asimetria <- function(x) {
x <- x[!is.na(x)]
mean(((x - mean(x)) / sd(x))^3)
}
n_atipicos <- function(x) {
q <- quantile(x, c(0.25, 0.75), na.rm = TRUE, names = FALSE)
ric <- q[2] - q[1]
sum(x < q[1] - 1.5 * ric | x > q[2] + 1.5 * ric, na.rm = TRUE)
}Carga y verificación
La base rotacion (paquete paqueteMODELOS)
tiene 1.470 filas y 24 columnas. Los nombres originales contienen
tildes, eñes y espacios (por ejemplo Viaje de Negocios,
Años_Experiencia, Satisfación_Laboral), lo que
suele generar problemas de codificación al escribir código. Por esa
razón se renombran por posición a nombres simples sin tildes; la
Tabla 1 conserva la equivalencia entre el nombre
original y el nombre usado en el análisis.
data("rotacion")
datos_raw <- as_tibble(rotacion)
stopifnot(ncol(datos_raw) == 24, nrow(datos_raw) > 0)
nombres_originales <- names(datos_raw)
nombres_nuevos <- c(
"rotacion", "edad", "viaje_negocios", "departamento", "distancia_casa",
"educacion", "campo_educacion", "satisfaccion_ambiental", "genero", "cargo",
"satisfaccion_laboral", "estado_civil", "ingreso_mensual", "trabajos_anteriores",
"horas_extra", "porcentaje_aumento_salarial", "rendimiento_laboral",
"anos_experiencia", "capacitaciones", "equilibrio_trabajo_vida", "antiguedad",
"antiguedad_cargo", "anos_ultima_promocion", "anos_mismo_jefe"
)
names(datos_raw) <- nombres_nuevos
# Verificaciones de consistencia de los valores esperados
stopifnot(
all(na.omit(datos_raw$rotacion) %in% c("Si", "No")),
all(na.omit(datos_raw$horas_extra) %in% c("Si", "No")),
"Casado" %in% datos_raw$estado_civil
)
# Categoría de referencia de 'viaje_negocios': el nivel que NO es "Raramente" ni "Frecuentemente"
# (se espera "No_Viaja"; se detecta automáticamente para no depender de cómo esté escrito)
nivel_viaje_ref <- setdiff(unique(as.character(datos_raw$viaje_negocios)),
c("Raramente", "Frecuentemente"))
stopifnot(length(nivel_viaje_ref) == 1)Para el análisis univariado se clasifican las 24 variables en tres grupos:
- Cualitativas nominales (8):
rotacion,viaje_negocios,departamento,campo_educacion,genero,cargo,estado_civil,horas_extra. - Cualitativas ordinales (5):
educacion,satisfaccion_ambiental,satisfaccion_laboral,rendimiento_laboral,equilibrio_trabajo_vida. Están codificadas numéricamente, pero son escalas de categorías ordenadas (tipo Likert), por eso se describen con frecuencias y no con media y desviación. - Cuantitativas (11):
edad,distancia_casa,ingreso_mensual,trabajos_anteriores,porcentaje_aumento_salarial,anos_experiencia,capacitaciones,antiguedad,antiguedad_cargo,anos_ultima_promocion,anos_mismo_jefe.
vars_nominales <- c("rotacion", "viaje_negocios", "departamento", "campo_educacion",
"genero", "cargo", "estado_civil", "horas_extra")
vars_ordinales <- c("educacion", "satisfaccion_ambiental", "satisfaccion_laboral",
"rendimiento_laboral", "equilibrio_trabajo_vida")
vars_cuant <- c("edad", "distancia_casa", "ingreso_mensual", "trabajos_anteriores",
"porcentaje_aumento_salarial", "anos_experiencia", "capacitaciones",
"antiguedad", "antiguedad_cargo", "anos_ultima_promocion",
"anos_mismo_jefe")
stopifnot(setequal(c(vars_nominales, vars_ordinales, vars_cuant), nombres_nuevos))
# Etiquetas legibles para tablas y gráficos
etiquetas <- c(
rotacion = "Rotación", edad = "Edad", viaje_negocios = "Viaje de negocios",
departamento = "Departamento", distancia_casa = "Distancia a casa",
educacion = "Educación (nivel)", campo_educacion = "Campo de educación",
satisfaccion_ambiental = "Satisfacción ambiental", genero = "Género",
cargo = "Cargo", satisfaccion_laboral = "Satisfacción laboral",
estado_civil = "Estado civil", ingreso_mensual = "Ingreso mensual",
trabajos_anteriores = "Trabajos anteriores", horas_extra = "Horas extra",
porcentaje_aumento_salarial = "% aumento salarial",
rendimiento_laboral = "Rendimiento laboral", anos_experiencia = "Años de experiencia",
capacitaciones = "Capacitaciones", equilibrio_trabajo_vida = "Equilibrio trabajo-vida",
antiguedad = "Antigüedad (empresa)", antiguedad_cargo = "Antigüedad en el cargo",
anos_ultima_promocion = "Años desde última promoción",
anos_mismo_jefe = "Años con el mismo jefe",
ingreso_miles = "Ingreso mensual (miles)"
)Se construye la base de trabajo datos con:
- la variable respuesta binaria
y(1 = rota, 0 = no rota) - las variables de texto convertidas a factor, fijando las categorías de referencia (Horas extra: “No”; Estado civil: “Casado”, la categoría más frecuente; Viaje de negocios: la categoría de no viaje)
ingreso_miles = ingreso_mensual / 1000, para que el coeficiente del ingreso se interprete por cada 1.000 unidades y no por una unidad (cuyo efecto sería numéricamente minúsculo).
datos <- datos_raw %>%
mutate(
across(where(is.character), as.factor),
rotacion = factor(rotacion, levels = c("No", "Si")),
y = as.integer(rotacion == "Si"),
horas_extra = factor(horas_extra, levels = c("No", "Si")),
estado_civil = relevel(estado_civil, ref = "Casado"),
viaje_negocios = factor(viaje_negocios,
levels = c(nivel_viaje_ref, "Raramente", "Frecuentemente")),
ingreso_miles = ingreso_mensual / 1000
)
n_duplicadas <- sum(duplicated(datos_raw))diccionario <- tibble(
original = nombres_originales,
nuevo = nombres_nuevos,
clasificacion = case_when(
nombres_nuevos %in% vars_nominales ~ "Cualitativa nominal",
nombres_nuevos %in% vars_ordinales ~ "Cualitativa ordinal (escala codificada)",
nombres_nuevos %in% vars_cuant ~ "Cuantitativa"
),
distintos = map_int(nombres_nuevos, ~ n_distinct(datos_raw[[.x]])),
faltantes = map_int(nombres_nuevos, ~ sum(is.na(datos_raw[[.x]])))
)
estilo_tabla(diccionario,
caption = "Tabla 1. Diccionario de variables: nombre original, nombre usado en el análisis, clasificación, valores distintos y valores faltantes.",
col.names = c("Nombre original", "Nombre en el análisis", "Clasificación",
"Valores distintos", "Valores faltantes"))| Nombre original | Nombre en el análisis | Clasificación | Valores distintos | Valores faltantes |
|---|---|---|---|---|
| Rotación | rotacion | Cualitativa nominal | 2 | 0 |
| Edad | edad | Cuantitativa | 43 | 0 |
| Viaje de Negocios | viaje_negocios | Cualitativa nominal | 3 | 0 |
| Departamento | departamento | Cualitativa nominal | 3 | 0 |
| Distancia_Casa | distancia_casa | Cuantitativa | 29 | 0 |
| Educación | educacion | Cualitativa ordinal (escala codificada) | 5 | 0 |
| Campo_Educación | campo_educacion | Cualitativa nominal | 6 | 0 |
| Satisfacción_Ambiental | satisfaccion_ambiental | Cualitativa ordinal (escala codificada) | 4 | 0 |
| Genero | genero | Cualitativa nominal | 2 | 0 |
| Cargo | cargo | Cualitativa nominal | 9 | 0 |
| Satisfación_Laboral | satisfaccion_laboral | Cualitativa ordinal (escala codificada) | 4 | 0 |
| Estado_Civil | estado_civil | Cualitativa nominal | 3 | 0 |
| Ingreso_Mensual | ingreso_mensual | Cuantitativa | 1349 | 0 |
| Trabajos_Anteriores | trabajos_anteriores | Cuantitativa | 10 | 0 |
| Horas_Extra | horas_extra | Cualitativa nominal | 2 | 0 |
| Porcentaje_aumento_salarial | porcentaje_aumento_salarial | Cuantitativa | 15 | 0 |
| Rendimiento_Laboral | rendimiento_laboral | Cualitativa ordinal (escala codificada) | 2 | 0 |
| Años_Experiencia | anos_experiencia | Cuantitativa | 40 | 0 |
| Capacitaciones | capacitaciones | Cuantitativa | 7 | 0 |
| Equilibrio_Trabajo_Vida | equilibrio_trabajo_vida | Cualitativa ordinal (escala codificada) | 4 | 0 |
| Antigüedad | antiguedad | Cuantitativa | 37 | 0 |
| Antigüedad_Cargo | antiguedad_cargo | Cuantitativa | 19 | 0 |
| Años_ultima_promoción | anos_ultima_promocion | Cuantitativa | 16 | 0 |
| Años_acargo_con_mismo_jefe | anos_mismo_jefe | Cuantitativa | 18 | 0 |
# Salida 1. Verificación de la codificación de la variable respuesta
with(datos, table(rotacion, y))## y
## rotacion 0 1
## No 1233 0
## Si 0 237Salida 1. Verificación de la codificación: todas las observaciones “Si” deben quedar con \(y = 1\) y todas las “No” con \(y = 0\).
La base tiene 1.470 registros y 24 variables, sin valores faltantes
en ninguna de ellas y con 0 filas duplicadas, por lo que no fue
necesario imputar ni eliminar observaciones. La clasificación de la
Tabla 1 es adecuada: ingreso_mensual (1.349 valores
distintos) y edad (43) se comportan como cuantitativas
continuas o casi continuas; trabajos_anteriores (10
valores) y capacitaciones (7 valores) son conteos discretos
que se tratan como cuantitativas; y las cinco escalas ordinales tienen
entre 2 y 5 niveles. Cabe destacar que rendimiento_laboral
solo toma 2 valores distintos (3 y 4), por lo que en la práctica es una
variable casi dicotómica. La Salida 1 confirma que la codificación es
correcta: los 1.233 empleados que no rotan quedan con \(y = 0\), los 237 que rotan quedan con \(y = 1\) y no hay cruces entre
categorías.
1 Selección de variables
1.1 Criterio de selección
Se seleccionan 3 variables categóricas (distintas de rotación) y 3 variables cuantitativas con los siguientes criterios:
- Sustento conceptual: que exista un mecanismo plausible que vincule la variable con la decisión de rotar.
- Interpretabilidad para la gerencia: que sean factores sobre los cuales la empresa pueda actuar (carga de trabajo, viajes, compensación) o que permitan segmentar (perfil personal).
- Baja redundancia entre sí: el principal supuesto del modelo logit es la ausencia (o baja) multicolinealidad; por ello se evita incluir varias variables que midan lo mismo (por ejemplo, edad, años de experiencia y antigüedad en la empresa son todas indicadores de “etapa de carrera” y suelen estar muy correlacionadas con el ingreso).
- Cobertura de dimensiones distintas: carga laboral, perfil personal y compensación/estabilidad.
Las variables elegidas y la dirección esperada del efecto se resumen en la Tabla 2.
seleccion <- tribble(
~H, ~Variable, ~Tipo, ~Contraste, ~Signo,
"H1", "Horas extra", "Categórica (binaria)", "Sí vs. No (referencia: No)", "Positivo (+)",
"H2", "Estado civil", "Categórica (3 niveles)", "Soltero vs. Casado (referencia: Casado)", "Positivo (+)",
"H3", "Viaje de negocios", "Categórica (3 niveles)", "Frecuentemente / Raramente vs. No viaja", "Positivo (+)",
"H4", "Ingreso mensual", "Cuantitativa", "Por cada 1.000 unidades adicionales", "Negativo (-)",
"H5", "Antigüedad en el cargo", "Cuantitativa", "Por cada año adicional", "Negativo (-)",
"H6", "Distancia a casa", "Cuantitativa", "Por cada unidad adicional de distancia", "Positivo (+)"
)
estilo_tabla(seleccion,
caption = "Tabla 2. Variables seleccionadas, tipo, contraste de interés y signo esperado del coeficiente en el modelo logit (+ aumenta la probabilidad de rotar; - la disminuye).",
col.names = c("Hipótesis", "Variable", "Tipo", "Contraste", "Signo esperado del coeficiente"))| Hipótesis | Variable | Tipo | Contraste | Signo esperado del coeficiente |
|---|---|---|---|---|
| H1 | Horas extra | Categórica (binaria) | Sí vs. No (referencia: No) | Positivo (+) |
| H2 | Estado civil | Categórica (3 niveles) | Soltero vs. Casado (referencia: Casado) | Positivo (+) |
| H3 | Viaje de negocios | Categórica (3 niveles) | Frecuentemente / Raramente vs. No viaja | Positivo (+) |
| H4 | Ingreso mensual | Cuantitativa | Por cada 1.000 unidades adicionales | Negativo (-) |
| H5 | Antigüedad en el cargo | Cuantitativa | Por cada año adicional | Negativo (-) |
| H6 | Distancia a casa | Cuantitativa | Por cada unidad adicional de distancia | Positivo (+) |
1.2 Justificación e hipótesis por variable
Variables categóricas
Hipótesis 1 - Horas extra. Trabajar horas extra incrementa la carga laboral, reduce el tiempo de descanso y de vida personal y familiar, y genera desgaste. Ese desgaste aumenta la intención de abandonar el cargo. Hipótesis: los empleados que trabajan horas extra tienen mayor probabilidad (mayores odds) de rotar que quienes no las trabajan; se espera \(\beta_{\text{Sí vs. No}} > 0\) (OR > 1).
Hipótesis 2 - Estado civil. Los empleados casados suelen tener compromisos económicos y familiares que aumentan el valor de la estabilidad laboral y el costo de cambiar de cargo. Los solteros tienen, en general, mayor movilidad y menos restricciones para cambiar. Hipótesis: los solteros tienen mayor probabilidad de rotar que los casados; se espera \(\beta_{\text{Soltero vs. Casado}} > 0\). Para los divorciados no se plantea una dirección a priori, su comparación se reporta de forma exploratoria.
Hipótesis 3 - Viaje de negocios. Viajar con frecuencia implica ausencias del hogar, fatiga y desorganización de las rutinas personales, lo que deteriora el equilibrio entre trabajo y vida. Hipótesis: la probabilidad de rotar aumenta con la frecuencia de viaje (Frecuentemente > Raramente > No viaja); se espera \(\beta > 0\) para ambas categorías frente a la de referencia, y mayor para “Frecuentemente”.
Variables cuantitativas
Hipótesis 4 - Ingreso mensual. Un mayor ingreso eleva la satisfacción con la compensación, aumenta el costo de oportunidad de irse y reduce el atractivo de ofertas externas. Hipótesis: a mayor ingreso, menor probabilidad de rotar; se espera \(\beta < 0\) (OR < 1) por cada 1.000 unidades adicionales.
Hipótesis 5 - Antigüedad en el cargo. Con más años en el cargo el empleado acumula conocimiento específico, relaciones y arraigo; en cambio, los primeros años concentran los desajustes entre la persona y el cargo y las salidas tempranas. Hipótesis: a mayor antigüedad en el cargo, menor probabilidad de rotar; se espera \(\beta < 0\). (Un posible efecto no lineal, por ejemplo estancamiento en antigüedades muy altas, se explora gráficamente en la Sección 3.2.)
Hipótesis 6 - Distancia a casa. Una mayor distancia implica más tiempo y costo de desplazamiento, más fatiga y menos tiempo personal. Hipótesis: a mayor distancia, mayor probabilidad de rotar; se espera \(\beta > 0\).
Para cada coeficiente \(\beta_j\) del modelo se contrasta \(H_0: \beta_j = 0\) (la variable no modifica los odds de rotar) frente a \(H_1: \beta_j \neq 0\), mediante el estadístico de Wald. La dirección esperada (signo de la Tabla 2) se compara luego con el signo estimado. Todas las relaciones planteadas son de asociación, no de causalidad.
1.3 Verificación de redundancia entre variables cuantitativas
Para respaldar el criterio 3, la Tabla 3 muestra las correlaciones de Spearman (apropiadas para variables asimétricas) entre las tres variables cuantitativas elegidas y otras candidatas que miden “etapa de carrera” (edad, años de experiencia y antigüedad en la empresa).
cand <- c("ingreso_mensual", "antiguedad_cargo", "distancia_casa",
"edad", "anos_experiencia", "antiguedad")
cor_cand <- cor(datos[, cand], method = "spearman", use = "pairwise.complete.obs")
dimnames(cor_cand) <- list(unname(etiquetas[cand]), unname(etiquetas[cand]))
estilo_tabla(round(cor_cand, 2), digits = 2,
caption = "Tabla 3. Matriz de correlaciones de Spearman entre variables cuantitativas candidatas (las tres primeras son las seleccionadas).")| Ingreso mensual | Antigüedad en el cargo | Distancia a casa | Edad | Años de experiencia | Antigüedad (empresa) | |
|---|---|---|---|---|---|---|
| Ingreso mensual | 1.00 | 0.39 | 0.00 | 0.47 | 0.71 | 0.46 |
| Antigüedad en el cargo | 0.39 | 1.00 | 0.01 | 0.20 | 0.49 | 0.85 |
| Distancia a casa | 0.00 | 0.01 | 1.00 | -0.02 | 0.00 | 0.01 |
| Edad | 0.47 | 0.20 | -0.02 | 1.00 | 0.66 | 0.25 |
| Años de experiencia | 0.71 | 0.49 | 0.00 | 0.66 | 1.00 | 0.59 |
| Antigüedad (empresa) | 0.46 | 0.85 | 0.01 | 0.25 | 0.59 | 1.00 |
Entre las tres variables seleccionadas las correlaciones de Spearman son bajas o moderadas: Ingreso mensual y Antigüedad en el cargo tienen ρ = 0.39, y la Distancia a casa es prácticamente independiente de las otras dos (ρ entre 0.00 y 0.01). Todas quedan por debajo del umbral de referencia de 0.7, por lo que no se espera multicolinealidad relevante. En cambio, las candidatas descartadas están más entrelazadas: los Años de experiencia se correlacionan con el ingreso (ρ = 0.71, por encima del umbral) y con la edad (ρ = 0.66), y la Antigüedad en la empresa se correlaciona con la Antigüedad en el cargo (ρ = 0.85), es decir, miden casi lo mismo. La Edad tiene una correlación moderada con el ingreso (ρ = 0.47) y, en rigor, podría haberse incluido sin generar colinealidad severa; se dejó fuera porque se solapa conceptualmente con el ingreso y la antigüedad (todas reflejan la “etapa de carrera”). Queda como candidata para una segunda versión del modelo.
2 Análisis univariado
Se caracteriza toda la información de la base rotacion,
usando indicadores y gráficos según el tipo de variable: frecuencias y
gráficos de barras para las cualitativas; medidas de tendencia central,
dispersión y forma, con histogramas y diagramas de caja, para las
cuantitativas.
2.1 Variable respuesta: rotación
n_total <- nrow(datos)
n_rot <- sum(datos$y)
p_rot <- mean(datos$y)
ic_rot <- binom.test(n_rot, n_total)$conf.int
odds_rot <- p_rot / (1 - p_rot)
tab_rot <- datos %>%
count(rotacion = as.character(rotacion), name = "n") %>%
mutate(pct = n / sum(n) * 100) %>%
bind_rows(tibble(rotacion = "Total", n = n_total, pct = 100))
estilo_tabla(tab_rot, digits = 1,
caption = "Tabla 4. Distribución de la variable rotación.",
col.names = c("Rotación", "Frecuencia", "Porcentaje (%)"))| Rotación | Frecuencia | Porcentaje (%) |
|---|---|---|
| No | 1233 | 83.9 |
| Si | 237 | 16.1 |
| Total | 1470 | 100.0 |
datos %>%
count(rotacion) %>%
mutate(pct = n / sum(n)) %>%
ggplot(aes(x = rotacion, y = pct, fill = rotacion)) +
geom_col(width = 0.6) +
geom_text(aes(label = paste0(percent(pct, accuracy = 0.1), " (n = ", n, ")")),
vjust = -0.4, size = 3.5) +
scale_fill_manual(values = c(No = col_azul, Si = col_naranja)) +
scale_y_continuous(labels = percent, expand = expansion(mult = c(0, 0.15))) +
labs(x = "Rotación", y = "Porcentaje") +
theme(legend.position = "none")
Figura 1. Distribución porcentual de la variable rotación.
De los 1470 empleados de la base, 237 rotaron, lo que equivale a una tasa de rotación de 16.1% (IC 95% exacto: 14.3% a 18.1%). Los odds de rotar son \(p/(1-p) =\) 0.192; es decir, por cada empleado que rota hay aproximadamente 5.2 que no rotan.
La variable respuesta está desbalanceada (la clase “rota” es minoritaria), lo que tiene tres implicaciones para el resto del análisis:
- La exactitud (accuracy) por sí sola es engañosa, un modelo que predijera “No rota” para todos los empleados alcanzaría una exactitud de 83.9% sin identificar a ningún empleado en riesgo. Por eso la evaluación se hará con la curva ROC, el AUC, la sensibilidad y la especificidad.
- El punto de corte de 0.5 para clasificar probablemente no sea adecuado, porque la probabilidad promedio de rotar es cercana a 16%; el corte se definirá en el punto 6.
- Aun así, hay 237 eventos, un tamaño suficiente para estimar un modelo con varias covariables (se verifica la regla de eventos por variable en la Sección 4.7).
2.2 Variables cuantitativas
desc_cuant <- datos %>%
select(all_of(vars_cuant)) %>%
pivot_longer(everything(), names_to = "variable", values_to = "valor") %>%
group_by(variable) %>%
summarise(
n = sum(!is.na(valor)),
Media = mean(valor, na.rm = TRUE),
DE = sd(valor, na.rm = TRUE),
CV = DE / Media * 100,
Min = min(valor, na.rm = TRUE),
Q1 = quantile(valor, 0.25, na.rm = TRUE, names = FALSE),
Mediana = median(valor, na.rm = TRUE),
Q3 = quantile(valor, 0.75, na.rm = TRUE, names = FALSE),
Max = max(valor, na.rm = TRUE),
Asimetria = asimetria(valor),
Atipicos = n_atipicos(valor),
.groups = "drop"
) %>%
arrange(match(variable, vars_cuant)) %>%
mutate(variable = unname(etiquetas[variable]))
estilo_tabla(desc_cuant, digits = 2,
caption = "Tabla 5. Estadísticos descriptivos de las variables cuantitativas (DE: desviación estándar; CV: coeficiente de variación en %; Atípicos: observaciones fuera de Q1 - 1.5·RIC o Q3 + 1.5·RIC).",
col.names = c("Variable", "n", "Media", "DE", "CV (%)", "Mín.", "Q1",
"Mediana", "Q3", "Máx.", "Asimetría", "Atípicos"))| Variable | n | Media | DE | CV (%) | Mín. | Q1 | Mediana | Q3 | Máx. | Asimetría | Atípicos |
|---|---|---|---|---|---|---|---|---|---|---|---|
| Edad | 1470 | 36.92 | 9.14 | 24.74 | 18 | 30 | 36 | 43 | 60 | 0.41 | 0 |
| Distancia a casa | 1470 | 9.19 | 8.11 | 88.19 | 1 | 2 | 7 | 14 | 29 | 0.96 | 0 |
| Ingreso mensual | 1470 | 6502.93 | 4707.96 | 72.40 | 1009 | 2911 | 4919 | 8379 | 19999 | 1.37 | 114 |
| Trabajos anteriores | 1470 | 2.69 | 2.50 | 92.75 | 0 | 1 | 2 | 4 | 9 | 1.02 | 52 |
| % aumento salarial | 1470 | 15.21 | 3.66 | 24.06 | 11 | 12 | 14 | 18 | 25 | 0.82 | 0 |
| Años de experiencia | 1470 | 11.28 | 7.78 | 68.98 | 0 | 6 | 10 | 15 | 40 | 1.11 | 63 |
| Capacitaciones | 1470 | 2.80 | 1.29 | 46.06 | 0 | 2 | 3 | 3 | 6 | 0.55 | 238 |
| Antigüedad (empresa) | 1470 | 7.01 | 6.13 | 87.42 | 0 | 3 | 5 | 9 | 40 | 1.76 | 104 |
| Antigüedad en el cargo | 1470 | 4.23 | 3.62 | 85.67 | 0 | 2 | 3 | 7 | 18 | 0.92 | 21 |
| Años desde última promoción | 1470 | 2.19 | 3.22 | 147.29 | 0 | 0 | 1 | 3 | 15 | 1.98 | 107 |
| Años con el mismo jefe | 1470 | 4.12 | 3.57 | 86.54 | 0 | 2 | 3 | 7 | 17 | 0.83 | 14 |
datos %>%
select(all_of(vars_cuant)) %>%
pivot_longer(everything(), names_to = "variable", values_to = "valor") %>%
mutate(variable = factor(variable, levels = vars_cuant,
labels = unname(etiquetas[vars_cuant]))) %>%
ggplot(aes(x = valor)) +
geom_histogram(aes(y = after_stat(density)), bins = 25,
fill = col_azul, color = "white", alpha = 0.85) +
geom_density(color = col_naranja, linewidth = 0.8) +
facet_wrap(~ variable, scales = "free", ncol = 3) +
labs(x = NULL, y = "Densidad")
Figura 2. Histogramas (con curva de densidad) de las variables cuantitativas.
datos %>%
select(all_of(vars_cuant)) %>%
pivot_longer(everything(), names_to = "variable", values_to = "valor") %>%
mutate(variable = factor(variable, levels = vars_cuant,
labels = unname(etiquetas[vars_cuant]))) %>%
ggplot(aes(x = "", y = valor)) +
geom_boxplot(fill = col_azul, alpha = 0.6, outlier.alpha = 0.5) +
facet_wrap(~ variable, scales = "free_y", ncol = 4) +
labs(x = NULL, y = NULL) +
theme(axis.text.x = element_blank())
Figura 3. Diagramas de caja de las variables cuantitativas.
A partir de los resultados de la Tabla 5, Figuras 2 y 3, se puede interpretar lo siguiente:
- Edad. Distribución aproximadamente simétrica (media 36.9 y mediana 36 años; asimetría 0.41), con rango de 18 a 60 años, CV de 24.7% y sin atípicos.
- Distancia a casa. La media (9.19) supera a la mediana (7), el CV es de 88% y la asimetría de 0.96: la mayoría de los empleados vive cerca (moda en 1 y 2 unidades) y una minoría llega hasta 29. No hay atípicos según la regla de 1.5·RIC.
- Ingreso mensual. Sesgo positivo marcado (asimetría 1.37): la media (6.503) supera con holgura a la mediana (4.919), el CV es de 72% y hay 114 atípicos en la cola superior (máximo 19.999). Son valores extremos plausibles (ingresos altos de una minoría) y no errores evidentes, pero pueden influir en el ajuste, por lo que se revisan con la distancia de Cook.
- Antigüedad en el cargo. Mediana de 3 años (media 4.23; CV 86%; asimetría 0.92) y 21 atípicos (antigüedades de 15 a 18 años). El histograma es multimodal, con concentraciones cerca de 0 a 2 años y alrededor de 7 años, un patrón que podría mostrar cohortes de ingreso o movimientos internos en momentos específicos (no es posible confirmarlo con estos datos).
- Trabajos anteriores. Mediana de 2 (media 2.69), CV de 93%, asimetría 1.02 y 52 atípicos (empleados con 9 trabajos previos).
- % de aumento salarial. Toma valores entre 11 y 25 (media 15.2; mediana 14), asimetría 0.82 y sin atípicos; se concentra en los porcentajes más bajos (11 a 14).
- Años de experiencia, antigüedad en la empresa, años desde la última promoción y años con el mismo jefe. Son las variables más asimétricas: asimetrías de 1.11, 1.76, 1.98 y 0.83, con 63, 104, 107 y 14 atípicos, respectivamente. Los Años desde la última promoción tienen la mayor dispersión relativa (CV de 147%): la mediana es de 1 año, pero hay empleados con hasta 15 años sin promoción.
- Capacitaciones. Variable de conteo de 0 a 6 con mediana 3 y rango intercuartílico muy estrecho (2 a 3), lo que hace que la regla de 1.5·RIC marque 238 observaciones como atípicas. En este caso el criterio es poco informativo: se trata de valores discretos legítimos (0, 1, 5 y 6 capacitaciones), no de errores.
En conjunto, casi todas las variables de ingreso y de tiempo tienen media mayor que la mediana y asimetría positiva (se considera marcada cuando |g₁| > 1), por lo que la mediana describe mejor al empleado típico. El modelo logit no exige normalidad de las covariables, así que esta asimetría no invalida el análisis; sí conviene vigilar el efecto de los valores extremos y, en el caso del ingreso, la forma de la relación con la rotación. Los mínimos y máximos son coherentes con la edad de la plantilla (por ejemplo, un máximo de 40 años de experiencia frente a una edad máxima de 60).
2.3 Variables cualitativas nominales
vars_nom_x <- setdiff(vars_nominales, "rotacion") # rotación ya se describió en la Sección 2.1
tab_nom <- map_dfr(vars_nom_x, function(v) {
datos %>%
count(Categoria = as.character(.data[[v]]), name = "Frecuencia") %>%
mutate(Variable = etiquetas[[v]],
Porcentaje = Frecuencia / sum(Frecuencia) * 100,
.before = 1) %>%
select(Variable, Categoria, Frecuencia, Porcentaje) %>%
arrange(desc(Frecuencia))
})
estilo_tabla(tab_nom, digits = 1,
caption = "Tabla 6. Distribución de frecuencias de las variables cualitativas nominales (excepto rotación).",
col.names = c("Variable", "Categoría", "Frecuencia", "Porcentaje (%)")) %>%
collapse_rows(columns = 1, valign = "top")| Variable | Categoría | Frecuencia | Porcentaje (%) |
|---|---|---|---|
| Viaje de negocios | Raramente | 1043 | 71.0 |
| Frecuentemente | 277 | 18.8 | |
| No_Viaja | 150 | 10.2 | |
| Departamento | IyD | 961 | 65.4 |
| Ventas | 446 | 30.3 | |
| RH | 63 | 4.3 | |
| Campo de educación | Ciencias | 606 | 41.2 |
| Salud | 464 | 31.6 | |
| Mercadeo | 159 | 10.8 | |
| Tecnicos | 132 | 9.0 | |
| Otra | 82 | 5.6 | |
| Humanidades | 27 | 1.8 | |
| Género | M | 882 | 60.0 |
| F | 588 | 40.0 | |
| Cargo | Ejecutivo_Ventas | 326 | 22.2 |
| Investigador_Cientifico | 292 | 19.9 | |
| Tecnico_Laboratorio | 259 | 17.6 | |
| Director_Manofactura | 145 | 9.9 | |
| Representante_Salud | 131 | 8.9 | |
| Gerente | 102 | 6.9 | |
| Representante_Ventas | 83 | 5.6 | |
| Director_Investigación | 80 | 5.4 | |
| Recursos_Humanos | 52 | 3.5 | |
| Estado civil | Casado | 673 | 45.8 |
| Soltero | 470 | 32.0 | |
| Divorciado | 327 | 22.2 | |
| Horas extra | No | 1054 | 71.7 |
| Si | 416 | 28.3 |
graficos_frecuencias <- function(df, var, ordinal = FALSE, color = col_azul) {
tab <- df %>%
count(categoria = as.character(.data[[var]])) %>%
mutate(pct = n / sum(n))
if (ordinal) {
tab <- tab %>%
mutate(categoria = factor(categoria, levels = as.character(sort(as.numeric(categoria)))))
} else {
tab <- tab %>% mutate(categoria = reorder(categoria, n))
}
g <- ggplot(tab, aes(x = categoria, y = pct)) +
geom_col(fill = color, width = 0.7) +
scale_y_continuous(labels = percent, expand = expansion(mult = c(0, 0.2))) +
labs(title = etiquetas[[var]], x = NULL, y = "Porcentaje")
if (ordinal) {
g + geom_text(aes(label = percent(pct, accuracy = 0.1)), vjust = -0.4, size = 3)
} else {
g + geom_text(aes(label = percent(pct, accuracy = 0.1)), hjust = -0.1, size = 3) +
coord_flip()
}
}
wrap_plots(map(vars_nom_x, ~ graficos_frecuencias(datos, .x)), ncol = 2)
Figura 4. Distribución porcentual de las variables cualitativas nominales.
A partir de los resultados de la Tabla 6 y Figura 4 se puede interpretar:
- Viaje de negocios: la gran mayoría viaja raramente (71.0%; 1.043 empleados), un 18.8% viaja con frecuencia (277) y solo un 10.2% no viaja (150).
- Departamento: predomina Investigación y Desarrollo (IyD, 65.4%), seguido de Ventas (30.3%); Recursos Humanos (RH) es muy pequeño (4.3%; 63 empleados).
- Campo de educación: Ciencias (41.2%) y Salud (31.6%) reúnen cerca de tres cuartas partes de la plantilla (72.8%); Otra (5.6%) y Humanidades (1.8%; 27 empleados) son categorías pequeñas.
- Género: 60.0% hombres (882) y 40.0% mujeres (588).
- Cargo: hay nueve cargos, los tres más frecuentes (Ejecutivo de Ventas, 22.2%; Investigador Científico, 19.9%; Técnico de Laboratorio, 17.6%) suman el 59.7%. Los menos frecuentes son Recursos Humanos (3.5%; 52 empleados), Director de Investigación (5.4%) y Representante de Ventas (5.6%).
- Estado civil: Casado 45.8% (673), Soltero 32.0% (470) y Divorciado 22.2% (327).
- Horas extra: el 28.3% (416 empleados) trabaja horas extra y el 71.7% (1.054) no.
Hay categorías con pocas observaciones (Humanidades, RH, Recursos Humanos) que darían estimaciones inestables si se incluyeran por separado en un modelo. Entre las variables seleccionadas, la categoría con menos casos es “No viaja” (150 empleados), que además es la categoría de referencia de Viaje de negocios y condiciona la precisión de sus coeficientes.
2.4 Variables cualitativas ordinales
tab_ord <- map_dfr(vars_ordinales, function(v) {
datos %>%
count(Nivel = as.character(.data[[v]]), name = "Frecuencia") %>%
mutate(Variable = etiquetas[[v]],
Porcentaje = Frecuencia / sum(Frecuencia) * 100,
.before = 1) %>%
select(Variable, Nivel, Frecuencia, Porcentaje) %>%
arrange(as.numeric(Nivel))
})
estilo_tabla(tab_ord, digits = 1,
caption = "Tabla 7. Distribución de frecuencias de las variables cualitativas ordinales (niveles según la codificación original de la base).",
col.names = c("Variable", "Nivel", "Frecuencia", "Porcentaje (%)")) %>%
collapse_rows(columns = 1, valign = "top")| Variable | Nivel | Frecuencia | Porcentaje (%) |
|---|---|---|---|
| Educación (nivel) | 1 | 170 | 11.6 |
| 2 | 282 | 19.2 | |
| 3 | 572 | 38.9 | |
| 4 | 398 | 27.1 | |
| 5 | 48 | 3.3 | |
| Satisfacción ambiental | 1 | 284 | 19.3 |
| 2 | 287 | 19.5 | |
| 3 | 453 | 30.8 | |
| 4 | 446 | 30.3 | |
| Satisfacción laboral | 1 | 289 | 19.7 |
| 2 | 280 | 19.0 | |
| 3 | 442 | 30.1 | |
| 4 | 459 | 31.2 | |
| Rendimiento laboral | 3 | 1244 | 84.6 |
| 4 | 226 | 15.4 | |
| Equilibrio trabajo-vida | 1 | 80 | 5.4 |
| 2 | 344 | 23.4 | |
| 3 | 893 | 60.7 | |
| 4 | 153 | 10.4 |
wrap_plots(map(vars_ordinales, ~ graficos_frecuencias(datos, .x, ordinal = TRUE, color = col_naranja)),
ncol = 3)
Figura 5. Distribución porcentual de las variables cualitativas ordinales.
Interpretación de la Tabla 7 y la Figura 5:
- Educación: el nivel 3 es el más frecuente (38.9%), seguido del 4 (27.1%) y el 2 (19.2%); el nivel 5 casi no aparece (3.3%; 48 empleados).
- Satisfacción ambiental: el 61.1% se ubica en los niveles 3 y 4, y el 38.8% en los niveles 1 y 2 (distribución casi uniforme dentro de cada par).
- Satisfacción laboral: patrón muy similar: 61.3% en los niveles 3 y 4 y 38.7% en los niveles 1 y 2. Es decir, cerca de 4 de cada 10 empleados declara una satisfacción laboral baja o intermedia-baja, un grupo que merece atención por su posible relación con la rotación.
- Rendimiento laboral: solo aparecen los niveles 3 (84.6%) y 4 (15.4%); casi no hay variabilidad, por lo que es poco probable que discrimine entre quienes rotan y quienes no (se confirma en el tamizaje de la Sección 3.5).
- Equilibrio trabajo-vida: el nivel 3 concentra el 60.7% (893 empleados); el 28.8% se ubica en los niveles 1 y 2 (el nivel 1 solo tiene 80 empleados, 5.4%) y el 10.4% en el nivel 4.
2.5 Síntesis del análisis univariado
El empleado típico de la base es hombre (60.0%), de unos 37 años, casado (45.8%), que trabaja en Investigación y Desarrollo (65.4%), viaja raramente (71.0%) y no hace horas extra (71.7%); tiene un ingreso mensual mediano de 4.919, una mediana de 3 años en el cargo y vive a una distancia mediana de 7 unidades. Las variables de ingreso y de tiempo (antigüedades, experiencia, años desde la última promoción) son asimétricas a la derecha y presentan valores extremos plausibles, mientras que la edad es aproximadamente simétrica. Algunas categorías tienen muy pocos casos (Humanidades, RH, Recursos Humanos) y Rendimiento laboral casi no varía (solo niveles 3 y 4). Finalmente, la variable respuesta está desbalanceada (16.1% de rotación), lo que obliga a evaluar el modelo con la curva ROC y el AUC, y a definir con cuidado el punto de corte.
3 Análisis bivariado
El objetivo es analizar la relación de la rotación (\(y = 1\) si rota, \(y = 0\) si no rota) con cada una de las 6 variables seleccionadas y compararla con las hipótesis de la Sección 1. El análisis tiene tres capas complementarias:
- Exploratoria descriptiva: tasas de rotación por categoría y comparación de las distribuciones de las variables cuantitativas entre quienes rotan y quienes no.
- Contrastes clásicos: chi-cuadrado de independencia (categóricas) y pruebas t de Welch y Wilcoxon (cuantitativas).
- Modelos logit simples \(\ln\left(\frac{p}{1-p}\right) = \beta_0 + \beta_1 x\), uno por variable: el signo de \(\hat\beta_1\) indica la dirección de la asociación y \(\exp(\hat\beta_1)\) su magnitud en términos de odds.
Al final se aplica el mismo modelo simple al resto de variables de la base (tamizaje), para identificar todas las posibles variables asociadas con la rotación.
vars_cat_sel <- c("horas_extra", "estado_civil", "viaje_negocios")
vars_cuant_sel <- c("ingreso_miles", "antiguedad_cargo", "distancia_casa")
vars_sel <- c(vars_cat_sel, vars_cuant_sel)3.1 Variables categóricas frente a la rotación
tab_cont <- map_dfr(vars_cat_sel, function(v) {
datos %>%
group_by(Categoria = .data[[v]]) %>%
summarise(N = n(), Rotan = sum(y), No_rotan = N - Rotan,
Tasa = Rotan / N * 100, Odds = Rotan / No_rotan,
.groups = "drop") %>%
mutate(Categoria = as.character(Categoria), Variable = etiquetas[[v]], .before = 1)
})
estilo_tabla(tab_cont, digits = 2,
caption = "Tabla 8. Tasa de rotación y odds de rotar según las categorías de las variables categóricas seleccionadas.",
col.names = c("Variable", "Categoría", "N", "Rotan", "No rotan",
"Tasa de rotación (%)", "Odds (rotan / no rotan)")) %>%
collapse_rows(columns = 1, valign = "top")| Variable | Categoría | N | Rotan | No rotan | Tasa de rotación (%) | Odds (rotan / no rotan) |
|---|---|---|---|---|---|---|
| Horas extra | No | 1054 | 110 | 944 | 10.44 | 0.12 |
| Si | 416 | 127 | 289 | 30.53 | 0.44 | |
| Estado civil | Casado | 673 | 84 | 589 | 12.48 | 0.14 |
| Divorciado | 327 | 33 | 294 | 10.09 | 0.11 | |
| Soltero | 470 | 120 | 350 | 25.53 | 0.34 | |
| Viaje de negocios | No_Viaja | 150 | 12 | 138 | 8.00 | 0.09 |
| Raramente | 1043 | 156 | 887 | 14.96 | 0.18 | |
| Frecuentemente | 277 | 69 | 208 | 24.91 | 0.33 |
prueba_chi <- map_dfr(vars_cat_sel, function(v) {
tab <- table(datos[[v]], datos$rotacion)
ct <- chisq.test(tab, correct = FALSE)
n <- sum(tab)
tibble(
Variable = etiquetas[[v]],
Chi2 = as.numeric(ct$statistic),
gl = as.numeric(ct$parameter),
p_valor = fmt_p(ct$p.value),
V_Cramer = sqrt(as.numeric(ct$statistic) / (n * (min(dim(tab)) - 1))),
Min_esperada = min(ct$expected)
)
})
estilo_tabla(prueba_chi, digits = 3,
caption = "Tabla 9. Prueba chi-cuadrado de independencia entre rotación y cada variable categórica (V de Cramér como medida del tamaño del efecto; se verifica que la frecuencia esperada mínima sea al menos 5).",
col.names = c("Variable", "χ²", "gl", "p-valor", "V de Cramér", "Frecuencia esperada mínima"))| Variable | χ² | gl | p-valor | V de Cramér | Frecuencia esperada mínima |
|---|---|---|---|---|---|
| Horas extra | 89.04 | 1 | <0.001 | 0.246 | 67.07 |
| Estado civil | 46.16 | 2 | <0.001 | 0.177 | 52.72 |
| Viaje de negocios | 24.18 | 2 | <0.001 | 0.128 | 24.18 |
grafico_tasa <- function(df, var, ymax = 0.35) {
tab <- df %>%
group_by(categoria = .data[[var]]) %>%
summarise(n = n(), tasa = mean(y), .groups = "drop")
ggplot(tab, aes(x = categoria, y = tasa)) +
geom_col(fill = col_naranja, width = 0.65) +
geom_hline(yintercept = mean(df$y), linetype = "dashed", color = "gray30") +
geom_text(aes(label = paste0(percent(tasa, accuracy = 0.1), "\n(n = ", n, ")")),
vjust = 1.15, size = 3.2, color = "white", fontface = "bold") +
scale_y_continuous(labels = percent, limits = c(0, ymax), expand = expansion(mult = c(0, 0))) +
labs(title = etiquetas[[var]], x = NULL, y = "Tasa de rotación")
}
wrap_plots(map(vars_cat_sel, ~ grafico_tasa(datos, .x)), nrow = 1)
Figura 6. Tasa de rotación por categoría de las variables categóricas seleccionadas, con la misma escala vertical en los tres paneles (la línea punteada es la tasa global de rotación, 16.1%).
En las Tablas 8, 9 y la Figura 6, la tasa global de rotación es de 16.1%. Por lo que de esto se puede interpretar:
- Horas extra. Rota el 30.5% de quienes hacen horas extra (127 de 416) frente al 10.4% de quienes no las hacen (110 de 1.054): casi tres veces más (2.9). Los odds pasan de 0.12 a 0.44. La asociación es significativa (χ² = 89.04; 1 gl; p < 0.001) y es la más fuerte de las tres (V de Cramér = 0.246, tamaño de efecto entre pequeño y moderado).
- Estado civil. Los solteros rotan el 25.5% (120 de 470), el doble que los casados (12.5%; 84 de 673); los divorciados rotan el 10.1% (33 de 327), similar o algo menos que los casados. Asociación significativa (χ² = 46.16; 2 gl; p < 0.001; V = 0.177).
- Viaje de negocios. Se observa un gradiente creciente con la frecuencia de viaje: 8.0% (No viaja; 12 de 150), 15.0% (Raramente; 156 de 1.043) y 24.9% (Frecuentemente; 69 de 277). Quienes viajan con frecuencia rotan unas 3.1 veces más que quienes no viajan. Asociación significativa pero de menor magnitud (χ² = 24.18; 2 gl; p < 0.001; V = 0.128).
- Supuesto de la prueba. Las frecuencias esperadas mínimas son 67.07, 52.72 y 24.18, todas muy superiores a 5, por lo que la aproximación chi-cuadrado es válida.
Por encima de la tasa global quedan quienes hacen horas extra, los solteros y quienes viajan frecuentemente; por debajo, quienes no hacen horas extra, los casados, los divorciados y quienes no viajan. Quienes viajan raramente (15.0%) están prácticamente en la media. Estos patrones son consistentes con las hipótesis H1, H2 y H3.
3.2 Variables cuantitativas frente a la rotación
desc_grupo <- datos %>%
select(rotacion, all_of(vars_cuant_sel)) %>%
pivot_longer(-rotacion, names_to = "variable", values_to = "valor") %>%
group_by(variable, rotacion) %>%
summarise(n = n(), Media = mean(valor), DE = sd(valor),
Q1 = quantile(valor, 0.25, names = FALSE),
Mediana = median(valor),
Q3 = quantile(valor, 0.75, names = FALSE),
.groups = "drop") %>%
arrange(match(variable, vars_cuant_sel), rotacion) %>%
mutate(variable = unname(etiquetas[variable]), rotacion = as.character(rotacion))
estilo_tabla(desc_grupo, digits = 2,
caption = "Tabla 10. Estadísticos descriptivos de las variables cuantitativas seleccionadas según rotación.",
col.names = c("Variable", "Rotación", "n", "Media", "DE", "Q1", "Mediana", "Q3")) %>%
collapse_rows(columns = 1, valign = "top")| Variable | Rotación | n | Media | DE | Q1 | Mediana | Q3 |
|---|---|---|---|---|---|---|---|
| Ingreso mensual (miles) | No | 1233 | 6.83 | 4.82 | 3.21 | 5.2 | 8.83 |
| Si | 237 | 4.79 | 3.64 | 2.37 | 3.2 | 5.92 | |
| Antigüedad en el cargo | No | 1233 | 4.48 | 3.65 | 2.00 | 3.0 | 7.00 |
| Si | 237 | 2.90 | 3.17 | 0.00 | 2.0 | 4.00 | |
| Distancia a casa | No | 1233 | 8.92 | 8.01 | 2.00 | 7.0 | 13.00 |
| Si | 237 | 10.63 | 8.45 | 3.00 | 9.0 | 17.00 |
pruebas_cuant <- map_dfr(vars_cuant_sel, function(v) {
x_si <- datos[[v]][datos$rotacion == "Si"]
x_no <- datos[[v]][datos$rotacion == "No"]
tt <- t.test(x_si, x_no) # t de Welch (varianzas desiguales)
wt <- wilcox.test(x_si, x_no, exact = FALSE) # Wilcoxon-Mann-Whitney
sp <- sqrt(((length(x_si) - 1) * var(x_si) + (length(x_no) - 1) * var(x_no)) /
(length(x_si) + length(x_no) - 2))
tibble(
Variable = etiquetas[[v]],
Dif_medias = mean(x_si) - mean(x_no),
d_Cohen = (mean(x_si) - mean(x_no)) / sp,
t_Welch = unname(tt$statistic),
p_t = fmt_p(tt$p.value),
p_Wilcoxon = fmt_p(wt$p.value)
)
})
estilo_tabla(pruebas_cuant, digits = 3,
caption = "Tabla 11. Comparación de las variables cuantitativas entre quienes rotan y quienes no (diferencia de medias = rotan - no rotan; d de Cohen como tamaño del efecto).",
col.names = c("Variable", "Diferencia de medias", "d de Cohen", "t (Welch)",
"p-valor (t)", "p-valor (Wilcoxon)"))| Variable | Diferencia de medias | d de Cohen | t (Welch) | p-valor (t) | p-valor (Wilcoxon) |
|---|---|---|---|---|---|
| Ingreso mensual (miles) | -2.046 | -0.440 | -7.483 | <0.001 | <0.001 |
| Antigüedad en el cargo | -1.581 | -0.442 | -6.847 | <0.001 | <0.001 |
| Distancia a casa | 1.717 | 0.212 | 2.888 | 0.004 | 0.002 |
datos %>%
select(rotacion, all_of(vars_cuant_sel)) %>%
pivot_longer(-rotacion, names_to = "variable", values_to = "valor") %>%
mutate(variable = factor(variable, levels = vars_cuant_sel,
labels = unname(etiquetas[vars_cuant_sel]))) %>%
ggplot(aes(x = rotacion, y = valor, fill = rotacion)) +
geom_boxplot(alpha = 0.8, outlier.alpha = 0.4) +
stat_summary(fun = mean, geom = "point", shape = 23, size = 2.5, fill = "white") +
scale_fill_manual(values = c(No = col_azul, Si = col_naranja)) +
facet_wrap(~ variable, scales = "free_y") +
labs(x = "Rotación", y = NULL) +
theme(legend.position = "none")
Figura 7. Distribución de las variables cuantitativas seleccionadas según rotación (el rombo blanco indica la media).
El modelo logit supone que la relación entre cada variable cuantitativa y el logaritmo de los odds es lineal. La Figura 8 compara el ajuste logístico simple (paramétrico, en naranja) con un suavizado no paramétrico loess (azul, discontinuo): si ambas curvas son similares, el supuesto de linealidad es razonable.
datos %>%
select(y, all_of(vars_cuant_sel)) %>%
pivot_longer(-y, names_to = "variable", values_to = "valor") %>%
mutate(variable = factor(variable, levels = vars_cuant_sel,
labels = unname(etiquetas[vars_cuant_sel]))) %>%
ggplot(aes(x = valor, y = y)) +
geom_jitter(height = 0.04, alpha = 0.12, color = "gray30") +
geom_smooth(method = "glm", formula = y ~ x, method.args = list(family = binomial),
color = col_naranja, fill = col_naranja, alpha = 0.2) +
geom_smooth(method = "loess", formula = y ~ x, se = FALSE,
color = col_azul, linetype = "dashed") +
facet_wrap(~ variable, scales = "free_x") +
labs(x = NULL, y = "P(rotación = 1)")
Figura 8. Probabilidad de rotar según cada variable cuantitativa seleccionada: ajuste logístico simple (naranja, con banda de confianza) frente a suavizado loess (azul, discontinuo). Los puntos son las observaciones (con dispersión vertical para evitar superposición).
Interpretación de las Tablas 10 y 11, al igual que las Figuras 7 y 8:
- Ingreso mensual. Quienes rotan ganan menos: media de 4.79 mil (mediana 3.2 mil) frente a 6.83 mil (mediana 5.2 mil) de quienes no rotan. La diferencia de medias es de -2.05 mil, con d de Cohen = -0.44 (efecto entre pequeño y moderado); t de Welch = -7.48 y p < 0.001 (Wilcoxon también p < 0.001).
- Antigüedad en el cargo. Los que rotan tienen media de 2.90 años (mediana 2) frente a 4.48 (mediana 3) de los que no rotan; diferencia de -1.58 años, d = -0.44, t = -6.85, p < 0.001 (Wilcoxon p < 0.001). El primer cuartil de quienes rotan es 0: al menos una cuarta parte de las rotaciones ocurrió en el primer año en el cargo.
- Distancia a casa. Quienes rotan viven más lejos: media de 10.63 (mediana 9) frente a 8.92 (mediana 7); diferencia de +1.72, d = 0.21 (efecto pequeño), t = 2.89, p = 0.004 (Wilcoxon p = 0.002). Es la diferencia más modesta de las tres, aunque significativa.
Las pruebas paramétrica (Welch) y no paramétrica (Wilcoxon) coinciden en todos los casos, por lo que la conclusión no depende de la asimetría de las variables.
Para la Distancia a casa el ajuste logístico y el suavizado loess casi coinciden: el supuesto de linealidad es razonable. Para la Antigüedad en el cargo hay una concordancia general, aunque el loess queda algo por encima del ajuste logístico en el primer año y se aplana alrededor de 7% en antigüedades altas (más de 10 años) en lugar de seguir bajando, pero sin señales de que la rotación vuelva a aumentar con mucha antigüedad (no se observa el “estancamiento” mencionado en H5). Para el Ingreso mensual la desviación es mayor: el loess muestra una probabilidad de rotar muy alta en los ingresos más bajos (cercana a 0.5 hacia 1 mil), que cae rápido hasta cerca de 0.10 en 5 a 6 mil y luego se mantiene en 0.12 a 0.14 entre 8 y 15 mil antes de bajar de nuevo; el ajuste logístico, más suave, subestima el riesgo en los ingresos más bajos y en el tramo de 8 a 15 mil. Esto sugiere que el efecto del ingreso se concentra en los salarios más bajos. Los extremos del loess se apoyan en pocos datos y deben leerse con cautela; por eso la linealidad de las tres variables se contrasta formalmente en la Sección 4.8 (spline natural y transformaciones logarítmicas).
3.3 Modelos logit simples para las 6 variables seleccionadas
Para cada variable se estima \(\ln\left(\frac{p}{1-p}\right) = \beta_0 + \beta_1
x\) con glm(..., family = binomial(link = "logit")).
Se reportan el coeficiente, su error estándar, el estadístico \(z\) de Wald, el p-valor, el OR y su
intervalo de confianza al 95% (por verosimilitud perfilada, con
confint). Para variables con más de dos categorías, la
significancia global de la variable se evalúa con la prueba de razón de
verosimilitud (Tabla 13).
ajustes_simples <- map(set_names(vars_sel), function(v) {
glm(reformulate(v, response = "y"), family = binomial(link = "logit"), data = datos)
})
coef_simples <- imap_dfr(ajustes_simples, function(m, v) {
tidy(m, conf.int = TRUE) %>%
filter(term != "(Intercept)") %>%
mutate(variable = v, .before = 1)
}) %>%
mutate(OR = exp(estimate), OR_inf = exp(conf.low), OR_sup = exp(conf.high))
tabla_biv <- coef_simples %>%
transmute(
Variable = unname(etiquetas[variable]),
Contraste = etq_termino(term),
Beta = estimate, EE = std.error, z = statistic,
p = fmt_p(p.value), Sig = signif_stars(p.value),
OR = OR, OR_inf = OR_inf, OR_sup = OR_sup
)
estilo_tabla(tabla_biv, digits = 3,
caption = "Tabla 12. Modelos logit simples (bivariados): coeficientes, significancia de Wald y odds ratio con IC 95% (Sig.: *** p<0.001, ** p<0.01, * p<0.05, . p<0.1).",
col.names = c("Variable", "Contraste", "β", "EE", "z", "p-valor",
"Sig.", "OR", "IC 95% inf.", "IC 95% sup.")) %>%
collapse_rows(columns = 1, valign = "top")| Variable | Contraste | β | EE | z | p-valor | Sig. | OR | IC 95% inf. | IC 95% sup. |
|---|---|---|---|---|---|---|---|---|---|
| Horas extra | Horas extra: Sí (vs. No) | 1.327 | 0.147 | 9.056 | <0.001 | *** | 3.771 | 2.832 | 5.032 |
| Estado civil | Estado civil: Divorciado (vs. Casado) | -0.239 | 0.217 | -1.101 | 0.271 | 0.787 | 0.508 | 1.194 | |
| Estado civil: Soltero (vs. Casado) | 0.877 | 0.157 | 5.571 | <0.001 | *** | 2.404 | 1.769 | 3.281 | |
| Viaje de negocios | Viaje: Raramente (vs. no viaja) | 0.704 | 0.313 | 2.249 | 0.025 | * | 2.023 | 1.139 | 3.931 |
| Viaje: Frecuentemente (vs. no viaja) | 1.339 | 0.331 | 4.039 | <0.001 | *** | 3.815 | 2.061 | 7.641 | |
| Ingreso mensual (miles) | Ingreso mensual (por cada 1.000) | -0.127 | 0.022 | -5.879 | <0.001 | *** | 0.881 | 0.842 | 0.917 |
| Antigüedad en el cargo | Antigüedad en el cargo (por año) | -0.146 | 0.024 | -6.033 | <0.001 | *** | 0.864 | 0.823 | 0.905 |
| Distancia a casa | Distancia a casa (por unidad) | 0.025 | 0.008 | 2.973 | 0.003 | ** | 1.025 | 1.008 | 1.042 |
lr_simples <- imap_dfr(ajustes_simples, function(m, v) {
a <- anova(m, test = "Chisq")
tibble(variable = unname(etiquetas[v]),
gl = a$Df[2], LR_chi2 = a$Deviance[2], p = fmt_p(a$`Pr(>Chi)`[2]))
})
estilo_tabla(lr_simples, digits = 3,
caption = "Tabla 13. Prueba de razón de verosimilitud (significancia global) de cada modelo logit simple frente al modelo nulo.",
col.names = c("Variable", "gl", "χ² (LR)", "p-valor"))| Variable | gl | χ² (LR) | p-valor |
|---|---|---|---|
| Horas extra | 1 | 81.40 | <0.001 |
| Estado civil | 2 | 44.00 | <0.001 |
| Viaje de negocios | 2 | 23.76 | <0.001 |
| Ingreso mensual (miles) | 1 | 45.49 | <0.001 |
| Antigüedad en el cargo | 1 | 42.70 | <0.001 |
| Distancia a casa | 1 | 8.58 | 0.003 |
Todos los coeficientes de este análisis son crudos (cada variable por separado). La interpretación de las Tablas 12 y 13 es la siguiente:
- Horas extra. β = 1.327 (signo positivo), OR = 3.77 (IC 95%: 2.83 a 5.03), z = 9.06, p < 0.001. Los odds de rotar de quienes hacen horas extra son 3.77 veces los de quienes no las hacen. LR = 81.40 (1 gl), p < 0.001.
- Estado civil. Soltero vs. Casado: β = 0.877, OR = 2.40 (1.77 a 3.28), p < 0.001; los odds de rotar de un soltero son 2.4 veces los de un casado. Divorciado vs. Casado: β = -0.239, OR = 0.787 (0.508 a 1.194), p = 0.271; el IC incluye 1, por lo que no hay evidencia de diferencia. La significancia global (LR = 44.00; 2 gl; p < 0.001) proviene de los solteros.
- Viaje de negocios. Raramente vs. No viaja: β = 0.704, OR = 2.02 (1.14 a 3.93), p = 0.025. Frecuentemente vs. No viaja: β = 1.339, OR = 3.82 (2.06 a 7.64), p < 0.001. Los OR crecen con la frecuencia de viaje, pero los IC son amplios porque la categoría de referencia solo tiene 150 empleados y 12 rotaciones. LR global = 23.76 (2 gl), p < 0.001.
- Ingreso mensual (miles). β = -0.127 (signo negativo), OR = 0.881 (0.842 a 0.917), p < 0.001: por cada 1.000 unidades adicionales de ingreso, los odds de rotar se reducen en 11.9%. LR = 45.49, p < 0.001.
- Antigüedad en el cargo. β = -0.146, OR = 0.864 (0.823 a 0.905), p < 0.001: cada año adicional en el cargo reduce los odds de rotar en 13.6%. LR = 42.70, p < 0.001.
- Distancia a casa. β = 0.025 (signo positivo), OR = 1.025 (1.008 a 1.042), p = 0.003: cada unidad adicional de distancia aumenta los odds en 2.5%. LR = 8.58, p = 0.003.
Ordenadas por el estadístico LR, la asociación bivariada más fuerte es la de Horas extra (81.4), seguida de Ingreso (45.5), Estado civil (44.0), Antigüedad en el cargo (42.7), Viaje de negocios (23.8) y Distancia a casa (8.6).
3.4 Comparación de los resultados con las hipótesis planteadas
La Tabla 14 contrasta el signo esperado (Sección 1) con el signo estimado en el modelo logit simple y su significancia al 5%.
hipotesis <- tribble(
~H, ~Variable, ~Contraste, ~term, ~signo_esp,
"H1", "Horas extra", "Sí vs. No", "horas_extraSi", "+",
"H2", "Estado civil", "Soltero vs. Casado", "estado_civilSoltero", "+",
"H3", "Viaje de negocios", "Raramente vs. No viaja", "viaje_negociosRaramente", "+",
"H3", "Viaje de negocios", "Frecuentemente vs. No viaja", "viaje_negociosFrecuentemente", "+",
"H4", "Ingreso mensual", "+1.000 unidades", "ingreso_miles", "-",
"H5", "Antigüedad en el cargo","+1 año", "antiguedad_cargo", "-",
"H6", "Distancia a casa", "+1 unidad", "distancia_casa", "+"
)
comparacion <- hipotesis %>%
left_join(coef_simples %>% select(term, estimate, OR, p.value), by = "term")
# Verificación: todos los términos de las hipótesis deben existir en los modelos
stopifnot(!anyNA(comparacion$estimate))
comparacion <- comparacion %>%
mutate(
signo_obs = ifelse(estimate > 0, "+", "-"),
Resultado = case_when(
p.value >= 0.05 ~ "No significativo: no se confirma la hipótesis",
signo_obs == signo_esp ~ "Se confirma la hipótesis",
TRUE ~ "Se contradice la hipótesis"
)
) %>%
transmute(H, Variable, Contraste,
Signo_esperado = ifelse(signo_esp == "+", "Positivo (+)", "Negativo (-)"),
Signo_observado = ifelse(signo_obs == "+", "Positivo (+)", "Negativo (-)"),
Beta = estimate, OR, p = fmt_p(p.value), Resultado)
estilo_tabla(comparacion, digits = 3,
caption = "Tabla 14. Hipótesis planteadas frente a los resultados del análisis bivariado (modelos logit simples).",
col.names = c("Hipótesis", "Variable", "Contraste", "Signo esperado",
"Signo observado", "β", "OR", "p-valor", "Resultado"))| Hipótesis | Variable | Contraste | Signo esperado | Signo observado | β | OR | p-valor | Resultado |
|---|---|---|---|---|---|---|---|---|
| H1 | Horas extra | Sí vs. No | Positivo (+) | Positivo (+) | 1.327 | 3.771 | <0.001 | Se confirma la hipótesis |
| H2 | Estado civil | Soltero vs. Casado | Positivo (+) | Positivo (+) | 0.877 | 2.404 | <0.001 | Se confirma la hipótesis |
| H3 | Viaje de negocios | Raramente vs. No viaja | Positivo (+) | Positivo (+) | 0.704 | 2.023 | 0.025 | Se confirma la hipótesis |
| H3 | Viaje de negocios | Frecuentemente vs. No viaja | Positivo (+) | Positivo (+) | 1.339 | 3.815 | <0.001 | Se confirma la hipótesis |
| H4 | Ingreso mensual | +1.000 unidades | Negativo (-) | Negativo (-) | -0.127 | 0.881 | <0.001 | Se confirma la hipótesis |
| H5 | Antigüedad en el cargo | +1 año | Negativo (-) | Negativo (-) | -0.146 | 0.864 | <0.001 | Se confirma la hipótesis |
| H6 | Distancia a casa | +1 unidad | Positivo (+) | Positivo (+) | 0.025 | 1.025 | 0.003 | Se confirma la hipótesis |
Los siete contrastes incluidos en la Tabla 14 resultan en “Se confirma la hipótesis”: en todos los casos el signo estimado coincide con el esperado y el p-valor es menor que 0.05. Ninguna hipótesis se contradice en estos contrastes. Además, H3 muestra descriptivamente el gradiente esperado (OR de 3.82 para viajar frecuentemente frente a 2.02 para viajar raramente), pero estos dos OR comparan cada categoría con “No viaja” y no constituyen por sí solos una prueba formal de que “Frecuentemente” tenga mayores odds que “Raramente”. Por ello, se realiza a continuación el contraste directo entre ambas categorías.
m_h3 <- ajustes_simples[["viaje_negocios"]]
b_h3 <- coef(m_h3)
V_h3 <- vcov(m_h3)
L_h3 <- setNames(rep(0, length(b_h3)), names(b_h3))
L_h3["viaje_negociosFrecuentemente"] <- 1
L_h3["viaje_negociosRaramente"] <- -1
est_h3 <- sum(L_h3 * b_h3)
se_h3 <- sqrt(as.numeric(t(L_h3) %*% V_h3 %*% L_h3))
z_h3 <- est_h3 / se_h3
p_h3 <- 2 * pnorm(-abs(z_h3))
ci_h3 <- est_h3 + c(-1, 1) * qnorm(0.975) * se_h3
contraste_h3 <- tibble(
Contraste = "Frecuentemente vs. Raramente",
Beta = est_h3,
OR = exp(est_h3),
IC_inf = exp(ci_h3[1]),
IC_sup = exp(ci_h3[2]),
p = p_h3
)
estilo_tabla(contraste_h3, digits = 3,
caption = "Tabla 14A. Contraste directo entre las categorías de frecuencia de viaje: Frecuentemente frente a Raramente.",
col.names = c("Contraste", "β", "OR", "IC 95% inf.", "IC 95% sup.", "p-valor"))| Contraste | β | OR | IC 95% inf. | IC 95% sup. | p-valor |
|---|---|---|---|---|---|
| Frecuentemente vs. Raramente | 0.635 | 1.886 | 1.368 | 2.6 | 0 |
El contraste directo evalúa si los odds de rotar difieren entre quienes viajan frecuentemente y quienes viajan raramente, sin utilizar “No viaja” como referencia. El resultado es estadísticamente significativo (OR = 1.89, p = <0.001), por lo que las dos categorías presentan odds de rotación diferentes. En consecuencia, el orden descriptivo de H3 se reporta como un patrón observado, pero no se presenta como un gradiente formalmente demostrado sin este contraste.
3.5 Tamizaje de todas las variables: identificación de las variables asociadas
Para identificar las variables asociadas con la rotación en toda la base (no solo en las 6 seleccionadas), se ajusta un modelo logit simple por cada una de las 23 variables restantes. Las variables ordinales se tratan como numéricas (efecto lineal por nivel) y las nominales como factores. Como se realizan 23 contrastes simultáneos, además del p-valor de la prueba de razón de verosimilitud se reporta el p-valor ajustado por Benjamini-Hochberg (BH) para controlar la tasa de falsos descubrimientos.
Además, cada variable se clasifica según su naturaleza (criterio propio): palanca de gestión (la empresa puede modificarla: horas extra, viaje de negocios, ingreso, porcentaje de aumento, capacitaciones, años desde la última promoción, satisfacciones y equilibrio trabajo-vida), contextual (condiciones de la organización o del puesto: cargo, departamento, distancia a casa, antigüedades y años con el mismo jefe) y perfil del empleado (características personales o de trayectoria).
vars_todas <- setdiff(nombres_nuevos, "rotacion")
sel_orig <- c("horas_extra", "estado_civil", "viaje_negocios",
"ingreso_mensual", "antiguedad_cargo", "distancia_casa")
# Naturaleza del factor
palancas <- c("horas_extra", "viaje_negocios", "ingreso_mensual", "porcentaje_aumento_salarial",
"capacitaciones", "anos_ultima_promocion", "satisfaccion_ambiental",
"satisfaccion_laboral", "equilibrio_trabajo_vida")
contextuales <- c("cargo", "departamento", "distancia_casa", "antiguedad_cargo",
"anos_mismo_jefe", "antiguedad")
stopifnot(all(c(palancas, contextuales) %in% vars_todas))
tamizaje <- map_dfr(vars_todas, function(v) {
m <- glm(reformulate(v, response = "y"), family = binomial(link = "logit"), data = datos)
a <- anova(m, test = "Chisq")
b <- coef(m)
gl <- a$Df[2]
f <- datos[[v]]
contraste <- if (is.factor(f)) {
if (gl == 1) paste0(levels(f)[2], " vs. ", levels(f)[1]) else paste0(gl + 1, " categorías")
} else if (v %in% vars_ordinales) {
"Por cada nivel"
} else {
"Por cada unidad"
}
tibble(
variable = v,
Tipo = case_when(v %in% vars_nominales ~ "Nominal",
v %in% vars_ordinales ~ "Ordinal (numérica)",
TRUE ~ "Cuantitativa"),
Factor = case_when(v %in% palancas ~ "Palanca de gestión",
v %in% contextuales ~ "Contextual",
TRUE ~ "Perfil del empleado"),
Contraste = contraste,
gl = gl,
LR_chi2 = a$Deviance[2],
p_LR = a$`Pr(>Chi)`[2],
Signo = ifelse(gl == 1, ifelse(b[2] > 0, "Positivo", "Negativo"), "varios niveles")
)
}) %>%
mutate(p_BH = p.adjust(p_LR, method = "BH"),
En_modelo = ifelse(variable %in% sel_orig, "Sí", "")) %>%
arrange(p_LR)
tabla_tamizaje <- tamizaje %>%
transmute(Variable = unname(etiquetas[variable]), Tipo, Contraste,
chi = paste0(f_num(LR_chi2, 2), " (", gl, ")"),
p = fmt_p(p_LR), p_ajustado = fmt_p(p_BH), Sig = signif_stars(p_BH),
Signo, Factor, En_modelo)
estilo_tabla(tabla_tamizaje, digits = 2,
caption = "Tabla 15. Tamizaje bivariado de todas las variables (modelos logit simples), ordenado por significancia. Sig. corresponde al p-valor ajustado (BH). Signo: sentido del coeficiente cuando la variable tiene un solo parámetro (según el contraste indicado). Tipo de factor: criterio propio para orientar la estrategia.",
col.names = c("Variable", "Tipo", "Contraste", "χ² (gl)", "p-valor (LR)",
"p ajustado (BH)", "Sig.", "Signo de β", "Tipo de factor",
"Seleccionada"))| Variable | Tipo | Contraste | χ² (gl) | p-valor (LR) | p ajustado (BH) | Sig. | Signo de β | Tipo de factor | Seleccionada |
|---|---|---|---|---|---|---|---|---|---|
| Horas extra | Nominal | Si vs. No | 81.40 (1) | <0.001 | <0.001 | *** | Positivo | Palanca de gestión | Sí |
| Cargo | Nominal | 9 categorías | 88.91 (8) | <0.001 | <0.001 | *** | varios niveles | Contextual | |
| Años de experiencia | Cuantitativa | Por cada unidad | 50.51 (1) | <0.001 | <0.001 | *** | Negativo | Perfil del empleado | |
| Ingreso mensual | Cuantitativa | Por cada unidad | 45.49 (1) | <0.001 | <0.001 | *** | Negativo | Palanca de gestión | Sí |
| Antigüedad en el cargo | Cuantitativa | Por cada unidad | 42.70 (1) | <0.001 | <0.001 | *** | Negativo | Contextual | Sí |
| Años con el mismo jefe | Cuantitativa | Por cada unidad | 39.87 (1) | <0.001 | <0.001 | *** | Negativo | Contextual | |
| Estado civil | Nominal | 3 categorías | 44.00 (2) | <0.001 | <0.001 | *** | varios niveles | Perfil del empleado | Sí |
| Edad | Cuantitativa | Por cada unidad | 39.52 (1) | <0.001 | <0.001 | *** | Negativo | Perfil del empleado | |
| Antigüedad (empresa) | Cuantitativa | Por cada unidad | 32.08 (1) | <0.001 | <0.001 | *** | Negativo | Contextual | |
| Viaje de negocios | Nominal | 3 categorías | 23.76 (2) | <0.001 | <0.001 | *** | varios niveles | Palanca de gestión | Sí |
| Satisfacción laboral | Ordinal (numérica) | Por cada nivel | 15.53 (1) | <0.001 | <0.001 | *** | Negativo | Palanca de gestión | |
| Satisfacción ambiental | Ordinal (numérica) | Por cada nivel | 15.50 (1) | <0.001 | <0.001 | *** | Negativo | Palanca de gestión | |
| Distancia a casa | Cuantitativa | Por cada unidad | 8.58 (1) | 0.003 | 0.006 | ** | Positivo | Contextual | Sí |
| Departamento | Nominal | 3 categorías | 10.49 (2) | 0.005 | 0.009 | ** | varios niveles | Contextual | |
| Campo de educación | Nominal | 6 categorías | 14.90 (5) | 0.011 | 0.017 | * | varios niveles | Perfil del empleado | |
| Equilibrio trabajo-vida | Ordinal (numérica) | Por cada nivel | 5.90 (1) | 0.015 | 0.022 | * | Negativo | Palanca de gestión | |
| Capacitaciones | Cuantitativa | Por cada unidad | 5.32 (1) | 0.021 | 0.029 | * | Negativo | Palanca de gestión | |
| Trabajos anteriores | Cuantitativa | Por cada unidad | 2.71 (1) | 0.100 | 0.127 | Positivo | Perfil del empleado | ||
| Años desde última promoción | Cuantitativa | Por cada unidad | 1.67 (1) | 0.196 | 0.237 | Negativo | Palanca de gestión | ||
| Educación (nivel) | Ordinal (numérica) | Por cada nivel | 1.44 (1) | 0.230 | 0.265 | Negativo | Perfil del empleado | ||
| Género | Nominal | M vs. F | 1.29 (1) | 0.257 | 0.281 | Positivo | Perfil del empleado | ||
| % aumento salarial | Cuantitativa | Por cada unidad | 0.27 (1) | 0.604 | 0.632 | Negativo | Palanca de gestión | ||
| Rendimiento laboral | Ordinal (numérica) | Por cada nivel | 0.01 (1) | 0.912 | 0.912 | Positivo | Perfil del empleado |
determinantes <- tamizaje %>% filter(p_BH < 0.05) %>% pull(variable)
det_pos <- tamizaje %>% filter(p_BH < 0.05, Signo == "Positivo") %>% pull(variable)
det_neg <- tamizaje %>% filter(p_BH < 0.05, Signo == "Negativo") %>% pull(variable)
det_multi <- tamizaje %>% filter(p_BH < 0.05, Signo == "varios niveles") %>% pull(variable)
no_sig <- tamizaje %>% filter(p_BH >= 0.05) %>% pull(variable)
det_gest <- tamizaje %>% filter(p_BH < 0.05, Factor == "Palanca de gestión") %>% pull(variable)
det_ctx <- tamizaje %>% filter(p_BH < 0.05, Factor == "Contextual") %>% pull(variable)
det_per <- tamizaje %>% filter(p_BH < 0.05, Factor == "Perfil del empleado") %>% pull(variable)Variables con p-valor ajustado (BH) < 0.05: Horas extra, Cargo, Años de experiencia, Ingreso mensual, Antigüedad en el cargo, Años con el mismo jefe, Estado civil, Edad, Antigüedad (empresa), Viaje de negocios, Satisfacción laboral, Satisfacción ambiental, Distancia a casa, Departamento, Campo de educación, Equilibrio trabajo-vida, Capacitaciones.
Para la interpretación de la Tabla 15, tras ajustar por comparaciones múltiples (BH), 17 de las 23 variables tienen asociación significativa con la rotación al 5%, y las 6 variables seleccionadas están entre ellas. Según el signo del coeficiente (columna “Signo de β”):
- Efecto positivo (aumentan la probabilidad de rotar): Horas extra, Distancia a casa.
- Efecto negativo (la disminuyen): Años de experiencia, Ingreso mensual, Antigüedad en el cargo, Años con el mismo jefe, Edad, Antigüedad (empresa), Satisfacción laboral, Satisfacción ambiental, Equilibrio trabajo-vida, Capacitaciones.
- Variables categóricas con varios niveles (significativas de forma global): Cargo, Estado civil, Viaje de negocios, Departamento, Campo de educación; el efecto por categoría se lee en los modelos simples o en el múltiple.
- Sin asociación significativa: Trabajos anteriores, Años desde última promoción, Educación (nivel), Género, % aumento salarial, Rendimiento laboral (p ajustado entre 0.127 y 0.912).
Según la naturaleza del factor (columna “Tipo de factor”), las palancas de gestión con asociación significativa son: Horas extra, Ingreso mensual, Viaje de negocios, Satisfacción laboral, Satisfacción ambiental, Equilibrio trabajo-vida, Capacitaciones; las variables contextuales significativas son: Cargo, Antigüedad en el cargo, Años con el mismo jefe, Antigüedad (empresa), Distancia a casa, Departamento; y las de perfil del empleado significativas, útiles para focalizar las acciones, son: Años de experiencia, Estado civil, Edad, Campo de educación. Esta lectura se retomará en el punto 7 (estrategia de retención).
Por el estadístico LR, las asociaciones más fuertes son Cargo (χ² = 88.9, pero con 8 gl), Horas extra (81.4), Años de experiencia (50.5), Ingreso mensual (45.5), Estado civil (44.0), Antigüedad en el cargo (42.7), Años con el mismo jefe (39.9) y Edad (39.5). Conviene tener presente tres cosas:
- Varias de estas variables miden aspectos parecidos (experiencia, antigüedad, edad e ingreso están correlacionadas entre sí, Tabla 3), por lo que en un modelo múltiple muchas perderían significancia: una asociación bivariada significativa no equivale a un efecto independiente.
- Entre las variables no seleccionadas hay candidatas para una segunda versión del modelo: Cargo, las satisfacciones laboral y ambiental (χ² ≈ 15.5 cada una) y el equilibrio trabajo-vida (χ² = 5.9; p ajustado = 0.022). Estas últimas son, además, factores sobre los cuales la empresa puede actuar.
- Rendimiento laboral, % de aumento salarial, género, educación, años desde la última promoción y trabajos anteriores no muestran asociación con la rotación.
3.6 Síntesis del análisis bivariado
Las seis variables seleccionadas se asocian con la rotación y en la dirección planteada (Tabla 14). Tienen más riesgo de rotar quienes hacen horas extra (OR crudo 3.77), los solteros frente a los casados (2.40), quienes viajan frecuentemente o raramente frente a quienes no viajan (3.82 y 2.02) y quienes viven más lejos (1.025 por unidad); tienen menos riesgo quienes tienen mayor ingreso (0.881 por cada 1.000) y mayor antigüedad en el cargo (0.864 por año). Por el estadístico LR, Horas extra presenta el mayor estadístico de contraste (81.4), seguida de Ingreso mensual (45.5), Estado civil (44.0), Antigüedad en el cargo (42.7), Viaje de negocios (23.8) y Distancia a casa (8.6). Estos valores describen la evidencia de asociación global de cada variable; no deben interpretarse como una medida directa de “importancia” o causalidad, y el estadístico LR tampoco es directamente comparable entre variables con distinto número de grados de libertad. Estos resultados son crudos, y el gráfico de linealidad sugiere que la relación con el ingreso no es estrictamente lineal en el logit (Sección 3.2). El tamizaje ubica además a Cargo, Años de experiencia, Edad y las satisfacciones entre las variables asociadas que podrían explorarse en una segunda versión del modelo.
4 Estimación del modelo de regresión logística
4.1 Especificación del modelo
Se estima un modelo logit múltiple con las 6 variables seleccionadas como covariables:
\[ \ln\left(\frac{p_i}{1-p_i}\right) = \beta_0 + \beta_1\,\text{HE}_i + \beta_2\,\text{SOLT}_i + \beta_3\,\text{DIV}_i + \beta_4\,\text{RAR}_i + \] \[ \beta_5\,\text{FREC}_i + \beta_6\,\text{ING}_i + \beta_7\,\text{ANT}_i + \beta_8\,\text{DIST}_i \]
donde \(p_i = P(Y_i = 1 \mid X)\) es la probabilidad de que el empleado \(i\) rote; HE = horas extra (Sí); SOLT y DIV = variables indicadoras de estado civil (Soltero, Divorciado; referencia: Casado); RAR y FREC = indicadoras de viaje de negocios (Raramente, Frecuentemente; referencia: no viaja); ING = ingreso mensual en miles; ANT = antigüedad en el cargo (años); DIST = distancia a casa.
niveles_ref <- tibble(
Variable = c("Horas extra", "Estado civil", "Viaje de negocios"),
Referencia = c(levels(datos$horas_extra)[1], levels(datos$estado_civil)[1],
levels(datos$viaje_negocios)[1]),
Otras = c(paste(levels(datos$horas_extra)[-1], collapse = ", "),
paste(levels(datos$estado_civil)[-1], collapse = ", "),
paste(levels(datos$viaje_negocios)[-1], collapse = ", "))
)
estilo_tabla(niveles_ref,
caption = "Tabla 16. Categorías de referencia y categorías comparadas para las variables categóricas del modelo.",
col.names = c("Variable", "Categoría de referencia", "Categorías comparadas"))| Variable | Categoría de referencia | Categorías comparadas |
|---|---|---|
| Horas extra | No | Si |
| Estado civil | Casado | Divorciado, Soltero |
| Viaje de negocios | No_Viaja | Raramente, Frecuentemente |
4.2 Estimación por máxima verosimilitud
modelo_multi <- glm(
y ~ horas_extra + estado_civil + viaje_negocios +
ingreso_miles + antiguedad_cargo + distancia_casa,
family = binomial(link = "logit"),
data = datos
)
# Salida 2. Resumen del modelo
summary(modelo_multi)##
## Call:
## glm(formula = y ~ horas_extra + estado_civil + viaje_negocios +
## ingreso_miles + antiguedad_cargo + distancia_casa, family = binomial(link = "logit"),
## data = datos)
##
## Coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.62669 0.38350 -6.85 7.4e-12 ***
## horas_extraSi 1.43850 0.15815 9.10 < 2e-16 ***
## estado_civilDivorciado -0.29060 0.23051 -1.26 0.20743
## estado_civilSoltero 0.87414 0.17143 5.10 3.4e-07 ***
## viaje_negociosRaramente 0.70857 0.33148 2.14 0.03255 *
## viaje_negociosFrecuentemente 1.34700 0.35340 3.81 0.00014 ***
## ingreso_miles -0.10258 0.02337 -4.39 1.1e-05 ***
## antiguedad_cargo -0.10810 0.02727 -3.96 7.4e-05 ***
## distancia_casa 0.03272 0.00931 3.51 0.00044 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## (Dispersion parameter for binomial family taken to be 1)
##
## Null deviance: 1298.6 on 1469 degrees of freedom
## Residual deviance: 1078.8 on 1461 degrees of freedom
## AIC: 1097
##
## Number of Fisher Scoring iterations: 5Salida 2. Resumen del modelo logit múltiple
(summary), con coeficientes, errores estándar, estadístico
\(z\) de Wald, devianzas y AIC.
4.3 Coeficientes, odds ratio y significancia de los parámetros
coef_multi_all <- tidy(modelo_multi, conf.int = TRUE) %>%
mutate(OR = exp(estimate), OR_inf = exp(conf.low), OR_sup = exp(conf.high))
coef_multi <- coef_multi_all %>% filter(term != "(Intercept)")
tabla_multi <- coef_multi_all %>%
transmute(
Termino = etq_termino(term), Beta = estimate, EE = std.error, z = statistic,
p = fmt_p(p.value), Sig = signif_stars(p.value),
OR = OR, OR_inf = OR_inf, OR_sup = OR_sup,
Cambio = ifelse(term == "(Intercept)", NA, round((OR - 1) * 100, 1))
)
estilo_tabla(tabla_multi, digits = 3,
caption = "Tabla 17. Modelo logit múltiple: coeficientes, prueba de Wald, odds ratio con IC 95% (verosimilitud perfilada) y cambio porcentual en los odds de rotar (Sig.: *** p<0.001, ** p<0.01, * p<0.05, . p<0.1).",
col.names = c("Término", "β", "EE", "z", "p-valor", "Sig.", "OR",
"IC 95% inf.", "IC 95% sup.", "Cambio en los odds (%)"))| Término | β | EE | z | p-valor | Sig. | OR | IC 95% inf. | IC 95% sup. | Cambio en los odds (%) |
|---|---|---|---|---|---|---|---|---|---|
| Intercepto | -2.627 | 0.383 | -6.849 | <0.001 | *** | 0.072 | 0.033 | 0.148 | |
| Horas extra: Sí (vs. No) | 1.439 | 0.158 | 9.096 | <0.001 | *** | 4.214 | 3.096 | 5.758 | 321.4 |
| Estado civil: Divorciado (vs. Casado) | -0.291 | 0.231 | -1.261 | 0.207 | 0.748 | 0.471 | 1.165 | -25.2 | |
| Estado civil: Soltero (vs. Casado) | 0.874 | 0.171 | 5.099 | <0.001 | *** | 2.397 | 1.716 | 3.363 | 139.7 |
| Viaje: Raramente (vs. no viaja) | 0.709 | 0.331 | 2.138 | 0.033 | * | 2.031 | 1.100 | 4.075 | 103.1 |
| Viaje: Frecuentemente (vs. no viaja) | 1.347 | 0.353 | 3.811 | <0.001 | *** | 3.846 | 1.985 | 8.010 | 284.6 |
| Ingreso mensual (por cada 1.000) | -0.103 | 0.023 | -4.390 | <0.001 | *** | 0.903 | 0.860 | 0.943 | -9.7 |
| Antigüedad en el cargo (por año) | -0.108 | 0.027 | -3.964 | <0.001 | *** | 0.898 | 0.850 | 0.946 | -10.2 |
| Distancia a casa (por unidad) | 0.033 | 0.009 | 3.513 | <0.001 | *** | 1.033 | 1.014 | 1.052 | 3.3 |
- Cada \(\beta_j\) es el cambio en el logaritmo de los odds de rotar por un aumento unitario de la variable (o por pasar de la categoría de referencia a la categoría indicada), manteniendo constantes las demás variables (ceteris paribus).
- \(\text{OR}_j = e^{\beta_j}\): multiplica los odds de rotar. Un OR > 1 indica más riesgo de rotar y un OR < 1 indica menos riesgo. El cambio porcentual en los odds es \((\text{OR} - 1) \times 100\).
- El intercepto corresponde al logit de un empleado de referencia (sin horas extra, casado, que no viaja) con ingreso, antigüedad y distancia iguales a cero; como esos valores quedan fuera del rango observado, el intercepto no tiene interpretación sustantiva y solo ancla el nivel de la curva.
- Significancia de cada parámetro: \(H_0: \beta_j = 0\) frente a \(H_1: \beta_j \neq 0\), con el estadístico de Wald \(z = \frac{\hat\beta_j} {\widehat{\text{EE}}(\hat\beta_j)}\).
coef_multi %>%
mutate(Significancia = ifelse(p.value < 0.05, "Significativo (5%)", "No significativo"),
etiqueta = etq_termino(term)) %>%
ggplot(aes(x = OR, y = reorder(etiqueta, OR), color = Significancia)) +
geom_vline(xintercept = 1, linetype = "dashed", color = "gray40") +
geom_pointrange(aes(xmin = OR_inf, xmax = OR_sup), size = 0.6) +
geom_text(aes(label = sprintf("%.2f", OR)), vjust = -1, size = 3, show.legend = FALSE) +
scale_x_log10() +
scale_color_manual(values = c("Significativo (5%)" = col_naranja,
"No significativo" = "gray55"), name = NULL) +
labs(x = "Odds ratio (escala logarítmica)", y = NULL) +
theme(legend.position = "bottom")
Figura 9. Odds ratio ajustados del modelo logit múltiple con IC 95% (escala logarítmica). La línea punteada (OR = 1) indica ausencia de efecto.
Para los coeficientes de la Tabla 17, la Figura 9 y la Salida 2, todos los efectos se leen ceteris paribus. El modelo usa los 1.470 empleados (237 con rotación). El intercepto (β₀ = -2.627; OR = 0.072) no tiene interpretación sustantiva.
- Horas extra (β = 1.439; OR = 4.21; IC 95%: 3.10 a 5.76; p < 0.001). Signo positivo: a igualdad de las demás variables, los odds de rotar de quien hace horas extra son 4.2 veces los de quien no las hace (un aumento de 321%).
- Estado civil. Soltero vs. Casado (β = 0.874; OR = 2.40; IC: 1.72 a 3.36; p < 0.001): los odds de rotar aumentan 139.7%. Divorciado vs. Casado (β = -0.291; OR = 0.748; IC: 0.471 a 1.165; p = 0.207): no significativo, el IC incluye 1.
- Viaje de negocios. Raramente vs. No viaja (β = 0.709; OR = 2.03; IC: 1.10 a 4.08; p = 0.033): los odds aumentan 103%, aunque el límite inferior del IC (1.10) está cerca de 1, por lo que es una significancia marginal. Frecuentemente vs. No viaja (β = 1.347; OR = 3.85; IC: 1.99 a 8.01; p < 0.001): los odds aumentan 284.6%. Los IC son amplios porque la categoría de referencia es pequeña (150 empleados).
- Ingreso mensual, por cada 1.000 unidades (β = -0.103; OR = 0.903; IC: 0.860 a 0.943; p < 0.001). Signo negativo: cada 1.000 unidades más de ingreso reducen los odds de rotar en 9.75%; para una diferencia de 5.000 unidades, el OR es 0.60 (una reducción de 40.1%).
- Antigüedad en el cargo, por año (β = -0.108; OR = 0.898; IC: 0.850 a 0.946; p < 0.001). Cada año adicional reduce los odds en 10.2%; para 5 años, el OR es 0.58 (una reducción de 41.8%).
- Distancia a casa, por unidad (β = 0.033; OR = 1.033; IC: 1.014 a 1.052; p < 0.001). Cada unidad adicional aumenta los odds en 3.3%; para 10 unidades, el OR es 1.39 (un aumento de 38.7%).
Todos los coeficientes son significativos al 5% excepto el de Divorciado. Por magnitud del estadístico z, los efectos más claros son Horas extra (9.10), Soltero (5.10), Ingreso (-4.39), Antigüedad en el cargo (-3.96), Viaje frecuente (3.81) y Distancia (3.51); Viaje raramente (2.14) es el más débil. Los OR por unidad de variables con escalas distintas no son directamente comparables entre sí.
Comparado con las hipótesis de la Tabla 2, los signos estimados coinciden con los esperados en H1 a H6, y las seis hipótesis se mantienen al controlar por las demás variables. Destaca la fuerza del efecto de las horas extra (OR = 4.21), que en el modelo múltiple es el más alto entre los coeficientes ligados a condiciones laborales.
4.4 Significancia global del modelo
Se contrasta \(H_0: \beta_1 = \beta_2 = \dots = \beta_8 = 0\) frente a \(H_a\): algún \(\beta_j \neq 0\), con la prueba de razón de verosimilitud: la diferencia entre la devianza del modelo nulo y la del modelo ajustado se distribuye \(\chi^2\) con grados de libertad igual a la diferencia de parámetros.
lr_global <- with(modelo_multi, null.deviance - deviance)
gl_global <- with(modelo_multi, df.null - df.residual)
p_global <- pchisq(lr_global, gl_global, lower.tail = FALSE)
n_obs <- nobs(modelo_multi)
mcfadden <- 1 - modelo_multi$deviance / modelo_multi$null.deviance
cox_snell <- 1 - exp((modelo_multi$deviance - modelo_multi$null.deviance) / n_obs)
nagelkerke <- cox_snell / (1 - exp(-modelo_multi$null.deviance / n_obs))
indicadores <- tibble(
Indicador = c("Observaciones utilizadas",
"Devianza nula (gl)", "Devianza residual (gl)",
"Estadístico LR: χ² (gl)", "p-valor de la prueba LR",
"AIC", "BIC", "R² de McFadden", "R² de Nagelkerke"),
Valor = c(
as.character(n_obs),
paste0(f_num(modelo_multi$null.deviance), " (", modelo_multi$df.null, ")"),
paste0(f_num(modelo_multi$deviance), " (", modelo_multi$df.residual, ")"),
paste0(f_num(lr_global), " (", gl_global, ")"),
fmt_p(p_global),
f_num(AIC(modelo_multi)), f_num(BIC(modelo_multi)),
f_num(mcfadden, 3), f_num(nagelkerke, 3)
)
)
estilo_tabla(indicadores,
caption = "Tabla 18. Significancia global e indicadores de ajuste del modelo logit múltiple.",
col.names = c("Indicador", "Valor"))| Indicador | Valor |
|---|---|
| Observaciones utilizadas | 1470 |
| Devianza nula (gl) | 1298.58 (1469) |
| Devianza residual (gl) | 1078.84 (1461) |
| Estadístico LR: χ² (gl) | 219.75 (8) |
| p-valor de la prueba LR | <0.001 |
| AIC | 1096.84 |
| BIC | 1144.47 |
| R² de McFadden | 0.169 |
| R² de Nagelkerke | 0.237 |
La prueba de razón de verosimilitud (χ² = 219.75 con 8 gl; p < 0.001) rechaza H₀: el modelo con las seis variables se ajusta significativamente mejor que el modelo nulo, es decir, al menos un coeficiente es distinto de cero. La inclusión de las covariables reduce la devianza de 1298.58 a 1078.84 (16.9%). El AIC (1096.84) y el BIC (1144.47) solo sirven para comparar con modelos alternativos. Los R² de McFadden (0.169) y de Nagelkerke (0.237) indican un ajuste moderado, esperable en un fenómeno multifactorial como la rotación, donde influyen muchas variables no observadas; no son comparables con el R² de una regresión lineal. Dado que en los modelos logit la bondad de ajuste pasa a un segundo plano frente al signo y la significancia de los coeficientes.
4.5 Significancia de cada variable (prueba de razón de verosimilitud)
Para variables categóricas con más de dos niveles (estado civil y
viaje de negocios) hay más de un coeficiente, por lo que se evalúa el
aporte global de cada variable al modelo comparando el modelo completo
con el modelo sin esa variable (drop1).
lr_por_variable <- drop1(modelo_multi, test = "Chisq") %>%
as.data.frame() %>%
rownames_to_column("variable") %>%
filter(variable != "<none>") %>%
transmute(Variable = unname(etiquetas[variable]), gl = Df, LRT = LRT,
p = fmt_p(`Pr(>Chi)`), Sig = signif_stars(`Pr(>Chi)`))
estilo_tabla(lr_por_variable, digits = 3,
caption = "Tabla 19. Aporte de cada variable al modelo múltiple: prueba de razón de verosimilitud al excluir la variable (Sig.: *** p<0.001, ** p<0.01, * p<0.05, . p<0.1).",
col.names = c("Variable", "gl", "χ² (LR)", "p-valor", "Sig."))| Variable | gl | χ² (LR) | p-valor | Sig. |
|---|---|---|---|---|
| Horas extra | 1 | 83.93 | <0.001 | *** |
| Estado civil | 2 | 38.99 | <0.001 | *** |
| Viaje de negocios | 2 | 20.32 | <0.001 | *** |
| Ingreso mensual (miles) | 1 | 23.09 | <0.001 | *** |
| Antigüedad en el cargo | 1 | 16.94 | <0.001 | *** |
| Distancia a casa | 1 | 12.09 | <0.001 | *** |
Una vez controlado el efecto de las demás variables, las seis aportan de forma significativa al modelo (p < 0.001 en todos los casos). Por el estadístico LR el aporte más grande es el de Horas extra (83.93; 1 gl), seguido de Estado civil (38.99; 2 gl), Ingreso mensual (23.09), Viaje de negocios (20.32; 2 gl), Antigüedad en el cargo (16.94) y Distancia a casa (12.09). Frente al análisis bivariado (Tabla 13), el aporte de Horas extra se mantiene (81.4 a 83.9), mientras que el de Ingreso (45.5 a 23.1) y el de Antigüedad en el cargo (42.7 a 16.9) se reduce a la mitad o menos, señal de que comparten información (ρ = 0.39, Tabla 3). Se concluye que las seis variables deben conservarse en el modelo.
4.6 Efecto crudo frente a efecto ajustado
La Tabla 20 compara el OR del análisis bivariado (crudo, Sección 3.3) con el OR del modelo múltiple (ajustado por las demás variables). Un cambio marcado en magnitud, o un cambio de signo, sugiere que la relación observada en el análisis bivariado estaba parcialmente explicada por otras variables (confusión).
comparacion_or <- coef_multi %>%
select(term, beta_aj = estimate, OR_aj = OR, p_aj = p.value) %>%
left_join(coef_simples %>% select(term, beta_cr = estimate, OR_cr = OR, p_cr = p.value),
by = "term") %>%
mutate(cambio_signo = ifelse(sign(beta_aj) != sign(beta_cr), "Sí", "No"),
cambio_OR = round((OR_aj - OR_cr) / OR_cr * 100, 1)) %>%
transmute(Termino = etq_termino(term), OR_cr, p_cr = fmt_p(p_cr), OR_aj, p_aj = fmt_p(p_aj), cambio_OR, cambio_signo)
estilo_tabla(comparacion_or, digits = 3,
caption = "Tabla 20. Odds ratio crudos (modelos simples) frente a odds ratio ajustados (modelo múltiple).",
col.names = c("Término", "OR crudo", "p-valor (crudo)", "OR ajustado",
"p-valor (ajustado)", "Cambio del OR (%)", "¿Cambia el signo?"))| Término | OR crudo | p-valor (crudo) | OR ajustado | p-valor (ajustado) | Cambio del OR (%) | ¿Cambia el signo? |
|---|---|---|---|---|---|---|
| Horas extra: Sí (vs. No) | 3.771 | <0.001 | 4.214 | <0.001 | 11.8 | No |
| Estado civil: Divorciado (vs. Casado) | 0.787 | 0.271 | 0.748 | 0.207 | -5.0 | No |
| Estado civil: Soltero (vs. Casado) | 2.404 | <0.001 | 2.397 | <0.001 | -0.3 | No |
| Viaje: Raramente (vs. no viaja) | 2.023 | 0.025 | 2.031 | 0.033 | 0.4 | No |
| Viaje: Frecuentemente (vs. no viaja) | 3.815 | <0.001 | 3.846 | <0.001 | 0.8 | No |
| Ingreso mensual (por cada 1.000) | 0.881 | <0.001 | 0.903 | <0.001 | 2.5 | No |
| Antigüedad en el cargo (por año) | 0.864 | <0.001 | 0.898 | <0.001 | 3.9 | No |
| Distancia a casa (por unidad) | 1.025 | 0.003 | 1.033 | <0.001 | 0.8 | No |
Ningún efecto cambia de signo al ajustar por las demás variables, y la mayoría de los OR casi no se mueve: Soltero (-0.3%), Viaje raramente (+0.4%), Viaje frecuentemente (+0.8%) y Distancia (+0.8%). Esto indica poca confusión entre las variables del modelo y respalda la robustez de los resultados bivariados. Hay dos matices. Primero, el efecto de las horas extra se amplía (OR de 3.77 a 4.21; +11.8%), es decir, al comparar empleados con las mismas características en las demás variables, la diferencia asociada a las horas extra es aún mayor. Segundo, los efectos de Ingreso (0.881 a 0.903) y de Antigüedad en el cargo (0.864 a 0.898) se atenúan hacia 1 (+2.5% y +3.9%), coherente con que ambas comparten información sobre la etapa de carrera, aunque siguen siendo significativos. Divorciado continúa sin ser significativo (p = 0.271 a 0.207) y Viaje raramente se mantiene significativo (p = 0.025 a 0.033). La distancia gana significancia (p = 0.003 a < 0.001).
4.7 Verificación de supuestos y diagnósticos
El modelo logit requiere:
- Variable respuesta binaria (Bernoulli).
- Observaciones independientes.
- Ausencia o baja multicolinealidad.
- Relación lineal entre las covariables cuantitativas y el logit. Se verifican como sigue.
(1) y (2) Respuesta binaria e independencia. La respuesta es binaria por construcción. La base tiene un registro por empleado, por lo que se asume independencia entre observaciones; esta condición no puede contrastarse con los datos disponibles.
(3) Eventos por variable y multicolinealidad.
n_param <- length(coef(modelo_multi)) - 1
epv <- min(sum(datos$y == 1), sum(datos$y == 0)) / n_param
vif_res <- car::vif(modelo_multi)
tabla_vif <- if (is.matrix(vif_res)) {
as.data.frame(vif_res) %>% rownames_to_column("Variable")
} else {
tibble(Variable = names(vif_res), GVIF = as.numeric(vif_res), Df = 1,
ajustado = sqrt(as.numeric(vif_res)))
}
names(tabla_vif) <- c("Variable", "GVIF", "gl", "GVIF_ajustado")
tabla_vif <- tabla_vif %>% mutate(Variable = unname(etiquetas[Variable]))
estilo_tabla(tabla_vif, digits = 3,
caption = "Tabla 21. Factor de inflación de la varianza generalizado (GVIF) del modelo múltiple. GVIF ajustado = GVIF^(1/(2·gl)); valores por debajo de 2 indican multicolinealidad baja.",
col.names = c("Variable", "GVIF", "gl", "GVIF ajustado"))| Variable | GVIF | gl | GVIF ajustado |
|---|---|---|---|
| Horas extra | 1.029 | 1 | 1.014 |
| Estado civil | 1.027 | 2 | 1.007 |
| Viaje de negocios | 1.013 | 2 | 1.003 |
| Ingreso mensual (miles) | 1.099 | 1 | 1.048 |
| Antigüedad en el cargo | 1.091 | 1 | 1.045 |
| Distancia a casa | 1.020 | 1 | 1.010 |
El modelo tiene 8 parámetros (sin el intercepto) y 237 eventos de la clase menos frecuente, lo que da 29.6 eventos por parámetro; la regla práctica recomienda al menos 10.
En la Tabla 21, el GVIF ajustado de las seis variables está entre 1.003 y 1.048 (el mayor es el de Ingreso mensual), muy por debajo del valor de referencia de 2. No hay multicolinealidad: las estimaciones y sus errores estándar no están distorsionados por redundancia entre covariables, lo que concuerda con las correlaciones bajas a moderadas de la Tabla 3. Además, con 29.6 eventos por parámetro, el tamaño de muestra es suficiente para estimar el modelo con 8 parámetros.
(4) Linealidad en el logit. Se revisó gráficamente en la Figura 8 (Sección 3.2), que la Distancia a casa sigue de cerca el ajuste logístico; la Antigüedad en el cargo lo sigue en general, con un riesgo algo mayor en el primer año y un aplanamiento en antigüedades altas; y el Ingreso mensual muestra la desviación más clara, en los salarios más bajos. Como la inspección visual es subjetiva, la Sección 4.8 contrasta formalmente la linealidad de las tres variables y compara especificaciones alternativas.
La distancia de Cook identifica observaciones con un peso desproporcionado en la estimación. Se usa como umbral de referencia \(4/n\).
cooks <- cooks.distance(modelo_multi)
umbral_cook <- 4 / length(cooks)
n_influyentes <- sum(cooks > umbral_cook)
tibble(obs = seq_along(cooks), cook = cooks) %>%
ggplot(aes(x = obs, y = cook)) +
geom_segment(aes(xend = obs, yend = 0), color = col_azul, alpha = 0.6) +
geom_hline(yintercept = umbral_cook, linetype = "dashed", color = col_naranja) +
labs(x = "Observación", y = "Distancia de Cook")
Figura 10. Distancia de Cook de cada observación en el modelo múltiple (línea punteada: umbral de referencia 4/n).
Observaciones con distancia de Cook superior a \(4/n\) (0.0027): 104 de 1470; el valor máximo es 0.0212.
Para la Figura 10, aunque 104 observaciones (7.1%) superan el umbral de referencia 4/n, el valor máximo de la distancia de Cook (0.0212) es muy pequeño en términos absolutos (los valores cercanos a 1 son los que indican influencia preocupante). Superar 4/n en una proporción de este orden es habitual con muestras grandes y respuesta binaria. El gráfico muestra un patrón homogéneo, sin observaciones que dominen la estimación, por lo que no hay razón para excluir registros.
4.8 Sensibilidad a la forma funcional de las covariables cuantitativas
La Figura 8 sugirió que la relación de algunas variables cuantitativas con el logit podría no ser lineal (en particular la del ingreso: riesgo alto en los salarios más bajos que luego se aplana). Para contrastarlo formalmente se compara, para cada variable cuantitativa, el modelo lineal con una versión en la que esa variable se modela con un spline natural de 3 grados de libertad (2 parámetros adicionales), dejando las demás covariables como están. Como el modelo lineal está anidado en el spline, se usa la prueba de razón de verosimilitud; un \(\Delta\)AIC negativo indica que la versión flexible ajusta mejor incluso penalizando su complejidad.
Se considera que una variable presenta evidencia de no linealidad relevante si la prueba LR es significativa al 5% y el \(\Delta\)AIC es menor que -4.
vars_lin <- c("ingreso_miles", "antiguedad_cargo", "distancia_casa")
terminos_base <- c("horas_extra", "estado_civil", "viaje_negocios", vars_lin)
prueba_lin <- map_dfr(vars_lin, function(v) {
terminos <- ifelse(terminos_base == v, paste0("splines::ns(", v, ", df = 3)"), terminos_base)
m_ns <- glm(reformulate(terminos, response = "y"),
family = binomial(link = "logit"), data = datos)
a <- anova(modelo_multi, m_ns, test = "Chisq")
tibble(var = v, Variable = etiquetas[[v]], LR = a$Deviance[2], gl = a$Df[2],
p = a$`Pr(>Chi)`[2], dAIC = AIC(m_ns) - AIC(modelo_multi))
})
estilo_tabla(prueba_lin %>% transmute(Variable, LR, gl, p = fmt_p(p), dAIC),
digits = 3,
caption = "Tabla 22. Contraste de linealidad en el logit de cada variable cuantitativa: modelo múltiple lineal frente al mismo modelo con esa variable modelada con un spline natural de 3 gl (ΔAIC = AIC del spline - AIC del modelo lineal).",
col.names = c("Variable", "χ² (LR)", "gl", "p-valor", "ΔAIC"))| Variable | χ² (LR) | gl | p-valor | ΔAIC |
|---|---|---|---|---|
| Ingreso mensual (miles) | 26.882 | 2 | <0.001 | -22.882 |
| Antigüedad en el cargo | 11.936 | 2 | 0.003 | -7.936 |
| Distancia a casa | 1.046 | 2 | 0.593 | 2.954 |
El spline usa parámetros adicionales y sus coeficientes no se leen como OR por unidad. Por eso se evalúan también transformaciones logarítmicas, que no agregan parámetros: el logaritmo del ingreso (comprime los valores altos, de modo que un mismo aumento absoluto pesa más en salarios bajos) y el logaritmo de la antigüedad en el cargo más 1 (el +1 permite el valor 0). La Tabla 23 compara las especificaciones y la Tabla 24 muestra si los coeficientes de las demás variables se ven afectados por la elección.
m_log <- update(modelo_multi, . ~ . - ingreso_miles + log(ingreso_mensual))
m_log2 <- update(modelo_multi, . ~ . - ingreso_miles - antiguedad_cargo +
log(ingreso_mensual) + log1p(antiguedad_cargo))
m_ns_ing <- update(modelo_multi, . ~ . - ingreso_miles + splines::ns(ingreso_miles, df = 3))
lista_mod <- list("Lineal (modelo principal)" = modelo_multi,
"Logaritmo del ingreso" = m_log,
"Logaritmo del ingreso y de (antigüedad en el cargo + 1)" = m_log2,
"Spline natural del ingreso (3 gl)" = m_ns_ing)
comp_modelos <- tibble(
Modelo = names(lista_mod),
Parametros = map_int(lista_mod, ~ length(coef(.x))),
Devianza = map_dbl(lista_mod, deviance),
AIC = map_dbl(lista_mod, AIC),
BIC = map_dbl(lista_mod, BIC)
) %>%
mutate(dAIC = AIC - AIC[1])
estilo_tabla(comp_modelos, digits = 2,
caption = "Tabla 23. Comparación de especificaciones del modelo múltiple (ΔAIC frente al modelo lineal; valores menores indican mejor ajuste penalizado).",
col.names = c("Especificación", "Parámetros", "Devianza", "AIC", "BIC", "ΔAIC"))| Especificación | Parámetros | Devianza | AIC | BIC | ΔAIC |
|---|---|---|---|---|---|
| Lineal (modelo principal) | 9 | 1079 | 1097 | 1144 | 0.00 |
| Logaritmo del ingreso | 9 | 1068 | 1086 | 1133 | -11.33 |
| Logaritmo del ingreso y de (antigüedad en el cargo + 1) | 9 | 1060 | 1078 | 1126 | -18.58 |
| Spline natural del ingreso (3 gl) | 11 | 1052 | 1074 | 1132 | -22.88 |
# Robustez de los coeficientes de las demás variables (OR) frente a la especificación
terminos_com <- c("horas_extraSi", "estado_civilDivorciado", "estado_civilSoltero",
"viaje_negociosRaramente", "viaje_negociosFrecuentemente", "distancia_casa")
or_de <- function(m) unname(exp(coef(m)[terminos_com]))
robustez <- tibble(
Termino = etq_termino(terminos_com),
Lineal = or_de(modelo_multi),
Log_ingreso = or_de(m_log),
Log_ing_antig = or_de(m_log2),
Spline_ingreso = or_de(m_ns_ing)
)
estilo_tabla(robustez, digits = 3,
caption = "Tabla 24. Odds ratio de las variables que no cambian de forma funcional, según la especificación del ingreso y la antigüedad en el cargo.",
col.names = c("Término", "Lineal (principal)", "Log(ingreso)",
"Log(ingreso) y log(antigüedad + 1)", "Spline del ingreso"))| Término | Lineal (principal) | Log(ingreso) | Log(ingreso) y log(antigüedad + 1) | Spline del ingreso |
|---|---|---|---|---|
| Horas extra: Sí (vs. No) | 4.214 | 4.325 | 4.324 | 4.368 |
| Estado civil: Divorciado (vs. Casado) | 0.748 | 0.749 | 0.753 | 0.747 |
| Estado civil: Soltero (vs. Casado) | 2.397 | 2.385 | 2.380 | 2.382 |
| Viaje: Raramente (vs. no viaja) | 2.031 | 2.067 | 2.089 | 2.185 |
| Viaje: Frecuentemente (vs. no viaja) | 3.846 | 3.920 | 4.064 | 4.207 |
| Distancia a casa (por unidad) | 1.033 | 1.035 | 1.035 | 1.035 |
# Valores usados en el texto y en la decisión
lin_ing <- prueba_lin %>% filter(var == "ingreso_miles")
d_aic_log <- AIC(m_log) - AIC(modelo_multi)
d_aic_log2 <- AIC(m_log2) - AIC(modelo_multi)
d_aic_ns <- AIC(m_ns_ing) - AIC(modelo_multi)
max_cambio <- 100 * max(abs(as.matrix(robustez[, -1]) / robustez$Lineal - 1))
or_ing10 <- exp(unname(coef(m_log2)["log(ingreso_mensual)"]) * log(1.1))
fmt_p_txt <- function(p) ifelse(p < 0.001, "p < 0.001", paste0("p = ", formatC(p, format = "f", digits = 3)))
res_lin <- prueba_lin %>%
mutate(texto = paste0(Variable, ": ", ifelse(p < 0.05, "se rechaza la linealidad", "no se rechaza la linealidad"),
" (χ² = ", f_num(LR), ", 2 gl, ", fmt_p_txt(p), "; ΔAIC = ", f_num(dAIC), ")"))
vars_mejora <- prueba_lin %>% filter(p < 0.05, dAIC < -4) %>% pull(Variable)
vars_no_lin <- prueba_lin %>% filter(p < 0.05) %>% pull(Variable)
txt_lin <- if (length(vars_no_lin) > 0) {
paste0("rechaza la linealidad en el logit de ", paste(vars_no_lin, collapse = " y "),
", por lo que se conservan especificaciones alternativas con logaritmos para compararlas en el punto 5")
} else {
"no rechaza la linealidad en el logit de ninguna de las tres variables"
}
decision <- if (length(vars_mejora) > 0) {
paste0("El criterio se cumple para ", paste(vars_mejora, collapse = " y "),
". La especificación lineal se conserva como modelo principal de interpretación, porque es la que corresponde a las hipótesis de la Sección 1 y permite leer los OR por unidad. Las especificaciones con logaritmos, que tienen el mismo número de parámetros, se conservan como modelos alternativos y se compararán con el lineal por AUC en la muestra de prueba y en validación cruzada (Sección 5.8); si la ganancia predictiva fuera material, según el criterio definido en esa sección, se adoptarían como modelo final. Como los coeficientes de las demás variables cambian como máximo ",
f_num(max_cambio, 1), "% (Tabla 24), las conclusiones sobre horas extra, estado civil, viaje de negocios y distancia no dependen de esta elección.")
} else {
"no se cumple el criterio para ninguna variable; se conserva la especificación lineal como modelo final."
}# Perfil base: se modifica una característica a la vez en esta sección y en la siguiente
perfil_base <- tibble(
horas_extra = factor("No", levels = levels(datos$horas_extra)),
estado_civil = factor("Casado", levels = levels(datos$estado_civil)),
viaje_negocios = factor("Raramente", levels = levels(datos$viaje_negocios)),
ingreso_miles = median(datos$ingreso_miles),
antiguedad_cargo = median(datos$antiguedad_cargo),
distancia_casa = median(datos$distancia_casa)
)
grid_ing <- perfil_base %>%
select(-ingreso_miles) %>%
slice(rep(1, 120)) %>%
mutate(ingreso_miles = seq(min(datos$ingreso_miles), max(datos$ingreso_miles), length.out = 120),
ingreso_mensual = ingreso_miles * 1000)
curvas <- tibble(
ingreso_miles = grid_ing$ingreso_miles,
`Lineal (modelo principal)` = predict(modelo_multi, newdata = grid_ing, type = "response"),
`Logaritmo del ingreso` = predict(m_log, newdata = grid_ing, type = "response"),
`Spline natural (3 gl)` = predict(m_ns_ing, newdata = grid_ing, type = "response")
) %>%
pivot_longer(-ingreso_miles, names_to = "Especificacion", values_to = "prob")
ggplot(curvas, aes(x = ingreso_miles, y = prob, color = Especificacion, linetype = Especificacion)) +
geom_line(linewidth = 0.9) +
geom_rug(data = datos, aes(x = ingreso_miles), inherit.aes = FALSE, sides = "b", alpha = 0.08) +
scale_color_manual(values = c("Lineal (modelo principal)" = col_naranja,
"Logaritmo del ingreso" = "gray40",
"Spline natural (3 gl)" = col_azul)) +
scale_linetype_manual(values = c("Lineal (modelo principal)" = "solid",
"Logaritmo del ingreso" = "dotdash",
"Spline natural (3 gl)" = "dashed")) +
scale_y_continuous(labels = percent) +
labs(x = "Ingreso mensual (miles)", y = "Probabilidad de rotar", color = NULL, linetype = NULL) +
theme(legend.position = "bottom")
Figura 11. Probabilidad estimada de rotar según el ingreso mensual para un perfil base (sin horas extra, casado, viaja raramente, antigüedad en el cargo y distancia en su mediana), con tres especificaciones del ingreso. Las marcas del eje horizontal muestran dónde hay observaciones.
Interpretación de las Tablas 22 a 24 y de la Figura 11:
Contraste de linealidad por variable: Ingreso mensual (miles): se rechaza la linealidad (χ² = 26.88, 2 gl, p < 0.001; ΔAIC = -22.88); Antigüedad en el cargo: se rechaza la linealidad (χ² = 11.94, 2 gl, p = 0.003; ΔAIC = -7.94); Distancia a casa: no se rechaza la linealidad (χ² = 1.05, 2 gl, p = 0.593; ΔAIC = 2.95). Para el ingreso, la especificación con logaritmo tiene un \(\Delta\)AIC de -11.33 y la del spline un \(\Delta\)AIC de -22.88; la que usa logaritmos del ingreso y de la antigüedad en el cargo, sin parámetros adicionales respecto al modelo lineal, tiene un \(\Delta\)AIC de -18.58 (Tabla 23). Los coeficientes de las demás variables cambian como máximo 9.4% entre especificaciones (Tabla 24). Con esta última especificación, un aumento de 10% en el ingreso multiplica los odds de rotar por 0.930. La Figura 11 permite ver dónde difieren las curvas: las marcas del eje horizontal se concentran en los ingresos bajos y medios.
Decisión: El criterio se cumple para Ingreso mensual (miles) y Antigüedad en el cargo. La especificación lineal se conserva como modelo principal de interpretación, porque es la que corresponde a las hipótesis de la Sección 1 y permite leer los OR por unidad. Las especificaciones con logaritmos, que tienen el mismo número de parámetros, se conservan como modelos alternativos y se compararán con el lineal por AUC en la muestra de prueba y en validación cruzada (Sección 5.8); si la ganancia predictiva fuera material, según el criterio definido en esa sección, se adoptarían como modelo final. Como los coeficientes de las demás variables cambian como máximo 9.4% (Tabla 24), las conclusiones sobre horas extra, estado civil, viaje de negocios y distancia no dependen de esta elección.
4.9 Lectura del modelo en probabilidades
Los OR describen efectos multiplicativos sobre los odds, pero para la
toma de decisiones suele ser más intuitivo ver cuánto cambia la
probabilidad. Se toma el perfil base de la sección anterior (sin horas
extra, casado, viaja raramente y con ingreso, antigüedad en el cargo y
distancia iguales a su mediana) y se modifica una sola característica a
la vez, dejando las demás fijas. Esta lectura es además el puente hacia
donde se estimará la probabilidad de rotación de un empleado hipotético
con el mismo mecanismo (predict con
newdata).
q <- function(x, p) unname(quantile(x, p))
nivel <- function(var, valor) factor(valor, levels = levels(datos[[var]]))
lista_esc <- list(
base = perfil_base,
he = mutate(perfil_base, horas_extra = nivel("horas_extra", "Si")),
sol = mutate(perfil_base, estado_civil = nivel("estado_civil", "Soltero")),
nov = mutate(perfil_base, viaje_negocios = nivel("viaje_negocios", nivel_viaje_ref)),
fre = mutate(perfil_base, viaje_negocios = nivel("viaje_negocios", "Frecuentemente")),
ing25 = mutate(perfil_base, ingreso_miles = q(datos$ingreso_miles, 0.25)),
ing75 = mutate(perfil_base, ingreso_miles = q(datos$ingreso_miles, 0.75)),
ant25 = mutate(perfil_base, antiguedad_cargo = q(datos$antiguedad_cargo, 0.25)),
ant75 = mutate(perfil_base, antiguedad_cargo = q(datos$antiguedad_cargo, 0.75)),
dis25 = mutate(perfil_base, distancia_casa = q(datos$distancia_casa, 0.25)),
dis75 = mutate(perfil_base, distancia_casa = q(datos$distancia_casa, 0.75))
)
etiq_esc <- c(
base = "Perfil base",
he = "Horas extra: Sí",
sol = "Estado civil: Soltero",
nov = "Viaje de negocios: No viaja",
fre = "Viaje de negocios: Frecuentemente",
ing25 = paste0("Ingreso en el percentil 25 (", f_num(q(datos$ingreso_miles, 0.25), 1), " mil)"),
ing75 = paste0("Ingreso en el percentil 75 (", f_num(q(datos$ingreso_miles, 0.75), 1), " mil)"),
ant25 = paste0("Antigüedad en el cargo en el percentil 25 (", q(datos$antiguedad_cargo, 0.25), " años)"),
ant75 = paste0("Antigüedad en el cargo en el percentil 75 (", q(datos$antiguedad_cargo, 0.75), " años)"),
dis25 = paste0("Distancia a casa en el percentil 25 (", q(datos$distancia_casa, 0.25), ")"),
dis75 = paste0("Distancia a casa en el percentil 75 (", q(datos$distancia_casa, 0.75), ")")
)
prob_esc <- tibble(
clave = names(lista_esc),
Escenario = unname(etiq_esc[names(lista_esc)]),
Prob = map_dbl(lista_esc, ~ unname(predict(modelo_multi, newdata = .x, type = "response")))
) %>%
mutate(Cambio_pp = (Prob - Prob[clave == "base"]) * 100)
pr <- setNames(prob_esc$Prob, prob_esc$clave)
pct <- function(x) percent(unname(x), accuracy = 0.1)
dpp <- function(a, b) sprintf("%+.1f", 100 * unname(a - b))
estilo_tabla(prob_esc %>% transmute(Escenario, Prob = Prob * 100, Cambio_pp),
digits = 1,
caption = "Tabla 25. Probabilidad estimada de rotar (modelo múltiple) para el perfil base y al modificar una sola característica a la vez (p.p.: puntos porcentuales).",
col.names = c("Escenario", "Probabilidad de rotar (%)", "Cambio frente al perfil base (p.p.)"))| Escenario | Probabilidad de rotar (%) | Cambio frente al perfil base (p.p.) |
|---|---|---|
| Perfil base | 7.5 | 0.0 |
| Horas extra: Sí | 25.4 | 17.9 |
| Estado civil: Soltero | 16.2 | 8.7 |
| Viaje de negocios: No viaja | 3.8 | -3.6 |
| Viaje de negocios: Frecuentemente | 13.2 | 5.8 |
| Ingreso en el percentil 25 (2.9 mil) | 9.0 | 1.6 |
| Ingreso en el percentil 75 (8.4 mil) | 5.4 | -2.1 |
| Antigüedad en el cargo en el percentil 25 (2 años) | 8.2 | 0.8 |
| Antigüedad en el cargo en el percentil 75 (7 años) | 5.0 | -2.5 |
| Distancia a casa en el percentil 25 (2) | 6.4 | -1.1 |
| Distancia a casa en el percentil 75 (14) | 9.2 | 1.7 |
prob_esc %>%
ggplot(aes(x = Prob, y = reorder(Escenario, Prob), fill = clave == "base")) +
geom_col(width = 0.7) +
geom_vline(xintercept = unname(pr["base"]), linetype = "dashed", color = "gray30") +
geom_text(aes(label = percent(Prob, accuracy = 0.1)), hjust = -0.1, size = 3) +
scale_fill_manual(values = c("TRUE" = "gray60", "FALSE" = col_naranja), guide = "none") +
scale_x_continuous(labels = percent, expand = expansion(mult = c(0, 0.15))) +
labs(x = "Probabilidad estimada de rotar", y = NULL)
Figura 12. Probabilidad estimada de rotar por escenario (el perfil base se muestra en gris; la línea punteada marca su probabilidad).
Al analizar la Tabla 25 y la Figura 12, el perfil base tiene una probabilidad estimada de rotar de 7.5%. Manteniendo todo lo demás igual:
- Horas extra: la probabilidad pasa a 25.4% (+17.9 p.p.).
- Estado civil: ser soltero la lleva a 16.2% (+8.7 p.p.).
- Viaje de negocios: viajar frecuentemente la lleva a 13.2% (+5.8 p.p.), y no viajar la reduce a 3.8% (-3.6 p.p.).
- Ingreso: del percentil 25 al 75, la probabilidad va de 9.0% a 5.4%.
- Antigüedad en el cargo: del percentil 25 al 75, va de 8.2% a 5.0%.
- Distancia a casa: del percentil 25 al 75, va de 6.4% a 9.2%.
Las diferencias en puntos porcentuales dependen del perfil de referencia: un mismo OR produce cambios pequeños cuando la probabilidad base es baja y mayores cuando se acerca a 0.5. Por eso se reportan para un perfil concreto y no como un efecto único.
4.10 Síntesis de la estimación
El modelo múltiple es globalmente significativo (LR χ² = 219.75 con 8 gl; p < 0.001) y, al controlar por las demás variables, las seis variables aportan de forma significativa (p < 0.001 en la prueba LR); la única categoría sin efecto significativo es Divorciado frente a Casado. Aumentan los odds de rotar: las horas extra (OR = 4.21), viajar frecuentemente (3.85) o raramente (2.03) frente a no viajar, ser soltero frente a casado (2.40) y cada unidad adicional de distancia a casa (1.033). Disminuyen los odds: cada 1.000 unidades de ingreso mensual (0.903) y cada año de antigüedad en el cargo (0.898). Los signos coinciden con las seis hipótesis de la Sección 1 y los efectos ajustados son muy parecidos a los crudos, lo que indica asociaciones robustas. Deben tenerse presentes cuatro limitaciones:
- El ajuste global es moderado (McFadden 0.169; Nagelkerke 0.237), normal en un fenómeno multifactorial.
- Los IC de Viaje de negocios son amplios porque la categoría de referencia solo tiene 150 empleados (12 rotaciones).
- El contraste formal de la Sección 4.8 rechaza la linealidad en el logit de Ingreso mensual (miles) y Antigüedad en el cargo, por lo que se conservan especificaciones alternativas con logaritmos para compararlas en el punto 5.
- Los resultados describen asociaciones y no causalidad. La capacidad predictiva del modelo se evalúa con la curva ROC y AUC, y de ella dependerán la definición del punto de corte y las conclusiones y la estrategia de retención.
# Modelo que se usará en los puntos 5 a 7 (evaluación, predicciones y conclusiones).
# Por ahora es el modelo múltiple de la Sección 4.2; si se decidiera adoptar otra especificación
# (Sección 4.8), basta con cambiar esta línea.
modelo_final <- modelo_multi
formula_final <- formula(modelo_final)
# Especificación alternativa (logaritmos, mismo número de parámetros) que se comparará por AUC en el punto 5
modelo_alt <- m_log2
formula_final## y ~ horas_extra + estado_civil + viaje_negocios + ingreso_miles +
## antiguedad_cargo + distancia_casa5 Evaluación del modelo
El poder predictivo del modelo se evalúa sobre observaciones que no intervienen en su estimación. Para ello se divide la base en una muestra de entrenamiento, con la que se reestima el modelo, y una muestra de prueba, sobre la que se calculan la matriz de confusión, la curva ROC y el AUC. El modelo interpretado en la Sección 4 (estimado con los 1.470 empleados) se mantiene como modelo de interpretación; la muestra de entrenamiento solo se usa para medir cuánto se degrada su capacidad de clasificación ante datos nuevos.
5.1 Partición de la base en entrenamiento y prueba
La partición es aleatoria y estratificada por la variable respuesta: se toma el 70% de los empleados para entrenamiento y el 30% restante para prueba, conservando en ambas muestras la misma proporción de rotación. Como solo 16.1% de los empleados de la base rotó (237 de 1470, Sección 2.1), un muestreo simple podría dejar en la muestra de prueba, por azar, una tasa de rotación distinta a la de la base, lo que distorsionaría la sensibilidad y la especificidad.
library(caret) # createDataPartition(); más adelante confusionMatrix() y validación cruzada
# La semilla se fija de nuevo aquí: las figuras con geom_jitter() consumen números
# aleatorios, y sin esta línea la partición cambiaría si se modifica alguna figura previa.
set.seed(123)
# createDataPartition() muestrea por separado dentro de cada nivel de 'rotacion'
# (estratificación), de modo que ambas muestras conserven la tasa de rotación.
# Se usa el factor 'rotacion' y no 'y': con una respuesta numérica la función
# estratifica por cuantiles, no por clase.
# list = FALSE devuelve una matriz de índices; as.vector() la convierte en vector,
# porque los tibbles no admiten subconjuntos de filas con una matriz.
idx_train <- as.vector(createDataPartition(datos$rotacion, p = 0.7, list = FALSE))
datos_train <- datos[idx_train, ] # muestra de entrenamiento (70%)
datos_test <- datos[-idx_train, ] # muestra de prueba (30%)
# Verificación: ninguna observación está en ambas muestras y entre las dos suman la base
stopifnot(length(intersect(idx_train, setdiff(seq_len(nrow(datos)), idx_train))) == 0,
nrow(datos_train) + nrow(datos_test) == nrow(datos))
# Tabla comparativa: tamaño y tasa de rotación en la base completa y en cada muestra
tab_particion <- bind_rows(Total = datos, Entrenamiento = datos_train,
Prueba = datos_test, .id = "Muestra") %>%
mutate(Muestra = factor(Muestra, levels = c("Total", "Entrenamiento", "Prueba"))) %>%
group_by(Muestra) %>%
summarise(N = n(), Rotan = sum(y), No_rotan = N - Rotan,
Tasa = mean(y) * 100, .groups = "drop")
estilo_tabla(tab_particion, digits = 1,
caption = "Tabla 26. Tamaño y tasa de rotación de la base completa y de las muestras de entrenamiento y prueba (partición estratificada 70/30).",
col.names = c("Muestra", "N", "Rotan", "No rotan", "Tasa de rotación (%)"))| Muestra | N | Rotan | No rotan | Tasa de rotación (%) |
|---|---|---|---|---|
| Total | 1470 | 237 | 1233 | 16.1 |
| Entrenamiento | 1030 | 166 | 864 | 16.1 |
| Prueba | 440 | 71 | 369 | 16.1 |
La Tabla 26 confirma que la partición conserva la estructura de la base: la muestra de entrenamiento tiene 1030 empleados (166 rotaciones) y la de prueba 440 (71 rotaciones), ambas con una tasa de rotación de 16.1% y 16.1%, respectivamente, igual a la de la base completa. De esto se desprenden dos consecuencias para la evaluación:
Estimación en entrenamiento. Con 166 eventos y 8 parámetros, el modelo reestimado cuenta con 20.8 eventos por parámetro, por encima de la regla práctica de 10, de modo que la reducción del tamaño de muestra no compromete la estabilidad de los coeficientes.
Precisión de la evaluación. Las métricas que dependen de los casos positivos (sensibilidad) se calculan sobre solo 71 rotaciones de la muestra de prueba, por lo que tienen un margen de error apreciable. Para no depender de una única partición, en la Sección 5.6 el AUC se complementa con validación cruzada.
5.2 Reestimación del modelo con la muestra de entrenamiento
El modelo de la Sección 4 se estimó con los 1470 empleados; para evaluarlo sobre datos no vistos se reestima la misma especificación solo con la muestra de entrenamiento. La evaluación de las secciones siguientes corresponde a este modelo reestimado, por lo que primero se verifica que conserve las conclusiones del modelo interpretado. El criterio de estabilidad es la dirección y la magnitud de los efectos, no su significancia: con el 70% de los datos los errores estándar aumentan cerca de 1.2 veces (\(\sqrt{1/0.7}\)), de modo que un coeficiente puede perder significancia sin que el efecto haya cambiado. Se considera que un efecto es inestable si cambia de signo o si su OR en entrenamiento queda fuera del intervalo de confianza del modelo completo.
# update() reutiliza la fórmula y las categorías de referencia de 'modelo_final'
# y solo cambia la base de datos: se garantiza que se estima la misma especificación.
modelo_train <- update(modelo_final, data = datos_train)
# Verificaciones: todas las categorías están presentes en entrenamiento y los
# dos modelos tienen exactamente los mismos parámetros.
stopifnot(all(table(datos_train$viaje_negocios) > 0),
all(table(datos_train$estado_civil) > 0),
all(table(datos_train$horas_extra) > 0),
identical(names(coef(modelo_train)), names(coef(modelo_final))))
# Función auxiliar: OR, IC 95% (verosimilitud perfilada) y p-valor de un modelo,
# con un sufijo en los nombres para poder unir las tablas de los dos modelos.
or_modelo <- function(m, sufijo) {
tidy(m, conf.int = TRUE) %>%
filter(term != "(Intercept)") %>%
transmute(term, OR = exp(estimate), inf = exp(conf.low),
sup = exp(conf.high), p = p.value) %>%
rename_with(~ paste0(.x, "_", sufijo), -term)
}
comp_train <- or_modelo(modelo_final, "c") %>%
left_join(or_modelo(modelo_train, "t"), by = "term") %>%
mutate(cambio = (OR_t / OR_c - 1) * 100,
mismo_signo = sign(log(OR_t)) == sign(log(OR_c)),
dentro_ic = OR_t >= inf_c & OR_t <= sup_c)
estilo_tabla(
comp_train %>%
transmute(Termino = etq_termino(term),
OR_c,
IC_c = paste0(f_num(inf_c, 3), " a ", f_num(sup_c, 3)),
p_c = fmt_p(p_c),
OR_t,
IC_t = paste0(f_num(inf_t, 3), " a ", f_num(sup_t, 3)),
p_t = fmt_p(p_t),
cambio = sprintf("%+.1f", cambio), # un decimal y signo explícito
Estable = ifelse(mismo_signo & dentro_ic, "Sí", "No")),
digits = 3,
caption = "Tabla 27. Odds ratio del modelo estimado con la base completa (Sección 4) frente al mismo modelo reestimado con la muestra de entrenamiento. Cambio del OR (%) = (OR entrenamiento / OR completo - 1) × 100. Estable: el OR de entrenamiento conserva el signo y queda dentro del IC 95% del modelo completo.",
col.names = c("Término", "OR completo", "IC 95% completo", "p-valor completo",
"OR entrenamiento", "IC 95% entrenamiento", "p-valor entrenamiento",
"Cambio del OR (%)", "Estable"))| Término | OR completo | IC 95% completo | p-valor completo | OR entrenamiento | IC 95% entrenamiento | p-valor entrenamiento | Cambio del OR (%) | Estable |
|---|---|---|---|---|---|---|---|---|
| Horas extra: Sí (vs. No) | 4.214 | 3.096 a 5.758 | <0.001 | 4.041 | 2.781 a 5.902 | <0.001 | -4.1 | Sí |
| Estado civil: Divorciado (vs. Casado) | 0.748 | 0.471 a 1.165 | 0.207 | 1.042 | 0.607 a 1.756 | 0.880 | +39.3 | No |
| Estado civil: Soltero (vs. Casado) | 2.397 | 1.716 a 3.363 | <0.001 | 2.917 | 1.936 a 4.438 | <0.001 | +21.7 | Sí |
| Viaje: Raramente (vs. no viaja) | 2.031 | 1.100 a 4.075 | 0.033 | 3.329 | 1.470 a 9.010 | 0.008 | +63.9 | Sí |
| Viaje: Frecuentemente (vs. no viaja) | 3.846 | 1.985 a 8.010 | <0.001 | 7.327 | 3.065 a 20.597 | <0.001 | +90.5 | Sí |
| Ingreso mensual (por cada 1.000) | 0.903 | 0.860 a 0.943 | <0.001 | 0.905 | 0.855 a 0.954 | <0.001 | +0.3 | Sí |
| Antigüedad en el cargo (por año) | 0.898 | 0.850 a 0.946 | <0.001 | 0.884 | 0.827 a 0.942 | <0.001 | -1.5 | Sí |
| Distancia a casa (por unidad) | 1.033 | 1.014 a 1.052 | <0.001 | 1.038 | 1.015 a 1.061 | <0.001 | +0.4 | Sí |
# ---- Valores para el texto ----
n_estables <- sum(comp_train$mismo_signo & comp_train$dentro_ic)
# Vectores con nombre = término del modelo
or_c_v <- setNames(comp_train$OR_c, comp_train$term) # OR modelo completo
or_t_v <- setNames(comp_train$OR_t, comp_train$term) # OR entrenamiento
cb_v <- setNames(comp_train$cambio, comp_train$term) # cambio del OR (%)
p_c_v <- setNames(comp_train$p_c, comp_train$term) # p-valor completo
p_t_v <- setNames(comp_train$p_t, comp_train$term) # p-valor entrenamiento
ic_t_v <- setNames(paste0(f_num(comp_train$inf_t), " a ", f_num(comp_train$sup_t)),
comp_train$term) # IC entrenamiento (texto, 2 decimales)
# Mayor cambio absoluto entre los efectos más estables (horas extra y cuantitativas)
cambio_nucleo <- max(abs(cb_v[c("horas_extraSi", "ingreso_miles",
"antiguedad_cargo", "distancia_casa")]))
# Tamaño, rotaciones y tasa de la categoría de referencia de viaje,
# en la base completa y en entrenamiento
ref_viaje <- function(df) {
df %>% filter(viaje_negocios == nivel_viaje_ref) %>%
summarise(n = n(), rot = sum(y), tasa = mean(y))
}
ref_viaje_c <- ref_viaje(datos)
ref_viaje_tr <- ref_viaje(datos_train)# Orden de los términos: según el OR del modelo completo
orden_terminos <- comp_train %>% arrange(OR_c) %>% pull(term) %>% etq_termino()
bind_rows(`Base completa (Sección 4)` = tidy(modelo_final, conf.int = TRUE),
`Entrenamiento (70%)` = tidy(modelo_train, conf.int = TRUE),
.id = "Modelo") %>%
filter(term != "(Intercept)") %>%
mutate(across(c(estimate, conf.low, conf.high), exp),
etiqueta = factor(etq_termino(term), levels = orden_terminos)) %>%
ggplot(aes(x = estimate, y = etiqueta, color = Modelo)) +
geom_vline(xintercept = 1, linetype = "dashed", color = "gray40") +
# position_dodge() separa verticalmente los dos modelos dentro de cada término
geom_pointrange(aes(xmin = conf.low, xmax = conf.high),
position = position_dodge(width = 0.6), size = 0.5) +
# Cortes explícitos para que la zona OR < 1 también tenga referencias numéricas
scale_x_log10(breaks = c(0.5, 0.75, 1, 1.5, 2, 3, 5, 10, 20),
labels = c("0.5", "0.75", "1", "1.5", "2", "3", "5", "10", "20")) +
scale_color_manual(values = c("Base completa (Sección 4)" = col_naranja,
"Entrenamiento (70%)" = col_azul), name = NULL) +
labs(x = "Odds ratio (escala logarítmica)", y = NULL) +
theme(legend.position = "bottom",
panel.grid.minor = element_blank()) # evita líneas de rejilla sin etiqueta
Figura 13. Odds ratio con IC 95% del modelo estimado con la base completa (naranja) y del reestimado con la muestra de entrenamiento (azul), en escala logarítmica. La línea punteada (OR = 1) indica ausencia de efecto. Los términos se ordenan según el OR del modelo completo; los OR de las variables cuantitativas se expresan por unidad (o por 1.000 en el ingreso), por lo que su cercanía a 1 no indica un efecto pequeño ni es comparable con los de las variables categóricas.
La comparación entre los dos modelos muestra que 7 de los 8 términos conservan el signo y quedan dentro del intervalo de confianza del modelo completo (Tabla 27 y Figura 13). Los efectos de las horas extra, el ingreso mensual, la antigüedad en el cargo y la distancia a casa son prácticamente idénticos en los dos modelos: ninguno de sus OR cambia más de 4.1%. Por su parte, el OR de los solteros frente a los casados aumenta 21.7%, pero se mantiene dentro del intervalo del modelo completo. El único término marcado como inestable es Divorciado frente a Casado: su OR pasa de 0.748 a 1.042, lo que implica un cambio de signo del coeficiente. No obstante, en ambos modelos su intervalo contiene el 1 (p = 0.207 y p = 0.880), de modo que el cambio corresponde a variación aleatoria alrededor de la ausencia de efecto y la conclusión se mantiene: en ninguno de los dos modelos hay evidencia de diferencia entre divorciados y casados.
Los mayores cambios se presentan en Viaje de negocios. El OR de viajar raramente pasa de 2.03 a 3.33 (+63.9%) y el de viajar frecuentemente de 3.85 a 7.33 (+90.5%); ambos quedan cerca del límite superior del intervalo del modelo completo, con intervalos muy amplios en entrenamiento (1.47 a 9.01 y 3.07 a 20.60). El aumento proviene de la categoría de referencia: “No viaja” rota 8.0% en la base completa (12 de 150) y 5.7% en entrenamiento (6 de 106). Como los dos OR de la variable se calculan frente a los odds de ese grupo, una menor rotación en la referencia los eleva a ambos al mismo tiempo, y con tan pocos eventos basta una o dos rotaciones para producir ese cambio. La dirección y el orden del efecto se mantienen (Frecuentemente > Raramente > No viaja), pero su magnitud es imprecisa, lo que confirma la advertencia de la Sección 4.3 sobre la amplitud de estos intervalos.
En este sentido, el modelo reestimado conserva las conclusiones del modelo interpretado en la Sección 4 y puede usarse para evaluar el poder predictivo. Para la estrategia de retención, el efecto de los viajes se interpreta por su dirección y no por la magnitud de sus OR.
5.3 Probabilidades predichas en la muestra de prueba
El modelo reestimado asigna a cada empleado de la muestra de prueba una probabilidad de rotar. Antes de resumir su desempeño en una sola medida, se examina cómo se distribuyen esas probabilidades entre quienes efectivamente rotaron y quienes no: un modelo con capacidad de discriminación asigna probabilidades sistemáticamente más altas al primer grupo. Como medida resumen se usa el coeficiente de discriminación de Tjur, igual a la diferencia entre la probabilidad media predicha de quienes rotan y la de quienes no rotan (0 indica ausencia de discriminación y 1, separación perfecta). Adicionalmente, se compara la probabilidad media predicha con la tasa de rotación observada en la muestra de prueba, como verificación de la calibración global del modelo.
# predict(type = "response") devuelve la probabilidad estimada P(y = 1 | x) en lugar
# del logit. Se usa 'modelo_train', que no vio la muestra de prueba; usar 'modelo_final'
# contaminaría la evaluación, porque se estimó con toda la base.
datos_test <- datos_test %>%
mutate(prob = predict(modelo_train, newdata = datos_test, type = "response"))
# Resumen de las probabilidades predichas por clase observada y en total
resumen_prob <- function(df) {
df %>% summarise(n = n(), Media = mean(prob), DE = sd(prob),
Min = min(prob), Q1 = quantile(prob, 0.25, names = FALSE),
Mediana = median(prob), Q3 = quantile(prob, 0.75, names = FALSE),
Max = max(prob), Sobre_05 = mean(prob > 0.5) * 100)
}
tab_prob <- bind_rows(
datos_test %>% group_by(Rotacion = as.character(rotacion)) %>% resumen_prob(),
datos_test %>% resumen_prob() %>% mutate(Rotacion = "Total")
) %>%
mutate(Rotacion = factor(Rotacion, levels = c("No", "Si", "Total"))) %>%
arrange(Rotacion)
estilo_tabla(tab_prob %>% mutate(across(Media:Max, ~ .x * 100)), digits = 1,
caption = "Tabla 28. Probabilidad estimada de rotar (%) en la muestra de prueba según la rotación observada (modelo reestimado con la muestra de entrenamiento). DE: desviación estándar; > 0.5: porcentaje de empleados con probabilidad estimada superior a 0.5.",
col.names = c("Rotación observada", "n", "Media", "DE", "Mín.", "Q1",
"Mediana", "Q3", "Máx.", "> 0.5 (%)"))| Rotación observada | n | Media | DE | Mín. | Q1 | Mediana | Q3 | Máx. | > 0.5 (%) |
|---|---|---|---|---|---|---|---|---|---|
| No | 369 | 14.2 | 13.6 | 0.3 | 4.0 | 9.5 | 20.3 | 68.5 | 2.7 |
| Si | 71 | 29.8 | 20.5 | 0.3 | 15.3 | 24.5 | 42.1 | 86.4 | 16.9 |
| Total | 440 | 16.7 | 16.0 | 0.3 | 4.8 | 11.1 | 23.9 | 86.4 | 5.0 |
# ---- Valores para el texto ----
media_si <- mean(datos_test$prob[datos_test$y == 1]) # prob. media de quienes rotan
media_no <- mean(datos_test$prob[datos_test$y == 0]) # prob. media de quienes no rotan
tjur <- media_si - media_no # coeficiente de discriminación
prob_med <- mean(datos_test$prob) # prob. media predicha (total)
tasa_test <- mean(datos_test$y) # tasa observada en prueba
n_sobre05 <- sum(datos_test$prob > 0.5) # empleados sobre el corte 0.5
n_si_sobre05 <- sum(datos_test$prob > 0.5 & datos_test$y == 1)
# Cuartiles y dispersión por clase
p_no <- datos_test$prob[datos_test$y == 0]
p_si <- datos_test$prob[datos_test$y == 1]
med_no <- median(p_no); q3_no <- quantile(p_no, 0.75, names = FALSE)
med_si <- median(p_si); q1_si <- quantile(p_si, 0.25, names = FALSE)
de_no <- sd(p_no); de_si <- sd(p_si)
max_si <- max(p_si)
n_si_test <- sum(datos_test$y)ggplot(datos_test, aes(x = prob, fill = rotacion)) +
geom_density(alpha = 0.55, color = NA) +
# geom_rug() marca cada observación sobre el eje horizontal
geom_rug(aes(color = rotacion), sides = "b", alpha = 0.4, show.legend = FALSE) +
geom_vline(xintercept = 0.5, linetype = "dashed", color = "gray30") +
geom_vline(xintercept = tasa_test, linetype = "dotted", color = "gray30") +
scale_fill_manual(values = c(No = col_azul, Si = col_naranja),
labels = c(No = "No rota", Si = "Rota"), name = NULL) +
scale_color_manual(values = c(No = col_azul, Si = col_naranja)) +
scale_x_continuous(labels = percent, limits = c(0, 1)) +
labs(x = "Probabilidad estimada de rotar", y = "Densidad") +
theme(legend.position = "bottom")
Figura 14. Densidad de la probabilidad estimada de rotar en la muestra de prueba, según la rotación observada (azul: no rota; naranja: rota). Cada curva integra 1, por lo que se comparan las formas de las distribuciones y no el tamaño de los grupos; por el suavizado, las curvas pueden extenderse algo más allá de los valores observados. Las marcas del eje horizontal son las observaciones; la línea discontinua indica el corte convencional de 0.5 y la punteada, la tasa de rotación observada en la muestra de prueba.
Al analizar la Tabla 28 y la Figura 14, se destacan cuatro resultados:
Capacidad de discriminación. Los empleados que rotaron reciben una probabilidad media de 29.8%, frente a 14.2% de quienes no rotaron, lo que da un coeficiente de discriminación de Tjur de 0.156 (15.6 p.p.). La diferencia se aprecia en toda la distribución: la mediana de quienes rotan (24.5%) supera el tercer cuartil de quienes no rotan (20.3%), es decir, más de la mitad de los empleados que rotaron recibió una probabilidad mayor que la del 75% de los que no rotaron; asimismo, el primer cuartil de quienes rotan (15.3%) supera la mediana de quienes no rotan (9.5%).
Superposición de las distribuciones. La probabilidad de quienes no rotan se concentra en valores bajos, mientras que la de quienes rotan es más dispersa (DE de 20.5 frente a 13.6 p.p.) y tiene una cola larga que llega hasta 86.4%. Las dos curvas se cruzan cerca de la tasa de rotación observada: por debajo de ese valor predominan los empleados que no rotan y por encima, los que rotan. No obstante, la superposición entre ambos grupos es amplia, de modo que el modelo ordena a los empleados según su riesgo, pero no los separa sin error. Este resultado es coherente con el ajuste moderado reportado en la Sección 4.5.
Calibración global. La probabilidad media predicha en la muestra de prueba (16.7%) es cercana a la tasa de rotación observada (16.1%), con una diferencia de +0.6 p.p. El modelo no sobrestima ni subestima el riesgo promedio, por lo que sus probabilidades pueden interpretarse como tales en el punto 6.
Corte convencional de 0.5. Solo 22 empleados (5.0%) superan una probabilidad de 0.5, y de ellos 12 rotaron. Con ese corte se identificaría a 12 de los 71 empleados que rotaron (16.9%), porque la mayoría de quienes rotan recibe probabilidades muy inferiores a 0.5.
En este sentido, el modelo tiene capacidad para ordenar a los empleados según su riesgo de rotación, aunque con una superposición considerable entre grupos. La Sección 5.4 cuantifica esta capacidad con la curva ROC y el AUC, y el punto 6 aborda la elección de un corte adecuado, dado que el valor convencional de 0.5 dejaría sin identificar a la gran mayoría de los empleados que rotan.
5.4 Curva ROC y área bajo la curva (AUC)
La curva ROC representa, para cada posible punto de corte de la probabilidad, la proporción de rotaciones identificadas correctamente (sensibilidad) frente a la proporción de empleados que no rotan y que el modelo clasificaría como en riesgo (1 - especificidad). El área bajo la curva (AUC) resume la capacidad de discriminación del modelo en todos los cortes a la vez y se interpreta como la probabilidad de que, al tomar al azar un empleado que rota y uno que no, el modelo asigne una probabilidad mayor al primero. Se contrasta \(H_0: \text{AUC} = 0.5\) (el modelo no discrimina mejor que el azar) frente a \(H_1: \text{AUC} > 0.5\), con el error estándar de DeLong, y el AUC se reporta con su intervalo de confianza al 95% por el mismo método. Para su lectura se adopta la escala de Hosmer, Lemeshow y Sturdivant (2013): un AUC de 0.5 indica ausencia de discriminación; entre 0.5 y 0.7, discriminación pobre; entre 0.7 y 0.8, aceptable; entre 0.8 y 0.9, excelente; y de 0.9 o más, sobresaliente. Finalmente, se compara el AUC de la muestra de prueba con el de la muestra de entrenamiento: una caída significativa indicaría sobreajuste, es decir, que el modelo reproduce particularidades de los datos con que se estimó y no se generaliza a empleados nuevos.
library(pROC) # curvas ROC, AUC, intervalos de confianza y comparación de curvas
# roc() construye la curva ROC a partir de la respuesta observada y de la probabilidad
# predicha. levels = c(0, 1) fija 0 como clase negativa y 1 como positiva, y
# direction = "<" indica que se esperan probabilidades mayores en la clase positiva.
# Sin estos argumentos pROC los detecta automáticamente y puede invertir la curva.
roc_test <- roc(response = datos_test$y, predictor = datos_test$prob,
levels = c(0, 1), direction = "<", quiet = TRUE)
# Curva en entrenamiento: fitted() devuelve las probabilidades ajustadas del modelo
# para las mismas filas con que se estimó (datos_train, sin valores faltantes).
roc_train <- roc(response = datos_train$y, predictor = fitted(modelo_train),
levels = c(0, 1), direction = "<", quiet = TRUE)
# ci.auc() calcula el IC 95% del AUC con el método de DeLong (no paramétrico)
ci_test <- ci.auc(roc_test, method = "delong")
ci_train <- ci.auc(roc_train, method = "delong")
# Contraste H0: AUC = 0.5 frente a H1: AUC > 0.5 (unilateral).
# var() aplicado a un objeto roc devuelve la varianza de DeLong del AUC.
auc_test <- as.numeric(auc(roc_test))
auc_train <- as.numeric(auc(roc_train))
z_auc <- (auc_test - 0.5) / sqrt(var(roc_test))
p_auc <- pnorm(z_auc, lower.tail = FALSE)
# roc.test() compara dos AUC; paired = FALSE porque las curvas provienen de muestras
# distintas (entrenamiento y prueba), no de los mismos empleados.
comp_auc <- roc.test(roc_train, roc_test, method = "delong", paired = FALSE)
# Categoría de la escala de Hosmer, Lemeshow y Sturdivant (2013)
escala_auc <- function(a) {
as.character(cut(a, breaks = c(-Inf, 0.5, 0.7, 0.8, 0.9, Inf), right = FALSE,
labels = c("Sin discriminación", "Pobre", "Aceptable",
"Excelente", "Sobresaliente")))
}
tab_auc <- tibble(
Muestra = c("Entrenamiento", "Prueba"),
n = c(nrow(datos_train), nrow(datos_test)),
Eventos = c(sum(datos_train$y), sum(datos_test$y)),
AUC = c(auc_train, auc_test),
IC_inf = c(ci_train[1], ci_test[1]),
IC_sup = c(ci_train[3], ci_test[3]),
Gini = 2 * AUC - 1,
Escala = escala_auc(AUC)
)
estilo_tabla(tab_auc, digits = 3,
caption = "Tabla 29. Área bajo la curva ROC (AUC) del modelo reestimado en las muestras de entrenamiento y prueba, con IC 95% (método de DeLong). Gini = 2·AUC - 1. Escala: Hosmer, Lemeshow y Sturdivant (2013).",
col.names = c("Muestra", "n", "Rotaciones", "AUC", "IC 95% inf.",
"IC 95% sup.", "Gini", "Discriminación"))| Muestra | n | Rotaciones | AUC | IC 95% inf. | IC 95% sup. | Gini | Discriminación |
|---|---|---|---|---|---|---|---|
| Entrenamiento | 1030 | 166 | 0.787 | 0.745 | 0.829 | 0.574 | Aceptable |
| Prueba | 440 | 71 | 0.750 | 0.688 | 0.812 | 0.500 | Aceptable |
# ---- Valores para el texto ----
ic_auc_test <- paste0(f_num(ci_test[1], 3), " a ", f_num(ci_test[3], 3))
ic_auc_train <- paste0(f_num(ci_train[1], 3), " a ", f_num(ci_train[3], 3))
dif_auc <- auc_train - auc_test
p_comp_auc <- comp_auc$p.value
# Punto de la curva de prueba que corresponde al corte convencional de 0.5
sens_05 <- mean(p_si > 0.5) # proporción de rotaciones con probabilidad > 0.5
esp_05 <- mean(p_no <= 0.5) # proporción de no rotaciones con probabilidad <= 0.5# coords() extrae los puntos de la curva (especificidad y sensibilidad) para todos
# los cortes; transpose = FALSE los devuelve como data frame.
roc_df <- bind_rows(
Prueba = coords(roc_test, "all", ret = c("specificity", "sensitivity"),
transpose = FALSE),
Entrenamiento = coords(roc_train, "all", ret = c("specificity", "sensitivity"),
transpose = FALSE),
.id = "Muestra") %>%
mutate(Muestra = factor(Muestra, levels = c("Prueba", "Entrenamiento"),
labels = c(paste0("Prueba (AUC = ", f_num(auc_test, 3), ")"),
paste0("Entrenamiento (AUC = ", f_num(auc_train, 3), ")"))))
ggplot(roc_df, aes(x = 1 - specificity, y = sensitivity, color = Muestra)) +
geom_abline(intercept = 0, slope = 1, linetype = "dotted", color = "gray40") +
# geom_path() une los puntos en el orden de los cortes (no los reordena por x)
geom_path(linewidth = 0.9) +
annotate("point", x = 1 - esp_05, y = sens_05, size = 3) +
# annotate("label") dibuja el texto sobre un recuadro blanco, para que no lo cruce la diagonal
annotate("label", x = 1 - esp_05 + 0.03, y = sens_05, hjust = 0, vjust = 0.5,
size = 3.2, label.size = 0, fill = "white",
label = paste0("Corte 0.5\nSensibilidad: ", percent(sens_05, accuracy = 0.1),
"\nEspecificidad: ", percent(esp_05, accuracy = 0.1))) +
scale_color_manual(values = setNames(c(col_naranja, col_azul), levels(roc_df$Muestra)),
name = NULL) +
scale_x_continuous(labels = percent, limits = c(0, 1)) +
scale_y_continuous(labels = percent, limits = c(0, 1)) +
coord_equal() + # misma escala en ambos ejes: la diagonal queda a 45 grados
labs(x = "1 - especificidad (tasa de falsos positivos)",
y = "Sensibilidad (tasa de verdaderos positivos)") +
theme(legend.position = "bottom")
Figura 15. Curva ROC del modelo reestimado en la muestra de prueba (naranja) y en la de entrenamiento (azul). La diagonal punteada corresponde a un clasificador al azar (AUC = 0.5). El punto negro señala la combinación de sensibilidad y especificidad en la muestra de prueba con el corte convencional de 0.5.
Al analizar la Tabla 29 y la Figura 15, se destacan cinco resultados:
- Capacidad de discriminación en la muestra de prueba. El AUC es 0.750 (IC 95%: 0.688 a 0.812), que corresponde a una discriminación aceptable en la escala adoptada, con un índice de Gini de 0.500. En términos prácticos, si se toma al azar un empleado que rota y uno que no, el modelo asigna una probabilidad mayor al que rota en el 75.0% de los casos.
- Significancia. El contraste frente a un clasificador al azar rechaza \(H_0\) (z = 7.88; p < 0.001): el modelo discrimina mejor que el azar entre los empleados que rotan y los que no.
- Precisión de la estimación. El intervalo de confianza es amplio: su límite inferior (0.688) se ubica en la franja de discriminación pobre y el superior (0.812) en la de discriminación excelente. Con solo 71 rotaciones en la muestra de prueba, el AUC no puede ubicarse con mayor precisión, por lo que en la Sección 5.6 se complementa con validación cruzada.
- Ausencia de sobreajuste. En la muestra de entrenamiento el AUC es 0.787 (IC 95%: 0.745 a 0.829), apenas 0.037 mayor que en la de prueba, y la diferencia no es significativa (prueba de DeLong para muestras independientes, p = 0.333). La Figura 15 muestra que esa diferencia se concentra en la zona de tasas bajas de falsos positivos, mientras que en el resto del recorrido las dos curvas prácticamente coinciden. Esto implica que el desempeño del modelo es estable en cortes intermedios y algo menor en datos nuevos cuando el corte es alto, es decir, cuando se busca intervenir a un grupo muy reducido de empleados.
- Corte convencional de 0.5. Se ubica en el extremo inferior izquierdo de la curva (sensibilidad de 16.9% y especificidad de 97.3%). A partir de ese punto la curva asciende con fuerza, lo que indica que reducir el corte permitiría identificar una proporción mucho mayor de rotaciones a costa de un aumento moderado de la tasa de falsos positivos.
En este sentido, el modelo tiene una capacidad de discriminación aceptable, estadísticamente significativa y sin evidencia de sobreajuste, lo que respalda su uso para priorizar a los empleados según su riesgo de rotación. No obstante, el AUC resume el desempeño en todos los cortes a la vez; la decisión práctica exige fijar uno, y la ubicación del corte de 0.5 en la curva confirma que ese valor no es adecuado para este problema. La elección del corte se aborda en el punto 6.
5.5 Matriz de confusión con el corte convencional de 0.5
Para usar el modelo en la toma de decisiones es necesario fijar un corte y clasificar a cada empleado como “rota” o “no rota”. La matriz de confusión cruza esa clasificación con la rotación observada; se construye primero con el corte convencional de 0.5, que sirve como línea base para el punto 6. Además de los indicadores habituales, se contrasta si la exactitud supera la tasa de no información (la que se obtendría clasificando a todos los empleados en la clase mayoritaria), se calcula el kappa de Cohen, que descuenta el acuerdo esperado por azar, y se aplica la prueba de McNemar, que evalúa si los dos tipos de error se presentan con la misma frecuencia.
# Clase predicha con el corte de 0.5 y clase observada, como factores con los mismos niveles
datos_test <- datos_test %>%
mutate(pred_05 = factor(ifelse(prob > 0.5, 1, 0), levels = c(0, 1)),
obs = factor(y, levels = c(0, 1)))
# confusionMatrix() (caret) cruza la predicción (data) con la realidad (reference) y
# calcula los indicadores. positive = "1" fija 'rota' como clase positiva; por defecto
# caret toma el primer nivel ("0") y reportaría los indicadores de la clase 'no rota'.
cm_05 <- confusionMatrix(data = datos_test$pred_05, reference = datos_test$obs,
positive = "1")
# Celdas de la matriz (cm_05$table tiene la predicción en filas y la realidad en columnas)
VP_05 <- cm_05$table["1", "1"]; FN_05 <- cm_05$table["0", "1"]
FP_05 <- cm_05$table["1", "0"]; VN_05 <- cm_05$table["0", "0"]
# Tabla 30: matriz en el formato del curso (realidad en filas, predicción en columnas)
mat_05 <- tibble(
Real = c("Rota", "No rota", "Total"),
Pred_rota = c(paste0(VP_05, " (VP)"), paste0(FP_05, " (FP)"), VP_05 + FP_05),
Pred_no = c(paste0(FN_05, " (FN)"), paste0(VN_05, " (VN)"), FN_05 + VN_05),
Total = c(VP_05 + FN_05, FP_05 + VN_05, nrow(datos_test))
)
estilo_tabla(mat_05,
caption = "Tabla 30. Matriz de confusión del modelo en la muestra de prueba con el corte de 0.5 (filas: rotación observada; columnas: clasificación del modelo). VP: verdaderos positivos; FN: falsos negativos; FP: falsos positivos; VN: verdaderos negativos.",
col.names = c("Rotación observada", "Predicción: rota", "Predicción: no rota", "Total"))| Rotación observada | Predicción: rota | Predicción: no rota | Total |
|---|---|---|---|
| Rota | 12 (VP) | 59 (FN) | 71 |
| No rota | 10 (FP) | 359 (VN) | 369 |
| Total | 22 | 418 | 440 |
# Tabla 31: indicadores derivados de la matriz
ov <- cm_05$overall # indicadores globales (exactitud, kappa, NIR, p-valores)
bc <- cm_05$byClass # indicadores por clase (sensibilidad, especificidad, etc.)
# Categoría de la escala de Landis y Koch (1977) para el kappa
escala_kappa <- function(k) {
as.character(cut(k, breaks = c(-Inf, 0, 0.2, 0.4, 0.6, 0.8, Inf),
labels = c("sin acuerdo", "leve", "aceptable", "moderado",
"considerable", "casi perfecto")))
}
ind_05 <- tibble(
Indicador = c("Exactitud", "Tasa de error", "Tasa de no información",
"p-valor (exactitud > no información)", "Kappa de Cohen",
"Sensibilidad", "Especificidad", "Valor predictivo positivo (precisión)",
"Valor predictivo negativo", "Exactitud balanceada", "p-valor de McNemar"),
Calculo = c("(VP + VN) / Total", "(FP + FN) / Total", "Proporción de la clase mayoritaria",
"Binomial unilateral", "Acuerdo corregido por azar",
"VP / (VP + FN)", "VN / (VN + FP)", "VP / (VP + FP)",
"VN / (VN + FN)", "(Sensibilidad + Especificidad) / 2", "FN frente a FP"),
Valor = c(
paste0(percent(ov["Accuracy"], accuracy = 0.1), " (IC 95%: ",
percent(ov["AccuracyLower"], accuracy = 0.1), " a ",
percent(ov["AccuracyUpper"], accuracy = 0.1), ")"),
percent(1 - ov["Accuracy"], accuracy = 0.1),
percent(ov["AccuracyNull"], accuracy = 0.1),
fmt_p(ov["AccuracyPValue"]),
paste0(f_num(ov["Kappa"], 3), " (", escala_kappa(ov["Kappa"]), ")"),
percent(bc["Sensitivity"], accuracy = 0.1),
percent(bc["Specificity"], accuracy = 0.1),
percent(bc["Pos Pred Value"], accuracy = 0.1),
percent(bc["Neg Pred Value"], accuracy = 0.1),
percent(bc["Balanced Accuracy"], accuracy = 0.1),
fmt_p(ov["McnemarPValue"]))
)
estilo_tabla(ind_05,
caption = "Tabla 31. Indicadores de clasificación del modelo en la muestra de prueba con el corte de 0.5 (clase positiva: rota). Escala del kappa (Landis y Koch, 1977): 0 a 0.20, leve; 0.21 a 0.40, aceptable; 0.41 a 0.60, moderado; 0.61 a 0.80, considerable; más de 0.80, casi perfecto.",
col.names = c("Indicador", "Cálculo", "Valor"))| Indicador | Cálculo | Valor |
|---|---|---|
| Exactitud | (VP + VN) / Total | 84.3% (IC 95%: 80.6% a 87.6%) |
| Tasa de error | (FP + FN) / Total | 15.7% |
| Tasa de no información | Proporción de la clase mayoritaria | 83.9% |
| p-valor (exactitud > no información) | Binomial unilateral | 0.428 |
| Kappa de Cohen | Acuerdo corregido por azar | 0.197 (leve) |
| Sensibilidad | VP / (VP + FN) | 16.9% |
| Especificidad | VN / (VN + FP) | 97.3% |
| Valor predictivo positivo (precisión) | VP / (VP + FP) | 54.5% |
| Valor predictivo negativo | VN / (VN + FN) | 85.9% |
| Exactitud balanceada | (Sensibilidad + Especificidad) / 2 | 57.1% |
| p-valor de McNemar | FN frente a FP | <0.001 |
# ---- Valores para el texto ----
aciertos_05 <- VP_05 + VN_05
fn_rate_neg <- FN_05 / (FN_05 + VN_05) # proporción que rota entre los clasificados "no rota"Al analizar las Tablas 30 y 31, se destacan cinco resultados:
- Exactitud frente a la tasa de no información. El modelo clasifica correctamente al 84.3% de los empleados (IC 95%: 80.6% a 87.6%), pero esa cifra no supera de forma significativa la tasa de no información (83.9%; p = 0.428). De los 371 aciertos, 359 corresponden a empleados que no rotaron; con este corte, el modelo clasifica prácticamente igual que una regla que asignara “no rota” a toda la plantilla.
- Detección de las rotaciones. El modelo identifica a 12 de los 71 empleados que rotaron (sensibilidad de 16.9%), de modo que 59 salidas pasan inadvertidas. A cambio, genera muy pocas falsas alarmas: 10 de 369 empleados que no rotaron (especificidad de 97.3%).
- Asimetría de los errores. Los falsos negativos (59) superan ampliamente a los falsos positivos (10), y la prueba de McNemar confirma que la diferencia es significativa (p < 0.001): con el corte de 0.5, el modelo se equivoca de manera sistemática hacia la clase “no rota”.
- Valores predictivos. De los 22 empleados clasificados en riesgo, 12 rotaron (valor predictivo positivo de 54.5%), es decir, una intervención sobre ellos alcanzaría a un empleado en riesgo real en algo más de la mitad de los casos. En cambio, entre los 418 clasificados como “no rota”, el 14.1% terminó rotando, una proporción cercana a la tasa de rotación de la muestra (16.1%): que el modelo clasifique a un empleado como “no rota” apenas modifica el riesgo que ya se conocía sin él.
- Acuerdo corregido por azar. El kappa de Cohen es 0.197, un acuerdo leve en la escala de Landis y Koch (1977), y la exactitud balanceada es 57.1%, apenas por encima del 50% que correspondería a una clasificación al azar.
En este sentido, el contraste entre estos indicadores y el AUC de 0.750 (Sección 5.4) muestra que el bajo desempeño de la clasificación no se debe a la capacidad del modelo para ordenar a los empleados según su riesgo, sino al corte utilizado: con una tasa de rotación cercana a 16%, el umbral de 0.5 es tan alto que muy pocos empleados lo superan. Un programa de retención basado en este corte llegaría a 22 empleados y, en el mejor de los casos, evitaría 12 de las 71 salidas. El punto 6 analiza cómo cambian estos indicadores al modificar el corte y propone uno acorde con los objetivos de la gerencia.
5.6 Validación cruzada
El AUC de la Sección 5.4 proviene de una sola partición aleatoria, con 71 rotaciones en la muestra de prueba; con otra partición, su valor habría sido distinto. Para evaluar si es representativo, se aplica una validación cruzada estratificada de 10 pliegues con 10 repeticiones sobre la base completa: cada pliegue se usa una vez como prueba mientras el modelo se estima con los nueve restantes. Se reportan el AUC de cada pliegue, que muestra cuánto varía el desempeño entre muestras de prueba, y el AUC agrupado de cada repetición, calculado con las probabilidades de los 1470 empleados a la vez, que ofrece una estimación más estable del desempeño en datos nuevos.
set.seed(123) # reproducibilidad de las particiones
n_rep <- 10 # repeticiones
k <- 10 # pliegues por repetición
cv_pliegues <- list() # AUC de cada pliegue
cv_agrupado <- list() # AUC agrupado de cada repetición
for (r in seq_len(n_rep)) {
# createFolds() (caret) divide las filas en k grupos estratificados por 'rotacion';
# returnTrain = FALSE devuelve, para cada pliegue, los índices de prueba.
pliegues <- createFolds(datos$rotacion, k = k, list = TRUE, returnTrain = FALSE)
prob_oof <- numeric(nrow(datos)) # probabilidades "fuera de pliegue" de la repetición
for (f in seq_along(pliegues)) {
idx <- pliegues[[f]]
# Misma especificación del modelo final, estimada sin el pliegue f
m_cv <- update(modelo_final, data = datos[-idx, ])
p_cv <- predict(m_cv, newdata = datos[idx, ], type = "response")
prob_oof[idx] <- p_cv
cv_pliegues[[length(cv_pliegues) + 1]] <- tibble(
repeticion = r, pliegue = f,
auc = as.numeric(auc(roc(datos$y[idx], p_cv, levels = c(0, 1),
direction = "<", quiet = TRUE))))
}
cv_agrupado[[r]] <- tibble(
repeticion = r,
auc = as.numeric(auc(roc(datos$y, prob_oof, levels = c(0, 1),
direction = "<", quiet = TRUE))))
}
cv_pliegues <- bind_rows(cv_pliegues)
cv_agrupado <- bind_rows(cv_agrupado)
# Resumen de una distribución de AUC
resumen_auc <- function(x) {
tibble(n = length(x), Media = mean(x), DE = sd(x), Min = min(x),
P2.5 = quantile(x, 0.025, names = FALSE), Mediana = median(x),
P97.5 = quantile(x, 0.975, names = FALSE), Max = max(x),
Sobre_07 = mean(x >= 0.7) * 100)
}
tab_cv <- bind_rows(
resumen_auc(cv_pliegues$auc) %>% mutate(Medida = "AUC por pliegue", .before = 1),
resumen_auc(cv_agrupado$auc) %>% mutate(Medida = "AUC agrupado por repetición", .before = 1)
)
estilo_tabla(tab_cv, digits = 3,
caption = "Tabla 32. Validación cruzada estratificada de 10 pliegues con 10 repeticiones sobre la base completa: distribución del AUC por pliegue (100 valores) y del AUC agrupado por repetición (10 valores). P2.5 y P97.5: percentiles 2.5 y 97.5; ≥ 0.7: porcentaje de valores en la franja de discriminación aceptable o superior.",
col.names = c("Medida", "n", "Media", "DE", "Mín.", "P2.5", "Mediana",
"P97.5", "Máx.", "≥ 0.7 (%)"))| Medida | n | Media | DE | Mín. | P2.5 | Mediana | P97.5 | Máx. | ≥ 0.7 (%) |
|---|---|---|---|---|---|---|---|---|---|
| AUC por pliegue | 100 | 0.770 | 0.056 | 0.614 | 0.672 | 0.776 | 0.865 | 0.883 | 89 |
| AUC agrupado por repetición | 10 | 0.769 | 0.002 | 0.766 | 0.766 | 0.770 | 0.771 | 0.771 | 100 |
# ---- Valores para el texto ----
auc_cv_media <- mean(cv_agrupado$auc) # estimación principal
auc_cv_de <- sd(cv_agrupado$auc)
auc_cv_rango <- range(cv_agrupado$auc)
auc_pl_media <- mean(cv_pliegues$auc)
auc_pl_de <- sd(cv_pliegues$auc)
auc_pl_rango <- range(cv_pliegues$auc)
auc_pl_p <- quantile(cv_pliegues$auc, c(0.025, 0.975), names = FALSE)
pct_pl_07 <- mean(cv_pliegues$auc >= 0.7) * 100
pos_test <- mean(cv_pliegues$auc <= auc_test) * 100 # percentil del AUC de prueba
optimismo <- auc_train - auc_cv_media # AUC entrenamiento - AUC validación cruzadaggplot(cv_pliegues, aes(x = auc)) +
geom_histogram(bins = 20, fill = col_azul, color = "white", alpha = 0.85) +
geom_vline(xintercept = 0.7, linetype = "dotted", color = "gray40") +
geom_vline(xintercept = auc_cv_media, linetype = "dashed", color = "black") +
geom_vline(xintercept = auc_test, color = col_naranja, linewidth = 0.9) +
labs(x = "AUC del pliegue", y = "Número de pliegues")
Figura 16. Distribución del AUC en los 100 pliegues de la validación cruzada (10 pliegues × 10 repeticiones). La línea discontinua negra indica la media del AUC agrupado por repetición; la línea naranja, el AUC de la muestra de prueba (Sección 5.4); y la punteada gris, el umbral de 0.7 a partir del cual la discriminación se considera aceptable.
Al analizar la Tabla 32 y la Figura 16, se destacan cuatro resultados:
- Estimación del desempeño en datos nuevos. El AUC agrupado es 0.769 en promedio y prácticamente no varía entre repeticiones (DE = 0.002; rango de 0.766 a 0.771), lo que corresponde a una discriminación aceptable. Esta es la estimación más estable de la capacidad del modelo para ordenar a empleados que no intervinieron en su estimación: al comparar al azar un empleado que rota y uno que no, el modelo asigna mayor probabilidad al primero en cerca del 77% de los casos.
- Representatividad de la muestra de prueba. El AUC de la Sección 5.4 (0.750) se ubica en el percentil 38 de la distribución de los pliegues. La partición de la Sección 5.1 fue, por tanto, típica y algo desfavorable: con otra partición, lo más probable habría sido un AUC ligeramente mayor, de modo que la evaluación de la Sección 5.4 no sobrestima el desempeño del modelo.
- Variabilidad entre muestras de prueba. El AUC por pliegue varía entre 0.614 y 0.883, y el 95% central de los valores se ubica entre 0.672 y 0.865 (DE = 0.056). Esta dispersión, mayor que la del intervalo de la Sección 5.4, se debe a que cada pliegue contiene alrededor de 24 rotaciones, e ilustra por qué una sola muestra de prueba ofrece una estimación imprecisa. Aun así, el 89% de los pliegues alcanza la franja de discriminación aceptable o superior (AUC ≥ 0.7).
- Optimismo de la evaluación en entrenamiento. El AUC de la muestra de entrenamiento (0.787) supera al de la validación cruzada en apenas 0.018. Esa diferencia mide el optimismo de evaluar el modelo con los mismos datos con que se estimó y, por su baja magnitud, confirma la ausencia de sobreajuste señalada en la Sección 5.4.
En este sentido, la validación cruzada respalda la conclusión de la Sección 5.4: el modelo tiene una capacidad de discriminación aceptable y estable, que no depende de la partición utilizada, con un AUC esperado en datos nuevos de 0.769.
5.7 Calibración
La discriminación indica si el modelo ordena bien a los empleados; la calibración, si sus probabilidades corresponden a las frecuencias observadas, condición necesaria para interpretarlas literalmente en el punto 6. Se agrupan los empleados de la muestra de prueba en quintiles de probabilidad predicha y se compara, en cada grupo, la probabilidad media con la tasa de rotación observada. Como medida resumen se estima la pendiente de calibración (regresión logística de la rotación observada sobre el logit de la probabilidad predicha) y el intercepto de calibración; un modelo perfectamente calibrado tiene pendiente 1 e intercepto 0.
# ntile() (dplyr) asigna a cada empleado el quintil de su probabilidad predicha (1 a 5)
datos_test <- datos_test %>% mutate(quintil = ntile(prob, 5))
tab_cal <- datos_test %>%
group_by(quintil) %>%
summarise(n = n(), Rotan = sum(y),
Pred = mean(prob) * 100, # probabilidad media predicha
Obs = mean(y) * 100, # tasa observada
Obs_inf = binom.test(sum(y), n())$conf.int[1] * 100, # IC 95% exacto
Obs_sup = binom.test(sum(y), n())$conf.int[2] * 100,
.groups = "drop") %>%
mutate(Dif = Obs - Pred)
estilo_tabla(tab_cal, digits = 1,
caption = "Tabla 33. Calibración del modelo en la muestra de prueba por quintiles de probabilidad predicha: probabilidad media predicha frente a tasa de rotación observada (%), con IC 95% exacto de la tasa observada. Diferencia = observada - predicha (p.p.).",
col.names = c("Quintil", "n", "Rotan", "Predicha (%)", "Observada (%)",
"IC 95% inf.", "IC 95% sup.", "Diferencia (p.p.)"))| Quintil | n | Rotan | Predicha (%) | Observada (%) | IC 95% inf. | IC 95% sup. | Diferencia (p.p.) |
|---|---|---|---|---|---|---|---|
| 1 | 88 | 5 | 2.0 | 5.7 | 1.9 | 12.8 | 3.7 |
| 2 | 88 | 5 | 5.8 | 5.7 | 1.9 | 12.8 | -0.1 |
| 3 | 88 | 8 | 11.5 | 9.1 | 4.0 | 17.1 | -2.4 |
| 4 | 88 | 20 | 21.2 | 22.7 | 14.5 | 32.9 | 1.6 |
| 5 | 88 | 33 | 43.2 | 37.5 | 27.4 | 48.5 | -5.7 |
# Pendiente de calibración: qlogis() transforma la probabilidad en logit
lp_test <- qlogis(datos_test$prob)
m_pend <- glm(y ~ lp_test, family = binomial, data = datos_test)
# Intercepto de calibración: offset() fija la pendiente en 1 y estima solo el intercepto
m_int <- glm(y ~ offset(lp_test), family = binomial, data = datos_test)
# ---- Valores para el texto ----
pend_cal <- unname(coef(m_pend)["lp_test"])
pend_cal_ic <- unname(confint.default(m_pend)["lp_test", ]) # IC 95% de Wald
int_cal <- unname(coef(m_int)[1])
int_cal_ic <- unname(confint.default(m_int)[1, ])
max_dif_cal <- max(abs(tab_cal$Dif))
diag_dentro <- all(tab_cal$Pred >= tab_cal$Obs_inf & tab_cal$Pred <= tab_cal$Obs_sup)
# Concentración de las rotaciones en el quintil de mayor riesgo
rot_q5 <- tab_cal$Rotan[tab_cal$quintil == 5]
captura_q5 <- rot_q5 / sum(datos_test$y)
razon_q5q1 <- tab_cal$Obs[tab_cal$quintil == 5] / tab_cal$Obs[tab_cal$quintil == 1]ggplot(tab_cal, aes(x = Pred / 100, y = Obs / 100)) +
geom_abline(intercept = 0, slope = 1, linetype = "dotted", color = "gray40") +
geom_errorbar(aes(ymin = Obs_inf / 100, ymax = Obs_sup / 100),
width = 0.015, color = col_naranja) +
geom_point(size = 3, color = col_naranja) +
scale_x_continuous(labels = percent, limits = c(0, 0.6)) +
scale_y_continuous(labels = percent, limits = c(0, 0.6)) +
coord_equal() +
labs(x = "Probabilidad media predicha", y = "Tasa de rotación observada")
Figura 17. Calibración del modelo en la muestra de prueba: tasa de rotación observada frente a probabilidad media predicha en cada quintil, con IC 95% exacto de la tasa observada. La diagonal punteada representa una calibración perfecta.
Al analizar la Tabla 33 y la Figura 17, se destacan tres resultados:
- Calibración en todo el rango. En los cinco quintiles la diferencia entre la tasa observada y la probabilidad predicha es a lo maximo 5.7 p.p., y la diagonal de calibración perfecta queda dentro de todos los intervalos de confianza. El intercepto de calibración (-0.053; IC 95%: -0.333 a 0.227) contiene el 0, por lo que el modelo no sobrestima ni subestima el riesgo de forma sistemática.
- Predicciones ligeramente extremas. En el quintil de menor riesgo la tasa observada supera a la predicha y en el de mayor riesgo ocurre lo contrario, de modo que el modelo exagera levemente hacia ambos extremos. La pendiente de calibración (0.773; IC 95%: 0.532 a 1.015) es coherente con ese patrón, aunque su intervalo incluye el 1 y la desviación no es significativa.
- Concentración del riesgo. La tasa de rotación observada aumenta de 5.7% en el primer quintil a 37.5% en el quinto, es decir, 6.6 veces más. El 20% de empleados con mayor probabilidad predicha concentra 33 de las 71 rotaciones de la muestra de prueba (46.5%).
En este sentido, las probabilidades del modelo pueden interpretarse como frecuencias esperadas de rotación en el punto 6, con la cautela de que en los perfiles de riesgo más alto pueden sobrestimar el riesgo en algunos puntos porcentuales.
5.8 Especificación alternativa del ingreso y la antigüedad en el cargo
La Sección 4.8 rechazó la linealidad en el logit del ingreso mensual y de la antigüedad en el cargo, y mostró que la especificación con el logaritmo de ambas variables mejora el ajuste sin agregar parámetros. El AIC mide ajuste, no capacidad de identificar a los empleados que rotan, por lo que ambas especificaciones se comparan aquí con el mismo protocolo de evaluación: AUC en la muestra de prueba, con la prueba de DeLong pareada (los dos modelos se evalúan sobre los mismos empleados), y AUC agrupado en la validación cruzada de la Sección 5.6, con los mismos pliegues para ambos modelos. Se adopta la especificación logarítmica como modelo final solo si mejora el AUC de validación cruzada en más de 0.01 y la diferencia en la muestra de prueba es significativa al 5%; de lo contrario, se conserva la especificación lineal por su interpretabilidad.
# Modelo alternativo (Sección 4.8) reestimado con la muestra de entrenamiento
modelo_alt_train <- update(modelo_alt, data = datos_train)
datos_test <- datos_test %>%
mutate(prob_alt = predict(modelo_alt_train, newdata = datos_test, type = "response"))
roc_alt_test <- roc(datos_test$y, datos_test$prob_alt, levels = c(0, 1),
direction = "<", quiet = TRUE)
auc_alt_test <- as.numeric(auc(roc_alt_test))
ci_alt_test <- ci.auc(roc_alt_test, method = "delong")
# Prueba de DeLong pareada: paired = TRUE porque ambos modelos se evalúan sobre
# los mismos empleados de la muestra de prueba
comp_alt <- roc.test(roc_test, roc_alt_test, method = "delong", paired = TRUE)
# Validación cruzada con los mismos pliegues para los dos modelos. Con la misma semilla
# y el mismo esquema que la Sección 5.6, createFolds() genera exactamente los mismos pliegues.
auc_de <- function(p) as.numeric(auc(roc(datos$y, p, levels = c(0, 1),
direction = "<", quiet = TRUE)))
set.seed(123)
cv_comp <- list()
for (r in seq_len(n_rep)) {
pliegues <- createFolds(datos$rotacion, k = k, list = TRUE, returnTrain = FALSE)
p_lin <- numeric(nrow(datos))
p_alt <- numeric(nrow(datos))
for (f in seq_along(pliegues)) {
idx <- pliegues[[f]]
p_lin[idx] <- predict(update(modelo_final, data = datos[-idx, ]),
newdata = datos[idx, ], type = "response")
p_alt[idx] <- predict(update(modelo_alt, data = datos[-idx, ]),
newdata = datos[idx, ], type = "response")
}
cv_comp[[r]] <- tibble(repeticion = r, auc_lin = auc_de(p_lin), auc_alt = auc_de(p_alt))
}
cv_comp <- bind_rows(cv_comp) %>% mutate(dif = auc_alt - auc_lin)
# Verificación: el AUC del modelo lineal coincide con el de la Sección 5.6 (mismos pliegues)
stopifnot(isTRUE(all.equal(cv_comp$auc_lin, cv_agrupado$auc)))
tab_alt <- tibble(
Modelo = c("Lineal (modelo principal)",
"Logaritmo del ingreso y de (antigüedad en el cargo + 1)",
"Diferencia (logarítmico - lineal)"),
AIC = c(AIC(modelo_final), AIC(modelo_alt), AIC(modelo_alt) - AIC(modelo_final)),
AUC_pr = c(auc_test, auc_alt_test, auc_alt_test - auc_test),
IC_pr = c(ic_auc_test,
paste0(f_num(ci_alt_test[1], 3), " a ", f_num(ci_alt_test[3], 3)),
fmt_p_txt(comp_alt$p.value)),
AUC_cv = c(mean(cv_comp$auc_lin), mean(cv_comp$auc_alt), mean(cv_comp$dif)),
DE_cv = c(sd(cv_comp$auc_lin), sd(cv_comp$auc_alt), sd(cv_comp$dif))
)
estilo_tabla(tab_alt, digits = 3,
caption = "Tabla 34. Comparación de la especificación lineal y la logarítmica: AIC en la base completa, AUC en la muestra de prueba (IC 95% de DeLong; en la fila de diferencia, p-valor de la prueba de DeLong pareada) y AUC agrupado en la validación cruzada de 10 pliegues con 10 repeticiones (media y DE entre repeticiones; en la fila de diferencia, media y DE de la diferencia pareada por repetición).",
col.names = c("Especificación", "AIC", "AUC prueba", "IC 95% / p-valor",
"AUC validación cruzada", "DE"))| Especificación | AIC | AUC prueba | IC 95% / p-valor | AUC validación cruzada | DE |
|---|---|---|---|---|---|
| Lineal (modelo principal) | 1096.84 | 0.750 | 0.688 a 0.812 | 0.769 | 0.002 |
| Logaritmo del ingreso y de (antigüedad en el cargo + 1) | 1078.26 | 0.761 | 0.700 a 0.822 | 0.774 | 0.002 |
| Diferencia (logarítmico - lineal) | -18.58 | 0.011 | p = 0.150 | 0.005 | 0.000 |
# ---- Regla de decisión y valores para el texto ----
umbral_auc <- 0.01
dif_auc_cv <- mean(cv_comp$dif)
dif_auc_pr <- auc_alt_test - auc_test
p_pareado <- comp_alt$p.value
pct_alt_mejor <- mean(cv_comp$dif > 0) * 100 # repeticiones en que el logarítmico supera al lineal
adopta_alt <- dif_auc_cv > umbral_auc & p_pareado < 0.05 & dif_auc_pr > 0Al analizar la Tabla 34, se destacan tres resultados:
- Ajuste. La especificación logarítmica reduce el AIC en 18.58 unidades con el mismo número de parámetros, lo que confirma que describe mejor la relación del ingreso y de la antigüedad en el cargo con la rotación.
- Discriminación en la muestra de prueba. El AUC pasa de 0.750 a 0.761 (+0.011), pero la diferencia no es significativa (prueba de DeLong pareada, p = 0.150) y los intervalos de confianza de ambos modelos se superponen casi por completo.
- Discriminación en la validación cruzada. El AUC agrupado pasa de 0.769 a 0.774, una mejora de 0.005 que se presenta en el 100% de las repeticiones. La ganancia es, por tanto, sistemática, pero de una magnitud que no modifica de forma perceptible la capacidad del modelo para identificar a los empleados en riesgo.
La regla de decisión no se cumple: la mejora en validación cruzada (0.005) es inferior al umbral de 0.01 y la diferencia en la muestra de prueba no es significativa, por lo que se conserva la especificación lineal como modelo final, por su interpretabilidad y porque la especificación logarítmica no mejora de forma material la identificación de los empleados en riesgo. No obstante, el mejor ajuste de la especificación logarítmica tiene una lectura sustantiva: el efecto del ingreso y de la antigüedad en el cargo no es constante, sino más intenso en los salarios más bajos y en los primeros años en el cargo (Sección 4.8, Figura 11). Los OR por unidad de la Sección 4.3 deben leerse, en consecuencia, como efectos promedio, y esta concentración del riesgo se retoma en la estrategia del punto 7.
5.9 Síntesis de la evaluación
El modelo tiene una capacidad de discriminación aceptable: en la muestra de prueba su AUC es 0.750 (IC 95%: 0.688 a 0.812), significativamente superior al de un clasificador al azar, y la validación cruzada estima un AUC de 0.769 en datos nuevos, estable entre particiones. No hay evidencia de sobreajuste y sus probabilidades están bien calibradas, por lo que pueden interpretarse como frecuencias esperadas de rotación. En la práctica, el modelo ordena adecuadamente a los empleados según su riesgo: el 20% con mayor probabilidad estimada concentra el 46.5% de las rotaciones. No obstante, con el corte convencional de 0.5 solo identifica al 16.9% de quienes rotan, lo que hace necesario definir otro corte (punto 6). Deben tenerse presentes tres limitaciones:
- La muestra de prueba contiene solo 71 rotaciones, lo que hace imprecisas las métricas que dependen de los casos positivos; en el caso del AUC, la validación cruzada atenúa esta limitación.
- La discriminación es aceptable, pero no excelente: las probabilidades de quienes rotan y quienes no se superponen ampliamente, de modo que cualquier corte producirá errores de ambos tipos.
- La validación es interna: los datos provienen de una sola organización y un solo período, por lo que el desempeño podría diferir si cambian las condiciones laborales o la composición de la plantilla.
6 Predicciones y punto de corte
Con el corte convencional de 0.5, el modelo identifica solo el 16.9% de las rotaciones (Sección 5.5). Para que sea útil en la decisión que plantea la gerencia, es decir, intervenir o no a un empleado para motivar su permanencia, es necesario elegir un corte acorde con ese propósito. En esta sección se analiza cómo cambian los indicadores de clasificación con el corte, se elige uno con un criterio explícito, se evalúa en la muestra de prueba y, finalmente, se aplica a un empleado hipotético.
6.1 Sensibilidad de la clasificación al punto de corte
El corte se elige con la muestra de entrenamiento y se evalúa en la de prueba, que no interviene en la elección. Para evitar el optimismo de las probabilidades dentro de muestra, se aplica la validación cruzada de la Sección 5.6 a la muestra de entrenamiento, con 10 pliegues y 5 repeticiones: cada empleado recibe el promedio de las probabilidades asignadas por modelos que no lo incluyeron. Con esas probabilidades se calculan los indicadores de clasificación para cada corte, junto con el porcentaje de empleados que quedarían clasificados en riesgo, que expresa la carga de intervenciones que implicaría cada opción para la gerencia.
set.seed(123)
n_rep_c <- 5 # repeticiones de la validación cruzada en entrenamiento
# Matriz de probabilidades fuera de pliegue: una fila por empleado, una columna por repetición
prob_oof_mat <- matrix(NA_real_, nrow = nrow(datos_train), ncol = n_rep_c)
for (r in seq_len(n_rep_c)) {
pliegues <- createFolds(datos_train$rotacion, k = 10, list = TRUE, returnTrain = FALSE)
for (f in seq_along(pliegues)) {
idx <- pliegues[[f]]
m_f <- update(modelo_final, data = datos_train[-idx, ])
prob_oof_mat[idx, r] <- predict(m_f, newdata = datos_train[idx, ], type = "response")
}
}
# rowMeans() promedia, para cada empleado, sus probabilidades de las 5 repeticiones
datos_train <- datos_train %>% mutate(prob_oof = rowMeans(prob_oof_mat))
stopifnot(!anyNA(datos_train$prob_oof))
# AUC de las probabilidades fuera de pliegue (control de coherencia con la Sección 5.6)
auc_oof_train <- as.numeric(auc(roc(datos_train$y, datos_train$prob_oof, levels = c(0, 1),
direction = "<", quiet = TRUE)))# Indicadores de clasificación para un corte dado
metricas_corte <- function(p, y, corte) {
pred <- as.integer(p > corte)
VP <- sum(pred == 1 & y == 1); FN <- sum(pred == 0 & y == 1)
FP <- sum(pred == 1 & y == 0); VN <- sum(pred == 0 & y == 0)
tibble(corte = corte,
intervenidos = (VP + FP) / length(y),
sens = VP / (VP + FN),
esp = VN / (VN + FP),
vpp = ifelse(VP + FP > 0, VP / (VP + FP), NA_real_),
vpn = ifelse(VN + FN > 0, VN / (VN + FN), NA_real_),
exact = (VP + VN) / length(y),
exact_bal = (sens + esp) / 2,
youden = sens + esp - 1)
}
# Rejilla fina (se usa en la figura y en la elección del corte) y rejilla de la tabla
cortes_fino <- map_dfr(seq(0.01, 0.80, by = 0.01),
~ metricas_corte(datos_train$prob_oof, datos_train$y, .x))
cortes_tab <- map_dfr(seq(0.05, 0.60, by = 0.05),
~ metricas_corte(datos_train$prob_oof, datos_train$y, .x))
estilo_tabla(cortes_tab %>%
mutate(across(intervenidos:exact_bal, ~ .x * 100),
corte = f_num(corte, 2), # cortes con dos decimales
youden = f_num(youden, 3)), # Youden con tres decimales
digits = 1,
caption = "Tabla 35. Indicadores de clasificación según el punto de corte, calculados con las probabilidades fuera de pliegue de la muestra de entrenamiento (validación cruzada de 10 pliegues con 5 repeticiones). Intervenidos: porcentaje de empleados clasificados en riesgo. VPP y VPN: valores predictivos positivo y negativo. Youden = sensibilidad + especificidad - 1.",
col.names = c("Corte", "Intervenidos (%)", "Sensibilidad (%)", "Especificidad (%)",
"VPP (%)", "VPN (%)", "Exactitud (%)", "Exactitud balanceada (%)",
"Youden"))| Corte | Intervenidos (%) | Sensibilidad (%) | Especificidad (%) | VPP (%) | VPN (%) | Exactitud (%) | Exactitud balanceada (%) | Youden |
|---|---|---|---|---|---|---|---|---|
| 0.05 | 75.9 | 89.8 | 26.7 | 19.1 | 93.1 | 36.9 | 58.2 | 0.165 |
| 0.10 | 51.7 | 80.7 | 53.9 | 25.2 | 93.6 | 58.3 | 67.3 | 0.347 |
| 0.15 | 37.6 | 74.7 | 69.6 | 32.0 | 93.5 | 70.4 | 72.1 | 0.443 |
| 0.20 | 28.9 | 66.9 | 78.4 | 37.2 | 92.5 | 76.5 | 72.6 | 0.452 |
| 0.25 | 21.2 | 59.0 | 86.1 | 45.0 | 91.6 | 81.7 | 72.6 | 0.451 |
| 0.30 | 15.7 | 48.2 | 90.5 | 49.4 | 90.1 | 83.7 | 69.4 | 0.387 |
| 0.35 | 11.2 | 39.2 | 94.2 | 56.5 | 89.0 | 85.3 | 66.7 | 0.334 |
| 0.40 | 9.2 | 31.9 | 95.1 | 55.8 | 87.9 | 85.0 | 63.5 | 0.271 |
| 0.45 | 7.1 | 23.5 | 96.1 | 53.4 | 86.7 | 84.4 | 59.8 | 0.196 |
| 0.50 | 4.8 | 17.5 | 97.7 | 59.2 | 86.0 | 84.8 | 57.6 | 0.152 |
| 0.55 | 3.5 | 15.7 | 98.8 | 72.2 | 85.9 | 85.4 | 57.3 | 0.145 |
| 0.60 | 2.9 | 14.5 | 99.3 | 80.0 | 85.8 | 85.6 | 56.9 | 0.138 |
# ---- Valores para el texto ----
# Fila de la rejilla fina correspondiente a un corte dado
en_corte <- function(c) cortes_fino[which.min(abs(cortes_fino$corte - c)), ]
m_05 <- en_corte(0.5)
m_03 <- en_corte(0.3)
m_06 <- en_corte(0.6)
# Corte en que se cruzan la sensibilidad y la especificidad
m_cruce <- cortes_fino[which.min(abs(cortes_fino$sens - cortes_fino$esp)), ]
tasa_train <- mean(datos_train$y)cortes_fino %>%
filter(corte <= 0.6) %>% # mismo rango que la Tabla 35
select(corte, Sensibilidad = sens, Especificidad = esp, VPP = vpp,
`Empleados intervenidos` = intervenidos) %>%
pivot_longer(-corte, names_to = "Indicador", values_to = "valor") %>%
mutate(Indicador = factor(Indicador, levels = c("Sensibilidad", "Especificidad",
"VPP", "Empleados intervenidos"))) %>%
ggplot(aes(x = corte, y = valor, color = Indicador, linetype = Indicador)) +
geom_line(linewidth = 0.9, na.rm = TRUE) +
geom_vline(xintercept = 0.5, linetype = "dashed", color = "gray30") +
geom_vline(xintercept = tasa_train, linetype = "dotted", color = "gray30") +
scale_color_manual(values = c("Sensibilidad" = col_naranja, "Especificidad" = col_azul,
"VPP" = "gray45", "Empleados intervenidos" = "#2c3e50"),
name = NULL) +
scale_linetype_manual(values = c("Sensibilidad" = "solid", "Especificidad" = "solid",
"VPP" = "solid", "Empleados intervenidos" = "longdash"),
name = NULL) +
scale_x_continuous(breaks = seq(0, 0.6, 0.1)) +
scale_y_continuous(labels = percent, limits = c(0, 1)) +
labs(x = "Punto de corte", y = NULL) +
theme(legend.position = "bottom")
Figura 18. Sensibilidad, especificidad, valor predictivo positivo (VPP) y porcentaje de empleados intervenidos según el punto de corte, con las probabilidades fuera de pliegue de la muestra de entrenamiento. La línea discontinua indica el corte convencional de 0.5 y la punteada, la tasa de rotación de la muestra de entrenamiento. Se muestran cortes hasta 0.6, porque por encima de ese valor quedan muy pocos empleados clasificados en riesgo y el VPP deja de ser informativo.
Al analizar la Tabla 35 y la Figura 18, se destacan cuatro resultados:
- Validación de las probabilidades. El AUC de las probabilidades fuera de pliegue es 0.775, similar al estimado por validación cruzada en la Sección 5.6 (0.769), por lo que la elección del corte no está afectada por el optimismo de la evaluación dentro de muestra.
- Corte de 0.16. Al reducir el corte, la sensibilidad aumenta y la especificidad disminuye; ambas se igualan en un corte de 0.16, cercano a la tasa de rotación de la muestra (16.1%), con valores alrededor de 74%. Con ese corte se clasificaría en riesgo al 35.4% de los empleados. Esto concuerda con lo observado en la Sección 5.3, donde las distribuciones de probabilidad de quienes rotan y quienes no se cruzan cerca de la tasa de rotación y no en 0.5.
- Cortes de 0.3 y 0.5. Al pasar del corte de 0.5 al de 0.3, la sensibilidad sube de 17.5% a 48.2%, el porcentaje de empleados intervenidos pasa de 4.8% a 15.7% y el VPP baja de 59.2% a 49.4%: se identifica una proporción mucho mayor de las rotaciones con una pérdida moderada de precisión.
- Cortes superiores a 0.5. Con estos cortes se clasificaría en riesgo a menos del 5% de la plantilla, y a partir de 0.6 (2.9% de intervenidos) el VPP se calcula sobre muy pocos casos, por lo que pierde estabilidad. El modelo es muy preciso en esa zona, pero deja sin identificar a la gran mayoría de quienes rotan.
6.2 Elección del punto de corte
La Tabla 35 y la Figura 18 muestran que ningún corte maximiza a la vez la sensibilidad y la especificidad, por lo que la elección requiere un criterio explícito. En el caso de la rotación, dejar sin intervenir a un empleado que termina saliendo (falso negativo) suele ser más costoso que intervenir a uno que no iba a rotar (falso positivo), ya que una acción de retención cuesta menos que reemplazar a un empleado; por esta razón se busca un corte inferior a 0.5, que aumente la sensibilidad. Se adopta el corte que maximiza el índice de Youden (sensibilidad + especificidad - 1), que corresponde al punto de la curva ROC más alejado de la diagonal del azar, como se anticipó en las Secciones 2.1 y 4.10.
Para leer el corte en términos de gestión se calcula su razón de costos implícita: bajo el supuesto de que la intervención evita la salida del empleado, un corte \(c\) minimiza el costo esperado cuando perder a un empleado cuesta \(k = (1 - c)/c\) veces lo que cuesta intervenirlo; si la intervención solo retiene a una parte de los empleados intervenidos, la razón necesaria para justificar el corte aumenta en la misma proporción. El corte elegido se compara con los que resultarían de otros supuestos de costo (\(k\) = 2, 4 y 9) y con el corte convencional de 0.5.
Para evaluar la estabilidad de la elección se identifica el rango de cortes cuyo índice de Youden se ubica a menos de 0.01 del máximo; dentro de ese rango, los cortes tienen un desempeño prácticamente equivalente.
# Corte que maximiza el índice de Youden en la rejilla fina (probabilidades fuera de pliegue)
corte_youden <- cortes_fino$corte[which.max(cortes_fino$youden)]
youden_max <- max(cortes_fino$youden)
k_youden <- (1 - corte_youden) / corte_youden
# Rango de cortes cuyo Youden está a menos de 0.01 del máximo (estabilidad de la elección)
meseta_youden <- range(cortes_fino$corte[cortes_fino$youden >= youden_max - 0.01])
# Cortes que se comparan: Youden, tres supuestos de costo (c = 1 / (1 + k)) y el convencional
cortes_crit <- c(corte_youden, 1 / 3, 1 / 5, 1 / 10, 0.5)
tab_crit <- map_dfr(cortes_crit, ~ metricas_corte(datos_train$prob_oof, datos_train$y, .x)) %>%
mutate(Criterio = c("Máximo índice de Youden", "Razón de costos k = 2",
"Razón de costos k = 4", "Razón de costos k = 9",
"Corte convencional"),
k = (1 - corte) / corte,
.before = 1)
estilo_tabla(tab_crit %>%
transmute(Criterio,
Corte = f_num(corte, 2),
k = f_num(k, 1),
across(c(intervenidos, sens, esp, vpp), ~ .x * 100),
Youden = f_num(youden, 3)),
digits = 1,
caption = "Tabla 36. Punto de corte según el criterio de elección, con las probabilidades fuera de pliegue de la muestra de entrenamiento. k: razón de costos implícita, es decir, cuántas veces más cuesta perder a un empleado que intervenirlo (k = (1 - corte) / corte). Intervenidos: porcentaje de empleados clasificados en riesgo. VPP: valor predictivo positivo.",
col.names = c("Criterio", "Corte", "k", "Intervenidos (%)", "Sensibilidad (%)",
"Especificidad (%)", "VPP (%)", "Youden"))| Criterio | Corte | k | Intervenidos (%) | Sensibilidad (%) | Especificidad (%) | VPP (%) | Youden |
|---|---|---|---|---|---|---|---|
| Máximo índice de Youden | 0.16 | 5.2 | 35.4 | 74.1 | 72.0 | 33.7 | 0.461 |
| Razón de costos k = 2 | 0.33 | 2.0 | 12.2 | 42.8 | 93.6 | 56.3 | 0.364 |
| Razón de costos k = 4 | 0.20 | 4.0 | 28.9 | 66.9 | 78.4 | 37.2 | 0.452 |
| Razón de costos k = 9 | 0.10 | 9.0 | 51.7 | 80.7 | 53.9 | 25.2 | 0.347 |
| Corte convencional | 0.50 | 1.0 | 4.8 | 17.5 | 97.7 | 59.2 | 0.152 |
# Corte adoptado para el resto de la sección
corte_elegido <- corte_youden
m_eleg <- en_corte(corte_elegido)
# ---- Valores para el texto ----
m_sup <- en_corte(meseta_youden[2]) # extremo superior de la meseta
crit <- function(i) tab_crit[i, ] # fila i de la Tabla 36Al analizar la Tabla 36, se destacan cuatro resultados:
- El índice de Youden alcanza su máximo (0.461) en un corte de 0.16, el mismo en que se igualan la sensibilidad y la especificidad (Sección 6.1) y prácticamente igual a la tasa de rotación de la muestra (16.1%), como se anticipó en la Sección 2.1. Con este corte se identifica al 74.1% de las rotaciones y se clasifica en riesgo al 35.4% de los empleados.
- El índice se mantiene a menos de 0.01 del máximo entre los cortes de 0.16 y 0.25. Dentro de ese rango la capacidad de separar a los dos grupos es prácticamente la misma, pero la carga de intervenciones cambia de forma apreciable: con un corte de 0.25 se intervendría al 21.2% de los empleados, a cambio de reducir la sensibilidad a 59.0%.
- El corte elegido equivale a suponer que perder a un empleado cuesta 5.2 veces lo que cuesta intervenirlo. Con un supuesto más conservador (k = 2), el corte sería 0.33 y la sensibilidad bajaría a 42.8%; con uno más exigente (k = 9), el corte sería 0.10, la sensibilidad subiría a 80.7% y se intervendría al 51.7% de la plantilla. En ambos casos el índice de Youden es menor que en el corte elegido.
- Frente al corte convencional de 0.5, el corte elegido multiplica la sensibilidad por 4.2 (de 17.5% a 74.1%). El VPP baja de 59.2% a 33.7%, pero sigue siendo cerca del doble de la tasa de rotación, es decir, una de cada tres intervenciones alcanza a un empleado que efectivamente iba a rotar.
En consecuencia, se adopta un corte de 0.16. Si la capacidad de intervención de la empresa fuera limitada, el corte podría elevarse hasta 0.25 sin una pérdida apreciable en la capacidad de discriminación, aceptando una menor proporción de rotaciones identificadas.
6.3 Evaluación del corte
El corte de 0.16 se eligió con la muestra de entrenamiento; su desempeño con empleados nuevos se evalúa en la muestra de prueba, que no intervino en la estimación del modelo ni en la elección del corte. Se construye la matriz de confusión con ese corte y sus indicadores se comparan con los obtenidos en entrenamiento (Tabla 36) y con los del corte de 0.5 en la misma muestra de prueba (Sección 5.5). Si los valores de prueba son cercanos a los de entrenamiento, el corte se generaliza a datos nuevos.
# Clase predicha en la muestra de prueba con el corte elegido
datos_test <- datos_test %>%
mutate(pred_eleg = factor(ifelse(prob > corte_elegido, 1, 0), levels = c(0, 1)))
cm_eleg <- confusionMatrix(data = datos_test$pred_eleg, reference = datos_test$obs,
positive = "1")
VP_e <- cm_eleg$table["1", "1"]; FN_e <- cm_eleg$table["0", "1"]
FP_e <- cm_eleg$table["1", "0"]; VN_e <- cm_eleg$table["0", "0"]
# Tabla 37: matriz de confusión con el corte elegido (mismo formato que la Tabla 30)
mat_eleg <- tibble(
Real = c("Rota", "No rota", "Total"),
Pred_rota = c(paste0(VP_e, " (VP)"), paste0(FP_e, " (FP)"), VP_e + FP_e),
Pred_no = c(paste0(FN_e, " (FN)"), paste0(VN_e, " (VN)"), FN_e + VN_e),
Total = c(VP_e + FN_e, FP_e + VN_e, nrow(datos_test))
)
estilo_tabla(mat_eleg,
caption = paste0("Tabla 37. Matriz de confusión del modelo en la muestra de prueba con el corte elegido (", f_num(corte_elegido, 2), "). Filas: rotación observada; columnas: clasificación del modelo."),
col.names = c("Rotación observada", "Predicción: rota", "Predicción: no rota", "Total"))| Rotación observada | Predicción: rota | Predicción: no rota | Total |
|---|---|---|---|
| Rota | 52 (VP) | 19 (FN) | 71 |
| No rota | 121 (FP) | 248 (VN) | 369 |
| Total | 173 | 267 | 440 |
# Indicadores de una matriz de confusión de caret
ind_cm <- function(cm) {
ov <- cm$overall; bc <- cm$byClass; tb <- cm$table
c(`Empleados intervenidos` = sum(tb["1", ]) / sum(tb),
Sensibilidad = bc[["Sensitivity"]],
Especificidad = bc[["Specificity"]],
VPP = bc[["Pos Pred Value"]],
VPN = bc[["Neg Pred Value"]],
`Exactitud balanceada` = bc[["Balanced Accuracy"]],
Exactitud = ov[["Accuracy"]],
Kappa = ov[["Kappa"]])
}
# Matriz del corte elegido en entrenamiento (probabilidades fuera de pliegue)
cm_eleg_tr <- confusionMatrix(
data = factor(ifelse(datos_train$prob_oof > corte_elegido, 1, 0), levels = c(0, 1)),
reference = factor(datos_train$y, levels = c(0, 1)),
positive = "1")
# Tabla 38: comparación de indicadores
fmt_ind <- function(x, ind) ifelse(ind == "Kappa", f_num(x, 3), percent(x, accuracy = 0.1))
indicadores <- names(ind_cm(cm_eleg))
tab_comp_corte <- tibble(
Indicador = indicadores,
c05_pr = fmt_ind(ind_cm(cm_05), indicadores),
eleg_tr = fmt_ind(ind_cm(cm_eleg_tr), indicadores),
eleg_pr = fmt_ind(ind_cm(cm_eleg), indicadores)
)
estilo_tabla(tab_comp_corte,
caption = paste0("Tabla 38. Indicadores de clasificación con el corte de 0.5 en la muestra de prueba y con el corte elegido (", f_num(corte_elegido, 2), ") en las muestras de entrenamiento (probabilidades fuera de pliegue) y de prueba (clase positiva: rota)."),
col.names = c("Indicador", "Corte 0.5 (prueba)",
paste0("Corte ", f_num(corte_elegido, 2), " (entrenamiento)"),
paste0("Corte ", f_num(corte_elegido, 2), " (prueba)")))| Indicador | Corte 0.5 (prueba) | Corte 0.16 (entrenamiento) | Corte 0.16 (prueba) |
|---|---|---|---|
| Empleados intervenidos | 5.0% | 35.4% | 39.3% |
| Sensibilidad | 16.9% | 74.1% | 73.2% |
| Especificidad | 97.3% | 72.0% | 67.2% |
| VPP | 54.5% | 33.7% | 30.1% |
| VPN | 85.9% | 93.5% | 92.9% |
| Exactitud balanceada | 57.1% | 73.0% | 70.2% |
| Exactitud | 84.3% | 72.3% | 68.2% |
| Kappa | 0.197 | 0.311 | 0.256 |
# ---- Valores para el texto ----
ie_pr <- ind_cm(cm_eleg) # indicadores del corte elegido en prueba
ie_tr <- ind_cm(cm_eleg_tr) # indicadores del corte elegido en entrenamiento
i05_pr <- ind_cm(cm_05) # indicadores del corte 0.5 en prueba
ic_sens_e <- binom.test(VP_e, VP_e + FN_e)$conf.int # IC 95% exacto de la sensibilidad
ic_esp_e <- binom.test(VN_e, VN_e + FP_e)$conf.int # IC 95% exacto de la especificidad
nir_pr <- cm_eleg$overall[["AccuracyNull"]]
# Costo del cambio de corte: intervenciones adicionales por rotación adicional detectada
extra_det <- VP_e - VP_05
extra_int <- (VP_e + FP_e) - (VP_05 + FP_05)
int_por_det <- extra_int / extra_detAl analizar las Tablas 37 y 38, se destacan cinco resultados:
- El corte se generaliza a la muestra de prueba. La sensibilidad es 73.2% (IC 95%: 61.4% a 83.1%), prácticamente igual a la de entrenamiento (74.1%). La especificidad baja de 72.0% a 67.2% (IC 95%: 62.2% a 72.0%), por lo que se clasifica en riesgo al 39.3% de los empleados, frente al 35.4% de entrenamiento; la diferencia es compatible con la variación muestral.
- Frente al corte de 0.5, el modelo identifica a 52 de las 71 rotaciones, en lugar de 12. La exactitud balanceada pasa de 57.1% a 70.2% y el kappa de 0.197 (acuerdo leve) a 0.256 (acuerdo aceptable).
- Para detectar 40 rotaciones adicionales se requieren 151 intervenciones adicionales (173 frente a 22), es decir, 3.8 intervenciones por cada rotación adicional identificada. Esta cifra es inferior a la razón de costos implícita en el corte (k = 5.2), de modo que el cambio de corte se justifica siempre que perder a un empleado cueste más de 3.8 veces lo que cuesta intervenirlo, suponiendo que la intervención evite la salida; si solo la evita en una fracción de los casos, ese umbral aumenta en la misma proporción.
- Entre los 267 empleados clasificados como “no rota”, solo el 7.1% terminó rotando, menos de la mitad de la tasa de rotación de la muestra (16.1%). Con el corte de 0.5, esa proporción era 14.1% (Sección 5.5); con el corte elegido, la clasificación “no rota” sí aporta información sobre el riesgo del empleado.
- La exactitud baja a 68.2%, por debajo de la tasa de no información (83.9%). Este resultado es esperable: el corte se eligió para identificar rotaciones y no para maximizar los aciertos globales, que con una clase minoritaria de 16% se maximizan clasificando a todos los empleados como “no rota”.
En consecuencia, el corte de 0.16 mantiene en datos nuevos el desempeño observado en entrenamiento: identifica cerca de tres de cada cuatro rotaciones y aproximadamente una de cada 3.3 intervenciones alcanza a un empleado que efectivamente iba a rotar.
6.4 Predicción para un empleado hipotético
La Sección 4.9 mostró cómo cambia la probabilidad de rotar al modificar una sola característica del perfil base. Aquí se aplica el modelo a un empleado hipotético que reúne varias características asociadas con mayor riesgo: hace horas extra, es soltero, viaja frecuentemente, tiene un ingreso mensual en el percentil 25 (2.9 mil), lleva un año en el cargo y vive a la distancia mediana (7 unidades). Su probabilidad de rotar se estima con el modelo de la Sección 4, con un intervalo de confianza al 95% calculado en la escala logit y transformado luego a probabilidad, y se compara con el corte de 0.16: si la supera, se recomienda intervenir.
Para orientar la intervención, se estima además cómo cambiaría la probabilidad si la empresa actuara sobre las palancas de gestión del perfil (Sección 3.5): eliminar las horas extra, reducir los viajes a la categoría “Raramente” y llevar el ingreso a la mediana de la plantilla, primero por separado y luego de forma conjunta. Estos escenarios describen asociaciones estimadas por el modelo y no efectos causales garantizados.
# Perfil del empleado hipotético: se parte de 'perfil_base' (Sección 4.9) y se modifican
# sus características con las funciones nivel() y q() definidas en esa sección
empleado <- perfil_base %>%
mutate(horas_extra = nivel("horas_extra", "Si"),
estado_civil = nivel("estado_civil", "Soltero"),
viaje_negocios = nivel("viaje_negocios", "Frecuentemente"),
ingreso_miles = q(datos$ingreso_miles, 0.25),
antiguedad_cargo = 1)
# Escenarios de intervención sobre las palancas de gestión
lista_int <- list(
actual = empleado,
sin_he = mutate(empleado, horas_extra = nivel("horas_extra", "No")),
viaje_rar = mutate(empleado, viaje_negocios = nivel("viaje_negocios", "Raramente")),
ing_med = mutate(empleado, ingreso_miles = median(datos$ingreso_miles)),
conjunto = mutate(empleado, horas_extra = nivel("horas_extra", "No"),
viaje_negocios = nivel("viaje_negocios", "Raramente"),
ingreso_miles = median(datos$ingreso_miles))
)
etiq_int <- c(
actual = "Empleado hipotético (situación actual)",
sin_he = "Sin horas extra",
viaje_rar = "Viaja raramente",
ing_med = paste0("Ingreso en la mediana (", f_num(median(datos$ingreso_miles), 1), " mil)"),
conjunto = "Las tres acciones combinadas"
)
# Probabilidad con IC 95%: predict(type = "link", se.fit = TRUE) devuelve el logit y su
# error estándar; el intervalo se construye en esa escala y plogis() lo lleva a probabilidad
pred_ic <- function(nd) {
p <- predict(modelo_final, newdata = nd, type = "link", se.fit = TRUE)
z <- qnorm(0.975)
tibble(prob = plogis(p$fit),
inf = plogis(p$fit - z * p$se.fit),
sup = plogis(p$fit + z * p$se.fit))
}
tab_emp <- imap_dfr(lista_int, ~ pred_ic(.x) %>% mutate(clave = .y, .before = 1)) %>%
mutate(Escenario = unname(etiq_int[clave]),
cambio = (prob - prob[clave == "actual"]) * 100,
decision = ifelse(prob > corte_elegido, "Intervenir", "No intervenir"))
estilo_tabla(tab_emp %>%
transmute(Escenario, prob = prob * 100, inf = inf * 100, sup = sup * 100,
cambio = ifelse(clave == "actual", NA, cambio), decision),
digits = 1,
caption = paste0("Tabla 39. Probabilidad estimada de rotar (%) del empleado hipotético en su situación actual y bajo escenarios de intervención sobre las palancas de gestión, con IC 95% (modelo de la Sección 4). Decisión: intervenir si la probabilidad supera el corte de ", f_num(corte_elegido, 2), ". p.p.: puntos porcentuales."),
col.names = c("Escenario", "Probabilidad (%)", "IC 95% inf.", "IC 95% sup.",
"Cambio frente a la situación actual (p.p.)", "Decisión"))| Escenario | Probabilidad (%) | IC 95% inf. | IC 95% sup. | Cambio frente a la situación actual (p.p.) | Decisión |
|---|---|---|---|---|---|
| Empleado hipotético (situación actual) | 70.2 | 60.6 | 78.3 | Intervenir | |
| Sin horas extra | 35.8 | 27.5 | 45.1 | -34.3 | Intervenir |
| Viaja raramente | 55.4 | 46.9 | 63.6 | -14.8 | Intervenir |
| Ingreso en la mediana (4.9 mil) | 65.7 | 55.7 | 74.4 | -4.5 | Intervenir |
| Las tres acciones combinadas | 19.3 | 15.1 | 24.5 | -50.8 | Intervenir |
# ---- Valores para el texto ----
te <- tab_emp %>% select(clave, prob, inf, sup, cambio, decision) %>%
split(.$clave) # te$actual, te$sin_he, te$viaje_rar, te$ing_med, te$conjunto
# Suma de las reducciones individuales (para contrastar con el escenario conjunto)
suma_ind <- -(te$sin_he$cambio + te$viaje_rar$cambio + te$ing_med$cambio)
# Probabilidad del empleado frente a la tasa de rotación de la plantilla
razon_tasa <- te$actual$prob / p_rotAl analizar la Tabla 39, se destacan cuatro resultados:
- La probabilidad estimada de que el empleado rote es 70.2% (IC 95%: 60.6% a 78.3%), 4.4 veces la tasa de rotación de la plantilla. Incluso el límite inferior del intervalo supera ampliamente el corte de 0.16, por lo que la recomendación de intervenir no depende de la incertidumbre de la estimación.
- Las palancas de gestión tienen efectos de magnitud muy distinta. Eliminar las horas extra reduce la probabilidad en 34.3 p.p., hasta 35.8%; reducir los viajes la baja en 14.8 p.p., y llevar el ingreso a la mediana, en solo 4.5 p.p. Para este perfil, la acción prioritaria es la redistribución de la carga de horas extra.
- Con las tres acciones combinadas, la probabilidad baja a 19.3% (IC 95%: 15.1% a 24.5%), un valor cercano a la tasa de rotación de la plantilla (16.1%) que, sin embargo, sigue por encima del corte. El riesgo remanente proviene de características que la empresa no puede modificar (ser soltero y llevar un año en el cargo), de modo que la intervención reduce el riesgo del empleado a un nivel promedio, pero no elimina la necesidad de seguimiento.
- La reducción conjunta (50.8 p.p.) es menor que la suma de las reducciones individuales (53.6 p.p.). Esto se debe a que el modelo es aditivo en el logit y no en la probabilidad: como se advirtió en la Sección 4.9, un mismo efecto sobre los odds produce cambios menores en puntos porcentuales cuando la probabilidad es más baja.
En este sentido, el modelo permite no solo identificar al empleado como candidato a intervención, sino también priorizar las acciones según su efecto estimado. Estos resultados describen asociaciones y no garantizan que la modificación de cada condición produzca el cambio estimado; su verificación requeriría el seguimiento de los empleados intervenidos.
6.5 Síntesis del punto de corte y las predicciones
Con una tasa de rotación de 16.1%, el corte convencional de 0.5 deja sin identificar a la mayoría de los empleados que rotan. Por ello se eligió, con el índice de Youden, un corte de 0.16, que coincide con la tasa de rotación de la plantilla y equivale a suponer que perder a un empleado cuesta unas 5 veces lo que cuesta intervenirlo. En la muestra de prueba, este corte identifica al 73.2% de las rotaciones, frente al 16.9% del corte de 0.5, a cambio de intervenir al 39.3% de la plantilla. Aplicado a un empleado hipotético que hace horas extra, es soltero, viaja frecuentemente, tiene un ingreso bajo y lleva un año en el cargo, el modelo estima una probabilidad de rotar de 70.2%, por lo que se recomienda intervenirlo; la acción con mayor efecto estimado es eliminar sus horas extra, que reduce esa probabilidad a la mitad. Deben tenerse presentes tres limitaciones:
- El corte refleja un supuesto implícito sobre los costos de la rotación; si la capacidad de intervención de la empresa es limitada, puede elevarse hasta 0.25 sin una pérdida apreciable de discriminación, aunque con menor sensibilidad.
- La sensibilidad en la muestra de prueba se estimó con solo 71 rotaciones, por lo que su intervalo de confianza es amplio (61.4% a 83.1%).
- Los escenarios de intervención describen asociaciones estimadas por el modelo, no efectos causales garantizados.
7 Conclusiones y estrategia de retención
7.1 Conclusiones
La rotación de empleados en la organización se asocia principalmente con la carga de trabajo, la frecuencia de viajes, la situación personal y las condiciones de compensación y permanencia en el cargo. En el modelo ajustado, las seis variables seleccionadas aportan de forma significativa y sus efectos tienen el signo planteado en las hipótesis de la Sección 1. Hacer horas extra es el factor con mayor asociación (OR = 4.21), seguido de viajar frecuentemente y de ser soltero; un mayor ingreso y una mayor antigüedad en el cargo se asocian con menor probabilidad de rotar, y la distancia a casa, con una probabilidad mayor. El efecto del ingreso y de la antigüedad no es constante: se concentra en los salarios más bajos y en los primeros años en el cargo (Secciones 4.8 y 5.8).
El modelo tiene una capacidad de discriminación aceptable (AUC de 0.769 en validación cruzada) y probabilidades bien calibradas. Con el corte de 0.16, identifica al 73.2% de los empleados que rotan en la muestra de prueba, frente al 16.9% del corte convencional de 0.5. En consecuencia, el modelo es útil para priorizar a los empleados según su riesgo y orientar las acciones de retención, aunque no permite anticipar con certeza qué empleado rotará.
7.2 Estrategia de retención
La estrategia se construye a partir de las variables asociadas significativamente con la rotación en el análisis bivariado (Sección 3.5). Entre ellas, las palancas de gestión, es decir, los factores sobre los que la empresa puede actuar, son: horas extra, ingreso mensual, viaje de negocios, satisfacción laboral, satisfacción ambiental, equilibrio trabajo-vida, capacitaciones. Las acciones se priorizan según la magnitud de su asociación en el modelo ajustado (Sección 4.3) y según su efecto estimado en el empleado hipotético (Sección 6.4). La Tabla 40 resume la estrategia.
or_aj <- function(t) coef_multi$OR[coef_multi$term == t] # OR ajustado de un término
signo_cap <- tamizaje$Signo[tamizaje$variable == "capacitaciones"] # signo en la Tabla 15
estrategia <- tribble(
~Linea, ~Variables, ~Evidencia, ~Accion, ~Prioridad,
"Carga de trabajo", "Horas extra",
paste0("OR ajustado = ", f_num(or_aj("horas_extraSi")),
"; sin horas extra, la probabilidad del empleado hipotético pasa de ",
pct(te$actual$prob), " a ", pct(te$sin_he$prob), " (Sección 6.4)"),
"Redistribuir la carga de horas extra: contratar personal en las áreas con más horas extra, fijar límites y rotar los turnos",
"Alta",
"Viajes de negocios", "Viaje de negocios",
paste0("OR ajustado = ", f_num(or_aj("viaje_negociosFrecuentemente")),
" (frecuentemente) y ", f_num(or_aj("viaje_negociosRaramente")),
" (raramente) frente a no viajar; dirección robusta, magnitud imprecisa (Sección 5.2)"),
"Rotar los viajes entre empleados, sustituir parte de ellos por reuniones virtuales y compensar el tiempo de viaje",
"Media",
"Compensación", "Ingreso mensual",
paste0("OR ajustado = ", f_num(or_aj("ingreso_miles"), 3),
" por cada 1.000; efecto más intenso en los salarios bajos (Sección 5.8)"),
"Ajustes salariales focalizados en los ingresos más bajos, en lugar de aumentos generales",
"Media",
"Primeros años en el cargo", "Antigüedad en el cargo",
paste0("OR ajustado = ", f_num(or_aj("antiguedad_cargo"), 3),
" por año; riesgo concentrado en los primeros años (Secciones 3.2 y 5.8)"),
"Programas de inducción, acompañamiento y plan de carrera durante los primeros años en el cargo",
"Media",
"Ambiente laboral", "Satisfacción laboral, satisfacción ambiental y equilibrio trabajo-vida",
"Asociación bivariada significativa (Sección 3.5); no incluidas en el modelo",
"Seguimiento del clima laboral, flexibilidad horaria y acciones de bienestar dirigidas a los empleados con baja satisfacción",
"Complementaria",
"Desarrollo", "Capacitaciones",
paste0("Asociación bivariada significativa, con signo ", tolower(signo_cap),
" (Sección 3.5); no incluida en el modelo"),
"Ampliar el acceso a capacitación, en especial para los empleados que han recibido menos formación",
"Complementaria"
)
estilo_tabla(estrategia,
caption = "Tabla 40. Estrategia de retención: líneas de acción, variables en que se sustentan, evidencia estadística, acciones propuestas y prioridad.",
col.names = c("Línea de acción", "Variables", "Evidencia", "Acción propuesta", "Prioridad")) %>%
column_spec(3, width = "16em") %>% # ancho de la columna de evidencia
column_spec(4, width = "18em") # ancho de la columna de acciones| Línea de acción | Variables | Evidencia | Acción propuesta | Prioridad |
|---|---|---|---|---|
| Carga de trabajo | Horas extra | OR ajustado = 4.21; sin horas extra, la probabilidad del empleado hipotético pasa de 70.2% a 35.8% (Sección 6.4) | Redistribuir la carga de horas extra: contratar personal en las áreas con más horas extra, fijar límites y rotar los turnos | Alta |
| Viajes de negocios | Viaje de negocios | OR ajustado = 3.85 (frecuentemente) y 2.03 (raramente) frente a no viajar; dirección robusta, magnitud imprecisa (Sección 5.2) | Rotar los viajes entre empleados, sustituir parte de ellos por reuniones virtuales y compensar el tiempo de viaje | Media |
| Compensación | Ingreso mensual | OR ajustado = 0.903 por cada 1.000; efecto más intenso en los salarios bajos (Sección 5.8) | Ajustes salariales focalizados en los ingresos más bajos, en lugar de aumentos generales | Media |
| Primeros años en el cargo | Antigüedad en el cargo | OR ajustado = 0.898 por año; riesgo concentrado en los primeros años (Secciones 3.2 y 5.8) | Programas de inducción, acompañamiento y plan de carrera durante los primeros años en el cargo | Media |
| Ambiente laboral | Satisfacción laboral, satisfacción ambiental y equilibrio trabajo-vida | Asociación bivariada significativa (Sección 3.5); no incluidas en el modelo | Seguimiento del clima laboral, flexibilidad horaria y acciones de bienestar dirigidas a los empleados con baja satisfacción | Complementaria |
| Desarrollo | Capacitaciones | Asociación bivariada significativa, con signo negativo (Sección 3.5); no incluida en el modelo | Ampliar el acceso a capacitación, en especial para los empleados que han recibido menos formación | Complementaria |
La prioridad alta de la carga de trabajo responde a que las horas extra son el factor con mayor asociación con la rotación y el que más reduce el riesgo estimado del empleado hipotético. Las líneas de prioridad media tienen efectos de menor magnitud o menos precisos, y el ambiente laboral y el desarrollo se consideran complementarios porque su evidencia proviene solo del análisis bivariado.
Para operar la estrategia, el modelo se aplica periódicamente a la plantilla y se prioriza a los empleados cuya probabilidad estimada supere el corte de 0.16. Las variables de perfil del empleado significativas en la Sección 3.5 (años de experiencia, estado civil, edad, campo de educación) sirven para focalizar las acciones, no como objeto de intervención. Dado que el modelo describe asociaciones y no efectos causales, se recomienda implementar las acciones primero como un piloto, con un grupo de comparación, para verificar su efecto real sobre la rotación antes de extenderlas a toda la organización.