Guía de clase: Modelo de conteo Poisson aplicado al conteo de meses sin cotizar al seguro social

1 Presentación

En el análisis del mercado laboral es frecuente estudiar fenómenos expresados como conteos. A diferencia de una variable dicotómica, que registra si ocurre o no un evento, una variable de conteo mide cuántas veces ocurre un evento en un periodo determinado.

En esta guía se desarrolla la aplicación del modelo de regresión de Poisson para analizar el número de meses sin cotizar al seguro social. El ejemplo se plantea desde una perspectiva laboral, considerando que el conteo de meses sin cotización puede interpretarse como una medida de discontinuidad contributiva o fragilidad en la trayectoria de formalidad laboral.

La unidad de análisis será el trabajador asalariado observado durante un periodo de referencia. La variable dependiente será el conteo de meses en los que no registra cotización al seguro social. El objetivo será explicar cómo dicho conteo varía en función de características sociodemográficas y laborales.

En términos didácticos, esta guía mantiene una estructura similar a la utilizada para la regresión logística multinomial: primero se presenta el fundamento estadístico del modelo, luego se prepara una base de datos, se ajusta el modelo, se interpretan sus parámetros, se realizan predicciones y se evalúa el ajuste.

2 Objetivos de aprendizaje

Al finalizar la guía, se espera que el estudiante sea capaz de:

1. Reconocer situaciones del mercado laboral en las que corresponde utilizar modelos de conteo.
2. Comprender la lógica probabilística del modelo de Poisson y su relación con la media esperada del conteo.
3. Ajustar un modelo de regresión de Poisson en R mediante la función glm().
4. Interpretar los parámetros estimados del modelo en términos del signo, de la razón de incidencias y de la variación esperada del conteo.
5. Generar predicciones para perfiles laborales específicos mediante el argumento newdata.
6. Evaluar supuestos básicos del modelo, especialmente la posible presencia de sobredispersión.
7. Identificar alternativas metodológicas cuando el modelo de Poisson no resulta adecuado.

3 Situación de análisis

Supóngase que se dispone de información sobre trabajadores asalariados del sector privado y se desea analizar cuántos meses no cotizaron al seguro social durante los últimos 24 meses. La variable dependiente puede tomar valores enteros no negativos:

\[ Y_i = 0, 1, 2, 3, \ldots \]

donde \(Y_i\) representa el número de meses sin cotización del trabajador \(i\).

Desde el punto de vista sustantivo, valores más altos de \(Y_i\) indican mayor discontinuidad contributiva. En cambio, valores cercanos a cero indican mayor estabilidad en la cotización al seguro social.

En este contexto, algunas preguntas de interés son:

• ¿Los trabajadores con contrato temporal presentan mayor cantidad esperada de meses sin cotizar?
• ¿La antigüedad laboral reduce la cantidad esperada de meses sin cotización?
• ¿Existen diferencias según nivel educativo, sexo o tamaño de empresa?
• ¿Qué perfil laboral presenta mayor riesgo de acumular meses sin cotización?

4 Fundamento teórico del modelo de Poisson

4.1 Variable de conteo

Una variable de conteo registra la cantidad de veces que ocurre un evento en una unidad de observación o en un intervalo de tiempo. En esta guía, el evento corresponde a un mes sin cotización al seguro social.

La variable dependiente debe cumplir tres características básicas:

• toma valores enteros;
• no admite valores negativos;
• representa una frecuencia o número de ocurrencias.

Por ejemplo, para el caso de los meses sin cotizar:

\[ Y_i \in \{0,1,2,3,\ldots,24\} \]

si el periodo de observación corresponde a 24 meses.

4.2 Distribución de Poisson

La distribución de Poisson se utiliza para modelar conteos. Si la variable aleatoria \(Y_i\) sigue una distribución de Poisson con parámetro \(\lambda_i\), se escribe:

\[ Y_i \sim Poisson(\lambda_i) \]

La función de probabilidad es:

\[ P(Y_i = y_i) = \frac{e^{-\lambda_i}\lambda_i^{y_i}}{y_i!}, \quad y_i = 0,1,2,\ldots \]

El parámetro \(\lambda_i\) representa tanto la media como la varianza de la distribución:

\[ E(Y_i) = \lambda_i \]

\[ Var(Y_i) = \lambda_i \]

Esta igualdad entre media y varianza es conocida como el supuesto de equidispersión. En datos laborales reales, este supuesto puede no cumplirse, por lo que será necesario evaluarlo.

4.3 Modelo de regresión de Poisson

El modelo de regresión de Poisson relaciona la media esperada del conteo con un conjunto de variables explicativas. Debido a que \(\lambda_i\) debe ser positiva, se utiliza una función de enlace logarítmica:

\[ \log(\lambda_i) = \beta_0 + \beta_1X_{1i} + \beta_2X_{2i} + \cdots + \beta_kX_{ki} \]

Por lo tanto:

\[ \lambda_i = E(Y_i|X_i) = exp(\beta_0 + \beta_1X_{1i} + \cdots + \beta_kX_{ki}) \]

La expresión anterior indica que el modelo no estima directamente el conteo observado, sino el valor esperado del conteo para un perfil determinado de trabajador.

4.4 Supuestos principales

Para aplicar el modelo de Poisson deben considerarse los siguientes supuestos:

