Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Incremental Area Under the Curve (iAUC) in R

Tags:

r

The area under the curve can be calculated using the trapz function of the pracma package in R. iAUC is much more accurate in many instances, particularly in biology. However there is no, to me known, R function to calculate this. The trap function follows:

Example:

a <- c(1, 4, 5, 6)
b <- c(2, 4, 6, 8)
plot(a, b)
trapz(a, b)

Answer: 21

Does anyone know how to use R to calculate the iAUC? A guide can be found here.

Thanks in advance.

Update In response to Cheng-Yow's answer below: Another example

# Recall that trapz-function is available in the pracma package
# Here is the function kindly provided by Cheng-Yow Timothy Lin
library(pracma) # should be loaded already
traps <- function(a, b) {
     for (i in 1:length(b)-1){
          iAUC[i] <-((b[i]+b[i+1])/2)*(a[i+1]-a[i])
     }
     return(sum(iAUC))
}

# Data example
a <- c(-15,  -5,   1,   2,   3,   5,   8,  10,  15,  20,  30,  45,  60,  90, 120)
b <- c(50.20604,  49.59338,  47.39944,  56.38831,  69.43493,  73.92512,  61.92072,  67.92632, 115.45669,
195.03242, 322.15894, 291.30094, 289.20284, 238.70562, 156.23798)

traps(a, b)

trapz(a, b)

They yield the same results- is this the incremental area under the curve? It is not how they have explained the trapz-function...

Thanks for any enlightenment!

like image 414
Adam Robinsson Avatar asked Oct 27 '25 05:10

Adam Robinsson


1 Answers

    a <- c(1,2,3,4,5)
b <- c(1,2,-1,-2,2)
plot(a,b)
lines(a, b)
abline(b[1],0)
iAUC <- vector()

for (i in 1:(length(a)-1)){

  if((b[i+1]-b[1] >= 0) && (b[i]-b[1] >= 0))
  {iAUC[i] <-((b[i]-b[1]+b[i+1]-b[1])/2)*(a[i+1]-a[i])} 

  else if((b[i+1]-b[1] < 0) && (b[i]-b[1] >= 0))
  {iAUC[i] <-(b[i]-b[1])*((b[i]-b[1])/(b[i]-b[i+1])*(a[i+1]-a[i])/2)} 

  else if((b[i+1]-b[1] >= 0) && (b[i]-b[1] < 0))
  {iAUC[i] <-(b[i+1]-b[1])*((b[i+1]-b[1])/(b[i+1]-b[i])*(a[i+1]-a[i])/2)} 

  else if((b[i]-b[1] < 0) && (b[i+1]-b[1] < 0))
  {iAUC[i] <- 0}
}
sum(iAUC)
iAUC    
like image 164
Cheng-Yow Timothy Lin Avatar answered Oct 29 '25 20:10

Cheng-Yow Timothy Lin