Diagnóstico del modelo

Author

Rubén Cool 

library(ggplot2)
library(mosaic) 
library(knitr)
library(kableExtra)
library(modeest) # Necesario para la función mfv()
library(e1071)
library(nortest)
library(lmtest)

1. Introducción y Planteamiento del Problema

El análisis de regresión lineal simple permite modelar la relación entre una variable de respuesta cuantitativa Y y una única variable explicativa X. En este reporte utilizaremos el conjunto de datos mtcars disponible en R, el cual contiene información sobre el consumo de combustible y las características de diseño de 32 automóviles de 1973–1974.

Nuestro objetivo es analizar cómo influye el peso del vehículo (X = \text{wt}, en miles de libras) sobre el rendimiento de combustible (Y = \text{mpg}, en millas por galón).

2. Análisis Exploratorio de Datos (AED)

Antes de ajustar un modelo matemático, es fundamental examinar la distribución univariada de las variables y su relación bivariada mediante gráficos de dispersión.

# Inspeccionamos las variables de interés
favstats(~ mpg, data = mtcars)
  min     Q1 median   Q3  max     mean       sd  n missing
 10.4 15.425   19.2 22.8 33.9 20.09062 6.026948 32       0
favstats(~ wt, data = mtcars)
   min      Q1 median   Q3   max    mean        sd  n missing
 1.513 2.58125  3.325 3.61 5.424 3.21725 0.9784574 32       0

Gráfico de dispersión

Visualizamos la relación entre ambas variables para comprobar si una tendencia lineal es adecuada.

ggplot(mtcars, aes(x = wt, y = mpg)) +
  geom_point(color = "#16a085", size = 3, alpha = 0.8) +
  geom_smooth(method = "lm", se = TRUE, color = "#2c3e50", fill = "#16a085", alpha = 0.2) +
  labs(
    title = "Relación entre el Peso y el Rendimiento del Automóvil",
    subtitle = "Dataset mtcars",
    x = "Peso del vehículo (miles de libras)",
    y = "Rendimiento (Millas por galón - mpg)"
  ) +
  theme_minimal()
`geom_smooth()` using formula = 'y ~ x'

Se observa una clara relación inversa (negativa): a medida que el peso del vehículo aumenta, el rendimiento en millas por galón disminuye de forma aproximadamente lineal.

3. Ajuste del Modelo de Regresión Lineal Simple

El modelo matemático poblacional se define como:

y_i=\beta_0 + \beta_1x_i +\varepsilon_i

Donde:

  • Y_i es el rendimiento (mpg).

  • X_i es el peso (wt).

  • \beta_0 es el intercepto.

  • \beta_1 es la pendiente.

  • \varepsilon_i es el término de error aleatorio.

Estimación de los coeficientes

Ajustamos el modelo en R utilizando la función lm():

modelo <- lm(mpg ~ wt, data = mtcars)
summary(modelo)

Call:
lm(formula = mpg ~ wt, data = mtcars)

Residuals:
    Min      1Q  Median      3Q     Max 
-4.5432 -2.3647 -0.1252  1.4096  6.8727 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  37.2851     1.8776  19.858  < 2e-16 ***
wt           -5.3445     0.5591  -9.559 1.29e-10 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 3.046 on 30 degrees of freedom
Multiple R-squared:  0.7528,    Adjusted R-squared:  0.7446 
F-statistic: 91.38 on 1 and 30 DF,  p-value: 1.294e-10

Interpretación de los parámetros

  • Intercepto (\hat{\beta}_0): 37.29 representa el promedio de rendimiento estimado para un vehículo hipotético con un peso de 0 miles de libras (interpretación matemática, fuera del rango observado).

  • Pendiente (\hat{\beta}_1): -5.34 indica que, por cada incremento de 1,000 libras en el peso del automóvil, se estima que el rendimiento disminuya en promedio 5.34 millas por galón.

4. Diagnóstico del Modelo