1. La variable dependiente es un conteo no negativo.
2. Las observaciones son independientes entre sí.
3. La media condicional del conteo está correctamente especificada.
4. La media y la varianza condicional son iguales o aproximadamente iguales.
5. No existe un exceso de ceros que el modelo no pueda representar adecuadamente.

Cuando la varianza es mayor que la media, se presenta sobredispersión. En ese caso, puede utilizarse un modelo cuasi-Poisson o un modelo binomial negativo.

5 Preparación del ambiente de trabajo

Para el desarrollo de la guía se utilizarán las siguientes librerías de R.

# Instalar previamente si fuera necesario:
# install.packages(c("dplyr", "ggplot2", "knitr", "broom", "MASS"))

library(dplyr)
library(ggplot2)
library(knitr)
library(broom)
library(MASS)

6 Construcción de una base de datos ficticia

En ausencia de una base administrativa real, se construye una base de datos ficticia con fines didácticos. La base representa trabajadores asalariados observados durante los últimos 24 meses.

La variable dependiente será meses_sin_cotizar, que registra el número de meses sin cotización al seguro social.

set.seed(123)

n <- 1200

datos <- data.frame(
  sexo = sample(c("Hombre", "Mujer"), n, replace = TRUE, prob = c(0.52, 0.48)),
  edad = round(pmin(pmax(rnorm(n, mean = 36, sd = 11), 18), 65)),
  nivel_educ = sample(
    c("Hasta primaria", "Secundaria", "Terciaria"),
    n, replace = TRUE, prob = c(0.25, 0.50, 0.25)
  ),
  tipo_contrato = sample(
    c("Permanente", "Temporal", "Sin contrato escrito"),
    n, replace = TRUE, prob = c(0.48, 0.30, 0.22)
  ),
  tamano_empresa = sample(
    c("Micro", "Pequeña", "Mediana/Grande"),
    n, replace = TRUE, prob = c(0.45, 0.35, 0.20)
  ),
  sector = sample(
    c("Industria", "Comercio", "Servicios", "Construcción"),
    n, replace = TRUE, prob = c(0.18, 0.32, 0.38, 0.12)
  ),
  antiguedad = round(pmax(rnorm(n, mean = 4, sd = 3), 0), 1)
)

# Definir categorías de referencia
datos <- datos %>%
  mutate(
    sexo = factor(sexo, levels = c("Hombre", "Mujer")),
    nivel_educ = factor(nivel_educ, levels = c("Terciaria", "Secundaria", "Hasta primaria")),
    tipo_contrato = factor(tipo_contrato, levels = c("Permanente", "Temporal", "Sin contrato escrito")),
    tamano_empresa = factor(tamano_empresa, levels = c("Mediana/Grande", "Pequeña", "Micro")),
    sector = factor(sector, levels = c("Industria", "Comercio", "Servicios", "Construcción"))
  )

# Generación del conteo esperado de meses sin cotizar
eta <- -0.50 +
  0.015 * (datos$edad - 35) +
  0.15 * (datos$sexo == "Mujer") +
  0.18 * (datos$nivel_educ == "Secundaria") +
  0.42 * (datos$nivel_educ == "Hasta primaria") +
  0.45 * (datos$tipo_contrato == "Temporal") +
  0.85 * (datos$tipo_contrato == "Sin contrato escrito") +
  0.18 * (datos$tamano_empresa == "Pequeña") +
  0.35 * (datos$tamano_empresa == "Micro") +
  0.15 * (datos$sector == "Comercio") +
  0.05 * (datos$sector == "Servicios") +
  0.25 * (datos$sector == "Construcción") -
  0.07 * datos$antiguedad

lambda <- exp(eta)

datos$meses_sin_cotizar <- rpois(n, lambda = lambda)

# Como se observa un periodo máximo de 24 meses, se limita el valor máximo.
datos$meses_sin_cotizar <- pmin(datos$meses_sin_cotizar, 24)

head(datos)

##     sexo edad     nivel_educ tipo_contrato tamano_empresa       sector
## 1 Hombre   48     Secundaria    Permanente          Micro     Comercio
## 2  Mujer   36      Terciaria    Permanente          Micro    Servicios
## 3 Hombre   36     Secundaria    Permanente          Micro Construcción
## 4  Mujer   19     Secundaria    Permanente        Pequeña Construcción
## 5  Mujer   45     Secundaria    Permanente          Micro     Comercio
## 6 Hombre   34 Hasta primaria    Permanente        Pequeña Construcción
##   antiguedad meses_sin_cotizar
## 1        6,0                 0
## 2        1,3                 2
## 3        7,9                 2
## 4        0,0                 3
## 5        0,7                 4
## 6        3,0                 0

La base contiene variables que permiten aproximar características laborales asociadas a la continuidad contributiva. En un análisis con datos reales, estas variables deberían provenir de registros administrativos, encuestas laborales o de una integración entre ambas fuentes.

7 Exploración descriptiva de la variable dependiente

Antes de ajustar el modelo, es necesario analizar el comportamiento de la variable de conteo.

resumen_y <- datos %>%
  summarise(
    observaciones = n(),
    media = mean(meses_sin_cotizar),
    varianza = var(meses_sin_cotizar),
    minimo = min(meses_sin_cotizar),
    maximo = max(meses_sin_cotizar),
    proporcion_ceros = mean(meses_sin_cotizar == 0)
  )

