Controlling covariates in linear regression in R

Question

I have made these two models:

(model1 <- summary(lm(mpg ~ drat + wt + cyl, mtcars)))

Call:
lm(formula = mpg ~ drat + wt + cyl, data = mtcars)

Residuals:
    Min      1Q  Median      3Q     Max 
-4.2944 -1.5576 -0.4667  1.5678  6.1014 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  39.7677     6.8729   5.786 3.26e-06 ***
drat         -0.0162     1.3231  -0.012 0.990317    
wt           -3.1947     0.8293  -3.852 0.000624 ***
cyl          -1.5096     0.4464  -3.382 0.002142 ** 
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.613 on 28 degrees of freedom
Multiple R-squared:  0.8302,    Adjusted R-squared:  0.812 
F-statistic: 45.64 on 3 and 28 DF,  p-value: 6.569e-11


(model2 <- summary(lm(mpg ~  wt + cyl + drat, mtcars)))

Call:
lm(formula = mpg ~ wt + cyl + drat, data = mtcars)

Residuals:
    Min      1Q  Median      3Q     Max 
-4.2944 -1.5576 -0.4667  1.5678  6.1014 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  39.7677     6.8729   5.786 3.26e-06 ***
wt           -3.1947     0.8293  -3.852 0.000624 ***
cyl          -1.5096     0.4464  -3.382 0.002142 ** 
drat         -0.0162     1.3231  -0.012 0.990317    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 2.613 on 28 degrees of freedom
Multiple R-squared:  0.8302,    Adjusted R-squared:  0.812 
F-statistic: 45.64 on 3 and 28 DF,  p-value: 6.569e-11

My understanding is that R uses "sequential" partitioning for the variance in mpg. So in model1, drat should be unadjusted, wt should be adjusted for drat and cyl should be adjusted for drat and wt. In model2, wt should be unadjusted, cyl should be adjusted for wt and drat should be adjusted for wt and cyl.

However, the coefficients in each model appear to be exactly the same, suggesting coefficients are not being adjusted at all. Are the coefficients not being adjusted at all?

Answer 1

The question, as phrased, is slightly ambiguous. It states that "the coefficients in each model appear to be exactly the same". There are two ways that statement could be interpreted, with respect to: (1) the Estimates of the coefficients, or (2) the tests of the coefficients.

1. Regarding the Estimates of the coefficients, they are being adjusted for the other variables in the model, but you don't see any difference because you have the same variables in both model1 and model2; the order they are listed in doesn't matter. The parameter estimates will only differ if the variables are correlated and the sets of variables that are included in the models differ. Consider:

   model1 <- lm(mpg ~ drat + wt + cyl,  mtcars)
   model2 <- lm(mpg ~  wt + cyl + drat, mtcars)
   model3 <- lm(mpg ~  wt       + drat, mtcars)
   cor(mtcars$wt, mtcars$cyl)
   # [1] 0.7824958
   summary(model2)
   # Coefficients:
   #               Estimate Std. Error t value Pr(>|t|)    
   #   (Intercept)  39.7677     6.8729   5.786 3.26e-06 ***
   #   wt           -3.1947     0.8293  -3.852 0.000624 ***
   #   cyl          -1.5096     0.4464  -3.382 0.002142 ** 
   #   drat         -0.0162     1.3231  -0.012 0.990317    
   summary(model3)
   # Coefficients:
   #               Estimate Std. Error t value Pr(>|t|)    
   #   (Intercept)   30.290      7.318   4.139 0.000274 ***
   #   wt            -4.783      0.797  -6.001 1.59e-06 ***
   #   drat           1.442      1.459   0.989 0.330854    
   

   Notice that the Estimate for wt is -3.1974 in model2 and -4.783 in model3. To learn more about how the parameter estimates for variables can change depending on whether a correlated variable is included or not, it may help you to read by answer here.

2. Regarding the tests of the coefficients, it depends on which function you use to get your output. You used summary(). What is reported then are the $t$-tests associated with the parameter estimates. Those are not computed by partitioning sums of squares. However, they are equivalent to $F$-tests using type III SS. The order variables are listed in doesn't matter for $t$-tests or $F$-tests that are based on type III SS. You can also use anova() to get significance tests of your parameter estimates. That is where R uses type I SS. And because type I SS are sequential, the order the variables are listed in does matter (although again, only if the variables are correlated). Consider:

   summary(model1)
   # Coefficients:
   #               Estimate Std. Error t value Pr(>|t|)    
   #   (Intercept)  39.7677     6.8729   5.786 3.26e-06 ***
   #   drat         -0.0162     1.3231  -0.012 0.990317    
   #   wt           -3.1947     0.8293  -3.852 0.000624 ***
   #   cyl          -1.5096     0.4464  -3.382 0.002142 ** 
   summary(model2)
   # Coefficients:
   #               Estimate Std. Error t value Pr(>|t|)    
   #   (Intercept)  39.7677     6.8729   5.786 3.26e-06 ***
   #   wt           -3.1947     0.8293  -3.852 0.000624 ***
   #   cyl          -1.5096     0.4464  -3.382 0.002142 ** 
   #   drat         -0.0162     1.3231  -0.012 0.990317    
   anova(model1)
   # Analysis of Variance Table
   #             Df Sum Sq Mean Sq F value    Pr(>F)    
   #   drat       1 522.48  522.48  76.525 1.691e-09 ***
   #   wt         1 334.33  334.33  48.967 1.308e-07 ***
   #   cyl        1  78.07   78.07  11.435  0.002142 ** 
   #   Residuals 28 191.17    6.83                      
   anova(model2)
   # Analysis of Variance Table
   #             Df Sum Sq Mean Sq  F value    Pr(>F)    
   #   wt         1 847.73  847.73 124.1627 8.382e-12 ***
   #   cyl        1  87.15   87.15  12.7645  0.001304 ** 
   #   drat       1   0.00    0.00   0.0001  0.990317    
   # Residuals 28 191.17    6.83                       
   

   Notice that the $p$-value for wt is 0.000624 in both summary() outputs, but is 1.308e-07 in anova(model1) and it 8.382e-12 in anova(model2). To learn more about sums of squares in general, it may help to read my answer here. Lastly, note that you can get an ANOVA table in R that uses other types of SS, such as II and III, by using Anova() in the car package.