Análisis de Regresión Lineal Simple

El análisis de regresión lineal simple es una técnica estadística que se utiliza para modelar la relación entre una variable dependiente y una variable independiente. En este documento se presenta un ejemplo de cómo realizar un análisis de regresión lineal simple en R.

La ecuación de regresión lineal simple es:

\[Y = \alpha + \beta X + \epsilon\]

Donde:

Objetivos del Análisis de Regresión

  • Estimar la ecuación de regresión mediante estimadores de \(\alpha\) (a) y \(\beta\) (b) a partir de una muestra y generar intervalos de confianza para los mismos.

  • Estimar que tanto está la variable dependiente controlada por la variable independiente, en otras palabras que tanto la variación en y se explica por la variación en x.

  • Usar la ecuación de regresión para predecir valores de y a partir de valores de x.

Supuestos del Análisis de Regresión Lineal

  1. Linealidad: La relación entre las variables es lineal.

  2. Homocedasticidad: La varianza de los residuales es constante para los valores de x. Esto es equivalente a la homogeneidad de varianza en los grupos del ANOVA.

  3. Normalidad: Los residuales siguen una distribución normal.

  4. Independencia: Los residuales son independientes entre sí.

Estimación de los parámetros de la regresión

Para estimar los parámetros de la regresión, se utiliza el método de mínimos cuadrados. El objetivo es encontrar los valores de \(\alpha\) y \(\beta\) que minimizan la suma de los cuadrados de los residuos:

\[\sum_{i=1}^{n} (Y_i - \hat{Y}_i)^2\]

Donde:

  • \(Y_i\) es el valor observado de la variable dependiente.
  • \(\hat{Y}_i\) es el valor predicho de la variable dependiente.

Cálculo manual de los coeficientes de regresión

Para calcular los coeficientes de regresión manualmente, se utilizan las siguientes fórmulas:

\[\beta = \frac{\sum{(X_i - \bar{X})(Y_i - \bar{Y})}}{\sum{(X_i - \bar{X})^2}}\]

\[\alpha = \bar{Y} - \beta \bar{X}\]

Donde:

  • \(\beta\) es la pendiente de la recta de regresión.
  • \(\alpha\) es la intersección de la recta de regresión con el eje Y.
  • \(X_i\) y \(Y_i\) son los valores de las variables independiente y dependiente, respectivamente.
  • \(\bar{X}\) y \(\bar{Y}\) son las medias de las variables independiente y dependiente, respectivamente.

Ejemplos de Análisis de Regresión Lineal Simple

Ejemplo 1

Vamos a usar los datos de la longitud (cm) del ala derecha de gorriones de diferentes edades (días). En este caso, podemos asumir que la variable dependiente es la longitud del ala y la independiente, la edad de las aves.

edad <- c(3.0,4.0,5.0,6.0,8.0,9.0,10.0,11.0,12.0,14.0,15.0,16.0,17.0)
l.ala <- c(1.4,1.5,2.2,2.4,3.1,3.2,3.2,3.9,4.1,4.7,4.5,5.2,5.0)
sparrow <- data.frame(edad,l.ala)
#tabla de datos
knitr::kable(sparrow, col.names = c("Edad (días)", "Longitud del ala (cm)"), caption = 'Tabla 1.  Longitud del ala derecha de gorriones de diferente edad')
Tabla 1. Longitud del ala derecha de gorriones de diferente edad
Edad (días) Longitud del ala (cm)
3 1.4
4 1.5
5 2.2
6 2.4
8 3.1
9 3.2
10 3.2
11 3.9
12 4.1
14 4.7
15 4.5
16 5.2
17 5.0

Gráfica de dispersión

# Gráfica de dispersión
library(ggplot2)
ggplot(sparrow, aes(x = edad, y = l.ala)) +
  geom_point() +
  labs(x = "Edad (días)", y = "Longitud del ala (cm)")

Figura 1. Gráfica de dispersión de la longitud del ala derecha de gorriones en función de la edad.

Cálculo de los coeficientes de regresión con R

Para calcular los coeficientes de regresión con R, se utiliza la función lm().

# Ajuste del modelo de regresión lineal
modelo <- lm(l.ala ~ edad, data = sparrow)
# Coeficientes de regresión
coeficientes <- coef(modelo)
summary(modelo)
## 
## Call:
## lm(formula = l.ala ~ edad, data = sparrow)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.30699 -0.21538  0.06553  0.16324  0.22507 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.71309    0.14790   4.821 0.000535 ***
## edad         0.27023    0.01349  20.027 5.27e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2184 on 11 degrees of freedom
## Multiple R-squared:  0.9733, Adjusted R-squared:  0.9709 
## F-statistic: 401.1 on 1 and 11 DF,  p-value: 5.267e-10