kable(resumen_y, digits = 3, caption = "Resumen del conteo de meses sin cotizar")

Resumen del conteo de meses sin cotizar

observaciones	media	varianza	minimo	maximo	proporcion_ceros
1200	1,304	1,693	0	8	0,322

La comparación entre media y varianza es importante para evaluar preliminarmente el supuesto de equidispersión. Si la varianza es considerablemente mayor que la media, puede existir sobredispersión.

ggplot(datos, aes(x = meses_sin_cotizar)) +
  geom_bar() +
  labs(
    title = "Distribución del número de meses sin cotizar",
    x = "Meses sin cotizar",
    y = "Cantidad de trabajadores"
  ) +
  theme_minimal()

8 Análisis descriptivo bivariado

Como primera aproximación, se puede comparar el promedio de meses sin cotizar según algunas características laborales.

tabla_contrato <- datos %>%
  group_by(tipo_contrato) %>%
  summarise(
    trabajadores = n(),
    promedio_meses_sin_cotizar = mean(meses_sin_cotizar),
    mediana = median(meses_sin_cotizar),
    proporcion_ceros = mean(meses_sin_cotizar == 0),
    .groups = "drop"
  )

kable(tabla_contrato, digits = 3, caption = "Meses sin cotizar según tipo de contrato")

Meses sin cotizar según tipo de contrato

tipo_contrato	trabajadores	promedio_meses_sin_cotizar	mediana	proporcion_ceros
Permanente	594	0,973	1	0,396
Temporal	343	1,376	1	0,292
Sin contrato escrito	263	1,958	2	0,198

ggplot(tabla_contrato, aes(x = tipo_contrato, y = promedio_meses_sin_cotizar)) +
  geom_col() +
  labs(
    title = "Promedio de meses sin cotizar según tipo de contrato",
    x = "Tipo de contrato",
    y = "Promedio de meses sin cotizar"
  ) +
  theme_minimal()

Este análisis permite identificar diferencias descriptivas. Sin embargo, no controla simultáneamente por otras características. Para ello se ajusta un modelo de regresión de Poisson.

9 Especificación del modelo

El modelo a estimar será:

\[ \log(\lambda_i) = \beta_0 + \beta_1 sexo_i + \beta_2 edad_i + \beta_3 nivel\_educ_i + \beta_4 tipo\_contrato_i + \beta_5 tamano\_empresa_i + \beta_6 sector_i + \beta_7 antiguedad_i \]

Donde:

\[ \lambda_i = E(Y_i|X_i) \]

representa la cantidad esperada de meses sin cotizar para el trabajador \(i\), dadas sus características.

10 Ajuste del modelo de Poisson en R

El modelo de Poisson se ajusta mediante la función glm(), especificando la familia poisson y el enlace logarítmico.

modelo_poisson <- glm(
  meses_sin_cotizar ~ sexo + edad + nivel_educ + tipo_contrato +
    tamano_empresa + sector + antiguedad,
  family = poisson(link = "log"),
  data = datos
)

summary(modelo_poisson)

## 
## Call:
## glm(formula = meses_sin_cotizar ~ sexo + edad + nivel_educ + 
##     tipo_contrato + tamano_empresa + sector + antiguedad, family = poisson(link = "log"), 
##     data = datos)
## 
## Coefficients:
##                                    Estimate Std. Error z value
## (Intercept)                       -0,862963   0,147713  -5,842
## sexoMujer                          0,108279   0,050827   2,130
## edad                               0,015950   0,002359   6,762
## nivel_educSecundaria               0,140296   0,064876   2,163
## nivel_educHasta primaria           0,307440   0,071760   4,284
## tipo_contratoTemporal              0,384224   0,062271   6,170
## tipo_contratoSin contrato escrito  0,765528   0,060926  12,565
## tamano_empresaPequeña              0,212522   0,078717   2,700
## tamano_empresaMicro                0,348614   0,074682   4,668
## sectorComercio                     0,102079   0,075949   1,344
## sectorServicios                   -0,067898   0,075028  -0,905
## sectorConstrucción                 0,277147   0,091112   3,042
## antiguedad                        -0,072722   0,009737  -7,468
##                                               Pr(>|z|)    
## (Intercept)                         0,0000000051528454 ***
## sexoMujer                                      0,03314 *  
## edad                                0,0000000000135694 ***
## nivel_educSecundaria                           0,03058 *  
## nivel_educHasta primaria            0,0000183310752937 ***
## tipo_contratoTemporal               0,0000000006821668 ***
## tipo_contratoSin contrato escrito < 0,0000000000000002 ***
## tamano_empresaPequeña                          0,00694 ** 
## tamano_empresaMicro                 0,0000030419690403 ***
## sectorComercio                                 0,17894    
## sectorServicios                                0,36548    
## sectorConstrucción                             0,00235 ** 
## antiguedad                          0,0000000000000811 ***
## ---
## Signif. codes:  0 '***' 0,001 '**' 0,01 '*' 0,05 '.' 0,1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 1695,3  on 1199  degrees of freedom
## Residual deviance: 1387,5  on 1187  degrees of freedom
## AIC: 3411,1
## 
## Number of Fisher Scoring iterations: 5

El resumen del modelo presenta los coeficientes estimados en escala logarítmica. Por lo tanto, la interpretación directa de los parámetros debe realizarse considerando que el modelo fue estimado sobre el logaritmo de la media esperada del conteo.

