1.- MUESTRA DE DATA

Los datos del fichero EdadPesoGrasas.txt corresponden a tres variables medidas en 25 individuos: edad, peso y cantidad 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"

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)

La correlación entre cada par de variables es directa o positiva, peor la única correlación fuerte y positiva que encontramos se da entre la edad y las grasas.

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.0000000

Peso vs edad: La correlación es 24%, indicando que la correlación es directa y débil. Nos muestra que cuando aumenta la edad en un 24% de las ocasiones, se incrementa el peso del individuo.

Peso vs grasas: La correlación es 26,53%, indicando que la correlación es directa y débil. Nos muestra que cuando aumenta las grasas en un 26,53% de las ocasiones, se incrementa el peso del individuo.

Edad vs grasas: La correlación es directa y fuerte. Indica que de las veces que aumenta la edad, aumenta en un 83.73% las grasas del individuo.

cor.test(grasas$peso, grasas$edad)
## 
##  Pearson's product-moment correlation
## 
## data:  grasas$peso and grasas$edad
## t = 1.1857, df = 23, p-value = 0.2478
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.1713697  0.5801269
## sample estimates:
##       cor 
## 0.2400133

Se obtiene un p valor de 0.2478 > 0.05. Por tanto, la correlación no es estadísticamente significativa al 5%.

cor.test(grasas$peso, grasas$grasas)
## 
##  Pearson's product-moment correlation
## 
## data:  grasas$peso and grasas$grasas
## t = 1.3196, df = 23, p-value = 0.2
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.1450415  0.5977634
## sample estimates:
##       cor 
## 0.2652935

Se obtiene un p valor de 0.20 siendo superior al nivel de significación del 5%. Esto indica que la correlación no es estadísticamente significativa.

cor.test(grasas$edad, grasas$grasas)
## 
##  Pearson's product-moment correlation
## 
## data:  grasas$edad and grasas$grasas
## t = 7.346, df = 23, p-value = 1.794e-07
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  0.6608860 0.9260782
## sample estimates:
##       cor 
## 0.8373534

El p valor es cero siendo inferior al 5%. En este caso, la correlación es estadísticamente significativa al 5%.

#3.- ANÁLISIS DESCRIPTIVO (EDA)

summary(grasas)
##       peso            edad           grasas     
##  Min.   :27.00   Min.   :20.00   Min.   :181.0  
##  1st Qu.:63.00   1st Qu.:30.00   1st Qu.:254.0  
##  Median :69.00   Median :37.00   Median :303.0  
##  Mean   :68.68   Mean   :39.12   Mean   :310.7  
##  3rd Qu.:76.00   3rd Qu.:50.00   3rd Qu.:374.0  
##  Max.   :89.00   Max.   :60.00   Max.   :451.0

#4.- REGRESIÓN SIMPLE

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-07

El intercepto asciende a 102.57. Esto indica que, independientemente de la edad, el nivel de grasas del individuo será un prmoedio de 102.57. El coeficiente es estadísticamente diferente de cero.

La edad: se estima que el incremento de la edad de un año, con el resto de factores constantes provoca un incremento de 5.32 unidades la grasa del individuo. La variable edad es estadísticamente significativa al 5%.

El coeficiente de correlación R2 es de un 70.12%, indicando que la variación de la edad permite explicar el 70.12% de las variaciones en las grasas (R2 > 70%, la fiabilidad es elevada).

Respecto al contraste de significación global, se obtiene un estadístico F de 53.96 con un p valor prácticamente nulo. Se rechaza la hipótesis nula al 5%; por lo tanto, las pendientes son conjuntamente significativas.

plot(grasas$edad, grasas$grasas, xlab = "Edad", ylab = "Grasas")
abline(regresion)

# Predicciones
#=========================
nuevas.edades <- data.frame(edad = seq(18, 65))
predict(regresion, nuevas.edades)
##        1        2        3        4        5        6        7        8 
## 198.3473 203.6680 208.9887 214.3093 219.6300 224.9507 230.2714 235.5921 
##        9       10       11       12       13       14       15       16 
## 240.9127 246.2334 251.5541 256.8748 262.1954 267.5161 272.8368 278.1575 
##       17       18       19       20       21       22       23       24 
## 283.4781 288.7988 294.1195 299.4402 304.7608 310.0815 315.4022 320.7229 
##       25       26       27       28       29       30       31       32 
## 326.0435 331.3642 336.6849 342.0056 347.3263 352.6469 357.9676 363.2883 
##       33       34       35       36       37       38       39       40 
## 368.6090 373.9296 379.2503 384.5710 389.8917 395.2123 400.5330 405.8537 
##       41       42       43       44       45       46       47       48 
## 411.1744 416.4950 421.8157 427.1364 432.4571 437.7778 443.0984 448.4191
nuevas.edades <- data.frame(edad = seq(18, 65))

