Parameters and AUC and IC50 of a dose response curve

Question

I have these dose response data:

df <- data.frame(viability=c(14,81,58,78,71,83,64,16,32,100,100,81,86,83,100,90,15,100,38,100,91,84,92,100),
                 dose=c(10,0.62,2.5,0.16,0.039,0.0024,0.0098,0.00061,10,0.62,2.5,0.16,0.039,0.0024,0.0098,0.00061,10,0.62,2.5,0.16,0.039,0.0024,0.0098,0.00061),
                 stringsAsFactors=F)

I then use the drc package's drm function to fit a log-logistic curve to these data:

library(drc)
fit <- drm(viability~dose,data=df,fct=LL.4(names=c("slope","low","high","ED50")),type="continuous")

> summary(fit)

Model fitted: Log-logistic (ED50 as parameter) (4 parms)

Parameter estimates:

                  Estimate Std. Error  t-value p-value
slope:(Intercept)  5.15328   18.07742  0.28507  0.7785
low:(Intercept)   20.19430   12.61122  1.60130  0.1250
high:(Intercept)  83.33181    4.96736 16.77586  0.0000
ED50:(Intercept)   2.98733    1.99685  1.49602  0.1503

Residual standard error:

 21.0743 (20 degrees of freedom)

I then generate predictions so I'll be able to plot the curve:

pred.df <- expand.grid(dose=exp(seq(log(max(df$dose)),log(min(df$dose)),length=100))) 
pred <- predict(fit,newdata=pred.df,interval="confidence") 
pred.df$viability <- pmax(pred[,1],0)
pred.df$viability <- pmin(pred.df$viability,100)
pred.df$viability.low <- pmax(pred[,2],0)
pred.df$viability.low <- pmin(pred.df$viability.low,100)
pred.df$viability.high <- pmax(pred[,3],0)
pred.df$viability.high <- pmin(pred.df$viability.high,100)

I also use the PharmacoGx Bioconductor package to compute AUC and IC50 for both the curve and its high and low bounds:

library(PharmacoGx)
auc.mid <- computeAUC(rev(pred.df$dose),rev(pred.df$viability))/((max(pred.df$viability)-min(pred.df$viability))*(max(pred.df$dose)-min(pred.df$dose)))
auc.low <- computeAUC(rev(pred.df$dose),rev(pred.df$viability.low))/((max(pred.df$viability.low)-min(pred.df$viability.low))*(max(pred.df)-min(pred.df$dose)))
auc.high <- computeAUC(rev(pred.df$dose),rev(pred.df$viability.high))/((max(pred.df$viability.high)-min(pred.df$viability.high))*(max(pred.df$dose)-min(pred.df$dose)))
ic50.mid <- computeIC50(rev(pred.df$dose),rev(pred.df$viability))
ic50.low <- computeIC50(rev(pred.df$dose),rev(pred.df$viability.low))
ic50.high <- computeIC50(rev(pred.df$dose),rev(pred.df$viability.high))

Ceating a table with all the parameters so I can plot everything together:

ann.df <- data.frame(param=c("slope","low","high","ED50","auc.mid","auc.high","auc.low","ic50.mid","ic50.high","ic50.low"),value=signif(c(summary(fit)$coefficient[,1],auc.mid,auc.high,auc.low,ic50.mid,ic50.high,ic50.low),2),stringsAsFactors=F)

And finally plotting it all:

library(ggplot2)
library(grid)
library(gridExtra)
pl <- ggplot(df,aes(x=dose,y=viability))+geom_point()+geom_ribbon(data=pred.df,aes(x=dose,y=viability,ymin=viability.low,ymax=viability.high),alpha=0.2)+labs(y="viability")+
  geom_line(data=pred.df,aes(x=dose,y=viability))+coord_trans(x="log")+theme_bw()+scale_x_continuous(name="dose",breaks=sort(unique(df$dose)),labels=format(signif(sort(unique(df$dose)),3),scientific=T))
ggdraw(pl)+draw_grob(tableGrob(ann.df,rows=NULL),x=0.1,y=0.175,width=0.3,height=0.4)

Which gives:

My questions are:

I thought that slope should be negative. How come it's 5.2?
the auc.mid, auc.high, and auc.lowcumputed as:

auc.mid <- computeAUC(rev(pred.df$dose),rev(pred.df$viability)) auc.low <- computeAUC(rev(pred.df$dose),rev(pred.df$viability.low)) auc.high <- computeAUC(rev(pred.df$dose),rev(pred.df$viability.high))

give 21.47818, 37.52389, and 2.678228, respectively.

Since these are not in the [0,1] range I thought that divinding them by the area under the highest corresponding viability will give what I'm looking for, i.e., relative AUC, but these values seem too low relative to what the figure shows. What are these AUCs then?

Also, how come auc.mid > auc.low > auc.high? I would think that it should be auc.high > auc.mid > auc.low

The IC50 values also seem a little low. Do they make sense?

Bonus question: how do I avoid the trailing zeros in slope, low, high, ED50, ic50.mid, and ic50.high in the figure?

As the slope you take summary(fit)$coefficient[,1]. In most models that would be the intercept. — Axeman

Petr Smirnov Petr Smirnov · Accepted Answer · 2017-04-19T20:41:32

The parameter you are pulling out is the hill slope parameter, or the coefficient in front of the concentration variable in the exponential, not the actual slope of the curve.
The AUC provided is in the [0-100] range, for the area above the curve. I ran the code and got the order as auc.low>auc.mid>auc.high. Traditionally the area under the response curve was reported, or 1-viability.
It is important to note that the PharmacoGx package uses a 3 parameter hill slope model, similar to LL.3 in drc. Therefore, the plot will not correspond to the function fit by PharmacoGx to calculate the IC50 or AUC.

Source: PharmacoGx dev.

Parameters and AUC and IC50 of a dose response curve

1 Answers