11 Interpretación de parámetros estimados

11.1 Interpretación en función del signo

El signo del coeficiente permite conocer la dirección de la relación entre una variable explicativa y el número esperado de meses sin cotizar.

• Si \(\beta_j > 0\), la variable se asocia con un aumento del número esperado de meses sin cotizar.
• Si \(\beta_j < 0\), la variable se asocia con una disminución del número esperado de meses sin cotizar.
• Si \(\beta_j = 0\), no se observa cambio en el número esperado de meses sin cotizar, manteniendo constantes las demás variables.

Esta interpretación permite identificar la orientación de la asociación, aunque no su intensidad.

11.2 Interpretación mediante razón de incidencias

Para interpretar la magnitud del efecto se utiliza la exponenciación del coeficiente:

\[ IRR_j = exp(\beta_j) \]

El valor \(IRR_j\) se conoce como razón de incidencias o razón de tasas. Su interpretación general es:

• Si \(IRR_j > 1\), aumenta el conteo esperado.
• Si \(IRR_j < 1\), disminuye el conteo esperado.
• Si \(IRR_j = 1\), no hay cambio en el conteo esperado.

La variación porcentual aproximada se calcula como:

\[ (exp(\beta_j)-1) \times 100 \]

tabla_modelo <- tidy(modelo_poisson, conf.int = TRUE, exponentiate = FALSE) %>%
  mutate(
    IRR = exp(estimate),
    variacion_porcentual = (IRR - 1) * 100,
    orientacion = case_when(
      estimate > 0 ~ "Aumenta el conteo esperado",
      estimate < 0 ~ "Disminuye el conteo esperado",
      TRUE ~ "Sin variación"
    )
  ) %>%
  dplyr::select(term, estimate, std.error, statistic, p.value, IRR, variacion_porcentual, orientacion)

kable(tabla_modelo, digits = 4, caption = "Coeficientes del modelo de Poisson e interpretación mediante IRR")

Coeficientes del modelo de Poisson e interpretación mediante IRR

term	estimate	std.error	statistic	p.value	IRR	variacion_porcentual	orientacion
(Intercept)	-0,8630	0,1477	-5,8422	0,0000	0,4219	-57,8090	Disminuye el conteo esperado
sexoMujer	0,1083	0,0508	2,1304	0,0331	1,1144	11,4359	Aumenta el conteo esperado
edad	0,0160	0,0024	6,7624	0,0000	1,0161	1,6078	Aumenta el conteo esperado
nivel_educSecundaria	0,1403	0,0649	2,1625	0,0306	1,1506	15,0614	Aumenta el conteo esperado
nivel_educHasta primaria	0,3074	0,0718	4,2843	0,0000	1,3599	35,9939	Aumenta el conteo esperado
tipo_contratoTemporal	0,3842	0,0623	6,1702	0,0000	1,4685	46,8474	Aumenta el conteo esperado
tipo_contratoSin contrato escrito	0,7655	0,0609	12,5650	0,0000	2,1501	115,0130	Aumenta el conteo esperado
tamano_empresaPequeña	0,2125	0,0787	2,6998	0,0069	1,2368	23,6794	Aumenta el conteo esperado
tamano_empresaMicro	0,3486	0,0747	4,6680	0,0000	1,4171	41,7102	Aumenta el conteo esperado
sectorComercio	0,1021	0,0759	1,3440	0,1789	1,1075	10,7470	Aumenta el conteo esperado
sectorServicios	-0,0679	0,0750	-0,9050	0,3655	0,9344	-6,5645	Disminuye el conteo esperado
sectorConstrucción	0,2771	0,0911	3,0418	0,0024	1,3194	31,9360	Aumenta el conteo esperado
antiguedad	-0,0727	0,0097	-7,4685	0,0000	0,9299	-7,0140	Disminuye el conteo esperado

11.3 Ejemplo de interpretación para una variable cualitativa

Se considera la variable tipo_contrato. La categoría de referencia es Permanente. Por lo tanto, los coeficientes de las categorías Temporal y Sin contrato escrito se interpretan en comparación con trabajadores con contrato permanente, manteniendo constantes las demás variables del modelo.

coef_contrato <- tabla_modelo %>%
  filter(grepl("tipo_contrato", term)) %>%
  dplyr::select(term, estimate, IRR, variacion_porcentual, orientacion, p.value)

kable(coef_contrato, digits = 4, caption = "Interpretación de la variable tipo de contrato")

Interpretación de la variable tipo de contrato

term	estimate	IRR	variacion_porcentual	orientacion	p.value
tipo_contratoTemporal	0,3842	1,4685	46,8474	Aumenta el conteo esperado	0
tipo_contratoSin contrato escrito	0,7655	2,1501	115,0130	Aumenta el conteo esperado	0

Si el coeficiente asociado a tipo_contratoTemporal es positivo, esto indica que los trabajadores con contrato temporal presentan una cantidad esperada de meses sin cotizar mayor que los trabajadores con contrato permanente, manteniendo constantes las demás variables.

En términos de razón de incidencias, si el IRR es mayor que 1, la cantidad esperada de meses sin cotizar es proporcionalmente mayor para trabajadores con contrato temporal respecto de la categoría de referencia.

