DoseResponse / drc

Fitting dose-response models in R
https://doseresponse.github.io/drc/
21 stars 16 forks source link

issues with boxcox and functional programmming #13

Open peperg opened 4 years ago

peperg commented 4 years ago

Hi!

First of all, thank you so much for this package, it is a must-have tool for all of us in the tox world.

I just came across an issue when trying to use drm() and boxcox() in a functional programming framework via PURRR. Here is a repex:

library(drc)
library(tidyverse)

alba_m1 <- S.alba %>% 
  group_by(Herbicide) %>% 
  nest() %>% 
  mutate(fit = map(.x = data, ~ drm(DryMatter~Dose, data = .x, fct = LL.4())
                   ))

alba_m2 <- S.alba %>% 
  group_by(Herbicide) %>% 
  nest() %>% 
  mutate(fit = map(.x = data, ~ drm(DryMatter~Dose, data = .x, fct = LL.4())),
         bc_fit = map(.x = fit, ~ boxcox(.x)))

alba_m3 <- S.alba %>% 
  group_by(Herbicide) %>% 
  nest() %>% 
  mutate(fit = map(.x = data, ~ boxcox(drm(DryMatter~Dose, data = .x, fct = LL.4()))
  ))

alba_m4a <- boxcox(drm(DryMatter~Dose, data = alba_m1$data[[1]], fct = LL.4()))
alba_m4b <- boxcox(drm(DryMatter~Dose, data = alba_m1$data[[2]], fct = LL.4()))

When used this way (a workflow that is becoming more and more common), alba_m1 works just fine, but alba_m2 and m3 fail at the boxcox step.

when run outside the map(), with the individual data subsets (as in alba_m4a and m4b) the boxcox works perfectly.

My current hypothesis as to what could be happening is that boxcox() fails to update the drc object. What i believe is happening is that the $call is stored in reference to the internal map() call, (e.g. "drm(formula = DryMatter ~ Dose, data = .x, fct = LL.4())") making it impossible for the update.drc() used by boxcox() to actually use it, even when nested inside the same map() (which should still give boxcox access to ".x"

Obviously, that update() part of boxplot is one of the main aspects, so i don't know if there is much room to change it. My request perhaps is to offer the possibility to run a boxcox = TRUE inside of the drm() call as it seems it was the case with the old multdrc(). perhaps that way the whole process can be run at once inside of a map().

I do know that drm() offers the option to run a boxcox transformation, but the way i understand it right now, you need to specify the lambda, instead of letting it choose the optimal. please correct me if i am wrong.

Also, i do know that you could run this simple example as it is presented in the help documents by telling drm() what the groups are instead of going the group_by - nest - map route. This is just a repex, the actual case where i came across this issue is a much more complex situation that would be hard to approach that way.

Thank you in advance!