Is there a way to obtain residual plots for all interaction terms?

Question

I am working on an exercise asking me "Plot the residuals against Y_hat, each predictor variable, and each two-factor interaction term on separate graphs." Here is a snippet of the data set I am using:

> dput(head(Commercial_Properties, 10))
structure(list(Rental_Rates = c(13.5, 12, 10.5, 15, 14, 10.5, 
14, 16.5, 17.5, 16.5), Age = c(1, 14, 16, 4, 11, 15, 2, 1, 1, 
8), Op_Expense_Tax = c(5.02, 8.19, 3, 10.7, 8.97, 9.45, 8, 6.62, 
6.2, 11.78), Vacancy_Rate = c(0.14, 0.27, 0, 0.05, 0.07, 0.24, 
0.19, 0.6, 0, 0.03), Total_Sq_Ft = c(123000, 104079, 39998, 57112, 
60000, 101385, 31300, 248172, 215000, 251015), residuals = c(`1` = -1.03567244005944, 
`2` = -1.51380641405037, `3` = -0.591053402133659, `4` = -0.133568082335235, 
`5` = 0.313283765150399, `6` = -3.18718522392237, `7` = -0.538356748944345, 
`8` = 0.236302385996349, `9` = 1.98922037248654, `10` = 0.105829602747806
)), row.names = c(NA, -10L), class = c("tbl_df", "tbl", "data.frame"
))

From here I created the proper linear model that includes two factor interaction terms:

commercial_properties_lm_two_degree_interaction <- 
  lm(data=Commercial_Properties, 
     formula=Rental_Rates ~ (Age + Op_Expense_Tax + Vacancy_Rate + Total_Sq_Ft)^2)

Next what I was hoping to accomplish was to plot the residuals not just of the linear terms, but also of the interaction terms. I attempted to do this using the residualPlots() function in the car package

library(car)
residualPlots(model=commercial_properties_lm_two_degree_interaction, 
              terms=~ (Age + Op_Expense_Tax + Vacancy_Rate + Total_Sq_Ft)^2)

When applied in this way the output only produced the residual plots against the linear terms, it didn't plot any interactions. So I then attempted to do it manually, but I got an error:

residualPlots(model=commercial_properties_lm_two_degree_interaction,
              terms=~ Age + Op_Expense_Tax + Vacancy_Rate + Tota_Sq_Ft + 
                Age:Op_Expense_Tax + Age:Vacancy_Rate)

Error in termsToMf(model, terms) : argument 'terms' not interpretable.

Now if I were to do things completely manually I was able to get an interaction plot for example:

with(data=Commercial_Properties, plot(x=Op_Expense_Tax * Vacancy_Rate, y=residuals))

plotted successfully. My issue is that sure I can do this completely manually for a reasonably small amount of variables, but it will get extremely tedious once the amount of variables begins to get larger.

So my question is if there is a way to use an already created function in R to make residual plots of the interaction terms or would I be left to doing it completely manually or most likely having to write some sort of loop ?

Note, I'm not asking about partial residuals. I haven't gotten to that point in my text I'm using. Just plain interaction terms against residuals.

Answer 1

You could do an eval(parse()) approach using the 'term.labels' attribute.

With gsub(':', '*', a[grep(':', a)]) pull out the interaction terms and replace : with * so it can be evaluated.

a <- attr(terms(commercial_properties_lm_two_degree_interaction), 'term.labels')

op <- par(mfrow=c(2, 3))
with(Commercial_Properties, 
     lapply(gsub(':', '*', a[grep(':', a)]), function(x)
            plot(eval(parse(text=x)), residuals, xlab=x)))
par(op)

Edit

This is how we would do this with a for loop in R (but see comments below):

as <- gsub(':', '*', a[grep(':', a)])
op <- par(mfrow=c(2, 3))
for (x in as) {
  with(Commercial_Properties, 
       plot(eval(parse(text=x)), residuals, xlab=x)
  )
}
par(op)