11.4 Ejemplo de interpretación para una variable cuantitativa

Se considera la variable antiguedad, medida en años. Su coeficiente indica cómo cambia el número esperado de meses sin cotizar ante un aumento de un año en la antigüedad laboral, manteniendo constantes las demás variables.

coef_antiguedad <- tabla_modelo %>%
  filter(term == "antiguedad") %>%
  dplyr::select(term, estimate, IRR, variacion_porcentual, orientacion, p.value)

kable(coef_antiguedad, digits = 4, caption = "Interpretación de la variable antigüedad")

Interpretación de la variable antigüedad

term	estimate	IRR	variacion_porcentual	orientacion	p.value
antiguedad	-0,0727	0,9299	-7,014	Disminuye el conteo esperado	0

Si el coeficiente de antiguedad es negativo, se interpreta que una mayor antigüedad laboral se asocia con una menor cantidad esperada de meses sin cotizar. En términos laborales, esto puede reflejar que la estabilidad en el empleo está asociada con una mayor continuidad contributiva.

12 Prueba global del modelo

Para evaluar si el conjunto de variables explicativas mejora el ajuste respecto de un modelo sin covariables, se puede utilizar una prueba de razón de verosimilitud.

El modelo nulo se especifica como:

\[ \log(\lambda_i) = \beta_0 \]

Luego se compara el modelo nulo con el modelo completo.

modelo_nulo <- glm(
  meses_sin_cotizar ~ 1,
  family = poisson(link = "log"),
  data = datos
)

anova(modelo_nulo, modelo_poisson, test = "Chisq")

## Analysis of Deviance Table
## 
## Model 1: meses_sin_cotizar ~ 1
## Model 2: meses_sin_cotizar ~ sexo + edad + nivel_educ + tipo_contrato + 
##     tamano_empresa + sector + antiguedad
##   Resid. Df Resid. Dev Df Deviance              Pr(>Chi)    
## 1      1199     1695,3                                      
## 2      1187     1387,5 12   307,86 < 0,00000000000000022 ***
## ---
## Signif. codes:  0 '***' 0,001 '**' 0,01 '*' 0,05 '.' 0,1 ' ' 1

La hipótesis nula de la prueba global sostiene que los coeficientes de las variables explicativas son conjuntamente iguales a cero. Un valor-p pequeño permite rechazar la hipótesis nula y concluir que el modelo con variables explicativas mejora el ajuste respecto del modelo nulo.

13 Predicción con el modelo ajustado

13.1 El argumento newdata

Uno de los usos más importantes del modelo ajustado es la predicción para perfiles laborales específicos. Para ello se utiliza la función predict() con el argumento newdata.

El argumento newdata debe contener un nuevo data.frame con las mismas variables explicativas utilizadas en el modelo. Además, las variables categóricas deben conservar los mismos niveles que en la base original.

13.2 Perfiles laborales hipotéticos

Se construyen tres perfiles de trabajadores:

1. Trabajador con contrato permanente, terciaria, empresa mediana o grande y mayor antigüedad.
2. Trabajador con contrato temporal, secundaria, empresa pequeña y menor antigüedad.
3. Trabajador sin contrato escrito, hasta primaria, microempresa y baja antigüedad.

perfiles <- data.frame(
  sexo = factor(c("Hombre", "Mujer", "Mujer"), levels = levels(datos$sexo)),
  edad = c(35, 35, 35),
  nivel_educ = factor(c("Terciaria", "Secundaria", "Hasta primaria"), levels = levels(datos$nivel_educ)),
  tipo_contrato = factor(c("Permanente", "Temporal", "Sin contrato escrito"), levels = levels(datos$tipo_contrato)),
  tamano_empresa = factor(c("Mediana/Grande", "Pequeña", "Micro"), levels = levels(datos$tamano_empresa)),
  sector = factor(c("Industria", "Comercio", "Construcción"), levels = levels(datos$sector)),
  antiguedad = c(6, 3, 1)
)

perfiles

##     sexo edad     nivel_educ        tipo_contrato tamano_empresa       sector
## 1 Hombre   35      Terciaria           Permanente Mediana/Grande    Industria
## 2  Mujer   35     Secundaria             Temporal        Pequeña     Comercio
## 3  Mujer   35 Hasta primaria Sin contrato escrito          Micro Construcción
##   antiguedad
## 1          6
## 2          3
## 3          1

La predicción se realiza en escala de respuesta, es decir, en cantidad esperada de meses sin cotizar.

perfiles$prediccion_meses <- predict(
  modelo_poisson,
  newdata = perfiles,
  type = "response"
)

kable(perfiles, digits = 3, caption = "Predicción del número esperado de meses sin cotizar")

Predicción del número esperado de meses sin cotizar

sexo	edad	nivel_educ	tipo_contrato	tamano_empresa	sector	antiguedad	prediccion_meses
Hombre	35	Terciaria	Permanente	Mediana/Grande	Industria	6	0,477
Mujer	35	Secundaria	Temporal	Pequeña	Comercio	3	1,529
Mujer	35	Hasta primaria	Sin contrato escrito	Micro	Construcción	1	4,177

La columna prediccion_meses representa la cantidad esperada de meses sin cotización para cada perfil laboral. Estos valores no son categorías, sino promedios esperados condicionales al conjunto de variables incluidas en el modelo.

14 Evaluación del ajuste del modelo

