Análisis de regresión lineal
setwd("~/Documents/pye1pm")library(pacman)
p_load("DT","prettydoc","xfun","readr")Tenemos 3 variables en nuestro conjunto de datos:
- Peso en Kg. de la persona
- Edades de las personas
- Contenido total de Grasas en sangre
grasas <- read.table("http://verso.mat.uam.es/~joser.berrendero/datos/EdadPesoGrasas.txt", header = TRUE)
names(grasas)## [1] "peso" "edad" "grasas"Otra manera de obtener los datos es desde un archivo, en este caso es un archivo de valores separados por comas (CSV)
grasas1 <- read_csv("grasas.csv")##
## ── Column specification ────────────────────────────────────────────────────────
## cols(
## peso = col_double(),
## edad = col_double(),
## grasas = col_double()
## )datatable(grasas1)Con el fin de conocer las relaciones existentes entre cada par de variables podemos representar una matriz de diagramas de dispersión. Al parecer existe una relación lineal bastante clara entre la edad y las grasas, pero no entre los otros dos pares de variables. Por otra parte el fichero contiene un dato atípico.
pairs(grasas)
## Cálculo y representación de la recta de mínimos cuadrados
El comando básico es lm (linear models).
El primer argumento de este comando es una fórmula y ~ x en la que se especifica cuál es la variable respuesta o dependiente (y) y cuál es la variable regresora o independiente (x).
El segundo argumento, llamado data especifica cuál es el fichero en el que se encuentran las variables. El resultado lo guardamos en un objeto llamado regresion. Este objeto es una lista que contiene toda la información relevante sobre el análisis. Mediante el comando summary obtenemos un resumen de los principales resultados:
regresion <- lm(grasas ~ edad, data = grasas )
summary(regresion)##
## Call:
## lm(formula = grasas ~ edad, data = grasas)
##
## Residuals:
## Min 1Q Median 3Q Max
## -63.478 -26.816 -3.854 28.315 90.881
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 102.5751 29.6376 3.461 0.00212 **
## edad 5.3207 0.7243 7.346 1.79e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 43.46 on 23 degrees of freedom
## Multiple R-squared: 0.7012, Adjusted R-squared: 0.6882
## F-statistic: 53.96 on 1 and 23 DF, p-value: 1.794e-07Inferencia en el modelo de regresión simple
Suponemos ahora que los datos proceden de un modelo de regresión simple de la forma
\[ y_i = \beta_0 + \beta_1 x_i + \epsilon_i, \ \ \ \ i=1,\ldots,n, \] donde los errores aleatorios ϵi son independientes con distribución normal de media 0 y varianza σ2
Bajo este modelo:
Los errores típicos de los estimadores de los parámetros β0 y β1 se encuentran en la columna Std Error de la salida anterior. En el ejemplo, sus valores son 29.638 y 0.724 respectivamente.
La columna t value contiene el estadístico t, es decir, cociente entre cada estimador y su error típico. Estos cocientes son la base para llevar a cabo los contrastes H0:β0=0 y H0:β1=0 . Los correspondientes p-valores aparecen en la columna Pr(>|t|). En este caso son muy pequeños por lo que se rechazan ambas hipótesis para los niveles de significación habituales.
El estimador de la desviación típica de los errores σ aparece como Residual standard error y su valor en el ejemplo es 43.5
Los intervalos de confianza para los parámetros se obtienen con el comando confint. El parámetro level permite elegir el nivel de confianza (por defecto es 0.95):
confint(regresion)## 2.5 % 97.5 %
## (Intercept) 41.265155 163.885130
## edad 3.822367 6.818986confint(regresion, level=0.90)## 5 % 95 %
## (Intercept) 51.780153 153.370132
## edad 4.079335 6.562018Indice de correlación
de tipo pearson
cor(grasas)## peso edad grasas
## peso 1.0000000 0.2400133 0.2652935
## edad 0.2400133 1.0000000 0.8373534
## grasas 0.2652935 0.8373534 1.0000000Los parámetros de la ecuación de la recta de mínimos cuadrados que relaciona la cantidad de grasas en la sangre en función del peso vienen dados por la columna ´Estimate´ de la tabla ´Coefficients´ de la salida anterior.
Por lo tanto, en este ejemplo la ecuación de la recta de mínimos cuadrados es:
\[ y = 102.575 + 5.321 x \]
y = 102.575 + 5.321*(50)
y## [1] 368.625Gráfica
Los siguientes comandos representan la nube de puntos (comando plot) y añaden la representación gráfica de la recta de mínimos cuadrados (comando abline aplicado al objeto generado por lm):
plot(grasas$edad, grasas$grasas, xlab="Edad de las personas", ylab="Contenido de grasa en sangre")
abline(regresion)
El coeficiente de determinación (es decir, el coeficiente de correlación al cuadrado) mide la bondad del ajuste de la recta a los datos. A partir de la salida anterior, vemos que su valor en este caso es Multiple R-squared: 0.701.
Estimación de predicciones
Supongamos que queremos utilizar la recta de mínimos cuadrados para predecir la cantidad de grasas para individuos de edades 31,31,32,…,50 . Basta crear un fichero de datos que contenga las nuevas variables regresoras y usar el comando predict:
nuevas.edades <- data.frame(edad = seq(20,60))
predict(regresion, nuevas.edades)## 1 2 3 4 5 6 7 8
## 208.9887 214.3093 219.6300 224.9507 230.2714 235.5921 240.9127 246.2334
## 9 10 11 12 13 14 15 16
## 251.5541 256.8748 262.1954 267.5161 272.8368 278.1575 283.4781 288.7988
## 17 18 19 20 21 22 23 24
## 294.1195 299.4402 304.7608 310.0815 315.4022 320.7229 326.0435 331.3642
## 25 26 27 28 29 30 31 32
## 336.6849 342.0056 347.3263 352.6469 357.9676 363.2883 368.6090 373.9296
## 33 34 35 36 37 38 39 40
## 379.2503 384.5710 389.8917 395.2123 400.5330 405.8537 411.1744 416.4950
## 41
## 421.8157Por ejemplo, para un individuo de 30 años, predecimos una cantidad de grasas de 262.2
- Los intervalos de confianza para la respuesta media y los intervalos de predicción para la respuesta se pueden obtener usando el comando predict. Por ejemplo, el siguiente código calcula y representa los dos tipos de intervalos para el rango de edades que va de 20 a 60 años (los de predicción en rojo):
nuevas.edades <- data.frame(edad = seq(20,60))
#Gráfico de dispersión y recta
plot(grasas$edad, grasas$grasas, xlab="Edad de las personas", ylab="Contenido de grasa en sangre")
abline(regresion)
# Intervalos de confianza de la respuesta media:
# ic es una matriz con tres columnas: la primera es la prediccion, las otras dos son los extremos del intervalo
ic <- predict(regresion, nuevas.edades, interval = "confidence")
lines(nuevas.edades$edad, ic[,2], lty=2)
lines(nuevas.edades$edad, ic[,3], lty=2)