How to reference R's automatically generated dummy variables in a linear model

Question

R is creating dummy variables for the color variable (color is either white or red) in my linear regression model. For example, color:pH returns the interaction term colorred:pH. Some questions:

Why is R only returning colorwhite for some interaction terms and colorred with other terms? Why doesn't it do colorwhite and colorred for each term interacting with color?

Next, if I want to rebuild a model using only some of the dummy variables R generated, how can I reference them in a linear model? If I retype the model manually and include colorwhite:pH as a term, R returns an error that the object colorwhite isn't found.

Answer 1

If you use the proper multiplicative notation, the model coefficients need to account for the interpretation of the intercept term. Assume WLOG that only color and ph are in the model (this is a superfluous example you've provided). In this case, if "color==red" is the default group. There's technically only 1 dummy in the model, 1 if color is white, 0 otherwise.

Then, fitting the pH interaction, the "colorwhite" parameter is interpreted as the expected difference in the outcome comparing white to red having a pH of exactly 0. There is also a pH parameter interpreted as an expected difference in the outcome comparing groups differing by 1 unit in pH having color red. Lastly, the "colorwhite:pH" parameter is interpreted as a difference in differences for those groups, i.e. the incremental change in the pH slope comparing whites to reds.

I think you should rewrite your color formula to remove ":" and replace them with "*"

> set.seed(1)
> a <- sample(letters[1:3], 100, replace=TRUE)
> b <- sample(LETTERS[1:3], 100, replace=TRUE)
> y <- rnorm(100)
> lm(y ~ a * b)

Call:
lm(formula = y ~ a * b)

Coefficients:
(Intercept)           ab           ac           bB           bC        ab:bB  
     0.1684      -0.3894      -0.2614      -0.2807      -0.3981       0.8720  
      ac:bB        ab:bC        ac:bC  
     0.2099       0.6215       0.4547  

> lm(y ~ a : b) ## wrong

Call:
lm(formula = y ~ a:b)

Coefficients:
(Intercept)        aa:bA        ab:bA        ac:bA        aa:bB        ab:bB  
   -0.03642      0.20484     -0.18456     -0.05654     -0.07587      0.40675  
      ac:bB        aa:bC        ab:bC        ac:bC  
   -0.12738     -0.19331      0.03883           NA  

> tapply(y, interaction(a, b), mean)
         a.A          b.A          c.A          a.B          b.B          c.B 
 0.168415779 -0.220978904 -0.092958625 -0.112286696  0.370329364 -0.163796519 
         a.C          b.C          c.C 
-0.229732632  0.002406992 -0.036419652