14.1 Comparación entre media y varianza

El modelo de Poisson supone que la media y la varianza condicional del conteo son iguales. Como aproximación inicial, se compara la media y la varianza observada de la variable dependiente.

media_y <- mean(datos$meses_sin_cotizar)
varianza_y <- var(datos$meses_sin_cotizar)
relacion_var_media <- varianza_y / media_y

data.frame(
  media = media_y,
  varianza = varianza_y,
  razon_varianza_media = relacion_var_media
) %>%
  kable(digits = 3, caption = "Relación entre varianza y media observada")

Relación entre varianza y media observada

media	varianza	razon_varianza_media
1,304	1,693	1,298

Una razón entre varianza y media considerablemente mayor que 1 puede ser una señal preliminar de sobredispersión.

14.2 Estadístico de dispersión de Pearson

Una forma más adecuada de evaluar sobredispersión consiste en calcular el estadístico de dispersión basado en los residuos de Pearson:

\[ \hat{\phi} = \frac{\sum r_i^2}{gl} \]

Donde \(r_i\) son los residuos de Pearson y \(gl\) son los grados de libertad residuales.

dispersion_pearson <- sum(residuals(modelo_poisson, type = "pearson")^2) /
  modelo_poisson$df.residual

dispersion_pearson

## [1] 1,01874

Como regla práctica, valores cercanos a 1 son compatibles con el supuesto de equidispersión. Valores claramente superiores a 1 sugieren sobredispersión.

14.3 Residuos del modelo

datos$residuos_pearson <- residuals(modelo_poisson, type = "pearson")
datos$ajustados <- fitted(modelo_poisson)

ggplot(datos, aes(x = ajustados, y = residuos_pearson)) +
  geom_point(alpha = 0.4) +
  geom_hline(yintercept = 0, linetype = "dashed") +
  labs(
    title = "Residuos de Pearson versus valores ajustados",
    x = "Valores ajustados",
    y = "Residuos de Pearson"
  ) +
  theme_minimal()

El gráfico permite observar si existen patrones sistemáticos en los residuos. Un patrón marcado podría sugerir problemas de especificación del modelo.

15 Alternativas ante sobredispersión

Si se detecta sobredispersión, las inferencias del modelo de Poisson pueden ser demasiado optimistas, especialmente porque los errores estándar podrían estar subestimados. En ese caso, pueden considerarse dos alternativas: el modelo cuasi-Poisson y el modelo binomial negativo.

15.1 Modelo cuasi-Poisson

El modelo cuasi-Poisson mantiene la misma media condicional que el modelo de Poisson, pero permite que la varianza sea proporcional a la media.

modelo_quasi <- glm(
  meses_sin_cotizar ~ sexo + edad + nivel_educ + tipo_contrato +
    tamano_empresa + sector + antiguedad,
  family = quasipoisson(link = "log"),
  data = datos
)

summary(modelo_quasi)

## 
## Call:
## glm(formula = meses_sin_cotizar ~ sexo + edad + nivel_educ + 
##     tipo_contrato + tamano_empresa + sector + antiguedad, family = quasipoisson(link = "log"), 
##     data = datos)
## 
## Coefficients:
##                                    Estimate Std. Error t value
## (Intercept)                       -0,862963   0,149091  -5,788
## sexoMujer                          0,108279   0,051301   2,111
## edad                               0,015950   0,002381   6,700
## nivel_educSecundaria               0,140296   0,065481   2,143
## nivel_educHasta primaria           0,307440   0,072429   4,245
## tipo_contratoTemporal              0,384224   0,062852   6,113
## tipo_contratoSin contrato escrito  0,765528   0,061494  12,449
## tamano_empresaPequeña              0,212522   0,079452   2,675
## tamano_empresaMicro                0,348614   0,075379   4,625
## sectorComercio                     0,102079   0,076658   1,332
## sectorServicios                   -0,067898   0,075727  -0,897
## sectorConstrucción                 0,277147   0,091962   3,014
## antiguedad                        -0,072722   0,009828  -7,399
##                                               Pr(>|t|)    
## (Intercept)                          0,000000009099072 ***
## sexoMujer                                      0,03501 *  
## edad                                 0,000000000032132 ***
## nivel_educSecundaria                           0,03235 *  
## nivel_educHasta primaria             0,000023594732249 ***
## tipo_contratoTemporal                0,000000001323229 ***
## tipo_contratoSin contrato escrito < 0,0000000000000002 ***
## tamano_empresaPequeña                          0,00758 ** 
## tamano_empresaMicro                  0,000004160154183 ***
## sectorComercio                                 0,18324    
## sectorServicios                                0,37011    
## sectorConstrucción                             0,00264 ** 
## antiguedad                           0,000000000000258 ***
## ---
## Signif. codes:  0 '***' 0,001 '**' 0,01 '*' 0,05 '.' 0,1 ' ' 1
## 
## (Dispersion parameter for quasipoisson family taken to be 1,018746)
## 
##     Null deviance: 1695,3  on 1199  degrees of freedom
## Residual deviance: 1387,5  on 1187  degrees of freedom
## AIC: NA
## 
## Number of Fisher Scoring iterations: 5

Este modelo ajusta los errores estándar para contemplar la sobredispersión, aunque no permite comparar modelos mediante AIC, porque no se basa en una verosimilitud completa.

