knitr::opts_chunk$set(
  echo = TRUE,
  warning = FALSE,
  message = FALSE,
  fig.align = "center",
  fig.width = 8,
  fig.height = 5
)

set.seed(2026)

# ============================================================
# PAQUETES
# ============================================================

paquetes <- c(
  "ggplot2",
  "dplyr",
  "knitr",
  "leaps",
  "lmtest"
)

faltantes <- paquetes[
  !sapply(
    paquetes,
    requireNamespace,
    quietly = TRUE
  )
]

if(length(faltantes) > 0){
  install.packages(faltantes)
}

library(ggplot2)
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(knitr)
library(leaps)
library(lmtest)
## Loading required package: zoo
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric

1 Introducción

En las clase anteriores hemos estudiado cómo estimar e interpretar modelos de regresión.

Sin embargo, en una aplicación real rara vez conocemos de antemano cuál es la ecuación correcta.

En ciencias actuariales podemos disponer de numerosas variables:

  • edad,
  • antigüedad,
  • suma asegurada,
  • número de siniestros anteriores,
  • kilometraje,
  • ubicación,
  • uso del vehículo,
  • experiencia del asegurado,
  • inflación,
  • tiempo.

La pregunta deja de ser únicamente:

¿Cómo estimamos una regresión?

y se convierte en:

¿Cómo construimos una especificación estadística adecuada para representar el fenómeno que estamos estudiando?

Construir un modelo significa decidir:

  1. qué variables deben incorporarse;
  2. qué forma debe tener su relación con la variable dependiente;
  3. si existen relaciones curvilíneas;
  4. si existen interacciones;
  5. si debemos transformar variables;
  6. cuándo agregar o eliminar predictores;
  7. cómo seleccionar entre modelos alternativos;
  8. si los residuos satisfacen los supuestos del modelo;
  9. si existe autocorrelación;
  10. si la especificación tiene sentido actuarial.

2 Objetivos de aprendizaje

Al finalizar esta clase deberías ser capaz de:

  • explicar qué es un modelo lineal general;
  • modelar relaciones curvilíneas;
  • formular e interpretar interacciones;
  • transformar la variable dependiente;
  • identificar modelos no lineales que pueden linealizarse;
  • determinar cuándo agregar o eliminar variables;
  • utilizar pruebas \(F\) parciales;
  • comprender el papel de los valores \(p\);
  • aplicar selección hacia adelante;
  • aplicar eliminación hacia atrás;
  • aplicar regresión por pasos;
  • aplicar regresión de mejores subconjuntos;
  • comparar modelos utilizando \(R^2\), \(R^2\) ajustado, AIC, BIC y \(C_p\);
  • representar diseños experimentales mediante regresión;
  • comprender la autocorrelación;
  • calcular e interpretar Durbin-Watson;
  • construir una estrategia integral de modelación actuarial.

3 De la regresión múltiple a la construcción de modelos

Consideremos el modelo clásico:

\[ Y_i = \beta_0+ \beta_1X_{1i}+ \beta_2X_{2i}+ \cdots+ \beta_pX_{pi}+ \varepsilon_i. \]

Los supuestos básicos incluyen:

\[ E(\varepsilon_i)=0, \]

\[ Var(\varepsilon_i)=\sigma^2, \]

y:

\[ Cov(\varepsilon_i,\varepsilon_j)=0, \qquad i\neq j. \]

Entonces:

\[ E(Y\mid X_1,\ldots,X_p) = \beta_0+ \beta_1X_1+ \cdots+ \beta_pX_p. \]

Una primera impresión podría ser que la regresión lineal únicamente permite representar relaciones rectas.

Esto no es correcto.

4 Modelo lineal general

El modelo lineal general puede escribirse como:

\[ \boxed{ Y= \beta_0+ \beta_1Z_1+ \beta_2Z_2+ \cdots+ \beta_pZ_p+ \varepsilon } \]

donde cada:

\[ Z_j \]

puede ser una función de las variables originales.

Por ejemplo:

\[ Z_1=X, \]

\[ Z_2=X^2, \]

\[ Z_3=X^3, \]

\[ Z_4=X_1X_2, \]

\[ Z_5=\log(X_3). \]

Por tanto, podemos construir modelos como:

\[ Y= \beta_0+ \beta_1X+ \beta_2X^2+ \varepsilon, \]

o:

\[ Y= \beta_0+ \beta_1X_1+ \beta_2X_2+ \beta_3X_1X_2+ \varepsilon. \]

5 ¿Por qué sigue siendo un modelo lineal?

Consideremos:

\[ Y= \beta_0+ \beta_1X+ \beta_2X^2+ \varepsilon. \]

Definamos:

\[ Z_1=X, \]

y:

\[ Z_2=X^2. \]

Entonces:

\[ Y= \beta_0+ \beta_1Z_1+ \beta_2Z_2+ \varepsilon. \]

La relación entre \(X\) y \(Y\) puede ser curva.

Pero los parámetros:

\[ \beta_0,\beta_1,\beta_2 \]

aparecen elevados únicamente a la primera potencia.

Por tanto:

\[ \boxed{ \text{linealidad se refiere a los parámetros} } \]

y no necesariamente a la forma gráfica de la relación entre \(X\) y \(Y\).

6 Representación matricial

Podemos escribir:

\[ \mathbf{Y} = \mathbf{X}\boldsymbol{\beta} + \boldsymbol{\varepsilon}. \]

Donde:

\[ \mathbf{Y} = \begin{pmatrix} Y_1\\ Y_2\\ \vdots\\ Y_n \end{pmatrix} \]

y:

\[ \boldsymbol{\beta} = \begin{pmatrix} \beta_0\\ \beta_1\\ \vdots\\ \beta_p \end{pmatrix}. \]

La matriz de diseño es:

\[ \mathbf{X} = \begin{pmatrix} 1 & X_{11} & X_{21} & \cdots & X_{p1}\\ 1 & X_{12} & X_{22} & \cdots & X_{p2}\\ \vdots & \vdots & \vdots & & \vdots\\ 1 & X_{1n} & X_{2n} & \cdots & X_{pn} \end{pmatrix}. \]

7 Demostración del estimador de mínimos cuadrados

El método de mínimos cuadrados busca:

\[ \hat{\boldsymbol{\beta}} = \arg\min_{\boldsymbol{\beta}} SSE. \]

Donde:

