Tidyverse way to get lm()-residuals using a subset of columns as predictors but keeping all in output

Viewed 231

I have some data like this:

dat = tibble(
      var1 = rep(c('A','B'),each=5)
    , var2 = rnorm(10)
    , var3 = rnorm(10)
    , var4 = rnorm(10)
    , var5 = rnorm(10)
)

I can get what I want by explicitly naming the columns to be used in the lm formula:

dat %>%
    #dat has columns: var1 through var5
    dplyr::group_by(var1) %>%
    dplyr::mutate(
        resids = resid(lm( var2 ~ var3 + var4 ))
    )

But I actually have many columns in my real data set, and the number and names of the ones I'll be using will vary. I do know the names of the ones I don't want, so I thought this would work:

dat %>%
    #dat has columns: var1 through var5
    dplyr::group_by(var1) %>%
    dplyr::mutate(
        resids = resid(lm( 
            formula = var2 ~ .
            , data = . %>% select(-var1, -var5)
     ))
    )

But that doesn't seem to work. Any suggestions?

3 Answers

Since you know the variables that aren't going to be on the right-hand side of the model, one option is to get those names into a vector and then build the formula for lm(). This means doing an extra step outside your pipe chain.

Building the formula is a relatively common approach for, e.g., making functions for fitting models. I wrote a blog post to show this approach here.

In your case, you can pull the names of the variables you want in the model as a character vector based on the variables you don't want.

modvars = dat %>%
     select(-var1, -var5, -var2) %>%
     names()
modvars
[1] "var3" "var4"

You can build the formula for model fitting using as.formula() after pasting the variables together. This is what that would look like:

as.formula(paste("var2 ~", paste(modvars, collapse = "+") ) )
var2 ~ var3 + var4

Even easier is the reformulate() approach from the comments thanks to @BenBolker.

reformulate(modvars, response = "var2")
var2 ~ var3 + var4

You could build this outside your pipe chain or put it directly within the chain. Here I do the latter.

dat %>%
     dplyr::group_by(var1) %>%
     dplyr::mutate(
          resids = resid(lm( 
               formula = reformulate(modvars, response = "var2") )
          )
     )
# A tibble: 10 x 6
# Groups:   var1 [2]
   var1     var2    var3    var4    var5 resids
   <chr>   <dbl>   <dbl>   <dbl>   <dbl>  <dbl>
 1 A      0.0792  0.265   0.637  -0.106   0.386
 2 A     -0.845   0.386   1.20    1.55   -0.232
 3 A      0.465   1.12   -0.750   0.726  -0.141
 4 A     -0.365  -1.19    0.174   0.347  -0.126
 5 A      0.395  -0.0515 -0.464  -0.0934  0.112
 6 B     -2.83   -0.0664 -0.0958  0.588  -1.99 
 7 B      0.383   1.16   -0.339   0.492   0.838
 8 B      1.35    0.270   2.40    0.626  -0.512
 9 B      0.620  -1.33    1.32   -0.148   0.688
10 B      0.664  -0.0487  0.426  -0.158   0.973

The residuals match your original approach where you wrote the formula out:

dat %>%
     #dat has columns: var1 through var5
     dplyr::group_by(var1) %>%
     dplyr::mutate(
          resids = resid(lm( var2 ~ var3 + var4 ))
     )

# A tibble: 10 x 6
# Groups:   var1 [2]
   var1     var2    var3    var4    var5 resids
   <chr>   <dbl>   <dbl>   <dbl>   <dbl>  <dbl>
 1 A      0.0792  0.265   0.637  -0.106   0.386
 2 A     -0.845   0.386   1.20    1.55   -0.232
 3 A      0.465   1.12   -0.750   0.726  -0.141
 4 A     -0.365  -1.19    0.174   0.347  -0.126
 5 A      0.395  -0.0515 -0.464  -0.0934  0.112
 6 B     -2.83   -0.0664 -0.0958  0.588  -1.99 
 7 B      0.383   1.16   -0.339   0.492   0.838
 8 B      1.35    0.270   2.40    0.626  -0.512
 9 B      0.620  -1.33    1.32   -0.148   0.688
10 B      0.664  -0.0487  0.426  -0.158   0.973

The issue here is not the subsettng but rather the confusion of the . in the formula. Is this referring to the variables or to the data.

There are other issues also with regard to this, ie the select. You should dodata = {.}%>%subset(-var1,var5) and not what you have, or simply do data = subset(.,-var1,-var5). HOw to solve this issue:

use nest + unnest

By nesting_by, the grouping variable is automatically removed from the data:

dat %>%
  nest_by(var1) %>%
  mutate(resid = list(resid(lm(formula = var2 ~ .-var5, data = data))))%>%
  unnest(c(data, resid))

use group_by + summarize

dat %>%
     group_by(var1) %>%
     summarise(resid = list(resid(lm(formula = var2~.-var1-var5, data =.))),.groups="drop")%>%
     unnest(resid)

Your code doesn't work because the special symbol "." in data = . %>% select(-var1, -var5) is actually the input data which has not been grouped. After dplyr 1.0.0 this issue has been solved by cur_data(), which gives the current data for the current group (exclusing grouping variables).

dat %>%
  group_by(var1) %>%
  mutate(
    resids = resid(lm(var2 ~ ., cur_data() %>% select(-var5)))
  )

Note that I use select(-var5) instead of select(ivar1, -var5) because cur_data() has excluded grouping variables, i.e. var1, so select(ivar1, -var5) will get an error:

Can't subset columns that don't exist.


As @IRTFM comments, you can also select variables before passing to lm. Remember to use cur_data() rather than ".".

dat %>%
  group_by(var1) %>%
  select(-var5) %>%
  mutate(
    resids = resid(lm(var2 ~ ., cur_data()))
  )

Here you still don't need to exclude var1 in select() because grouping variables cannot be excluded.

Related