15.2 Modelo binomial negativo

El modelo binomial negativo es una alternativa cuando la varianza es mayor que la media y se desea modelar explícitamente la sobredispersión.

modelo_nb <- glm.nb(
  meses_sin_cotizar ~ sexo + edad + nivel_educ + tipo_contrato +
    tamano_empresa + sector + antiguedad,
  data = datos
)

summary(modelo_nb)

## 
## Call:
## glm.nb(formula = meses_sin_cotizar ~ sexo + edad + nivel_educ + 
##     tipo_contrato + tamano_empresa + sector + antiguedad, data = datos, 
##     init.theta = 87.95256438, link = log)
## 
## Coefficients:
##                                    Estimate Std. Error z value
## (Intercept)                       -0,864204   0,148883  -5,805
## sexoMujer                          0,109115   0,051280   2,128
## edad                               0,015947   0,002381   6,699
## nivel_educSecundaria               0,140246   0,065394   2,145
## nivel_educHasta primaria           0,307012   0,072396   4,241
## tipo_contratoTemporal              0,384648   0,062741   6,131
## tipo_contratoSin contrato escrito  0,765484   0,061511  12,445
## tamano_empresaPequeña              0,214021   0,079330   2,698
## tamano_empresaMicro                0,349678   0,075281   4,645
## sectorComercio                     0,101726   0,076593   1,328
## sectorServicios                   -0,068551   0,075632  -0,906
## sectorConstrucción                 0,277327   0,092009   3,014
## antiguedad                        -0,072650   0,009820  -7,398
##                                               Pr(>|z|)    
## (Intercept)                          0,000000006452092 ***
## sexoMujer                                      0,03335 *  
## edad                                 0,000000000020997 ***
## nivel_educSecundaria                           0,03198 *  
## nivel_educHasta primaria             0,000022281563177 ***
## tipo_contratoTemporal                0,000000000875028 ***
## tipo_contratoSin contrato escrito < 0,0000000000000002 ***
## tamano_empresaPequeña                          0,00698 ** 
## tamano_empresaMicro                  0,000003400906816 ***
## sectorComercio                                 0,18414    
## sectorServicios                                0,36474    
## sectorConstrucción                             0,00258 ** 
## antiguedad                           0,000000000000138 ***
## ---
## Signif. codes:  0 '***' 0,001 '**' 0,01 '*' 0,05 '.' 0,1 ' ' 1
## 
## (Dispersion parameter for Negative Binomial(87,9526) family taken to be 1)
## 
##     Null deviance: 1672,7  on 1199  degrees of freedom
## Residual deviance: 1369,8  on 1187  degrees of freedom
## AIC: 3413
## 
## Number of Fisher Scoring iterations: 1
## 
## 
##               Theta:  88 
##           Std. Err.:  234 
## 
##  2 x log-likelihood:  -3384,973

Se puede comparar el AIC del modelo de Poisson y del modelo binomial negativo.

AIC(modelo_poisson, modelo_nb)

##                df      AIC
## modelo_poisson 13 3411,118
## modelo_nb      14 3412,973

Un menor AIC indica mejor ajuste relativo, aunque la selección del modelo debe considerar también la coherencia teórica y la interpretación sustantiva.

16 Predicción comparada entre modelos

En caso de sobredispersión, puede ser útil comparar las predicciones del modelo de Poisson y del modelo binomial negativo para los mismos perfiles laborales.

perfiles$pred_poisson <- predict(modelo_poisson, newdata = perfiles, type = "response")
perfiles$pred_binomial_negativo <- predict(modelo_nb, newdata = perfiles, type = "response")

kable(perfiles, digits = 3, caption = "Predicción comparada entre Poisson y binomial negativo")

Predicción comparada entre Poisson y binomial negativo

sexo	edad	nivel_educ	tipo_contrato	tamano_empresa	sector	antiguedad	prediccion_meses	pred_poisson	pred_binomial_negativo
Hombre	35	Terciaria	Permanente	Mediana/Grande	Industria	6	0,477	0,477	0,476
Mujer	35	Secundaria	Temporal	Pequeña	Comercio	3	1,529	1,529	1,531
Mujer	35	Hasta primaria	Sin contrato escrito	Micro	Construcción	1	4,177	4,177	4,178

17 Interpretación aplicada de los resultados

La lectura sustantiva del modelo debe vincular los resultados estadísticos con el fenómeno laboral analizado. En este caso, el conteo de meses sin cotizar al seguro social puede entenderse como un indicador de discontinuidad contributiva.

Desde esta perspectiva:

• un coeficiente positivo en tipo_contratoTemporal indicaría que los trabajadores con contrato temporal acumulan más meses esperados sin cotización que los trabajadores con contrato permanente;
• un coeficiente positivo en tipo_contratoSin contrato escrito indicaría una mayor discontinuidad contributiva respecto de trabajadores con contrato permanente;
• un coeficiente negativo en antiguedad indicaría que la permanencia en el empleo se asocia con menor cantidad esperada de meses sin cotización;
• un coeficiente positivo en tamano_empresaMicro indicaría que los trabajadores de microempresas presentan mayor cantidad esperada de meses sin cotizar en comparación con trabajadores de empresas medianas o grandes.

La interpretación debe realizarse siempre manteniendo constantes las demás variables incluidas en el modelo.

18 Consideraciones metodológicas para datos reales

