Построение кривых доза-ответ с помощью ggplot2 и drc

Включена недавняя статья авторов пакета drc инструкции по извлечению параметров для использования ggplot2. Они не работают в ggplot2, но извлекают данные из модели. Это их решение применительно к вашим данным.

demo1 <- reshape2::melt(demo,id.vars = "X") # get numbers ready for use.
demo.LL.4 <- drm(data = demo1,value~X,fct=LL.4(),na.action = na.omit) # run model.

Функция predict может извлекать параметры из drm моделей. Он несовместим с несколькими кривыми, которые были подогнаны с помощью curveid.

# predictions and confidence intervals.
demo.fits <- expand.grid(conc=exp(seq(log(1.00e-04), log(1.00e-09), length=100))) 
# new data with predictions
pm <- predict(demo.LL.4, newdata=demo.fits, interval="confidence") 
    demo.fits$p <- pm[,1]
    demo.fits$pmin <- pm[,2]
    demo.fits$pmax <- pm[,3]

Они советуют сместить нулевую концентрацию, чтобы избежать проблем с corre_trans.

demo1$XX <- demo1$X
demo1$XX[demo1$XX == 0] <- 1.00e-09

Затем идет построение кривой, если пропустить geom_ribbon, ошибки не будут отображаться.

ggplot(demo1, aes(x = XX, y = value)) +
  geom_point() +
  geom_ribbon(data=demo.fits, aes(x=conc, y=p, ymin=pmin, ymax=pmax), alpha=0.2) +
  geom_line(data=demo.fits, aes(x=conc, y=p)) +
  coord_trans(x="log") 

Чтобы построить график нескольких кривых вместе, процесс можно повторить. Добавьте идентификаторы к каждому набору.

demo.fits_1 <- data.frame(label = "curve1", demo.fits)

Затем используйте rbind, чтобы объединить все извлеченные параметры. Оттуда ggplot может обрабатывать цвета.

Я собираюсь ответить на свой вопрос и, надеюсь, это поможет другим, столкнувшимся с той же проблемой.

Конечно, можно построить кривые доза-реакция с помощью ggplot2 и пакета drc с простым добавлением либо geom_, либо stat_smooth (method=drm, fct=LL.4(),se=FALSE) при построении графика в линейном масштабе, либо geom_ или stat_smooth (method=drm, fct=L.4(),se=FALSE), если добавлено scale_x_log10().

Чтобы использовать шкалу log10, я преобразовал свои данные в:

demo <- demo %>% 
      mutate(X = 
       ifelse(X == 0, 
              yes = (sort(demo$X[which.min(sort(demo$X)) + 1]/100)),
              no = X
              )
            )         #looks for the pre-lowest value in X and divides it by 100

В этом случае я заменил значение X = 0 на X = 1/100 от предпоследнего значения X (в данном случае 1e-10). Однако вы можете легко отбросить значение 0, которое испортит логарифмическое построение, полностью исключив его из набора данных, как это делает Prism. Одна вещь, которую следует отметить, как я выяснил, заключается в том, что ggplot сначала масштабирует оси, а затем добавляет данные, поэтому код ломается при попытке log10 (0).

Другая тонкость заключается в том, что функция stat_smooth прекрасно справляется с моделями DRM с использованием method = drm, но она не знает, как соответствовать доверительным интервалам «SE». Таким образом, выбор se = FALSE позволяет строить график и, по моему скромному мнению, в любом случае делает график менее беспорядочным - просто добавьте полосы ошибок.

И, наконец, изменение fct = LL.4() на fct = L.4() позволяет строить график в масштабе log10, потому что снова сначала выбирается масштаб, а затем выполняется подгонка. Таким образом, даже несмотря на то, что значения оси не логарифмические, ggplot фактически преобразовал набор данных в log10, поэтому функция подгонки теперь должна быть просто logit-4P (т.е. L.4 ()) вместо log-logit-4P (LL .4 ()).

Функции geom_smooth () и stat_smooth () естественным образом принимают тот же цвет, что и набор данных, устраняя необходимость настраивать цвет подобранной функции, чтобы он соответствовал цвету точек данных.

В итоге:

demo <- demo %>% 
      mutate(X = 
       ifelse(X == 0, 
              yes = (sort(demo$X[which.min(sort(demo$X)) + 1]/100)),
              no = X
              )
            )
demo.long <- reshape2::melt(demo,id.vars = "X") #reshapes the demo dataset to long format
ggplot(data = demo.long,
       aes(x = X, y = value, col = variable)
      ) + 
   geom_point() + 
   geom_smooth(method = drm, fct = L.4(), se = FALSE) +
   scale_x_log10() #plots out the dataset with the corresponding 4-parameter log-logit dose response curves