## 7.2 O modelo mais simples possível library(MASS) data(Animals) anim.m2 <- lm(log(brain)~log(body),data=Animals, subset=!(log(Animals$body)>8&log(Animals$brain)<6)) anim.m0 <- lm(log(brain)~1, data=Animals, subset=!(log(Animals$body)>8&log(Animals$brain)<6)) anim.m0 <- update(anim.m2, .~. -log(body)) anova(anim.m0,anim.m2) anova(anim.m2) #1) Os comandos retornam a mesma tabela summary(anim.m0) mean(log(Animals$brain[!(log(Animals$body)>8&log(Animals$brain)<6)])) sd(log(Animals$brain[!(log(Animals$body)>8&log(Animals$brain)<6)])) #2) o valor do comando 'mean' retorna o mesmo valor que o intercepto estimado para o modelo anim.m0, e o valor obtido através do comando sd é o mesmo do erro parão residual do modelo. ## 7B Exercícios de regressão múltipla ##1 babies <- read.csv("babies.txt", sep="") babies <- babies[babies$gestation!="999",] #Retirando NAs babies$bwt <- babies$bwt/35.274 #Convertendo o peso para Kg m1 <- lm(babies$bwt~babies$gestation) summary(m1) prevv <- seq(min(babies$gestation),max(babies$gestation), length.out= 100) s2 <- var(babies$bwt, na.rm = T) n1 <- 1/length(babies$gestation) xx2 <- (prevv-mean(babies$gestation))^2 ssx <- sum(xx2) sey <- sqrt(s2*(n1+xx2/ssx)) sey t <- -qt(0.975, 1221-2) coef(m1) ic <- -0.28531452 + 0.01316161*prevv ic sup <- ic+(t*sey) inf <- ic-(t*sey) par(las=1) plot(babies$bwt~babies$gestation, xlab= "Tempo de gestação (dias)", ylab = "Peso do bebê (Kg)") abline(m1, col="red", lwd= 2) lines(x=prevv, y=sup, col="blue") lines(x=prevv, y=inf, col="blue") ##2 init.h = c(600, 700, 800, 950, 1100, 1300, 1500) h.d = c(253, 337, 395, 451, 495, 534, 573) mod1 <- lm(h.d~init.h) mod2 <- lm(h.d~init.h+I(init.h^2)) mod3 <- lm(h.d~init.h+I(init.h^2)+I(init.h^3)) anova(mod1,mod2) anova(mod2,mod3) plot(h.d~init.h) abline(mod1) cf.m2 <- coef(mod2) cf.m3 <- coef(mod3) curve(cf.m2[1]+cf.m2[2]*x+cf.m2[3]*x^2, add=T, lty=2) curve(cf.m3[1]+cf.m3[2]*x+cf.m3[3]*x^2 +cf.m3[4]*x^3, add=T, lwd=2, lty=1, col= "blue") #Sim, o polinômio de terceiro grau descreve ainda melhor os dados de Galileu que o de segundo grau. ##3 babies <- read.csv("babies.txt", sep="") babies <- babies[babies$gestation!="999"&babies$age!="99"&babies$height!="99"&babies$weight!="999"&babies$smoke!="9",] mcheio <- lm(bwt~gestation+parity+age+height+weight+smoke, data = babies) library(MASS) Scope=list(upper=~gestation+parity+age+height+weight+smoke,lower=~1) melhormodelo <- stepAIC(mcheio,Scope,direction="both") melhormodelo$anova #resultados