Los coeficientes de regresión son \(\alpha = 0.7130945\) y \(\beta = 0.270229\).

Gráfica de la recta de regresión con la ecuación y el intervalo de confianza

# Gráfica de la recta de regresión con ecuación
eq <- substitute(
  italic(y) == a + b %.% italic(x),
  list(
    a = round(as.numeric(coeficientes[1]), 3),
    b = round(as.numeric(coeficientes[2]), 3)
  )
)
# plot
ggplot(sparrow, aes(x = edad, y = l.ala)) +
  geom_point() +
  geom_smooth(method = "lm", se = TRUE, color = "red") +
  annotate("text", x = 5.5, y = 5, label = as.character(as.expression(eq)), parse = TRUE) +
  labs(x = "Edad (días)", y = "Longitud del ala (cm)")
## `geom_smooth()` using formula = 'y ~ x'

Figura 2. Gráfica de dispersión de la longitud del ala derecha de gorriones en función de la edad con la recta de regresión, intervalo de confianza (95%) y la ecuación de la recta.

Intervalo de confianza para la pendiente (\(\beta\))

# Intervalo de confianza para la pendiente
confint(modelo, "edad", level = 0.95)
##          2.5 %    97.5 %
## edad 0.2405308 0.2999272

Prueba de hipótesis para la pendiente

Para realizar una prueba de hipótesis para la pendiente, se plantean las siguientes hipótesis:

  • Hipótesis nula (\(H_0\)): La pendiente de la recta de regresión es igual a cero (\(\beta = 0\)).
  • Hipótesis alternativa (\(H_1\)): La pendiente de la recta de regresión es diferente de cero (\(\beta \neq 0\)).
# Prueba de hipótesis para la pendiente
summary(modelo)$coefficients
##              Estimate Std. Error   t value     Pr(>|t|)
## (Intercept) 0.7130945 0.14790445  4.821319 5.347485e-04
## edad        0.2702290 0.01349312 20.027169 5.267053e-10

El valor p para la pendiente es 5.2670535^{-10}.

Prueba del supuesto de heterocedasticidad

Para probar el supuesto de homocedasticidad, se puede realizar un análisis de los residuales.

# Análisis de los residuales
residuos <- residuals(modelo)
plot(modelo, which = 1)

# Prueba de Breusch-Pagan
library(lmtest)
bptest(modelo)
## 
##  studentized Breusch-Pagan test
## 
## data:  modelo
## BP = 1.6349, df = 1, p-value = 0.201

Aunque el gráfico de residuos vs. valores ajustados muestra cierta dispersión, la prueba de Breusch-Pagan no rechaza la hipótesis nula de homocedasticidad.

Predicción de la longitud del ala para una edad dada con intervalo de confianza

Para predecir la longitud del ala para una edad dada, se utiliza la función predict().

# Predicción de la longitud del ala para una edad dada
nueva_edad <- 7
prediccion <- predict(modelo, newdata = data.frame(edad = nueva_edad), interval = "confidence", level = 0.95)
prediccion
##        fit      lwr      upr
## 1 2.604698 2.444344 2.765051

Predicción de la edad para una longitud del ala dada con intervalo de confianza

Para predecir la edad para una longitud del ala dada, se utiliza la función predict() con el modelo invertido.

# Predicción de la edad para una longitud del ala dada
nueva_l_ala <- c(4.5, 3)
prediccion_inversa <- predict(lm(edad ~ l.ala, data = sparrow), newdata = data.frame(l.ala = nueva_l_ala), interval = "confidence", level = 0.95)
prediccion_inversa
##         fit       lwr      upr
## 1 13.906551 13.257517 14.55558
## 2  8.503874  7.990108  9.01764

Alternativa

Si los supuestos de normalidad y homocedasticidad no se cumplen, se puede utilizar otras alternativas que incluyen:

  • Transformaciones de las variables (logarítmica, cuadrática, etc.).

  • Regresión polinómica.

  • Modelos no-paramétricos como la regresión de Loess.

Regresión polinómica

# Regresión polinómica
modelo_polinomico <- lm(l.ala ~ poly(edad, 2), data = sparrow)
# Gráfica de la regresión polinómica
ggplot(sparrow, aes(x = edad, y = l.ala)) +
  geom_point() +
  geom_smooth(method = "lm", formula = y ~ poly(x, 2), se = TRUE, color = "green") +
  labs(x = "Edad (días)", y = "Longitud del ala (cm)")

Figura 3. Gráfica de dispersión de la longitud del ala derecha de gorriones en función de la edad con la regresión polinómica y el intervalo de confianza (95%).

Coeficientes de la regresión polinómica:

