Análisis de regresión lineal
setwd("~/pye1pm")library(pacman)
p_load("DT","prettydoc","xfun","readr")Se tienen 3 variables entre nuestros datos
la edad en años cumplidos el peso en kilogramos El total de las 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 objetivo de conocer las relaciones existentes entre las distintas variables, se elabora una matriz de diagramas de dispersion, y se puede observar como no hay una relacion lineal entre el peso y las grasas.
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 es la formula y ~ x donde se deja clara la variable dependiente y la independiente
El siguiente argumento es data, nos dice cual es el fichero donde se pueden encontrar las variables. Los resultados se guardan en un objeto llamado regresion, donde se contienen todos los datos relevantes del analisis. Con summary se obtiene:
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
Veremos ahora que los datos vienen de un sistema de regresión simple con forma
\[ y_i = \beta_0 + \beta_1 x_i + \epsilon_i, \ \ \ \ i=1,\ldots,n, \] donde errores aleatorios ϵi son independientes a distribución norma media 0 y varianza σ2
Bajo este modelo:
Los errores de los estimadores de parametros β0 y β1 se encuentran en la columna std Error de la salida anterior. En el ejemplo, Sus valores son 21.769 y -0.9209324.
La columna t value contiene al estaidstico t, cociente entre los estimadores y su error tipico, que son la base para los contrastes H0:β0=0 y H0:β1=0, los p-valores aparecen en la columna Pr(>|t|)
Los intervalos de confianza se obtienen con el comando confint. Mientras que el parametro Level es para elegir el nivel de confianza que por defecto esta en .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´ de la salida anterior.
Los parametros de la ecuacion de recta minimos cuadrados que relaciona la cantidad de grasas en funcion del peso se da por la columna ‘estimate’ de la tabla ‘coefficients’
En este ejemplo la ecuación de la recta mínimos cuadrados es:
\[ y = 199.298 + 1.622 x \]
y = 199.298 + 1.622*(44)
y## [1] 270.666Gráfica
Los comandos representan la nube de puntos y añaden representacion grafica de la recta minimos cuadrados siendo un comando aplicado al objeto generado por lm
Los comandos siguientes
plot(grasas$peso, grasas$grasas, xlab="peso en kilos de las personas", ylab="Contenido de grasa en sangre")
abline(regresion)
El coeficiente de determinacion mide la bondad del ajuste de recta a los datos.
Estimación de predicciones
Si queremos utilizar la recta de minimos cuadrados. Basta crear un fichero que contenga las nuevas variables y usar el comando predict:
nuevos.pesos <- data.frame(peso = seq(20,90))
predict(regresion, nuevos.pesos)## 1 2 3 4 5 6 7 8
## 231.7444 233.3667 234.9890 236.6114 238.2337 239.8561 241.4784 243.1008
## 9 10 11 12 13 14 15 16
## 244.7231 246.3454 247.9678 249.5901 251.2125 252.8348 254.4572 256.0795
## 17 18 19 20 21 22 23 24
## 257.7018 259.3242 260.9465 262.5689 264.1912 265.8136 267.4359 269.0582
## 25 26 27 28 29 30 31 32
## 270.6806 272.3029 273.9253 275.5476 277.1700 278.7923 280.4146 282.0370
## 33 34 35 36 37 38 39 40
## 283.6593 285.2817 286.9040 288.5264 290.1487 291.7710 293.3934 295.0157
## 41 42 43 44 45 46 47 48
## 296.6381 298.2604 299.8828 301.5051 303.1274 304.7498 306.3721 307.9945
## 49 50 51 52 53 54 55 56
## 309.6168 311.2391 312.8615 314.4838 316.1062 317.7285 319.3509 320.9732
## 57 58 59 60 61 62 63 64
## 322.5955 324.2179 325.8402 327.4626 329.0849 330.7073 332.3296 333.9519
## 65 66 67 68 69 70 71
## 335.5743 337.1966 338.8190 340.4413 342.0637 343.6860 345.3083Por ejemplo, para un individuo de 30 kilos, predecimos una cantidad de grasas de 278.7923
Los intervalos de confianza para respuesta media y de prediccion se pueden obtener con el comando predict. El siguiente codigo calcula y representa los dos tipos de intervalos para el rango de peso de 20 a 90:
nuevos.pesos <- data.frame(peso = seq(20,90))
#Gráfico de dispersión y recta
plot(grasas$peso, grasas$grasas, xlab="peso en kilos 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)