#7b. Exercícios de Regressão Múltipla #------------------------------------------------------------------------------ #Uma estimativa da incerteza na previsão do modelo #------------------------------------------------------------------------------ #Upload dos dados bebes <- read.table("babies.csv", header= TRUE, as.is = TRUE, sep= "\t") str(bebes) #verificando a estrutura dos dados bebes #conversão de onça para quilograma em bwt babies$bwt <- babies$bwt/35.274 #retirando os dados faltantes babies <- babies[babies$gestation!=999,] mod1 <- lm(babies$bwt~babies$gestation) summary(mod1) #x : valor da variável preditora para o qual se quer estimar o erro padrão; #criando um vetor para representar a variável preditora como uma sequência de 100 valores #compreendidos entre os valores máximo e mínimo da variável preditora observada x <- seq(min(babies$gestation), max(babies$gestation), length.out= 100) #s2: variância dos valores de y observados vari <- var(babies$bwt, na.rm = T) #n: número de observações n <- 1/length(babies$gestation) #X menos a média dos valores de x observados elevado ao quadrado m <- (x-mean(babies$gestation))^2 #soma dos desvios quadráticos dos valores de x observados SSX <- sum(m) #------------------------------------------------------------------------------- #Galileu estava Certo? #------------------------------------------------------------------------------- ############## Ajuste de Polinômios ########################### #Galileu Galilei investigou a distância percorrida por um corpo lançado de uma rampa #a diferentes alturas. Para isto ele usou um plano inclinado a uma certa altura do chão, #de onde soltou bolas de metal. Com este experimento ele verificou que a relação entre #a altura de lançamento e a distância percorrida não é linear. #As alturas de lançamentos e as distâncias percorridas obtidas por Galileu são: init.h = c(600, 700, 800, 950, 1100, 1300, 1500)#alturas de lançamentos h.d = c(253, 337, 395, 451, 495, 534, 573)#distâncias percorridas #Podemos ajustar os dois modelos e usar o comando anova para testar se o acréscimo do termo #quadrático implicou em melhora. O ajuste da reta é feito do modo usual: mod1 <- lm(h.d~init.h) #Os comando para ajustar a equação da parábola são: mod2 <- lm(h.d~init.h+I(init.h^2)) # E mod3 <- lm(h.d~init.h+I(init.h^2)+I(init.h^3)) #E agora podemos comparar os modelos dois a dois: anova(mod1,mod2) #Analysis of Variance Table #Model 1: h.d ~ init.h #Model 2: h.d ~ init.h + I(init.h^2) # Res.Df RSS Df Sum of Sq F Pr(>F) #1 5 5671.2 #2 4 744.1 1 4927.1 26.487 0.00676 ** #--- #Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 anova(mod2,mod3) #Analysis of Variance Table #Model 1: h.d ~ init.h + I(init.h^2) #Model 2: h.d ~ init.h + I(init.h^2) + I(init.h^3) # Res.Df RSS Df Sum of Sq F Pr(>F) #1 4 744.08 #2 3 48.25 1 695.82 43.26 0.00715 ** #--- #Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 #Pelos dados fornecidos pelas anovas, conseguimos perceber que o #modelo de parábola mod3 é superior aos outros dois. RSS, que indica #a variação dos dados que nãoé explicada pelos dados é a menor entre #os modelos. #### Plotando plot(h.d~init.h) #para mod1 abline(mod1)#para mod1 #coeficientes mod2 cf.m2 <- coef(mod2) # (Intercept) init.h I(init.h^2) #-2.401719e+02 1.052016e+00 -3.436937e-04 #coeficientes mod3 cf.m3 <- coef(mod3) # (Intercept) init.h I(init.h^2) I(init.h^3) #-7.815729e+02 2.771023e+00 -2.066508e-03 5.477104e-07 #curvas de mod2 e mod3 considerando equação y=a+b*X+c*X^2+d*X^3 curve(cf.m2[1]+cf.m2[2]*x+cf.m2[3]*x^2, add=T, lty=2)#mod2 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")#mod3 Podemos observar que houver um melhor ajuste dos dados em mod3, o polinômio com maior grau. #------------------------------------------------------------------------------- #Massa de recém-Nascidos #------------------------------------------------------------------------------- #Upload dos dados bebes <- read.table("babies.csv", header= TRUE, as.is = TRUE, sep= "\t") str(bebes) #verificando a estrutura dos dados bebes #Retirando todos os dados faltantes babies <- babies[babies$bwt!=999,] #Retirando os 999 babies <- babies[babies$gestation!=999,] #Retirando os 999 babies <- babies[babies$parity!=9,] #Retirando os 9 babies <- babies[babies$age!=99,] #Retirando os 99 babies <- babies[babies$height!=99,] #Retirando os 99 babies <- babies[babies$weight!=999,] #Retirando os 999 babies <- babies[babies$smoke!=9,] #Retirando os 9 bebes #verificando dados novamente