Using purrr and modelr to make predictions from dose-response models (drm)

458 Views Asked by At

I have a data set that looks, in abbreviate form, like this:

library(tidyverse)    
dat_s<-tibble(
      type=c(rep("A", 9), rep("B", 8), rep("C", 10)),
      ref=c("ref3", "ref3", "ref1",  "ref2", "ref2", "ref1", "ref2", "ref2", "ref2", "ref2", 
            "ref1", "ref2", "ref2", "ref3", "ref2", "ref3", "ref1", "ref3", 
            "ref2", "ref3", "ref1", "ref1", "ref3", "ref1", "ref1", "ref2", "ref2"),
      info=as.character(sample(100, 27)),
      liv=c(3.0e-05, 2.9e-07, 2.2e-07, 2.7e-07, 2.6e-06, 4.8e-07, 1.4e-05, 2.6e-06, 7.7e-06, 2.2e-06, 
            1.5e-07, 1.6e-07, 1.8e-06, 6.1e-08, 4.9e-06, 4.9e-06, 1.8e-06, 1.5e-07,
            4.3e-08, 1.8e-06, 1.0e-07, 1.6e-07, 9.7e-07, 1.0e-06, 6.4e-07, 1.2e-07, 5.7e-06),
      prod=c(0.00, 2, 3, 4.80, 2.10, 5.10, 0.00, 0.13, 2.00, 0.13, 0.00, 4.10, 4.60, 2.10, 0.26, 0.00, 
             4.60, 0.00, 4.60, 2.10, 4.80, 0.00, 0.00, 1.80, 3.60, 4.10, 0.00)
    )%>%
  mutate(livp1=liv+1)

I want to calculate dose response relationships for each combination of type and ref, make predictions to plot a curve, and calculate residuals. The info column is to reflect that I have additional columns in this data frame which I need to preserve, but are not important in the dose-response analysis.

I start by creating the models using a function and a nested data frame:

dr_s<-function(df){drc::drm(data=df, prod~log(livp1), fct=drc::LL.3())}

dat_mods<-
  dat_s%>%
  group_by(type, ref)%>%
  nest() %>%
  mutate(dr_mod=map(data, dr_s))

Which works to create the models and put them in the data frame. To use add_predictions with models of the type drm, the input has to be a data.frame (rather than a tibble). When I try to add predictions for each livp1 variable (according to the comments below):

dat_mods%>%

mutate(mod_preds=map2(data, dr_mod, ~add_predictions(data=as.data.frame(.x), model=.y))))

I get a non-numeric argument to binary operator error message. This code works fine when the info column is not character class. However, I need to retain this information with the predicted data, and would like to avoid pulling it from the data frame if possible.

Any guidance is appreciated!

1

There are 1 best solutions below

6
On BEST ANSWER

This one was pretty silly, the predict method for models of class drm doesn't work on objects of class tibble. So you have to convert .x to data.frame.

dat_s%>%
  group_by(type, ref)%>%
  nest() %>%
  mutate(dr_mod=map(data, dr_s),
         mod_preds=map2(data, dr_mod, 
                        ~modelr::add_predictions(data = as.data.frame(.x[,"livp1"]),
                                                 model = .y)))
## A tibble: 9 x 5
## Groups:   type, ref [9]
#  type  ref   data             dr_mod mod_preds       
#  <chr> <chr> <list>           <list> <list>          
#1 A     ref3  <tibble [2 × 4]> <drc>  <df[,2] [2 × 2]>
#2 A     ref1  <tibble [2 × 4]> <drc>  <df[,2] [2 × 2]>
#3 A     ref2  <tibble [5 × 4]> <drc>  <df[,2] [5 × 2]>
#4 B     ref2  <tibble [4 × 4]> <drc>  <df[,2] [4 × 2]>
#5 B     ref1  <tibble [2 × 4]> <drc>  <df[,2] [2 × 2]>
#6 B     ref3  <tibble [2 × 4]> <drc>  <df[,2] [2 × 2]>
#7 C     ref3  <tibble [3 × 4]> <drc>  <df[,2] [3 × 2]>
#8 C     ref2  <tibble [3 × 4]> <drc>  <df[,2] [3 × 2]>
#9 C     ref1  <tibble [4 × 4]> <drc>  <df[,2] [4 × 2]>