Para validar que las inferencias estadísticas del modelo sean confiables, debemos evaluar los supuestos sobre los residuos (\varepsilon_i): normalidad, homocedasticidad (varianza constante) e independencia.

Gráficos de Diagnóstico Estándar

R permite generar los cuatro gráficos clásicos de diagnóstico de una manera muy limpia:

par(mfrow = c(2, 2))
plot(modelo)

par(mfrow = c(1, 1)) # Restaurar configuración gráfica

Análisis detallado de los supuestos:

  1. Normalidad de los errores (Q-Q Normal): En el gráfico Normal Q-Q, la mayoría de los puntos se alinean razonablemente bien sobre la línea diagonal punteada, lo que indica que el supuesto de normalidad de los residuos se cumple de forma aceptable.

  2. Puntos influyentes (Residuales vs Leverage): Identifica observaciones atípicas que podrían ejercer una influencia excesiva sobre la pendiente de la recta de regresión (por ejemplo, vehículos con pesos o rendimientos muy extremos dentro de la muestra).

  3. Homocedasticidad (Scale-Location): La dispersión de los residuos estandarizados se mantiene relativamente constante a lo largo de los valores ajustados, sugiriendo homocedasticidad (varianza de los errores constante).

  4. Independencia (Residuales vs Valores Ajustados): En el primer gráfico (Residuals vs Fitted), los puntos se distribuyen de manera relativamente aleatoria alrededor de la línea curva de color rojo, lo que nos permite sospechar que no se cumple la suposición de independencia.

5. Análisis de la variable respuesta

En este apartado se pretende realizar un análisis exploratorio de la variables rendimiento, en el cual se realizan diversos gráficos e indicadores para observar el comportamiento y la distribución de cada variable.

grafico_1<-ggplot(data=mtcars)+
  geom_histogram(aes(x=mpg),binwidth = 1.2,fill = "#4D6291", color = "black")+
  theme(axis.text.x = element_text(size = 10))+
  labs(title = "Histograma variable Rendimiento", x = "Millas por galón", y = "Frecuencia")+
  theme(plot.title = element_text(hjust = 0.5))

grafico_1

En relación con los cuartiles, encontramos que el 25% de los autos tienen un rendimiento menor o igual a 15.425 millas por galón, el 50% de los autos tienen un rendimiento menor o igual a 19.200 millas por galón y por último el 75% de los autos tienen un rendimiento menor o igual a 22.80 millas por galón.

kable(data.frame(Cuartiles =quantile(mtcars$mpg,prob=c(0.25,0.5,0.75))),
      caption= "Indicadores de posición",aling= "cc")
Indicadores de posición
Cuartiles
25% 15.425
50% 19.200
75% 22.800

Esta representación de lo cuartiles se evidencia de una mejor forma en el siguiente gráfico, donde podemos afirmar lo anteriormente dicho y anexándole que hay una pequeña cantidad de autos con rendimientos muy altos (puntos rojos):

grafico_2<-ggplot(data=mtcars)+
  geom_boxplot(aes(x=mpg), col= "blue", outlier.color = "red")+
  labs(title = "Diagrama de cajas - Variable Rendimiento", x = "Millas por galón", y="Porcentaje")+
  theme(plot.title = element_text(hjust = 0.5))
grafico_2

Podemos observar que la desviación de los rendimientos varía por encima o por debajo de la media en 6.03 Millas por galón. Por otro lado, el coeficiente de variación (CV = 30%) indica que el rendimiento presenta una variabilidad moderada alrededor de su media (20.09 millas por galón). Esto refleja que los datos muestran cierta heterogeneidad en su distribución.”