En una aplicación con registros administrativos o encuestas laborales, es necesario considerar los siguientes aspectos:

1. Definición del periodo de observación: por ejemplo, últimos 12 o 24 meses.
2. Construcción del conteo: identificar para cada trabajador cuántos meses no registra cotización.
3. Exposición al riesgo: si no todos los trabajadores son observados durante la misma cantidad de meses, debe utilizarse un offset.
4. Ceros estructurales: algunos trabajadores pueden tener cero meses sin cotizar por condiciones institucionales o contractuales muy estables.
5. Sobredispersión: en datos laborales suele ser frecuente que la varianza supere la media.
6. Dependencia temporal: si se dispone de datos panel, puede existir correlación entre observaciones del mismo trabajador.

18.1 Uso de offset cuando cambia el periodo observado

Si algunos trabajadores fueron observados durante distintos números de meses, el modelo debe incorporar un offset. Por ejemplo, si meses_observados indica la exposición de cada trabajador, el modelo sería:

modelo_poisson_offset <- glm(
  meses_sin_cotizar ~ sexo + edad + nivel_educ + tipo_contrato +
    tamano_empresa + sector + antiguedad + offset(log(meses_observados)),
  family = poisson(link = "log"),
  data = datos_reales
)

En este caso, el modelo estima una tasa de meses sin cotizar por unidad de exposición, en lugar de un conteo bruto.

19 Ejercicios de fijación

19.1 Ejercicio 1: análisis descriptivo del conteo

Utilizando la base datos, realice lo siguiente:

1. Calcule la media, varianza y proporción de ceros de la variable meses_sin_cotizar.
2. Compare el promedio de meses sin cotizar según nivel_educ.
3. Elabore un gráfico de barras con la distribución de la variable dependiente.
4. Comente si el uso de un modelo de conteo parece apropiado.

19.2 Ejercicio 2: ajuste e interpretación del modelo de Poisson

Ajuste un modelo de Poisson considerando como variable dependiente meses_sin_cotizar y como variables explicativas:

• sexo
• edad
• nivel_educ
• tipo_contrato
• tamano_empresa
• sector
• antiguedad

Luego responda:

1. ¿Qué variables presentan coeficientes positivos?
2. ¿Qué variables presentan coeficientes negativos?
3. Interprete el resultado de una variable cualitativa.
4. Interprete el resultado de una variable cuantitativa.
5. Calcule la razón de incidencias para cada parámetro estimado.

19.3 Ejercicio 3: predicción para perfiles laborales

Construya tres perfiles laborales hipotéticos y estime la cantidad esperada de meses sin cotizar para cada uno. Luego responda:

1. ¿Qué perfil presenta mayor cantidad esperada de meses sin cotización?
2. ¿Qué características laborales podrían explicar esta diferencia?
3. ¿Cómo podría utilizarse esta información para focalizar acciones de formalización laboral?

19.4 Ejercicio 4: evaluación de sobredispersión

A partir del modelo ajustado:

1. Calcule el estadístico de dispersión de Pearson.
2. Evalúe si existe evidencia de sobredispersión.
3. Ajuste un modelo cuasi-Poisson.
4. Ajuste un modelo binomial negativo.
5. Compare los resultados y discuta cuál modelo resulta más apropiado.

20 Código integrado para replicación rápida

El siguiente bloque resume el flujo básico de trabajo para ajustar el modelo de Poisson.

# 1. Cargar librerías
library(dplyr)
library(broom)
library(MASS)

# 2. Ajustar modelo de Poisson
modelo_poisson <- glm(
  meses_sin_cotizar ~ sexo + edad + nivel_educ + tipo_contrato +
    tamano_empresa + sector + antiguedad,
  family = poisson(link = "log"),
  data = datos
)

# 3. Revisar resultados
summary(modelo_poisson)

# 4. Calcular razón de incidencias
tidy(modelo_poisson) %>%
  mutate(
    IRR = exp(estimate),
    variacion_porcentual = (IRR - 1) * 100
  )

# 5. Evaluar sobredispersión
dispersion_pearson <- sum(residuals(modelo_poisson, type = "pearson")^2) /
  modelo_poisson$df.residual

dispersion_pearson

# 6. Modelo alternativo si existe sobredispersión
modelo_nb <- glm.nb(
  meses_sin_cotizar ~ sexo + edad + nivel_educ + tipo_contrato +
    tamano_empresa + sector + antiguedad,
  data = datos
)

summary(modelo_nb)

21 Bibliografía sugerida

Agresti, A. (2013). Categorical Data Analysis. Wiley.

Cameron, A. C., & Trivedi, P. K. (2013). Regression Analysis of Count Data. Cambridge University Press.

Hilbe, J. M. (2014). Modeling Count Data. Cambridge University Press.

Long, J. S. (1997). Regression Models for Categorical and Limited Dependent Variables. SAGE Publications.

Venables, W. N., & Ripley, B. D. (2002). Modern Applied Statistics with S. Springer.

R Core Team. R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing.

Wickham, H., François, R., Henry, L., Müller, K., & Vaughan, D. dplyr: A Grammar of Data Manipulation.

Wickham, H. ggplot2: Elegant Graphics for Data Analysis. Springer.

Robinson, D., Hayes, A., & Couch, S. broom: Convert Statistical Objects into Tidy Tibbles.