# install.packages("lmtest")
library(lmtest)
# install.packages("olsrr")
library(olsrr)
# install.packages("corrplot")
library(corrplot)
library(tidyverse)
Como exemplo, vamos usar uma base de treino nativa do R.
mtcars <- mtcars
Relembrando o comando da regressão linear para investigar a relação entre o consumo de combustível e a potência de alguns modelos de carro:
mtcars %>%
lm(mpg ~ hp, data = .) %>%
summary()
##
## Call:
## lm(formula = mpg ~ hp, data = .)
##
## Residuals:
## Min 1Q Median 3Q Max
## -5.7121 -2.1122 -0.8854 1.5819 8.2360
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 30.09886 1.63392 18.421 < 2e-16 ***
## hp -0.06823 0.01012 -6.742 1.79e-07 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 3.863 on 30 degrees of freedom
## Multiple R-squared: 0.6024, Adjusted R-squared: 0.5892
## F-statistic: 45.46 on 1 and 30 DF, p-value: 1.788e-07
Para extrair os resíduos do modelo, podemos armazenar o modelo em um objeto e chamar o atributo “residuals” do objeto com o modelo:
reg <- mtcars %>%
lm(mpg ~ hp, data = .)
reg$residuals
## Mazda RX4 Mazda RX4 Wag Datsun 710 Hornet 4 Drive
## -1.59374995 -1.59374995 -0.95363068 -1.19374995
## Hornet Sportabout Valiant Duster 360 Merc 240D
## 0.54108812 -4.83489134 0.91706759 -1.46870730
## Merc 230 Merc 280 Merc 280C Merc 450SE
## -0.81717412 -2.50678234 -3.90678234 -1.41777049
## Merc 450SL Merc 450SLC Cadillac Fleetwood Lincoln Continental
## -0.51777049 -2.61777049 -5.71206353 -5.02978075
## Chrysler Imperial Fiat 128 Honda Civic Toyota Corolla
## 0.29364342 6.80420581 3.84900992 8.23597754
## Toyota Corona Dodge Challenger AMC Javelin Camaro Z28
## -1.98071757 -4.36461883 -4.66461883 -0.08293241
## Pontiac Firebird Fiat X1-9 Porsche 914-2 Lotus Europa
## 1.04108812 1.70420581 2.10991276 8.01093488
## Ford Pantera L Ferrari Dino Maserati Bora Volvo 142E
## 3.71340487 1.54108812 7.75761261 -1.26197823
Com essa informação, podemos contruir um gráfico para analisar nossos
resíduos em comparação aos valores previstos pelo modelo usando o
ggplot
.
residuos <- reg$residuals
valores_previstos <- predict(reg)
residuos_fitted <- cbind(residuos, valores_previstos) %>% as.data.frame()
residuos_fitted <- as.data.frame(residuos)
ggplot(residuos_fitted, aes(x = valores_previstos, y = residuos))+
geom_point() +
stat_smooth(method='lm')
## `geom_smooth()` using formula 'y ~ x'
Também podemos tirar um gráfico de dispersão entre resíduos e os valores previstos pelo modelo apenas com a função plot(). (Considere apenas o gráfico Residuals vs Fitted)
plot(reg)
Para verificar a matriz de correlação entre as variáveis de uma tabela, usamos:
cor(mtcars)
## mpg cyl disp hp drat wt
## mpg 1.0000000 -0.8521620 -0.8475514 -0.7761684 0.68117191 -0.8676594
## cyl -0.8521620 1.0000000 0.9020329 0.8324475 -0.69993811 0.7824958
## disp -0.8475514 0.9020329 1.0000000 0.7909486 -0.71021393 0.8879799
## hp -0.7761684 0.8324475 0.7909486 1.0000000 -0.44875912 0.6587479
## drat 0.6811719 -0.6999381 -0.7102139 -0.4487591 1.00000000 -0.7124406
## wt -0.8676594 0.7824958 0.8879799 0.6587479 -0.71244065 1.0000000
## qsec 0.4186840 -0.5912421 -0.4336979 -0.7082234 0.09120476 -0.1747159
## vs 0.6640389 -0.8108118 -0.7104159 -0.7230967 0.44027846 -0.5549157
## am 0.5998324 -0.5226070 -0.5912270 -0.2432043 0.71271113 -0.6924953
## gear 0.4802848 -0.4926866 -0.5555692 -0.1257043 0.69961013 -0.5832870
## carb -0.5509251 0.5269883 0.3949769 0.7498125 -0.09078980 0.4276059
## qsec vs am gear carb
## mpg 0.41868403 0.6640389 0.59983243 0.4802848 -0.55092507
## cyl -0.59124207 -0.8108118 -0.52260705 -0.4926866 0.52698829
## disp -0.43369788 -0.7104159 -0.59122704 -0.5555692 0.39497686
## hp -0.70822339 -0.7230967 -0.24320426 -0.1257043 0.74981247
## drat 0.09120476 0.4402785 0.71271113 0.6996101 -0.09078980
## wt -0.17471588 -0.5549157 -0.69249526 -0.5832870 0.42760594
## qsec 1.00000000 0.7445354 -0.22986086 -0.2126822 -0.65624923
## vs 0.74453544 1.0000000 0.16834512 0.2060233 -0.56960714
## am -0.22986086 0.1683451 1.00000000 0.7940588 0.05753435
## gear -0.21268223 0.2060233 0.79405876 1.0000000 0.27407284
## carb -0.65624923 -0.5696071 0.05753435 0.2740728 1.00000000
Podemos visualizar o resultado de maneira mais agradável com:
corrplot(cor(mtcars), method = "number")
A função bptest
do pacote lmplot
nos ajuda
a verificar a heterocedasticidade do modelo:
bptest(reg)
##
## studentized Breusch-Pagan test
##
## data: reg
## BP = 0.049298, df = 1, p-value = 0.8243
Nesta função, aplica-se um teste de hipótese no qual a hipótese nula corresponde aos resíduos estarem distribuídos com variância igual, e a hipótese alternativa aos resíduos estarem distribuídos com variância desigual.
Por fim, para testar a multicolinearidade entre os dados disponíveis,
podemos usar a função ols_vif_tol
do pacote
olsrr
:
reg2 <- mtcars %>% lm(mpg ~., data=.)
ols_vif_tol(reg2)
## Variables Tolerance VIF
## 1 cyl 0.06504559 15.373833
## 2 disp 0.04625295 21.620241
## 3 hp 0.10170833 9.832037
## 4 drat 0.29632966 3.374620
## 5 wt 0.06594180 15.164887
## 6 qsec 0.13283814 7.527958
## 7 vs 0.20137444 4.965873
## 8 am 0.21512374 4.648487
## 9 gear 0.18665589 5.357452
## 10 carb 0.12644228 7.908747