# Coeficientes de la regresión polinómica
coeficientes_polinomico <- coef(modelo_polinomico)
coeficientes_polinomico
##    (Intercept) poly(edad, 2)1 poly(edad, 2)2 
##      3.4153846      4.3740386     -0.3683301
# p-value
p_value_polinomico <- summary(modelo_polinomico)$coefficients[2,4]
p_value_polinomico
## [1] 7.801604e-10
# R^2
r2_polinomico <- summary(modelo_polinomico)$r.squared
r2_polinomico
## [1] 0.9802084

Modelos no-lineales Regresión de Loess

# Regresión de Loess
modelo_loess <- loess(l.ala ~ edad, data = sparrow)
# Gráfica de la regresión de Loess
ggplot(sparrow, aes(x = edad, y = l.ala)) +
  geom_point() +
  geom_smooth(method = "loess", se = TRUE, color = "blue") +
  labs(x = "Edad (días)", y = "Longitud del ala (cm)")
## `geom_smooth()` using formula = 'y ~ x'

Figura 4. Gráfica de dispersión de la longitud del ala derecha de gorriones en función de la edad con la regresión de Loess y el intervalo de confianza (95%).

Regresión múltiple

La regresión múltiple es una extensión de la regresión lineal simple que permite modelar la relación entre una variable dependiente y múltiples variables independientes. La ecuación de regresión múltiple es: \[Y = \alpha + \beta_1 X_1 + \beta_2 X_2 + ... + \beta_n X_n + \epsilon\] Donde: - \(Y\) es la variable dependiente. - \(X_1, X_2, ..., X_n\) son las variables independientes. - \(\alpha\) es la intersección de la recta de regresión con el eje Y. - \(\beta_1, \beta_2, ..., \beta_n\) son las pendientes de la recta de regresión. - \(\epsilon\) es el error aleatorio o residual.

Ejemplo de Regresión Múltiple con datos de mortalidad y salud en ciudades pequeñas

# Cargar los datos
library(MASS)
library(readxl)
small_cities <- read_excel("death_small_cities.xlsx", sheet = "data")

# modelo
modelo_multiplo <- lm(death1K ~ income1K + doctor100K + hosp100K + density, data = small_cities)
# Resumen del modelo
summary(modelo_multiplo)
## 
## Call:
## lm(formula = death1K ~ income1K + doctor100K + hosp100K + density, 
##     data = small_cities)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.6404 -0.7904  0.3053  0.9164  2.7906 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 12.2662552  2.0201467   6.072 1.95e-07 ***
## income1K    -0.3302302  0.2345518  -1.408   0.1656    
## doctor100K   0.0073916  0.0069336   1.066   0.2917    
## hosp100K     0.0005837  0.0007219   0.809   0.4228    
## density     -0.0094629  0.0048868  -1.936   0.0587 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.601 on 48 degrees of freedom
## Multiple R-squared:  0.1437, Adjusted R-squared:  0.07235 
## F-statistic: 2.014 on 4 and 48 DF,  p-value: 0.1075

Selección de mejor predictor

# Selección de mejor predictor
library(MASS)
modelo_seleccion <- stepAIC(modelo_multiplo, direction = "both")
## Start:  AIC=54.65
## death1K ~ income1K + doctor100K + hosp100K + density
## 
##              Df Sum of Sq    RSS    AIC
## - hosp100K    1    1.6763 124.75 53.369
## - doctor100K  1    2.9139 125.99 53.892
## <none>                    123.07 54.652
## - income1K    1    5.0825 128.16 54.797
## - density     1    9.6144 132.69 56.639
## 
## Step:  AIC=53.37
## death1K ~ income1K + doctor100K + density
## 
##              Df Sum of Sq    RSS    AIC
## <none>                    124.75 53.369
## - doctor100K  1    5.1882 129.94 53.529
## - income1K    1    6.1544 130.91 53.921
## + hosp100K    1    1.6763 123.07 54.652
## - density     1    8.3192 133.07 54.791
# Resumen del modelo
summary(modelo_seleccion)
## 
## Call:
## lm(formula = death1K ~ income1K + doctor100K + density, data = small_cities)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.7852 -0.6675  0.4010  0.9348  2.7495 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 12.565900   1.978836   6.350 6.74e-08 ***
## income1K    -0.359140   0.230990  -1.555   0.1264    
## doctor100K   0.009284   0.006504   1.428   0.1598    
## density     -0.008580   0.004746  -1.808   0.0768 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.596 on 49 degrees of freedom
## Multiple R-squared:  0.132,  Adjusted R-squared:  0.0789 
## F-statistic: 2.485 on 3 and 49 DF,  p-value: 0.07163