\[ SSE= (\mathbf{Y}-\mathbf{X}\boldsymbol{\beta})' (\mathbf{Y}-\mathbf{X}\boldsymbol{\beta}). \]

Desarrollamos:

\[ SSE= \mathbf{Y}'\mathbf{Y} - 2\boldsymbol{\beta}'\mathbf{X}'\mathbf{Y} + \boldsymbol{\beta}'\mathbf{X}'\mathbf{X}\boldsymbol{\beta}. \]

Derivamos respecto de:

\[ \boldsymbol{\beta}. \]

Entonces:

\[ \frac{\partial SSE} {\partial\boldsymbol{\beta}} = -2\mathbf{X}'\mathbf{Y} + 2\mathbf{X}'\mathbf{X}\boldsymbol{\beta}. \]

Igualando a cero:

\[ -2\mathbf{X}'\mathbf{Y} + 2\mathbf{X}'\mathbf{X} \hat{\boldsymbol{\beta}} =0. \]

Entonces:

\[ \mathbf{X}'\mathbf{X} \hat{\boldsymbol{\beta}} = \mathbf{X}'\mathbf{Y}. \]

Finalmente:

\[ \boxed{ \hat{\boldsymbol{\beta}} = (\mathbf{X}'\mathbf{X})^{-1} \mathbf{X}'\mathbf{Y} } \]

si:

\[ \mathbf{X}'\mathbf{X} \]

es invertible.

Esta misma expresión funciona cuando una columna de la matriz contiene:

\[ X^2, \]

\[ X_1X_2, \]

o cualquier transformación previamente definida.

8 Modelado de relaciones curvilíneas

Supongamos que queremos estudiar el costo esperado anual de una póliza en función de la edad del asegurado.

Podríamos sospechar una relación en forma de U:

  • mayor riesgo en edades jóvenes;
  • menor riesgo en edades intermedias;
  • incremento posterior del riesgo.

Un modelo estrictamente lineal puede resultar insuficiente.

9 Caso actuarial 1: edad y costo anual

Generaremos una cartera simulada.

set.seed(1234)

n <- 500

cartera_edad <- data.frame(
  edad = sample(
    18:80,
    n,
    replace = TRUE
  )
)

cartera_edad$costo <-
  900 -
  28 * cartera_edad$edad +
  0.36 * cartera_edad$edad^2 +
  rnorm(
    n,
    mean = 0,
    sd = 180
  )

head(cartera_edad)

10 Exploración gráfica

ggplot(
  cartera_edad,
  aes(
    x = edad,
    y = costo
  )
) +
  geom_point(
    alpha = 0.35
  ) +
  labs(
    title = "Costo anual de siniestros y edad",
    x = "Edad",
    y = "Costo anual"
  )

11 Modelo de primer orden

Comenzamos con:

\[ Y_i= \beta_0+ \beta_1Edad_i+ \varepsilon_i. \]

m1 <- lm(
  costo ~ edad,
  data = cartera_edad
)

summary(m1)
## 
## Call:
## lm(formula = costo ~ edad, data = cartera_edad)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -595.78 -136.75  -18.66  140.24  774.89 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 134.6010    27.1255   4.962 9.59e-07 ***
## edad          7.4885     0.5189  14.430  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 208.1 on 498 degrees of freedom
## Multiple R-squared:  0.2948, Adjusted R-squared:  0.2934 
## F-statistic: 208.2 on 1 and 498 DF,  p-value: < 2.2e-16

La ecuación estimada tiene la forma:

\[ \hat Y= b_0+b_1Edad. \]

En este modelo:

\[ b_1 \]

representa el cambio esperado en el costo por cada año adicional de edad.

El problema es que supone que ese cambio es constante.

12 Recta estimada

ggplot(
  cartera_edad,
  aes(
    x = edad,
    y = costo
  )
) +
  geom_point(
    alpha = 0.25
  ) +
  geom_smooth(
    method = "lm",
    se = FALSE
  ) +
  labs(
    title = "Modelo lineal: edad y costo",
    x = "Edad",
    y = "Costo anual"
  )

13 Residuos y forma funcional

Recordemos:

\[ e_i= y_i-\hat y_i. \]

Si el modelo describe apropiadamente la estructura sistemática, los residuos deberían aparecer aproximadamente distribuidos alrededor de cero sin un patrón claro.

datos_residuales_m1 <- data.frame(
  ajustado = fitted(m1),
  residuo = rstandard(m1)
)

ggplot(
  datos_residuales_m1,
  aes(
    x = ajustado,
    y = residuo
  )
) +
  geom_point(
    alpha = 0.5
  ) +
  geom_hline(
    yintercept = 0,
    linetype = 2
  ) +
  labs(
    title = "Residuos del modelo lineal",
    x = "Valor ajustado",
    y = "Residual estandarizado"
  )

Un patrón curvo sugiere que el modelo lineal de primer orden no está capturando completamente la relación.

14 Modelo de segundo orden

Introducimos:

\[ Edad^2. \]

Entonces:

\[ \boxed{ Y_i= \beta_0+ \beta_1Edad_i+ \beta_2Edad_i^2+ \varepsilon_i } \]

m2 <- lm(
  costo ~
    edad +
    I(edad^2),
  data = cartera_edad
)

summary(m2)
## 
## Call:
## lm(formula = costo ~ edad + I(edad^2), data = cartera_edad)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -611.70 -117.46    3.86  116.12  575.19 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 907.91838   63.67436   14.26   <2e-16 ***
## edad        -28.44624    2.78755  -10.21   <2e-16 ***
## I(edad^2)     0.36269    0.02777   13.06   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 179.8 on 497 degrees of freedom
## Multiple R-squared:  0.475,  Adjusted R-squared:  0.4729 
## F-statistic: 224.9 on 2 and 497 DF,  p-value: < 2.2e-16

En R:

I(edad^2)

significa que deseamos incorporar literalmente:

\[ Edad^2. \]

15 Modelo cuadrático gráficamente

ggplot(
  cartera_edad,
  aes(
    x = edad,
    y = costo
  )
) +
  geom_point(
    alpha = 0.25
  ) +
  geom_smooth(
    method = "lm",
    formula = y ~ poly(x, 2, raw = TRUE),
    se = FALSE
  ) +
  labs(
    title = "Modelo cuadrático del costo de siniestros",
    x = "Edad",
    y = "Costo anual"
  )

16 Interpretación del modelo cuadrático

Tenemos:

\[ E(Y|X) = \beta_0+ \beta_1X+ \beta_2X^2. \]

Para conocer el efecto de \(X\), derivamos:

\[ \frac{\partial E(Y|X)} {\partial X} = \beta_1+ 2\beta_2X. \]

Por tanto:

\[ \boxed{ \text{el efecto de }X\text{ depende del nivel de }X } \]

Ya no existe una única pendiente para todas las edades.

17 Efectos marginales

b <- coef(m2)

b
## (Intercept)        edad   I(edad^2) 
## 907.9183794 -28.4462391   0.3626929

Calculemos el efecto marginal para diferentes edades.

efecto_25 <-
  b["edad"] +
  2 *
  b["I(edad^2)"] *
  25

efecto_50 <-
  b["edad"] +
  2 *
  b["I(edad^2)"] *
  50

efecto_70 <-
  b["edad"] +
  2 *
  b["I(edad^2)"] *
  70

efectos_edad <- data.frame(
  Edad = c(
    25,
    50,
    70
  ),
  Efecto_marginal = c(
    efecto_25,
    efecto_50,
    efecto_70
  )
)

kable(
  efectos_edad,
  digits = 2,
  caption = "Efecto marginal de la edad"
)
Efecto marginal de la edad
Edad Efecto_marginal
25 -10.31
50 7.82
70 22.33

18 Punto mínimo de la curva

Tenemos:

\[ E(Y|X) = \beta_0+ \beta_1X+ \beta_2X^2. \]

Derivamos:

\[ \frac{dE(Y|X)} {dX} = \beta_1+ 2\beta_2X. \]

En un mínimo o máximo:

\[ \beta_1+ 2\beta_2X = 0. \]

Despejamos:

\[ \boxed{ X^* = -\frac{\beta_1} {2\beta_2} } \]

edad_minima <-
  -b["edad"] /
  (
    2 *
    b["I(edad^2)"]
  )

edad_minima
##     edad 
## 39.21533

Interpretación:

El modelo estima la edad aproximada en la que el costo esperado alcanza su nivel mínimo.

19 ¿Debemos agregar el término cuadrático?

Tenemos dos modelos anidados.

Modelo reducido:

\[ M_R: Y= \beta_0+ \beta_1X+ \varepsilon. \]

Modelo completo:

\[ M_C: Y= \beta_0+ \beta_1X+ \beta_2X^2+ \varepsilon. \]

La hipótesis es:

\[ H_0: \beta_2=0 \]

contra:

\[ H_a: \beta_2\neq0. \]

anova(
  m1,
  m2
)

Si el valor \(p\) es suficientemente pequeño, tenemos evidencia de que la incorporación del término cuadrático mejora el modelo.

20 Comparación de ajuste

comparacion_curvatura <- data.frame(
  Modelo = c(
    "Lineal",
    "Cuadrático"
  ),
  R2 = c(
    summary(m1)$r.squared,
    summary(m2)$r.squared
  ),
  R2_ajustado = c(
    summary(m1)$adj.r.squared,
    summary(m2)$adj.r.squared
  ),
  AIC = c(
    AIC(m1),
    AIC(m2)
  )
)

kable(
  comparacion_curvatura,
  digits = 3,
  caption = "Comparación del modelo lineal y cuadrático"
)
Comparación del modelo lineal y cuadrático
Modelo R2 R2_ajustado AIC
Lineal 0.295 0.293 6761.164
Cuadrático 0.475 0.473 6615.620

21 Interacción

Una interacción aparece cuando:

\[ \boxed{ \text{el efecto de una variable depende del nivel de otra} } \]

Un modelo sin interacción sería:

\[ Y= \beta_0+ \beta_1X_1+ \beta_2X_2+ \varepsilon. \]

Un modelo con interacción es:

\[ \boxed{ Y= \beta_0+ \beta_1X_1+ \beta_2X_2+ \beta_3X_1X_2+ \varepsilon } \]

El término:

\[ X_1X_2 \]

es la interacción.

22 Caso actuarial: antigüedad y uso comercial

Supongamos que queremos explicar el costo anual utilizando:

  • antigüedad del vehículo;
  • uso particular o comercial.

Es razonable pensar que el efecto de un año adicional de antigüedad sea distinto entre ambos tipos de uso.

set.seed(2026)

n_inter <- 700

cartera_inter <- data.frame(
  antiguedad = runif(
    n_inter,
    0,
    15
  ),
  comercial = rbinom(
    n_inter,
    1,
    0.30
  )
)

cartera_inter$costo <-
  500 +
  35 *
  cartera_inter$antiguedad +
  300 *
  cartera_inter$comercial +
  45 *
  cartera_inter$antiguedad *
  cartera_inter$comercial +
  rnorm(
    n_inter,
    0,
    180
  )

head(cartera_inter)

23 Modelo sin interacción

\[ E(Y) = \beta_0+ \beta_1Antiguedad+ \beta_2Comercial. \]

m_sin_inter <- lm(
  costo ~
    antiguedad +
    comercial,
  data = cartera_inter
)

summary(m_sin_inter)
## 
## Call:
## lm(formula = costo ~ antiguedad + comercial, data = cartera_inter)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -729.45 -125.35   -1.42  133.69  512.72 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)   409.68      15.53   26.38   <2e-16 ***
## antiguedad     48.35       1.74   27.79   <2e-16 ***
## comercial     602.28      16.25   37.05   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 195.9 on 697 degrees of freedom
## Multiple R-squared:  0.7605, Adjusted R-squared:  0.7598 
## F-statistic:  1106 on 2 and 697 DF,  p-value: < 2.2e-16

24 Modelo con interacción

\[ E(Y) = \beta_0+ \beta_1Antiguedad+ \beta_2Comercial+ \beta_3Antiguedad\times Comercial. \]

m_inter <- lm(
  costo ~
    antiguedad *
    comercial,
  data = cartera_inter
)

summary(m_inter)
## 
## Call:
## lm(formula = costo ~ antiguedad * comercial, data = cartera_inter)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -548.87 -109.03   -0.82  114.60  628.99 
## 
## Coefficients:
##                      Estimate Std. Error t value Pr(>|t|)    
## (Intercept)           499.782     16.077  31.087   <2e-16 ***
## antiguedad             36.094      1.894  19.056   <2e-16 ***
## comercial             290.494     30.159   9.632   <2e-16 ***
## antiguedad:comercial   41.255      3.475  11.874   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 178.8 on 696 degrees of freedom
## Multiple R-squared:  0.8008, Adjusted R-squared:  0.7999 
## F-statistic: 932.7 on 3 and 696 DF,  p-value: < 2.2e-16

En R:

antiguedad * comercial

equivale a:

antiguedad +
comercial +
antiguedad:comercial

25 Demostración de la interacción

Tenemos:

\[ E(Y) = \beta_0+ \beta_1X+ \beta_2D+ \beta_3XD. \]

Donde:

\[ D= \begin{cases} 0 & \text{vehículo particular}\\ 1 & \text{vehículo comercial} \end{cases} \]

Para un vehículo particular:

\[ D=0. \]

Entonces:

\[ E(Y|D=0) = \beta_0+ \beta_1X. \]

Para un vehículo comercial:

\[ D=1. \]

Entonces:

\[ E(Y|D=1) = \beta_0+ \beta_1X+ \beta_2+ \beta_3X. \]

Agrupamos:

\[ \boxed{ E(Y|D=1) = (\beta_0+\beta_2) + (\beta_1+\beta_3)X } \]

Por tanto:

\[ \beta_2 \]

representa la diferencia en interceptos y:

\[ \beta_3 \]

representa la diferencia en pendientes.

26 Efectos marginales con interacción

Derivamos:

\[ \frac{\partial E(Y)} {\partial X} = \beta_1+ \beta_3D. \]

Para particulares:

\[ D=0, \]

entonces:

\[ \frac{\partial E(Y)} {\partial X} = \beta_1. \]

Para comerciales:

\[ D=1, \]

entonces:

\[ \boxed{ \frac{\partial E(Y)} {\partial X} = \beta_1+\beta_3 } \]

27 Cálculo en R

coef_inter <- coef(m_inter)

efecto_particular <-
  coef_inter["antiguedad"]

efecto_comercial <-
  coef_inter["antiguedad"] +
  coef_inter["antiguedad:comercial"]

efectos_interaccion <- data.frame(
  Tipo = c(
    "Particular",
    "Comercial"
  ),
  Efecto_antiguedad = c(
    efecto_particular,
    efecto_comercial
  )
)

kable(
  efectos_interaccion,
  digits = 2,
  caption = "Efecto marginal de la antigüedad"
)
Efecto marginal de la antigüedad
Tipo Efecto_antiguedad
Particular 36.09
Comercial 77.35

28 Visualización de la interacción

cartera_inter$tipo <-
  factor(
    cartera_inter$comercial,
    levels = c(
      0,
      1
    ),
    labels = c(
      "Particular",
      "Comercial"
    )
  )

ggplot(
  cartera_inter,
  aes(
    x = antiguedad,
    y = costo,
    linetype = tipo
  )
) +
  geom_point(
    alpha = 0.20
  ) +
  geom_smooth(
    method = "lm",
    se = FALSE
  ) +
  labs(
    title = "Interacción entre antigüedad y uso",
    x = "Antigüedad del vehículo",
    y = "Costo anual",
    linetype = "Uso"
  )

29 Significancia de la interacción

Podemos comparar:

anova(
  m_sin_inter,
  m_inter
)

La prueba es:

\[ H_0: \beta_3=0 \]

contra:

\[ H_a: \beta_3\neq0. \]

Si rechazamos \(H_0\), concluimos que el efecto de la antigüedad depende del tipo de uso.

30 Transformaciones de la variable dependiente

En ciencias actuariales muchas variables monetarias presentan:

  • asimetría positiva;
  • valores extremos;
  • varianza creciente con la media.

Por ejemplo:

\[ Y= \text{severidad del siniestro}. \]

En estos casos puede ser conveniente transformar la variable dependiente.

31 Caso actuarial: severidad

set.seed(111)

n_sev <- 500

severidad <- data.frame(
  suma_asegurada = runif(
    n_sev,
    20,
    500
  )
)

mu_log <-
  6.5 +
  0.003 *
  severidad$suma_asegurada

severidad$indemnizacion <-
  exp(
    mu_log +
    rnorm(
      n_sev,
      0,
      0.45
    )
  )

head(severidad)

32 Visualización

ggplot(
  severidad,
  aes(
    x = suma_asegurada,
    y = indemnizacion
  )
) +
  geom_point(
    alpha = 0.35
  ) +
  labs(
    title = "Indemnización y suma asegurada",
    x = "Suma asegurada",
    y = "Indemnización"
  )

La dispersión aumenta a medida que aumenta la indemnización esperada.

33 Modelo en escala original

sev_original <- lm(
  indemnizacion ~
    suma_asegurada,
  data = severidad
)

summary(sev_original)
## 
## Call:
## lm(formula = indemnizacion ~ suma_asegurada, data = severidad)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -1793.4  -514.4  -140.2   361.3  6080.4 
## 
## Coefficients:
##                Estimate Std. Error t value Pr(>|t|)    
## (Intercept)    379.6001    87.8573   4.321 1.88e-05 ***
## suma_asegurada   5.1941     0.2958  17.558  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 883.4 on 498 degrees of freedom
## Multiple R-squared:  0.3824, Adjusted R-squared:  0.3811 
## F-statistic: 308.3 on 1 and 498 DF,  p-value: < 2.2e-16

34 Residuos en escala original

resid_original <- data.frame(
  fitted = fitted(
    sev_original
  ),
  resid = rstandard(
    sev_original
  )
)

ggplot(
  resid_original,
  aes(
    x = fitted,
    y = resid
  )
) +
  geom_point(
    alpha = 0.45
  ) +
  geom_hline(
    yintercept = 0,
    linetype = 2
  ) +
  labs(
    title = "Residuos en escala original",
    x = "Valor ajustado",
    y = "Residual estandarizado"
  )

Una forma de abanico puede sugerir:

\[ Var(\varepsilon_i) \neq \sigma^2. \]

Esto se denomina:

\[ \boxed{ \text{heterocedasticidad} } \]

35 Transformación logarítmica

Definimos:

\[ Y^*= \log(Y). \]

Entonces:

\[ \boxed{ \log(Y) = \beta_0+ \beta_1X+ \varepsilon } \]

sev_log <- lm(
  log(indemnizacion) ~
    suma_asegurada,
  data = severidad
)

summary(sev_log)
## 
## Call:
## lm(formula = log(indemnizacion) ~ suma_asegurada, data = severidad)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.44098 -0.30063 -0.00506  0.29787  1.30889 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)    6.4510497  0.0440034  146.60   <2e-16 ***
## suma_asegurada 0.0031637  0.0001482   21.35   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4425 on 498 degrees of freedom
## Multiple R-squared:  0.478,  Adjusted R-squared:  0.4769 
## F-statistic: 455.9 on 1 and 498 DF,  p-value: < 2.2e-16

36 Residuos después de transformar

resid_log <- data.frame(
  fitted = fitted(
    sev_log
  ),
  resid = rstandard(
    sev_log
  )
)

ggplot(
  resid_log,
  aes(
    x = fitted,
    y = resid
  )
) +
  geom_point(
    alpha = 0.45
  ) +
  geom_hline(
    yintercept = 0,
    linetype = 2
  ) +
  labs(
    title = "Residuos después de la transformación logarítmica",
    x = "Valor ajustado",
    y = "Residual estandarizado"
  )

La transformación logarítmica puede:

  • reducir asimetría;
  • estabilizar la varianza;
  • reducir el impacto de valores extremadamente grandes;
  • transformar relaciones multiplicativas en aditivas.

37 Interpretación del modelo log-lineal

Tenemos:

\[ \log(Y) = \beta_0+ \beta_1X. \]

Aplicamos exponencial:

\[ Y= e^{\beta_0+\beta_1X}. \]

Entonces:

\[ Y= e^{\beta_0} e^{\beta_1X}. \]

Si \(X\) aumenta una unidad:

\[ \frac{Y(X+1)} {Y(X)} = e^{\beta_1}. \]

Por tanto, el cambio porcentual exacto es:

\[ \boxed{ 100 \left( e^{\beta_1}-1 \right)\% } \]

Para valores pequeños de \(\beta_1\):

\[ e^{\beta_1}-1 \approx \beta_1. \]

Entonces podemos utilizar la aproximación:

\[ 100\beta_1\%. \]

38 Cambio porcentual en nuestro ejemplo

beta_suma <-
  coef(
    sev_log
  )["suma_asegurada"]

cambio_porcentual <-
  100 *
  (
    exp(beta_suma) -
    1
  )

cambio_porcentual
## suma_asegurada 
##      0.3168739

39 37. Problema de retransfomación

Supongamos:

\[ \log(Y) = X\beta+ \varepsilon. \]

Entonces:

\[ Y= e^{X\beta+\varepsilon}. \]

Separando:

\[ Y= e^{X\beta} e^\varepsilon. \]

Por tanto:

\[ E(Y|X) = e^{X\beta} E(e^\varepsilon). \]

Si:

\[ \varepsilon \sim N(0,\sigma^2), \]

entonces:

\[ E(e^\varepsilon) = e^{\sigma^2/2}. \]

Por consiguiente:

\[ \boxed{ E(Y|X) = \exp \left( X\beta+ \frac{\sigma^2}{2} \right) } \]

Esto es importante actuarialmente porque normalmente estamos interesados en la media o costo esperado.

40 Predicción corregida

Supongamos:

\[ SumaAsegurada=250. \]

nuevo <- data.frame(
  suma_asegurada = 250
)

pred_log <-
  predict(
    sev_log,
    newdata = nuevo
  )

sigma2 <-
  summary(
    sev_log
  )$sigma^2

mediana_pred <-
  exp(
    pred_log
  )

media_pred <-
  exp(
    pred_log +
    sigma2 / 2
  )

predicciones_log <- data.frame(
  Mediana_predicha = mediana_pred,
  Media_predicha_corregida = media_pred
)

kable(
  predicciones_log,
  digits = 2
)
Mediana_predicha Media_predicha_corregida
1396.86 1540.5

41 Transformación recíproca

Otra posibilidad es:

\[ Y^* = \frac{1}{Y}. \]

Entonces:

\[ \frac{1}{Y} = \beta_0+ \beta_1X+ \varepsilon. \]

sev_rec <- lm(
  I(1 / indemnizacion) ~
    suma_asegurada,
  data = severidad
)

summary(sev_rec)
## 
## Call:
## lm(formula = I(1/indemnizacion) ~ suma_asegurada, data = severidad)
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -9.810e-04 -2.701e-04 -5.617e-05  1.583e-04  2.338e-03 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     1.523e-03  4.358e-05   34.95   <2e-16 ***
## suma_asegurada -2.629e-06  1.467e-07  -17.92   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.0004382 on 498 degrees of freedom
## Multiple R-squared:  0.3921, Adjusted R-squared:  0.3908 
## F-statistic: 321.2 on 1 and 498 DF,  p-value: < 2.2e-16

No existe una transformación universalmente superior.

Debemos considerar:

  • comportamiento residual;
  • interpretación;
  • estabilidad de la varianza;
  • valores extremos;
  • desempeño predictivo;
  • sentido actuarial.

42 Modelos no lineales intrínsecamente lineales

Consideremos:

\[ E(Y) = \beta_0 \beta_1^X. \]

Este modelo es exponencial.

Tomamos logaritmos:

\[ \log E(Y) = \log \left( \beta_0 \beta_1^X \right). \]

Por propiedades de los logaritmos:

\[ \log E(Y) = \log(\beta_0) + X\log(\beta_1). \]

Definamos:

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

y:

\[ \alpha_1= \log(\beta_1). \]

Entonces:

\[ \boxed{ \log E(Y) = \alpha_0+ \alpha_1X } \]

que es un modelo lineal.

43 Caso actuarial: inflación de siniestros

Supongamos que el costo medio aumenta alrededor de 6% por año.

\[ E(Y_t) = \beta_0 (1.06)^t. \]

set.seed(908)

t <- 0:15

inflacion <- data.frame(
  anio = t
)

inflacion$costo_medio <-
  10000 *
  1.06^t *
  exp(
    rnorm(
      length(t),
      0,
      0.04
    )
  )

inflacion

44 Visualización de la inflación

ggplot(
  inflacion,
  aes(
    x = anio,
    y = costo_medio
  )
) +
  geom_point(
    size = 2
  ) +
  geom_line() +
  labs(
    title = "Evolución del costo medio de siniestros",
    x = "Años desde el período base",
    y = "Costo medio"
  )

45 Linealización

Estimamos:

\[ \log(Y_t) = \alpha_0+ \alpha_1t+ \varepsilon_t. \]

m_inflacion <- lm(
  log(costo_medio) ~
    anio,
  data = inflacion
)

summary(m_inflacion)
## 
## Call:
## lm(formula = log(costo_medio) ~ anio, data = inflacion)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.046514 -0.031054  0.002517  0.014915  0.064803 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 9.226616   0.016080  573.79  < 2e-16 ***
## anio        0.056576   0.001827   30.97 2.68e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.03368 on 14 degrees of freedom
## Multiple R-squared:  0.9856, Adjusted R-squared:  0.9846 
## F-statistic: 959.4 on 1 and 14 DF,  p-value: 2.684e-14

46 Recuperación de parámetros

Sabemos:

\[ \beta_0= e^{\alpha_0} \]

y:

\[ \beta_1= e^{\alpha_1}. \]

alpha <- coef(
  m_inflacion
)

beta0 <- exp(
  alpha[1]
)

beta1 <- exp(
  alpha[2]
)

parametros_inflacion <- data.frame(
  beta0 = beta0,
  beta1 = beta1,
  Crecimiento_porcentual =
    100 *
    (
      beta1 - 1
    )
)

kable(
  parametros_inflacion,
  digits = 3
)
beta0 beta1 Crecimiento_porcentual
(Intercept) 10164.09 1.058 5.821

Si:

\[ \hat\beta_1 \approx 1.06, \]

entonces:

El costo medio aumenta aproximadamente un 6% por período.

47 Determinación de cuándo agregar o eliminar variables

Supongamos un modelo reducido:

\[ M_R: Y= \beta_0+ \beta_1X_1+ \cdots+ \beta_qX_q+ \varepsilon. \]

Y un modelo completo:

\[ M_C: Y= \beta_0+ \beta_1X_1+ \cdots+ \beta_qX_q+ \beta_{q+1}X_{q+1} +\cdots+ \beta_pX_p+ \varepsilon. \]

Queremos evaluar si las variables adicionales contribuyen significativamente.

48 Hipótesis

La hipótesis nula es:

\[ H_0: \beta_{q+1} = \beta_{q+2} = \cdots = \beta_p = 0. \]

La hipótesis alternativa:

\[ H_a: \text{al menos uno de estos coeficientes es diferente de cero}. \]

49 Prueba F parcial

Definamos:

\[ SSE_R \]

como el error del modelo reducido y:

\[ SSE_C \]

como el error del modelo completo.

Como el modelo completo contiene más predictores:

\[ SSE_C \leq SSE_R. \]

La reducción en error es:

\[ SSE_R-SSE_C. \]

El estadístico es:

\[ \boxed{ F= \frac{ (SSE_R-SSE_C)/(p-q) }{ SSE_C/(n-p-1) } } \]

El numerador mide la mejora promedio generada por las nuevas variables.

El denominador corresponde al error cuadrático medio del modelo completo.

50 Cartera para el problema mayor

set.seed(444)

n_cartera <- 600

cartera <- data.frame(
  edad = runif(
    n_cartera,
    18,
    75
  ),
  antiguedad = runif(
    n_cartera,
    0,
    20
  ),
  siniestros_previos = rpois(
    n_cartera,
    0.7
  ),
  suma_asegurada = runif(
    n_cartera,
    30,
    400
  ),
  kilometraje = runif(
    n_cartera,
    5,
    50
  ),
  comercial = rbinom(
    n_cartera,
    1,
    0.22
  ),
  zona_alta = rbinom(
    n_cartera,
    1,
    0.30
  ),
  score = rnorm(
    n_cartera,
    650,
    60
  )
)

cartera$costo <-
  300 +
  5 *
  cartera$edad +
  95 *
  cartera$siniestros_previos +
  1.8 *
  cartera$suma_asegurada +
  180 *
  cartera$comercial +
  140 *
  cartera$zona_alta +
  rnorm(
    n_cartera,
    0,
    250
  )

head(cartera)

51 Modelo reducido

modelo_R <- lm(
  costo ~
    edad +
    siniestros_previos +
    suma_asegurada,
  data = cartera
)

summary(modelo_R)
## 
## Call:
## lm(formula = costo ~ edad + siniestros_previos + suma_asegurada, 
##     data = cartera)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -812.53 -191.80   -3.96  191.30  825.89 
## 
## Coefficients:
##                    Estimate Std. Error t value Pr(>|t|)    
## (Intercept)        412.1931    40.4148  10.199  < 2e-16 ***
## edad                 4.2387     0.6875   6.165  1.3e-09 ***
## siniestros_previos 111.5635    13.0401   8.555  < 2e-16 ***
## suma_asegurada       1.7207     0.1054  16.318  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 276.2 on 596 degrees of freedom
## Multiple R-squared:  0.3923, Adjusted R-squared:  0.3892 
## F-statistic: 128.2 on 3 and 596 DF,  p-value: < 2.2e-16

52 Modelo completo

Agregaremos:

  • uso comercial;
  • zona de riesgo alto.
modelo_C <- lm(
  costo ~
    edad +
    siniestros_previos +
    suma_asegurada +
    comercial +
    zona_alta,
  data = cartera
)

summary(modelo_C)
## 
## Call:
## lm(formula = costo ~ edad + siniestros_previos + suma_asegurada + 
##     comercial + zona_alta, data = cartera)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -733.00 -169.30   11.71  170.96  770.68 
## 
## Coefficients:
##                     Estimate Std. Error t value Pr(>|t|)    
## (Intercept)        325.46976   38.00368   8.564  < 2e-16 ***
## edad                 4.34013    0.63184   6.869 1.63e-11 ***
## siniestros_previos 104.36451   11.99225   8.703  < 2e-16 ***
## suma_asegurada       1.71992    0.09682  17.764  < 2e-16 ***
## comercial          219.65109   25.38472   8.653  < 2e-16 ***
## zona_alta          147.24839   23.12854   6.367 3.87e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 253.6 on 594 degrees of freedom
## Multiple R-squared:  0.4894, Adjusted R-squared:  0.4851 
## F-statistic: 113.9 on 5 and 594 DF,  p-value: < 2.2e-16

53 Prueba F con R

anova(
  modelo_R,
  modelo_C
)

Estamos contrastando:

\[ H_0: \beta_{Comercial} = \beta_{ZonaAlta} = 0. \]

54 Demostración manual de la prueba F

Calculamos:

\[ SSE_R. \]

SSE_R <-
  sum(
    residuals(
      modelo_R
    )^2
  )

SSE_R
## [1] 45455169

Calculamos:

\[ SSE_C. \]

SSE_C <-
  sum(
    residuals(
      modelo_C
    )^2
  )

SSE_C
## [1] 38188769

Número de observaciones:

n_f <-
  nrow(
    cartera
  )

n_f
## [1] 600

Número de predictores del modelo completo:

p_f <-
  length(
    coef(
      modelo_C
    )
  ) -
  1

p_f
## [1] 5

Número de predictores del modelo reducido:

q_f <-
  length(
    coef(
      modelo_R
    )
  ) -
  1

q_f
## [1] 3

Estadístico:

F_manual <-
  (
    (SSE_R - SSE_C) /
    (p_f - q_f)
  ) /
  (
    SSE_C /
    (n_f - p_f - 1)
  )

F_manual
## [1] 56.51192

Valor \(p\):

p_value_F <-
  pf(
    F_manual,
    df1 = p_f - q_f,
    df2 = n_f - p_f - 1,
    lower.tail = FALSE
  )

p_value_F
## [1] 3.409253e-23

55 Relación entre la prueba F y la prueba t

Cuando agregamos una sola variable:

\[ H_0: \beta_j=0. \]

La prueba \(F\) tiene un grado de libertad en el numerador.

Entonces:

\[ \boxed{ F=t^2 } \]

Construimos dos modelos.

m_base <- lm(
  costo ~
    edad +
    siniestros_previos,
  data = cartera
)

m_add <- lm(
  costo ~
    edad +
    siniestros_previos +
    suma_asegurada,
  data = cartera
)

Prueba F:

resultado_F <-
  anova(
    m_base,
    m_add
  )

resultado_F

Prueba t:

resultado_t <-
  summary(
    m_add
  )$coefficients[
    "suma_asegurada",
    "t value"
  ]

resultado_t
## [1] 16.31843

Verificación:

equivalencia_ft <- data.frame(
  F = resultado_F$F[2],
  t_cuadrado = resultado_t^2
)

kable(
  equivalencia_ft,
  digits = 5
)
F t_cuadrado
266.2911 266.2911

56 Uso de valores p

Una regla sencilla podría ser:

  • incorporar una variable si \(p<0.05\);
  • eliminarla si \(p>0.05\).

Sin embargo, esta regla debe utilizarse cuidadosamente.

El valor \(p\) de un coeficiente responde a:

\[ H_0: \beta_j=0 \]

condicionado a que el resto de variables actualmente en el modelo permanezcan incluidas.

Por tanto:

\[ \boxed{ \text{el valor }p\text{ depende de la especificación} } \]

57 Significancia condicional

Supongamos:

\[ Y= \beta_0+ \beta_1X_1+ \beta_2X_2+ \varepsilon. \]

La prueba:

\[ H_0: \beta_1=0 \]

pregunta realmente:

¿Existe evidencia de una relación entre \(X_1\) y \(Y\) después de controlar por \(X_2\)?

No es equivalente a preguntar simplemente si:

\[ Cor(X_1,Y)\neq0. \]

58 Significancia estadística versus importancia actuarial

Una variable puede presentar:

\[ p>0.05 \]

pero seguir siendo relevante por:

  • teoría actuarial;
  • regulación;
  • estructura tarifaria;
  • control de confusión;
  • estabilidad del modelo;
  • comparabilidad temporal;
  • políticas de negocio.

También puede ocurrir:

\[ p<0.001 \]

pero que el efecto sea económicamente muy pequeño.

Por tanto:

\[ \boxed{ \text{significancia estadística} \neq \text{importancia actuarial} } \]

59 Análisis de un problema mayor

Tenemos ahora ocho predictores potenciales:

  • edad;
  • antigüedad;
  • siniestros previos;
  • suma asegurada;
  • kilometraje;
  • uso comercial;
  • zona alta;
  • score.

Queremos explicar:

\[ Costo. \]

modelo_grande <- lm(
  costo ~
    edad +
    antiguedad +
    siniestros_previos +
    suma_asegurada +
    kilometraje +
    comercial +
    zona_alta +
    score,
  data = cartera
)

summary(modelo_grande)
## 
## Call:
## lm(formula = costo ~ edad + antiguedad + siniestros_previos + 
##     suma_asegurada + kilometraje + comercial + zona_alta + score, 
##     data = cartera)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -737.23 -172.37    8.95  171.70  767.07 
## 
## Coefficients:
##                     Estimate Std. Error t value Pr(>|t|)    
## (Intercept)        346.80504  122.22225   2.837   0.0047 ** 
## edad                 4.31678    0.63377   6.811 2.39e-11 ***
## antiguedad          -1.26746    1.84669  -0.686   0.4928    
## siniestros_previos 104.84603   12.04997   8.701  < 2e-16 ***
## suma_asegurada       1.72112    0.09709  17.728  < 2e-16 ***
## kilometraje         -0.52889    0.78787  -0.671   0.5023    
## comercial          219.43783   25.44561   8.624  < 2e-16 ***
## zona_alta          147.66846   23.17496   6.372 3.76e-10 ***
## score                0.00982    0.17222   0.057   0.9545    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 254 on 591 degrees of freedom
## Multiple R-squared:  0.4902, Adjusted R-squared:  0.4833 
## F-statistic: 71.04 on 8 and 591 DF,  p-value: < 2.2e-16

60 Análisis descriptivo

summary(
  cartera
)
##       edad         antiguedad       siniestros_previos suma_asegurada  
##  Min.   :18.01   Min.   : 0.08497   Min.   :0.0000     Min.   : 30.15  
##  1st Qu.:31.24   1st Qu.: 5.36738   1st Qu.:0.0000     1st Qu.:110.59  
##  Median :44.36   Median :10.67298   Median :1.0000     Median :205.84  
##  Mean   :45.96   Mean   :10.35698   Mean   :0.7283     Mean   :207.25  
##  3rd Qu.:60.75   3rd Qu.:15.22699   3rd Qu.:1.0000     3rd Qu.:299.70  
##  Max.   :74.96   Max.   :19.96092   Max.   :5.0000     Max.   :399.58  
##   kilometraje       comercial        zona_alta          score      
##  Min.   : 5.151   Min.   :0.0000   Min.   :0.0000   Min.   :467.2  
##  1st Qu.:15.076   1st Qu.:0.0000   1st Qu.:0.0000   1st Qu.:614.7  
##  Median :26.352   Median :0.0000   Median :0.0000   Median :651.2  
##  Mean   :26.881   Mean   :0.2117   Mean   :0.2783   Mean   :652.8  
##  3rd Qu.:38.183   3rd Qu.:0.0000   3rd Qu.:1.0000   3rd Qu.:695.1  
##  Max.   :49.972   Max.   :1.0000   Max.   :1.0000   Max.   :855.8  
##      costo        
##  Min.   : -82.21  
##  1st Qu.: 799.15  
##  Median :1037.55  
##  Mean   :1044.88  
##  3rd Qu.:1279.68  
##  Max.   :2078.68

61 Matriz de correlaciones

Este es el bloque que generaba el error en la versión anterior.

Ahora utilizamos explícitamente:

dplyr::select()

para evitar cualquier conflicto.

variables_numericas <-
  cartera %>%
  dplyr::select(
    costo,
    edad,
    antiguedad,
    siniestros_previos,
    suma_asegurada,
    kilometraje,
    comercial,
    zona_alta,
    score
  )

matriz_correlaciones <-
  cor(
    variables_numericas,
    use = "complete.obs"
  )

round(
  matriz_correlaciones,
  2
)
##                    costo  edad antiguedad siniestros_previos suma_asegurada
## costo               1.00  0.22      -0.03               0.26           0.53
## edad                0.22  1.00      -0.04              -0.02           0.05
## antiguedad         -0.03 -0.04       1.00              -0.02           0.01
## siniestros_previos  0.26 -0.02      -0.02               1.00          -0.01
## suma_asegurada      0.53  0.05       0.01              -0.01           1.00
## kilometraje         0.00 -0.02      -0.03               0.08           0.00
## comercial           0.25 -0.04      -0.02               0.04           0.00
## zona_alta           0.19  0.02       0.02               0.04          -0.01
## score              -0.03 -0.03       0.04               0.01          -0.04
##                    kilometraje comercial zona_alta score
## costo                     0.00      0.25      0.19 -0.03
## edad                     -0.02     -0.04      0.02 -0.03
## antiguedad               -0.03     -0.02      0.02  0.04
## siniestros_previos        0.08      0.04      0.04  0.01
## suma_asegurada            0.00      0.00     -0.01 -0.04
## kilometraje               1.00      0.01      0.01  0.05
## comercial                 0.01      1.00     -0.02 -0.02
## zona_alta                 0.01     -0.02      1.00 -0.01
## score                     0.05     -0.02     -0.01  1.00

62 ¿Por qué importa la correlación entre predictores?

Supongamos que:

\[ X_1 \]

y:

\[ X_2 \]

están altamente correlacionadas.

Puede ocurrir que \(X_1\) sea significativa cuando \(X_2\) no está incluida.

Pero al incorporar \(X_2\), \(X_1\) puede dejar de ser significativa.

Esto ocurre porque ambas variables contienen información parcialmente redundante.

Por tanto:

\[ \boxed{ \text{la significancia de una variable depende del resto del modelo} } \]

63 Procedimientos de selección de variables

Estudiaremos:

  1. selección hacia adelante;
  2. eliminación hacia atrás;
  3. regresión por pasos;
  4. mejores subconjuntos.

64 Principio de parsimonia

Supongamos:

\[ M_1: Y= \beta_0+ \beta_1X_1+ \beta_2X_2+ \varepsilon \]

y:

\[ M_2: Y= \beta_0+ \beta_1X_1+ \cdots+ \beta_{10}X_{10} + \varepsilon. \]

Si ambos presentan comportamiento predictivo similar, normalmente preferimos el modelo más simple.

Un modelo parsimonioso suele ser:

  • más interpretable;
  • más estable;
  • más sencillo de implementar;
  • menos susceptible al sobreajuste.

Pero:

\[ \boxed{ \text{más simple no significa automáticamente mejor} } \]

65 Selección hacia adelante

Forward selection comienza con:

\[ M_0: Y= \beta_0+ \varepsilon. \]

Luego:

  1. evalúa todos los predictores;
  2. agrega el que produzca la mayor mejora;
  3. vuelve a evaluar las variables restantes;
  4. continúa hasta que ninguna mejora suficientemente el modelo.

Esquemáticamente:

\[ M_0 \rightarrow M_1 \rightarrow M_2 \rightarrow \cdots \]

66 Modelo nulo

modelo_nulo <- lm(
  costo ~ 1,
  data = cartera
)

summary(modelo_nulo)
## 
## Call:
## lm(formula = costo ~ 1, data = cartera)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1127.09  -245.73    -7.33   234.80  1033.80 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  1044.88      14.43   72.43   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 353.4 on 599 degrees of freedom

67 Modelo máximo

modelo_maximo <- lm(
  costo ~
    edad +
    antiguedad +
    siniestros_previos +
    suma_asegurada +
    kilometraje +
    comercial +
    zona_alta +
    score,
  data = cartera
)

summary(modelo_maximo)
## 
## Call:
## lm(formula = costo ~ edad + antiguedad + siniestros_previos + 
##     suma_asegurada + kilometraje + comercial + zona_alta + score, 
##     data = cartera)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -737.23 -172.37    8.95  171.70  767.07 
## 
## Coefficients:
##                     Estimate Std. Error t value Pr(>|t|)    
## (Intercept)        346.80504  122.22225   2.837   0.0047 ** 
## edad                 4.31678    0.63377   6.811 2.39e-11 ***
## antiguedad          -1.26746    1.84669  -0.686   0.4928    
## siniestros_previos 104.84603   12.04997   8.701  < 2e-16 ***
## suma_asegurada       1.72112    0.09709  17.728  < 2e-16 ***
## kilometraje         -0.52889    0.78787  -0.671   0.5023    
## comercial          219.43783   25.44561   8.624  < 2e-16 ***
## zona_alta          147.66846   23.17496   6.372 3.76e-10 ***
## score                0.00982    0.17222   0.057   0.9545    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 254 on 591 degrees of freedom
## Multiple R-squared:  0.4902, Adjusted R-squared:  0.4833 
## F-statistic: 71.04 on 8 and 591 DF,  p-value: < 2.2e-16

68 Selección hacia adelante en R

Utilizaremos AIC como criterio.

modelo_forward <- step(
  modelo_nulo,
  scope = list(
    lower = formula(
      modelo_nulo
    ),
    upper = formula(
      modelo_maximo
    )
  ),
  direction = "forward",
  trace = 0
)

summary(modelo_forward)
## 
## Call:
## lm(formula = costo ~ suma_asegurada + siniestros_previos + comercial + 
##     edad + zona_alta, data = cartera)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -733.00 -169.30   11.71  170.96  770.68 
## 
## Coefficients:
##                     Estimate Std. Error t value Pr(>|t|)    
## (Intercept)        325.46976   38.00368   8.564  < 2e-16 ***
## suma_asegurada       1.71992    0.09682  17.764  < 2e-16 ***
## siniestros_previos 104.36451   11.99225   8.703  < 2e-16 ***
## comercial          219.65109   25.38472   8.653  < 2e-16 ***
## edad                 4.34013    0.63184   6.869 1.63e-11 ***
## zona_alta          147.24839   23.12854   6.367 3.87e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 253.6 on 594 degrees of freedom
## Multiple R-squared:  0.4894, Adjusted R-squared:  0.4851 
## F-statistic: 113.9 on 5 and 594 DF,  p-value: < 2.2e-16

Modelo seleccionado:

formula(
  modelo_forward
)
## costo ~ suma_asegurada + siniestros_previos + comercial + edad + 
##     zona_alta

69 Eliminación hacia atrás

Backward selection comienza con todas las variables.

\[ M_p \rightarrow M_{p-1} \rightarrow M_{p-2} \rightarrow \cdots \]

Una variable eliminada no vuelve a entrar.

modelo_backward <- step(
  modelo_maximo,
  direction = "backward",
  trace = 0
)

summary(modelo_backward)
## 
## Call:
## lm(formula = costo ~ edad + siniestros_previos + suma_asegurada + 
##     comercial + zona_alta, data = cartera)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -733.00 -169.30   11.71  170.96  770.68 
## 
## Coefficients:
##                     Estimate Std. Error t value Pr(>|t|)    
## (Intercept)        325.46976   38.00368   8.564  < 2e-16 ***
## edad                 4.34013    0.63184   6.869 1.63e-11 ***
## siniestros_previos 104.36451   11.99225   8.703  < 2e-16 ***
## suma_asegurada       1.71992    0.09682  17.764  < 2e-16 ***
## comercial          219.65109   25.38472   8.653  < 2e-16 ***
## zona_alta          147.24839   23.12854   6.367 3.87e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 253.6 on 594 degrees of freedom
## Multiple R-squared:  0.4894, Adjusted R-squared:  0.4851 
## F-statistic: 113.9 on 5 and 594 DF,  p-value: < 2.2e-16

Fórmula final:

formula(
  modelo_backward
)
## costo ~ edad + siniestros_previos + suma_asegurada + comercial + 
##     zona_alta

70 Regresión por pasos

La regresión stepwise combina ambas estrategias.

Una variable puede:

  • entrar al modelo;
  • permanecer;
  • salir posteriormente.

Esto permite reconsiderar decisiones tomadas en etapas anteriores.

modelo_step <- step(
  modelo_nulo,
  scope = list(
    lower = formula(
      modelo_nulo
    ),
    upper = formula(
      modelo_maximo
    )
  ),
  direction = "both",
  trace = 0
)

summary(modelo_step)
## 
## Call:
## lm(formula = costo ~ suma_asegurada + siniestros_previos + comercial + 
##     edad + zona_alta, data = cartera)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -733.00 -169.30   11.71  170.96  770.68 
## 
## Coefficients:
##                     Estimate Std. Error t value Pr(>|t|)    
## (Intercept)        325.46976   38.00368   8.564  < 2e-16 ***
## suma_asegurada       1.71992    0.09682  17.764  < 2e-16 ***
## siniestros_previos 104.36451   11.99225   8.703  < 2e-16 ***
## comercial          219.65109   25.38472   8.653  < 2e-16 ***
## edad                 4.34013    0.63184   6.869 1.63e-11 ***
## zona_alta          147.24839   23.12854   6.367 3.87e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 253.6 on 594 degrees of freedom
## Multiple R-squared:  0.4894, Adjusted R-squared:  0.4851 
## F-statistic: 113.9 on 5 and 594 DF,  p-value: < 2.2e-16

Fórmula:

formula(
  modelo_step
)
## costo ~ suma_asegurada + siniestros_previos + comercial + edad + 
##     zona_alta

71 Comparación de procedimientos

comparacion_step <- data.frame(
  Procedimiento = c(
    "Forward",
    "Backward",
    "Stepwise"
  ),
  AIC = c(
    AIC(
      modelo_forward
    ),
    AIC(
      modelo_backward
    ),
    AIC(
      modelo_step
    )
  ),
  R2_ajustado = c(
    summary(
      modelo_forward
    )$adj.r.squared,
    summary(
      modelo_backward
    )$adj.r.squared,
    summary(
      modelo_step
    )$adj.r.squared
  ),
  Numero_predictores = c(
    length(
      coef(
        modelo_forward
      )
    ) - 1,
    length(
      coef(
        modelo_backward
      )
    ) - 1,
    length(
      coef(
        modelo_step
      )
    ) - 1
  )
)

kable(
  comparacion_step,
  digits = 3,
  caption = "Comparación de procedimientos de selección"
)
Comparación de procedimientos de selección
Procedimiento AIC R2_ajustado Numero_predictores
Forward 8353.4 0.485 5
Backward 8353.4 0.485 5
Stepwise 8353.4 0.485 5

72 Advertencia sobre los procedimientos automáticos

Los métodos automáticos son herramientas exploratorias.

No deben interpretarse como una máquina que descubre el modelo verdadero.

Entre sus limitaciones:

  • prueban múltiples especificaciones;
  • pueden seleccionar variables por azar;
  • los resultados pueden cambiar con pequeñas variaciones de los datos;
  • pueden ignorar relaciones teóricas;
  • pueden eliminar variables sustantivamente importantes;
  • pueden favorecer el sobreajuste.

Por ello:

\[ \boxed{ \text{selección automática} + \text{criterio estadístico} + \text{criterio actuarial} + \text{validación} } \]

es preferible a confiar únicamente en un algoritmo.

73 Regresión de los mejores subconjuntos

Con:

\[ p \]

predictores existen:

\[ \boxed{ 2^p } \]

subconjuntos posibles.

Para ocho variables:

2^8
## [1] 256

Tenemos 256 combinaciones posibles.

74 Best subsets en R

best <- regsubsets(
  costo ~
    edad +
    antiguedad +
    siniestros_previos +
    suma_asegurada +
    kilometraje +
    comercial +
    zona_alta +
    score,
  data = cartera,
  nvmax = 8
)

best_res <- summary(
  best
)

75 Tabla de mejores subconjuntos

tabla_best <- data.frame(
  Variables = 1:8,
  R2 = best_res$rsq,
  R2_ajustado = best_res$adjr2,
  Cp = best_res$cp,
  BIC = best_res$bic
)

kable(
  tabla_best,
  digits = 3,
  caption = "Comparación de mejores subconjuntos"
)
Comparación de mejores subconjuntos
Variables R2 R2_ajustado Cp BIC
1 0.281 0.280 237.635 -185.067
2 0.354 0.351 155.452 -242.542
3 0.412 0.409 89.451 -293.216
4 0.455 0.451 42.286 -331.749
5 0.489 0.485 3.897 -364.958
6 0.490 0.485 5.451 -359.012
7 0.490 0.484 7.003 -353.070
8 0.490 0.483 9.000 -346.676

76 ¿Por qué no seleccionar utilizando únicamente \(R^2\)?

Tenemos:

\[ R^2 = 1- \frac{SSE}{SST}. \]

Cuando agregamos variables:

\[ SSE \]

no puede aumentar.

Por tanto:

\[ R^2 \]

no disminuye.

Eso significa que maximizar únicamente \(R^2\) favorece modelos grandes.

77 \(R^2\) ajustado

El coeficiente ajustado es:

\[ \boxed{ R^2_{adj} = 1- \frac{ SSE/(n-p-1) }{ SST/(n-1) } } \]

La nueva variable debe producir suficiente mejora como para compensar la pérdida de un grado de libertad.

Por eso:

\[ R^2_{adj} \]

sí puede disminuir.

78 Visualización del \(R^2\) ajustado

ggplot(
  tabla_best,
  aes(
    x = Variables,
    y = R2_ajustado
  )
) +
  geom_line() +
  geom_point(
    size = 2
  ) +
  labs(
    title = "R² ajustado según el tamaño del modelo",
    x = "Número de variables",
    y = "R² ajustado"
  )

Mejor tamaño:

which.max(
  best_res$adjr2
)
## [1] 5

79 Mallows \(C_p\)

Una expresión habitual es:

\[ \boxed{ C_p= \frac{SSE_p} {MSE_{full}} - (n-2p') } \]

donde:

\[ p' \]

incluye el intercepto.

Buscamos modelos con:

  • \(C_p\) relativamente bajo;
  • \(C_p\) aproximadamente cercano al número de parámetros.

Debe utilizarse como criterio de comparación y no como una regla absoluta.

80 BIC

El BIC puede expresarse aproximadamente como:

\[ \boxed{ BIC= n \log \left( \frac{SSE}{n} \right) + k\log(n) } \]

donde:

\[ k \]

es el número de parámetros.

Buscamos:

\[ \boxed{ BIC\text{ pequeño} } \]

porque el criterio combina ajuste y penalización por complejidad.

81 BIC según número de variables

ggplot(
  tabla_best,
  aes(
    x = Variables,
    y = BIC
  )
) +
  geom_line() +
  geom_point(
    size = 2
  ) +
  labs(
    title = "BIC según el tamaño del modelo",
    x = "Número de variables",
    y = "BIC"
  )

82 Modelo seleccionado mediante BIC

modelo_bic_id <-
  which.min(
    best_res$bic
  )

modelo_bic_id
## [1] 5

Variables:

coef(
  best,
  id = modelo_bic_id
)
##        (Intercept)               edad siniestros_previos     suma_asegurada 
##         325.469763           4.340129         104.364505           1.719923 
##          comercial          zona_alta 
##         219.651086         147.248394

83 La elección final

La elección final no debería depender únicamente de:

\[ p, \]

\[ R^2, \]

\[ R^2_{adj}, \]

\[ AIC, \]

\[ BIC. \]

También debemos analizar:

83.0.1 Plausibilidad

¿Los signos de los coeficientes tienen sentido?

83.0.2 Estabilidad

¿Cambian demasiado al utilizar otra muestra?

83.0.3 Predictividad

¿Cómo funciona el modelo con observaciones nuevas?

83.0.4 Parsimonia

¿Estamos incorporando variables que realmente aportan?

83.0.5 Interpretabilidad

¿Podemos explicar cada parámetro?

83.0.6 Implementación

¿Los predictores estarán disponibles cuando el modelo sea utilizado?

83.0.7 Gobernanza

¿Podemos justificar el uso de cada variable?

84 Construcción versus validación

Un modelo puede funcionar muy bien con los datos utilizados para estimarlo y mal con nuevos datos.

Esto se denomina:

\[ \boxed{ \text{sobreajuste} } \]

Por ello debemos separar:

\[ \text{entrenamiento} \]

y:

\[ \text{evaluación}. \]

85 División entrenamiento-prueba

Utilizaremos:

\[ 70\% \]

para entrenamiento y:

\[ 30\% \]

para prueba.

set.seed(77)

id_train <- sample(
  seq_len(
    nrow(
      cartera
    )
  ),
  size = floor(
    0.70 *
    nrow(
      cartera
    )
  )
)

train <-
  cartera[
    id_train,
    ,
    drop = FALSE
  ]

test <-
  cartera[
    -id_train,
    ,
    drop = FALSE
  ]

dim(train)
## [1] 420   9
dim(test)
## [1] 180   9

86 Modelo de entrenamiento

modelo_train <- lm(
  costo ~
    edad +
    siniestros_previos +
    suma_asegurada +
    comercial +
    zona_alta,
  data = train
)

summary(modelo_train)
## 
## Call:
## lm(formula = costo ~ edad + siniestros_previos + suma_asegurada + 
##     comercial + zona_alta, data = train)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -723.89 -177.85   14.44  168.15  751.02 
## 
## Coefficients:
##                    Estimate Std. Error t value Pr(>|t|)    
## (Intercept)        296.3370    47.2881   6.267 9.27e-10 ***
## edad                 4.5542     0.7902   5.763 1.61e-08 ***
## siniestros_previos 118.9314    14.0679   8.454 4.86e-16 ***
## suma_asegurada       1.7491     0.1176  14.873  < 2e-16 ***
## comercial          217.1281    30.6428   7.086 5.99e-12 ***
## zona_alta          139.0675    27.8438   4.995 8.71e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 259.3 on 414 degrees of freedom
## Multiple R-squared:  0.497,  Adjusted R-squared:  0.4909 
## F-statistic: 81.81 on 5 and 414 DF,  p-value: < 2.2e-16

Predicciones:

pred_test <-
  predict(
    modelo_train,
    newdata = test
  )

87 RMSE

\[ \boxed{ RMSE= \sqrt{ \frac{1}{n} \sum_{i=1}^{n} (y_i-\hat y_i)^2 } } \]

RMSE <- sqrt(
  mean(
    (
      test$costo -
      pred_test
    )^2
  )
)

RMSE
## [1] 241.2524

88 MAE

\[ \boxed{ MAE= \frac{1}{n} \sum_{i=1}^{n} |y_i-\hat y_i| } \]

MAE <- mean(
  abs(
    test$costo -
    pred_test
  )
)

MAE
## [1] 190.9633

El RMSE penaliza con mayor intensidad los errores grandes.

El MAE representa el error absoluto promedio.

89 Caso actuarial: siniestralidad mensual

Generaremos una serie con:

\[ \rho=0.70. \]

set.seed(707)

n_mes <- 60

rho <- 0.70

epsilon <- numeric(
  n_mes
)

z <- rnorm(
  n_mes,
  mean = 0,
  sd = 0.025
)

for(t in 2:n_mes){
  
  epsilon[t] <-
    rho *
    epsilon[t-1] +
    z[t]
}

serie <- data.frame(
  mes = 1:n_mes
)

serie$siniestralidad <-
  0.55 +
  0.0015 *
  serie$mes +
  epsilon

head(serie)

90 Evolución temporal

ggplot(
  serie,
  aes(
    x = mes,
    y = siniestralidad
  )
) +
  geom_line() +
  geom_point() +
  labs(
    title = "Razón de siniestralidad mensual",
    x = "Mes",
    y = "Razón de siniestralidad"
  )

91 Modelo ignorando la autocorrelación

Estimamos:

\[ Y_t= \beta_0+ \beta_1t+ \varepsilon_t. \]

m_serie <- lm(
  siniestralidad ~
    mes,
  data = serie
)

summary(m_serie)
## 
## Call:
## lm(formula = siniestralidad ~ mes, data = serie)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.056605 -0.026015 -0.000575  0.020024  0.054636 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 0.5768980  0.0075013  76.906  < 2e-16 ***
## mes         0.0009926  0.0002139   4.641 2.03e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.02869 on 58 degrees of freedom
## Multiple R-squared:  0.2708, Adjusted R-squared:  0.2582 
## F-statistic: 21.54 on 1 and 58 DF,  p-value: 2.028e-05

92 Residuos temporalmente ordenados

serie$residuo <-
  residuals(
    m_serie
  )

ggplot(
  serie,
  aes(
    x = mes,
    y = residuo
  )
) +
  geom_line() +
  geom_point() +
  geom_hline(
    yintercept = 0,
    linetype = 2
  ) +
  labs(
    title = "Residuos ordenados temporalmente",
    x = "Mes",
    y = "Residual"
  )

Rachas prolongadas de residuos del mismo signo pueden indicar autocorrelación positiva.

93 Prueba de Durbin-Watson

Sean:

\[ e_t= y_t-\hat y_t. \]

El estadístico es:

\[ \boxed{ d= \frac{ \sum_{t=2}^{n} (e_t-e_{t-1})^2 }{ \sum_{t=1}^{n} e_t^2 } } \]

El valor de \(d\) está aproximadamente entre:

\[ 0 \]

y:

\[ 4. \]

94 Interpretación de Durbin-Watson

Una guía aproximada es:

\[ d\approx2 \]

indica ausencia de autocorrelación.

\[ d<2 \]

sugiere autocorrelación positiva.

\[ d>2 \]

sugiere autocorrelación negativa.

En los extremos:

\[ d\rightarrow0 \]

indica fuerte autocorrelación positiva.

Mientras:

\[ d\rightarrow4 \]

indica fuerte autocorrelación negativa.

95 Demostración de la relación con \(\rho\)

Partimos del numerador:

\[ \sum_{t=2}^{n} (e_t-e_{t-1})^2. \]

Desarrollamos:

\[ (e_t-e_{t-1})^2 = e_t^2+ e_{t-1}^2- 2e_te_{t-1}. \]

Por tanto:

\[ \sum (e_t-e_{t-1})^2 = \sum e_t^2 + \sum e_{t-1}^2 - 2 \sum e_te_{t-1}. \]

Para muestras razonablemente grandes:

\[ \sum e_t^2 \approx \sum e_{t-1}^2. \]

Entonces:

\[ \sum (e_t-e_{t-1})^2 \approx 2 \sum e_t^2 - 2 \sum e_te_{t-1}. \]

Dividiendo por:

\[ \sum e_t^2, \]

obtenemos:

\[ d \approx 2 - 2 \frac{ \sum e_te_{t-1} }{ \sum e_t^2 }. \]

El cociente es aproximadamente:

\[ \hat\rho. \]

Por tanto:

\[ \boxed{ d \approx 2(1-\hat\rho) } \]

96 Ejemplos

Si:

\[ \hat\rho=0, \]

entonces:

\[ d \approx 2. \]

Si:

\[ \hat\rho=0.70, \]

entonces:

\[ d \approx 2(1-0.70) = 0.60. \]

Si:

\[ \hat\rho=-0.50, \]

entonces:

\[ d \approx 2(1+0.50) = 3. \]

97 Cálculo manual de Durbin-Watson

e <-
  residuals(
    m_serie
  )

DW_manual <-
  sum(
    diff(e)^2
  ) /
  sum(
    e^2
  )

DW_manual
## [1] 0.7798607

98 Prueba con lmtest

lmtest::dwtest(
  m_serie,
  alternative = "greater"
)
## 
##  Durbin-Watson test
## 
## data:  m_serie
## DW = 0.77986, p-value = 2.383e-08
## alternative hypothesis: true autocorrelation is greater than 0

Para autocorrelación positiva:

\[ H_0: \rho=0 \]

contra:

\[ H_a: \rho>0. \]

99 Función de autocorrelación

acf(
  residuals(
    m_serie
  ),
  main =
    "Función de autocorrelación de los residuos"
)

El rezago 1 estudia aproximadamente:

\[ Cor( e_t, e_{t-1} ). \]

100 ¿Por qué importa la autocorrelación?

Si los errores están correlacionados, los errores estándar de los coeficientes pueden estar incorrectamente estimados.

Recordemos:

\[ t= \frac{ \hat\beta_j }{ SE(\hat\beta_j) }. \]

Si:

\[ SE(\hat\beta_j) \]

es incorrecto, también pueden resultar incorrectos:

  • los estadísticos \(t\);
  • los valores \(p\);
  • los intervalos de confianza.

Por tanto, podríamos declarar una tendencia estadísticamente significativa cuando la evidencia es más débil de lo que parece.

101 Flujo completo de construcción de modelos

Una estrategia razonable puede estructurarse así:

101.1 Paso 1. Definir el problema

Preguntar:

¿Qué variable queremos explicar o predecir?

Por ejemplo:

\[ Y= \text{costo esperado}. \]

101.2 Paso 2. Identificar predictores

Basándonos en:

  • teoría;
  • conocimiento actuarial;
  • experiencia;
  • regulación;
  • evidencia previa.

101.3 Paso 3. Explorar los datos

Analizar:

  • distribuciones;
  • valores faltantes;
  • relaciones bivariadas;
  • valores extremos;
  • correlaciones.

101.4 Paso 4. Analizar forma funcional

Preguntar:

  • ¿la relación es lineal?;
  • ¿existe curvatura?;
  • ¿existe interacción?;
  • ¿debemos transformar \(Y\)?;
  • ¿la relación parece multiplicativa?

101.5 Paso 5. Estimar modelos candidatos

Utilizar modelos apropiados.

101.6 Paso 6. Comparar modelos anidados

Aplicar:

\[ \boxed{ \text{prueba F parcial} } \]

101.7 Paso 7. Analizar selección de variables

Considerar:

  • forward;
  • backward;
  • stepwise;
  • best subsets.

101.8 Paso 8. Evaluar parsimonia

Considerar:

\[ R^2_{adj}, \]

\[ AIC, \]

\[ BIC, \]

\[ C_p. \]

101.9 Paso 9. Analizar residuos

Evaluar:

  • forma funcional;
  • heterocedasticidad;
  • valores atípicos;
  • influencia;
  • independencia.

101.10 Paso 10. Evaluar autocorrelación

Cuando exista orden temporal, analizar:

\[ Durbin-Watson \]

y:

\[ ACF. \]

101.11 Paso 11. Validar

Comparar desempeño fuera de muestra.

101.12 Paso 12. Aplicar criterio actuarial

Finalmente preguntar:

¿El modelo tiene sentido para el fenómeno que estamos estudiando?