Data Exploration

library(ggplot2)
## Warning: package 'ggplot2' was built under R version 4.6.1
library(corrplot)
## Warning: package 'corrplot' was built under R version 4.6.1
## corrplot 0.95 loaded
library(car)
## Warning: package 'car' was built under R version 4.6.1
## Loading required package: carData
## Warning: package 'carData' was built under R version 4.6.1
data(mtcars)

Linear Regression

model <- lm(mpg ~ wt + hp, data = mtcars)
summary(model)
## 
## Call:
## lm(formula = mpg ~ wt + hp, data = mtcars)
## 
## Residuals:
##    Min     1Q Median     3Q    Max 
## -3.941 -1.600 -0.182  1.050  5.854 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 37.22727    1.59879  23.285  < 2e-16 ***
## wt          -3.87783    0.63273  -6.129 1.12e-06 ***
## hp          -0.03177    0.00903  -3.519  0.00145 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.593 on 29 degrees of freedom
## Multiple R-squared:  0.8268, Adjusted R-squared:  0.8148 
## F-statistic: 69.21 on 2 and 29 DF,  p-value: 9.109e-12

If everything else is constant for every increase in one unit of wt mpg ‘increases’ by -3.88. If everything else is constant for every increase in one unit of hp mpg ‘increases’ by -0.032.

par(mfrow = c(2,2))
plot(model)

The Linearity assumption is not met. There is a clear non linear relationship on the plot. The Normality assumption is mostly met. There is a skewness at the top. The Equal Variance is also mostly met. The trend line is quite flat and only rises at the high fitted values.

MSE

mse1 <- mean(residuals(model)^2)
mse1
## [1] 6.095242

MSE is reported to be 6.10. if we take the square root of we get roughly 2.47mpg

R²

summary(model)$r.squared
## [1] 0.8267855

R² is 0.8268. This means that weight and horsepower explain almost 83% of varioant in mpg.

Model with Interaction Term

int_model <- model2 <- lm(mpg ~ wt * hp, data = mtcars)
summary(int_model)
## 
## Call:
## lm(formula = mpg ~ wt * hp, data = mtcars)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.0632 -1.6491 -0.7362  1.4211  4.5513 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 49.80842    3.60516  13.816 5.01e-14 ***
## wt          -8.21662    1.26971  -6.471 5.20e-07 ***
## hp          -0.12010    0.02470  -4.863 4.04e-05 ***
## wt:hp        0.02785    0.00742   3.753 0.000811 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.153 on 28 degrees of freedom
## Multiple R-squared:  0.8848, Adjusted R-squared:  0.8724 
## F-statistic: 71.66 on 3 and 28 DF,  p-value: 2.981e-13

The interaction term wt:hp has a p-value = 0.0008. That is below the alpha of 0.05 which means it is significant.

Compare R² for both Models

summary(model)$r.squared
## [1] 0.8267855
summary(int_model)$r.squared
## [1] 0.8847637

R² as you can see from the above outputs the R² increase with the interaction variable. The equation for the new model is: mpg = 49.81 − 8.22·wt − 0.120·hp + 0.0278·(wt × hp) The effect of weight and horsepower now depends on the other. If they are both negative the interaction makes it so the effect is less on mpg.

Winsorizing

limits <- quantile(mtcars$hp, c(0.05, 0.95))
limits
##     5%    95% 
##  63.65 253.55
mtcars$hp_w <- mtcars$hp
mtcars$hp_w[mtcars$hp_w < limits[1]] <- limits[1]
mtcars$hp_w[mtcars$hp_w > limits[2]] <- limits[2]

win_model <- lm(mpg ~ wt + hp_w, data = mtcars)
summary(win_model)
## 
## Call:
## lm(formula = mpg ~ wt + hp_w, data = mtcars)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.8825 -1.6545 -0.0968  0.8367  5.7259 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 37.31722    1.56964  23.774  < 2e-16 ***
## wt          -3.58279    0.66427  -5.394  8.5e-06 ***
## hp_w        -0.03952    0.01059  -3.732 0.000824 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.546 on 29 degrees of freedom
## Multiple R-squared:  0.833,  Adjusted R-squared:  0.8215 
## F-statistic: 72.34 on 2 and 29 DF,  p-value: 5.348e-12
cbind(coef(model), coef(win_model))
##                    [,1]        [,2]
## (Intercept) 37.22727012 37.31722375
## wt          -3.87783074 -3.58278620
## hp          -0.03177295 -0.03951904
summary(model)$r.squared
## [1] 0.8267855
summary(win_model)$r.squared
## [1] 0.8330309

R² increases slightly, from 0.8268 to 0.8330, so the fit is a little better but nearly the same. The coefficients for hp become slightly more negative. The coefficients for qt get a little larger.

Check for Multicollinearity

vif(model)
##       wt       hp 
## 1.766625 1.766625

Both of the the VIFs are around 1.77. This means that multicollinearity is not an issue for this model.