plot(grasas$edad, grasas$grasas, xlab = "Edad", ylab = "Grasas")
abline(regresion)

ic <- predict(regresion, nuevas.edades, interval = 'confidence')
lines(nuevas.edades$edad, ic[,2], lty=2, col = 'lightgreen')
lines(nuevas.edades$edad, ic[,3], lty=2, col = 'lightgreen')


#Intervalos de predicción
ic <- predict(regresion, nuevas.edades, interval = 'prediction')
lines(nuevas.edades$edad, ic[,2], lty = 2, col = 'blue')
lines(nuevas.edades$edad, ic[,3], lty = 2, col = 'blue')

El intervalo de la predicción puntual es más estrecho y por ello más preciso que el intervalo de la esperanza de la predicción.

#4.- DIAGNOSIS DEL MODELO

residuos <- regresion$residuals
valores.ajustados <- fitted(regresion)
head(residuos)
##           1           2           3           4           5           6 
##   6.6737469 -18.9886687  25.7496889   0.8045681  45.1463073  66.4079497
head(valores.ajustados)
##        1        2        3        4        5        6 
## 347.3263 208.9887 379.2503 262.1954 405.8537 235.5921
plot(valores.ajustados, residuos)

Los errores son independientes.

plot(residuos)

mean(residuos)
## [1] -5.143802e-16

La esperanza matemática de los residuos es cero, aceptándose la hipótesis de esperanza condicional nula. Por tanto, el estimador de MCO es insesgado.

qqnorm(residuos) 
qqline(residuos)

La mayoría de los puntos se ajustan adecuadamente a los cuantiles de la normal, se acepta la hipótesis de normalidad para los residuos.

shapiro.test(residuos)
## 
##  Shapiro-Wilk normality test
## 
## data:  residuos
## W = 0.9634, p-value = 0.4863

El p valor es 0.48 siendo superior al nivel de significación del 5%. Indicando que la distribución de los residuos es aproximadamente normal.

# TEST DE HOMOCEDASTICIDAD DE WHITE
# ====================================

u2 <- residuos**2
head(u2)
##            1            2            3            4            5            6 
##   44.5388975  360.5695383  663.0464807    0.6473298 2038.1890649 4410.0157831
# Regresión auxiliar de White
#============================
modelo_white <- lm(u2 ~edad + edad**2, data = grasas)
summary(modelo_white)
## 
## Call:
## lm(formula = u2 ~ edad + edad^2, data = grasas)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -1817.7 -1487.3  -892.9   458.4  6493.9 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)
## (Intercept) 2083.408   1473.200   1.414    0.171
## edad          -8.835     36.002  -0.245    0.808
## 
## Residual standard error: 2160 on 23 degrees of freedom
## Multiple R-squared:  0.002612,   Adjusted R-squared:  -0.04075 
## F-statistic: 0.06022 on 1 and 23 DF,  p-value: 0.8083
white_stat <- length(u2)*summary(modelo_white)$r.squared
white_stat
## [1] 0.06528991
p_v <- pchisq(white_stat, df = 2)
p_v
## [1] 0.03211786
p_v_1 <- 1-p_v
p_v_1
## [1] 0.9678821

El p valor es del 96.7% siendo superior al nivel de significación del 5%, indicando homocedasticidad.

5.- CONCLUSIÓN

Se observa que a medida que aumenta la edad del individuo, su índice de grasas aumenta de forma estadísticamente significativa. La regresión obtenida tiene una fiabilidad elevada, puesto que el coeficiente de deteerminación o R2 es superior al 70%. En cuanto a la diagnosis del modelo, se verifica la independencia, pues los datos proceden de un MAS (srs), con una media de los residuales igual a cero, distribución normal y varianza homocedástica. Por lo tanto, la inferencia realizada en el modelo (contrastes e intervalos de confianza) es correcta.