library(lmtest)
## Carregando pacotes exigidos: zoo
##
## Anexando pacote: 'zoo'
## Os seguintes objetos são mascarados por 'package:base':
##
## as.Date, as.Date.numeric
dados <- read.csv("BGSgirls.csv", header = T, sep = ",")[-c(42,68),]
modelo <- lm(BMI18~HT2+WT2+HT9+WT9+ST9,data=dados)
summary(modelo)
##
## Call:
## lm(formula = BMI18 ~ HT2 + WT2 + HT9 + WT9 + ST9, data = dados)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.1855 -1.0815 -0.2294 0.9603 3.2969
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 26.627571 6.107523 4.360 5.0e-05 ***
## HT2 -0.124733 0.091238 -1.367 0.177
## WT2 -0.302688 0.192597 -1.572 0.121
## HT9 0.007765 0.067743 0.115 0.909
## WT9 0.290263 0.054356 5.340 1.4e-06 ***
## ST9 -0.013584 0.015807 -0.859 0.393
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.477 on 62 degrees of freedom
## Multiple R-squared: 0.4222, Adjusted R-squared: 0.3755
## F-statistic: 9.059 on 5 and 62 DF, p-value: 1.664e-06
hist(rstandard(modelo))
shapiro.test(rstandard(modelo))
##
## Shapiro-Wilk normality test
##
## data: rstandard(modelo)
## W = 0.98664, p-value = 0.6829
qqnorm(residuals(modelo))
qqline(residuals(modelo), col = "red")
plot(rstandard(modelo),fitted(modelo))
bptest(modelo)
##
## studentized Breusch-Pagan test
##
## data: modelo
## BP = 6.4566, df = 5, p-value = 0.2643
gqtest(modelo)
##
## Goldfeld-Quandt test
##
## data: modelo
## GQ = 0.68006, df1 = 28, df2 = 28, p-value = 0.8433
## alternative hypothesis: variance increases from segment 1 to 2
plot(rstandard(modelo))
dwtest(modelo)
##
## Durbin-Watson test
##
## data: modelo
## DW = 2.441, p-value = 0.9679
## alternative hypothesis: true autocorrelation is greater than 0
modelo1 <- lm(BMI18~WT2+HT9+WT9+ST9,data=dados)
summary(modelo1)
##
## Call:
## lm(formula = BMI18 ~ WT2 + HT9 + WT9 + ST9, data = dados)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.2235 -0.8877 -0.1233 0.9475 3.2538
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 24.210689 5.886245 4.113 0.000115 ***
## WT2 -0.421217 0.173161 -2.433 0.017847 *
## HT9 -0.050143 0.053230 -0.942 0.349786
## WT9 0.307244 0.053281 5.766 2.64e-07 ***
## ST9 -0.007989 0.015373 -0.520 0.605114
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.488 on 63 degrees of freedom
## Multiple R-squared: 0.4047, Adjusted R-squared: 0.3669
## F-statistic: 10.71 on 4 and 63 DF, p-value: 1.101e-06
modelo2 <- lm(BMI18~WT2+WT9+ST9,data=dados)
summary(modelo2)
##
## Call:
## lm(formula = BMI18 ~ WT2 + WT9 + ST9, data = dados)
##
## Residuals:
## Min 1Q Median 3Q Max
## -3.0286 -0.9935 -0.1336 0.8735 3.5814
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 18.86520 1.56299 12.070 < 2e-16 ***
## WT2 -0.44367 0.17136 -2.589 0.0119 *
## WT9 0.28203 0.04603 6.127 6.13e-08 ***
## ST9 -0.01375 0.01409 -0.976 0.3328
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.486 on 64 degrees of freedom
## Multiple R-squared: 0.3963, Adjusted R-squared: 0.368
## F-statistic: 14.01 on 3 and 64 DF, p-value: 4.018e-07
modelo3 <- lm(BMI18~WT2+WT9,data=dados)
summary(modelo3)
##
## Call:
## lm(formula = BMI18 ~ WT2 + WT9, data = dados)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.8964 -1.0368 -0.1033 0.9647 3.6818
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 18.83382 1.56209 12.057 < 2e-16 ***
## WT2 -0.47531 0.16820 -2.826 0.00626 **
## WT9 0.26929 0.04412 6.103 6.43e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 1.486 on 65 degrees of freedom
## Multiple R-squared: 0.3874, Adjusted R-squared: 0.3685
## F-statistic: 20.55 on 2 and 65 DF, p-value: 1.214e-07
library(rgl)
plot3d(dados$WT2, dados$WT9, dados$BMI18, type = "s", radius = 0.8, col = "blue", xlab = "WT2", ylab = "WT9", zlab = "BMI18")
WT2_seq <- seq(min(dados$WT2), max(dados$WT2), length.out = 40)
WT9_seq <- seq(min(dados$WT9), max(dados$WT9), length.out = 40)
grid <- expand.grid(WT2 = WT2_seq, WT9 = WT9_seq)
grid$y_hat <- predict(modelo3, grid)
surface3d(x = WT2_seq, y = WT9_seq, z = matrix(grid$y_hat, nrow = length(WT2_seq)), alpha=0.6, front="lines", back="lines", col="red")