Análisis de regresión lineal
setwd("~/pye1pm")library(pacman)
p_load("DT", "prettydoc", "xfun")Tenemos 3 variables en nuestro conjunto de datos: Edades de las personas Peso en Kg. de la persona *Grasas total de grasas en la sangre
grasas <- read.table("http://verso.mat.uam.es/~joser.berrendero/datos/EdadPesoGrasas.txt", header = TRUE)
names(grasas)## [1] "peso" "edad" "grasas"#Con el proposito de saber cuáles son las relaciones que existen entre cada variable, se representa a continuación una matriz de diagramas de dispersión.
pairs(grasas)
Cálculo y representación de la recta de mínimos cuadrados
El comando básico es lm (linear models).
El primer comando que se utiliza es una fórmula y ~ x, en donde esto nos dice cuál es la variable dependiente o independiente, en este caso la dependiente es la “y” grasas y la independiente es la “x” peso. El segundo argumento, “data”, es en dónde se especifica de dónde es que se sacan estas variables anteriormente mencionadas. Los dos argumentos anteriores son guardados en la variable “regresion”, y con el comando summary se obtiene un resumen de los principales resultados ##Aqui lo que se realiza es
regresion <- lm(grasas ~ peso, data = grasas )
summary(regresion)##
## Call:
## lm(formula = grasas ~ peso, data = grasas)
##
## Residuals:
## Min 1Q Median 3Q Max
## -127.729 -53.686 -9.239 46.537 128.404
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 199.298 85.818 2.322 0.0294 *
## peso 1.622 1.229 1.320 0.2000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 76.65 on 23 degrees of freedom
## Multiple R-squared: 0.07038, Adjusted R-squared: 0.02996
## F-statistic: 1.741 on 1 and 23 DF, p-value: 0.2Inferencia 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 21.7696 y -0.9209 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|).
El estimador de la desviación típica de los errores σ aparece como Residual standard error y su valor en el ejemplo es 76.65
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) 21.7696109 376.825392
## peso -0.9209324 4.165618confint(regresion, level=0.90)## 5 % 95 %
## (Intercept) 52.2166142 346.378389
## peso -0.4847468 3.729432Indice 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´ que se mostró anteriormente
En este ejemplo la ecuación de la recta de mínimos cuadrados es:
Aquí se puede observar la ecuación de la pendiente, en dónde el dato inicial(variable dependiente) es el error estimado interceptado y la variable independiente es el error estimado de peso en este caso, acompañado de la “x”, que en este caso será la variable en la que se le introduce el peso para poder calcular cuáles serán las grasas según el modelo que se ha planteado
\[ y = 199.298 + 1.622 x \]
##En este ejemplo se realizó utilizando el peso de 69 y nos arroja un resultado de 288.508 de grasas, el cual en este caso está dentro del error aceptado.
y = 199.298 + 1.622*(55)
y## [1] 288.508Grá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): En la siguiente gráfica se representa lo que son los datos de nuestra tabla proporcionada, junto a la recta que hemos generado con ayuda del comando lm
plot(grasas$peso, grasas$grasas, xlab="Peso de las personas", ylab="Contenido de grasa en sangre")
abline(regresion)
El coeficiente de determinación 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.07038.
Estimación de predicciones
Supongamos que queremos utilizar la recta de mínimos cuadrados para predecir la cantidad de grasas para individuos de pesos 50, 51, 52… 80 Basta crear un fichero de datos que contenga las nuevas variables regresoras y usar el comando predict:
nuevos.pesos <- data.frame(peso = seq(50,80))
predict(regresion, nuevos.pesos)## 1 2 3 4 5 6 7 8
## 280.4146 282.0370 283.6593 285.2817 286.9040 288.5264 290.1487 291.7710
## 9 10 11 12 13 14 15 16
## 293.3934 295.0157 296.6381 298.2604 299.8828 301.5051 303.1274 304.7498
## 17 18 19 20 21 22 23 24
## 306.3721 307.9945 309.6168 311.2391 312.8615 314.4838 316.1062 317.7285
## 25 26 27 28 29 30 31
## 319.3509 320.9732 322.5955 324.2179 325.8402 327.4626 329.0849Por ejemplo, para un individuo de peso de 69, predecimos una cantidad de grasas de 311.2391
- 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 pesos que va de 50 a 80 kg
nuevos.pesos <- data.frame(peso = seq(50,80))
#Gráfico de dispersión y recta
plot(grasas$peso, grasas$grasas, xlab="Peso 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, nuevos.pesos, interval = "confidence")
lines(nuevos.pesos$peso, ic[,2], lty=2)
lines(nuevos.pesos$peso, ic[,3], lty=2)