En cuanto a la forma de la distribución, el coeficiente de curtosis es de -0.37. Al presentar un valor negativo, la distribución se clasifica como levemente platicúrtica, lo que indica que la concentración de los datos alrededor de la media es ligeramente menor a la de una distribución normal (presentando una forma relativamente más aplanada.

# 1. Obtenemos los valores 
promedio <- mean(mtcars$mpg)
desviacion <- sd(mtcars$mpg)

# 2. Creamos la estructura Vertical (Indicador vs Valor)
tabla_vertical <- data.frame(
  Indicador = c("Moda", "Varianza", "Desviación Estándar", "Coeficiente de Variación (%)", "Curtosis"),
  Valor = c(
    mfv(mtcars$mpg)[1],
    var(mtcars$mpg),
    desviacion,
    (desviacion / promedio) * 100,
    kurtosis(mtcars$mpg)
  )
)

# 3. Presentación con Kable
kable(
  tabla_vertical,
  digits = 2, # Redondea automáticamente a 2 decimales
  caption = "Indicadores Descriptivos de Rendimiento (mpg)",
  align = c("l", "c"),
  col.names = c("Estadístico", "Valor")
) %>%
  kable_styling(bootstrap_options = c("striped", "hover"), full_width = FALSE, position = "center")
Indicadores Descriptivos de Rendimiento (mpg)
Estadístico Valor
Moda 10.40
Varianza 36.32
Desviación Estándar 6.03
Coeficiente de Variación (%) 30.00
Curtosis -0.37

“Al evaluar la presencia de valores atípicos mediante el método del Rango Intercuartílico (1.5 \times IQR), no se detectaron valores atípicos (outliers) en la variable rendimiento (mpg). Los datos oscilan entre un mínimo de 10.4 mpg y un máximo de 33.9 mpg, manteniéndose todos dentro de las barreras aceptables de la distribución.”

# Cálculo de cuartiles y Rango Intercuartílico (IQR)
Q1 <- quantile(mtcars$mpg, 0.25)
Q3 <- quantile(mtcars$mpg, 0.75)
IQR_val <- IQR(mtcars$mpg)

# Límites para detectar atípicos
limite_inferior <- Q1 - 1.5 * IQR_val
limite_superior <- Q3 + 1.5 * IQR_val

# Identificar valores fuera de los límites
atipicos <- mtcars$mpg[mtcars$mpg < limite_inferior | mtcars$mpg > limite_superior]
atipicos # Devuelve numeric(0) -> No hay valores atípicos
[1] 33.9

6. Supuesto de Normalidad de los residuales

Planteamiento de las Hipótesis

Para cualquier prueba de normalidad sobre los residuales del modelo (\varepsilon_i):

  • \text{H}_0 (Hipótesis Nula): Los residuales del modelo siguen una distribución normal (\varepsilon_i \sim N(0, \sigma^2)).

  • \text{H}_1 (Hipótesis Alternativa): Los residuales del modelo no siguen una distribución normal.

Regla de Decisión (\alpha = 0.05):

  • Si p\text{-valor} > 0.05: No se rechaza \text{H}_0. Se concluye que existe evidencia suficiente de que se cumple el supuesto de normalidad en los residuales.

  • Si p\text{-valor} \le 0.05: Se rechaza \text{H}_0. Existen evidencias de que los residuales no se distribuyen de forma normal.

# 3. Prueba Principal: Shapiro-Wilk (Recomendada para n < 50)
shapiro_res <- shapiro.test(modelo$residuals)
print(shapiro_res)

    Shapiro-Wilk normality test

data:  modelo$residuals
W = 0.94508, p-value = 0.1044

Conclusión:

Como el p\text{-valor} = 0.1044 es mayor que el nivel de significancia \alpha = 0.05, no se rechaza la hipótesis nula (\text{H}_0). Por lo tanto, con un nivel de confianza del 95%, existe suficiente evidencia estadística para afirmar que los residuales del modelo satisfacen el supuesto de normalidad.

7. Supuesto de Homocedasticidad

Planteamiento de las Hipótesis

  • \text{H}_0 (Hipótesis Nula): Los residuales son homocedásticos (la varianza de los errores es constante, \sigma_i^2 = \sigma^2).

  • \text{H}_1 (Hipótesis Alternativa): Los residuales son heterocedásticos (la varianza de los errores no es constante).

Regla de Decisión (\alpha = 0.05):

  • Si p\text{-valor} > 0.05: No se rechaza \text{H}_0 \rightarrow Se cumple el supuesto de homocedasticidad.

  • Si p\text{-valor} \le 0.05: Se rechaza \text{H}_0 \rightarrow Existe evidencia de heterocedasticidad.

# 2. Prueba formal de Breusch-Pagan
prueba_bp <- bptest(modelo)
print(prueba_bp)

    studentized Breusch-Pagan test

data:  modelo
BP = 0.040438, df = 1, p-value = 0.8406

Mediante la prueba de Breusch-Pagan se obtuvo un estadístico \chi^2 = 0.0404 con un p\text{-valor} = 0.8406. Al ser el p\text{-valor} > 0.05, no existe evidencia estadística para rechazar la hipótesis nula con un nivel de confianza del 95%.

8. Supuesto de Independencia

Planteamiento de las Hipótesis

Para la prueba de Durbin-Watson:

  • \text{H}_0 (Hipótesis Nula): \rho = 0 (Los errores son independientes / No existe autocorrelación en los residuales).

  • \text{H}_1 (Hipótesis Alternativa): \rho \neq 0 (Existe autocorrelación en los residuales).

Regla de Decisión (\alpha = 0.05):

  • Si p\text{-valor} > 0.05: No se rechaza \text{H}_0 \rightarrow Se cumple el supuesto de independencia.

  • Si p\text{-valor} \le 0.05: Se rechaza \text{H}_0 \rightarrow Existe autocorrelación entre los residuales.

# 2. Prueba de Durbin-Watson
prueba_dw <- dwtest(modelo)
print(prueba_dw)

    Durbin-Watson test

data:  modelo
DW = 1.2517, p-value = 0.0102
alternative hypothesis: true autocorrelation is greater than 0
# 3. Gráfico de autocorrelación (ACF) y residuales en función del orden
datos_indep <- data.frame(
  Orden = 1:nrow(mtcars),
  Residuales = residuals(modelo)
)

# Gráfico de residuales seguidos por orden de observación
ggplot(datos_indep, aes(x = Orden, y = Residuales)) +
  geom_point(color = "#2b5c8f", size = 2.5) +
  geom_line(color = "#2b5c8f", linetype = "dashed", alpha = 0.5) +
  geom_hline(yintercept = 0, linetype = "dashed", color = "#d9534f", linewidth = 1) +
  labs(
    title = "Residuales por Orden de Observación",
    subtitle = "Evaluación de Independencia de Errores",
    x = "Orden de la Observación",
    y = "Residuales"
  ) +
  theme_minimal()

Al ejecutar la prueba sobre el modelo ajustado, se obtuvo un estadístico DW = 1.2517 con un p\text{-valor} = 0.0102. Al ser el p\text{-valor} \le 0.05, se rechaza la hipótesis nula de independencia bajo la ordenación original de los datos. Cabe destacar que, al tratarse de un conjunto de datos transversal sin estructura temporal, este resultado refleja que el orden de registro de los automóviles captura cierta dependencia no modelada o la omisión de covariables relevantes en la especificación lineal.

Conclusiones

  • El modelo satisface satisfactoriamente los supuestos de normalidad en los residuales y homocedasticidad (varianza constante).

  • No obstante, se detectó un desvío en el supuesto de independencia de los errores (Durbin-Watson p < \alpha). Dado que el conjunto de datos mtcars es de corte transversal, este hallazgo evidencia la omisión de covariables relevantes (como la potencia del motor hp o la cilindrada disp) que agrupan a vehículos con características similares.

Recomendación Final:

Aunque el modelo lineal simple con el peso como único predictor posee una alta capacidad explicativa (R^2 \approx 75\%), para corregir la estructura atrapada en los residuales y realizar inferencias totalmente insesgadas, se recomienda transicionar hacia un modelo de regresión lineal múltiple que incorpore covariables adicionales como la potencia (hp).