Model syntax in R: how to input dynamic variables?

Viewed 38

I want to run linear models (in this case, multivariate models with two response variables) within a for loop in which a new data frame called bc_applied is created at each iteration, as well as the vector targets. In my code, the column names "target1" and "target2" change at every iteration, which means I can't explicitly write variable names, instead I want to extract them from the vector targets.

Here is an example:

targets <- c("target1","target2")
bc_applied <- data.frame("dsRNA" = c(rep("gene1",5),rep("gene2",5),rep("gene3",5)),
"target1" = runif(15), "target2" = runif(15))

But when running

lm(bc_applied[,targets] ~ dsRNA, data = bc_applied)

The following error is returned:

Error in model.frame.default(formula = bc_applied[, targets] ~ dsRNA,  : 
  invalid type (list) for variable 'bc_applied[, targets]'

The desired output is given by

lm(cbind(target1, target2) ~ dsRNA, data = bc_applied)
1 Answers

According to ?lm

If response is a matrix a linear model is fitted separately by least-squares to each column of the matrix.

With cbind, it is creating a matrix. So, we need an option that takes a matrix. After subsetting the dataset with the columns, convert it to a matrix with as.matrix and it should work

lm(as.matrix(bc_applied[,targets]) ~ dsRNA, data = bc_applied)

-output

#Call:
#lm(formula = as.matrix(bc_applied[, targets]) ~ dsRNA, data = bc_applied)

#Coefficients:
#             target1   target2 
#(Intercept)   0.45161   0.47457
#dsRNAgene2    0.36341   0.29226
#dsRNAgene3   -0.07115  -0.03003

Or another option is to create a formula with paste

lm(paste0('cbind(', toString(targets),') ~ dsRNA'), data = bc_applied)

-output

#Call:
#lm(formula = paste0("cbind(", toString(targets), ") ~ dsRNA"), 
#    data = bc_applied)

#Coefficients:
#             target1   target2 
#(Intercept)   0.45161   0.47457
#dsRNAgene2    0.36341   0.29226
#dsRNAgene3   -0.07115  -0.03003

or create the formula with glue

lm(glue::glue('cbind({toString(targets)}) ~ dsRNA'), bc_applied)

or another option is

lm(do.call(cbind, asplit(bc_applied[, targets], 2)) ~ dsRNA, bc_applied)

Crosschecking with cbind

lm(cbind(target1, target2)~ dsRNA, data = bc_applied)

-output

#Call:
#lm(formula = cbind(target1, target2) ~ dsRNA, data = bc_applied)

#Coefficients:
#             target1   target2 
#(Intercept)   0.45161   0.47457
#dsRNAgene2    0.36341   0.29226
#dsRNAgene3   -0.07115  -0.03003
Related