Dose response curve fitting examples

Alexander Betz

2019-06-04

library(bendr)
library(tibble)
library(dplyr)
library(ggplot2)

You can read your data from csv,tsv or excel:

fname = "/path/to/file"
require(readr)
fdata = readr::read_csv(fname)
fdata = readr::read_tsv(fname)
require(readxl)
fdata = readxl::read_excel(fname)
names(fdata)
View(fdata)

But for a quick test, we generate some dummy data:

sdata = as_tibble(
  data.frame(
    conc=c(1,10,100,1000,5000),
    effect=c(12,27,69,86,98)
  ))

then we set the formula for the dose response model

drc.formula = effect ~ 100 / (1 + 10^((logEC50-logconc) * slope))

if you want to adapt the dose response curve to the maximum and minimum effect values, use this formula instead

loe = min(sdata$effect)
hie = max(sdata$effect)
drc.formula.dyn = substitute(
  effect ~ bottom + ((top-bottom)/(1+10^((logEC50-logconc)*slope))),
  list(top = hie, bottom=loe))

Here is an example for fitting a single dose-response curve


data.logged = sdata %>% mutate(logconc = log10(conc))

fitObject = bendr::fitdr(drc.formula,data.logged)
summary(fitObject$curve.fit)
#> 
#> Formula: effect ~ 100/(1 + 10^((logEC50 - logconc) * slope))
#> 
#> Parameters:
#>         Estimate Std. Error t value Pr(>|t|)    
#> logEC50  1.55890    0.07911  19.705 0.000286 ***
#> slope    0.64660    0.06727   9.612 0.002390 ** 
#> ---
#> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
#> 
#> Residual standard error: 3.941 on 3 degrees of freedom
#> 
#> Number of iterations to convergence: 2 
#> Achieved convergence tolerance: NA
pSingle = plot_single_drc(fitObject) + xlab("log(concentration)") + ylab("effect")
pSingle 

#ggsave("/path/to/image.png",pSingle,width="8",height="6",units="in")

Reporting goodness of fit

R^2 is considered inadequate for nonlinear fits (https://www.ncbi.nlm.nih.gov/pmc/articles/PMC2892436/) Instead, use akaike’s An Information Criterion(AIC).

Fitting multiple dose response curves simulatneously


mdata = as_tibble(
  data.frame(
    conc=c(1,10,100,1000,5000),
    effectA=c(12,27,69,86,98),
    effectB=c(15,23,64,82,91)
  ))

fo = fitdr_multi(drc.formula,mdata,2:3,conc)
pMulti = plot_multi_drc(fo,plot.layout='multi')
pMulti