Section 2: Linear Regression using lm

0. Load necessary packages and dataset

# Load the mtcars dataset
data(mtcars)

1. Create a linear regression model to predict mpg (miles per gallon) with just wt and hp variables in the dataset. 10pts

# Q1
data(mtcars)
m1 <- lm(mpg ~ wt + hp, data = mtcars)

summary(m1)
## 
## 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

2. Interpret both coefficients, like what we did in class. 10pts

Interpretation of the wt coefficient (−3.87783)

Holding all other variables (horsepower, hp) constant, a 1-unit increase in wt is associated with an average decrease of 3.87783 mpg.

Or more practically, an additional 1000 lbs of car weight predicts about 3.88 fewer miles per gallon, on average, for cars with the same horsepower.

Interpretation of the hp coefficient (−0.03177)

Holding weight (wt) constant, a 1-unit increase in horsepower (hp) is associated with an average decrease of 0.03177 mpg.
An additional, 10 horsepower corresponds to about 0.32 fewer mpg on average, for cars with the same weight.

3. What assumptions are being made when we use linear regression? Are they met in this dataset? Use diagnostic plots and describe what you observe from the plots. 10pts

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

par(mfrow = c(1, 1))

Linear Regression Assumptions: (1) a linear relationship between predictors and the response, (2) independent errors, (3) constant variance of errors (homoscedasticity), (4) approximately normal errors

Diagnostic plots (plot(m1))

  • Residuals vs Fitted: The red smooth line shows a slight curved (U-shaped) pattern rather than staying flat around 0, suggesting moderate nonlinearity. The variances are similar across different fitted values.

  • Normal Q–Q: Points follow the line in the middle, but deviate in the tails—especially the upper tail—indicating residuals are not perfectly normal

  • Scale–Location: Compared to Residuals vs Fitted, the nonlinearity problem is less obvious.

  • Residuals vs Leverage (Cook’s distance): All points have low cook’s distance, indicating no influential points.

Conclusion: The assumptions are approximately reasonable, but there is evidence of moderate nonlinearity, mild heteroskedasticity, mild non-normality, and no obvious influential outliers.

4. Evaluate the model by reporting the MSE (Mean Square Error) 10pts

m1_mse <- mean((m1$fitted.values - mtcars$mpg)^2)

print(paste("Mean Squared Error (MSE) for the model:", round(m1_mse, 4)))
## [1] "Mean Squared Error (MSE) for the model: 6.0952"

The MSE for mpg ~ wt + hp is 6.0952. A smaller MSE indicates better in-sample fit (lower average squared error).

5. Evaluate the model by reporting the R^2 10pts

r2 <- summary(m1)$r.squared
r2
## [1] 0.8267855

The model’s R^2 is 0.8268.
This means the model using wt and hp explains about 82.68% of the variation in mpg in the mtcars dataset.

6. Try adding interaction term between wt and hp to your linear regression model fitted in question 1.

I fit a new model with an interaction term: mpg ~ wt * hp, which includes wt, hp, and wt:hp.

m2 <- lm(mpg ~ wt * hp, data = mtcars)   # wt * hp = wt + hp + wt:hp
summary_m2 <- summary(m2)
print(summary_m2)
## 
## 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

• Report the significance of the interaction term at 0.05 level. 10pts

  • Significance at 0.05: The interaction term (wt:hp) is statistically significant if its p-value < 0.05.
    From the model output, the p-value for wt:hp is 0.000811 so the interaction term is significant at the 0.05 level.

• Check how this interaction term influences the model’s performance in terms of R^2. 10pts

  • Effect on model performance (R^2):
    The original model (Q1) has R^2 = 0.8268.
    The interaction model has R^2 = 0.8848.
    So R^2 increases when adding the interaction.

• How do you interpret your new model? 10pts

  • Interpretation of the new model:
    With an interaction term, the effect of wt on mpg depends on hp, and the effect of hp on mpg depends on wt.
    Specifically, the predicted mpg is: \[ \widehat{mpg} = \beta_0 + \beta_1 wt + \beta_2 hp + \beta_3 (wt \cdot hp). \] The marginal effect of wt is \(\beta_1 + \beta_3 hp\), and the marginal effect of hp is \(\beta_2 + \beta_3 wt\).

7. And let’s assume hp has outliers, apply 5%/95% winsorization technique to fix it. Then fit a new linear regression model on the winsorized hp and original wt. No interaction term is needed. Compare the performance of the newly fitted model with the model in question 1. What differences do you observe in R^2 and the coefficients? 10pts

# 5%/95% winsorization for hp
mt2 <- mtcars

lo_hp <- quantile(mt2$hp, 0.05, na.rm = TRUE)
hi_hp <- quantile(mt2$hp, 0.95, na.rm = TRUE)

mt2$hp_win <- mt2$hp
mt2$hp_win[mt2$hp_win < lo_hp] <- lo_hp
mt2$hp_win[mt2$hp_win > hi_hp] <- hi_hp

# Fit new model (no interaction)
m3 <- lm(mpg ~ wt + hp_win, data = mt2)

# Summaries (like class)
summary(m1)
## 
## 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
summary(m3)
## 
## Call:
## lm(formula = mpg ~ wt + hp_win, data = mt2)
## 
## 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_win      -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

What differences do you observe in R^2 and the coefficients?

Original model (Q1): mpg ~ wt + hp
- R^2 = 0.8268 (Adj. R^2 = 0.8148)
- Coefficients: Intercept = 37.2273, wt = −3.8778, hp = −0.03177 (p = 0.00145)

Winsorized model (Q7): mpg ~ wt + hp_win
- R^2 = 0.8330 (Adj. R^2 = 0.8215)
- Coefficients: Intercept = 37.3172, wt = −3.5828, hp_win = −0.03952 (p = 0.000824)

8. Check multicollinearity using vif() on the model fitted in question 1. What do you find? 10pts

library(car)
## Loading required package: carData
vif_vals <- car::vif(m1)
vif_vals
##       wt       hp 
## 1.766625 1.766625

I checked multicollinearity on the model from question 1 (mpg ~ wt + hp) using vif().

Since these VIFs are well below common cutoffs (5 or 10), there is no serious multicollinearity between wt and hp in this model.

Q9 (no credit). Does an improved R^2 really improve model predictability?

Not necessarily. A higher in-sample R^2 only means the model fits the current data better; it can increase simply by adding more variables (even if they do not generalize), which is why adjusted R^2 is useful because it penalizes unnecessary predictors.