title: Parameters and AUC and IC50 of a dose response curve
body: <p>I have these dose response data:</p>
  
  <pre><code>df &lt;- 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)
  </code></pre>
  
  <p>I then use the <code>drc</code> package's <code>drm</code> function to fit a log-logistic curve to these data:</p>
  
  <pre><code>library(drc)
  fit &lt;- drm(viability~dose,data=df,fct=LL.4(names=c("slope","low","high","ED50")),type="continuous")
  
  &gt; 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)
  </code></pre>
  
  <p>I then generate predictions so I'll be able to plot the curve:</p>
  
  <pre><code>pred.df &lt;- expand.grid(dose=exp(seq(log(max(df$dose)),log(min(df$dose)),length=100))) 
  pred &lt;- predict(fit,newdata=pred.df,interval="confidence") 
  pred.df$viability &lt;- pmax(pred[,1],0)
  pred.df$viability &lt;- pmin(pred.df$viability,100)
  pred.df$viability.low &lt;- pmax(pred[,2],0)
  pred.df$viability.low &lt;- pmin(pred.df$viability.low,100)
  pred.df$viability.high &lt;- pmax(pred[,3],0)
  pred.df$viability.high &lt;- pmin(pred.df$viability.high,100)
  </code></pre>
  
  <p>I also use the <code>PharmacoGx</code> <code>Bioconductor</code> package to compute AUC and IC50 for both the curve and its high and low bounds:</p>
  
  <pre><code>library(PharmacoGx)
  auc.mid &lt;- 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 &lt;- 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 &lt;- 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 &lt;- computeIC50(rev(pred.df$dose),rev(pred.df$viability))
  ic50.low &lt;- computeIC50(rev(pred.df$dose),rev(pred.df$viability.low))
  ic50.high &lt;- computeIC50(rev(pred.df$dose),rev(pred.df$viability.high))
  </code></pre>
  
  <p>Ceating a table with all the parameters so I can plot everything together:</p>
  
  <pre><code>ann.df &lt;- 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)
  </code></pre>
  
  <p>And finally plotting it all:</p>
  
  <pre><code>library(ggplot2)
  library(grid)
  library(gridExtra)
  pl &lt;- 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)
  </code></pre>
  
  <p>Which gives:
  <a href="https://i.stack.imgur.com/Am38t.png" rel="nofollow noreferrer"><img src="https://i.stack.imgur.com/Am38t.png" alt="enter image description here"></a></p>
  
  <p>My questions are:</p>
  
  <ol>
  <li><p>I thought that slope should be negative. How come it's 5.2?</p></li>
  <li><p>the <code>auc.mid</code>, <code>auc.high</code>, and <code>auc.low</code>cumputed as:</p>
  
  <p>auc.mid &lt;- computeAUC(rev(pred.df$dose),rev(pred.df$viability))
  auc.low &lt;- computeAUC(rev(pred.df$dose),rev(pred.df$viability.low))
  auc.high &lt;- computeAUC(rev(pred.df$dose),rev(pred.df$viability.high))</p></li>
  </ol>
  
  <p>give 21.47818, 37.52389, and 2.678228, respectively.</p>
  
  <p>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?</p>
  
  <p>Also, how come <code>auc.mid</code> > <code>auc.low</code> > <code>auc.high</code>? I would think that it should be <code>auc.high</code> > <code>auc.mid</code> > <code>auc.low</code></p>
  
  <ol start="3">
  <li>The IC50 values also seem a little low. Do they make sense?</li>
  </ol>
  
  <p>Bonus question: how do I avoid the trailing zeros in <code>slope</code>, <code>low</code>, <code>high</code>, <code>ED50</code>, <code>ic50.mid</code>, and <code>ic50.high</code> in the figure?</p>