Bioequivalence and Bioavailability Forum

Hi John,

❝ I modified the PKNCA codes from the example above to calculate PK parameters for the following dataset. PKNCA's AUClast value was 5997.7885480, but my excel calculated value was 6825.422.???? I am not sure what went wrong.

In PKNCA the linear-up/log-down rule is the default (see this article). Let’s go for a second opinion:

calc.AUC <- function(t, C, rule = "linlog", digits = 5) {
  # if nonnumeric codes like BQL or NR are used, they
  # must be enclosed in 'single' or "double" quotes
  if (!length(C) == length(t))
    stop ("Both vectors must have the same number of elements.")
  if (!length(unique(t)) == length(t))
    stop ("Values in t-vector are not unique.")
  if (!rule %in% c("lin", "linlog"))
    stop ("rule must be either \"lin\" or \"linlog\".")
  # convert eventual nonnumerics to NAs
  C    <- unlist(suppressWarnings(lapply(C, as.numeric)))
  tmax <- min(t[C == max(C, na.rm = TRUE)], na.rm = TRUE)
  C[t  < tmax & is.na(C)] <- 0   # set NAs prior to tmax to zero
  x    <- NULL
  y    <- data.frame(t = t, C = C, section = "", pAUC = 0)
  y    <- y[with(y, order(t)), ] # belt plus suspenders
  NAs  <- which(is.na(C))        # location of eventual NAs
  if (sum(NAs) > 0) {            # remove NAs temporarily
    x      <- y[!complete.cases(y), ]
    x$pAUC <- x$AUC <- NA
    y      <- y[complete.cases(y), ]
  }
  for (i in 1:(nrow(y) - 1)) {
    if (rule == "linlog") {
      y$section[1]     <- "up"
      if (y$C[i+1]     < y$C[i]) { # decreasing
        y$pAUC[i+1]    <- (y$t[i+1] - y$t[i]) * (y$C[i+1] - y$C[i]) / log(y$C[i+1] / y$C[i])
        y$section[i+1] <- "down"
      } else {                     # increasing or equal
        y$pAUC[i+1]    <- 0.5 * (y$t[i+1] - y$t[i]) * (y$C[i+1] + y$C[i])
        y$section[i+1] <- "up"
      }
    } else {
        y$pAUC[i+1]  <- 0.5 * (y$t[i+1] - y$t[i]) * (y$C[i+1] + y$C[i])
    }
  }
  y$AUC    <- cumsum(y$pAUC)         # sum partial AUCs
  y        <- rbind(x, y)            # get the NAs back
  y        <- y[with(y, order(t)), ] # sort by time
  y[, 4:5] <- round(y[, 4:5], digits)
  if (rule == "lin") y <- y[, -3]    # not needed
  if (rule == "linlog") {
    attr(y, "rule") <- "linear-up/log-down rule"
  } else {
    attr(y, "rule") <- "linear rule"
  }
  return(y)
}

Now your data (IMHO, zero at the end doesn’t make sense):

t  <- c(0, 0.5, 1, 2, 4, 8, 24, 48, 72)
C  <- c(0, 307.838, 467.58, 439.886, 489.621, 306.011, 54.424, 3.628, 'BQL')


With my homebrew:

df <- calc.AUC(t, C, rule = "linlog", digits = 6)
cat(attr(df, "rule"), "\n"); print(df, row.names = FALSE)
linear-up/log-down rule
    t       C section      pAUC       AUC
  0.0   0.000      up    0.0000    0.0000
  0.5 307.838      up   76.9595   76.9595
  1.0 467.580      up  193.8545  270.8140
  2.0 439.886    down  453.5921  724.4061
  4.0 489.621      up  929.5070 1653.9131
  8.0 306.011    down 1562.6031 3216.5162
 24.0  54.424    down 2331.1067 5547.6230
 48.0   3.628    down  450.1656 5997.7885
 72.0      NA                NA        NA

df <- calc.AUC(t, C, rule = "lin", digits = 6)
cat(attr(df, "rule"), "\n"); print(df, row.names = FALSE)
linear rule
    t       C      pAUC       AUC
  0.0   0.000    0.0000    0.0000
  0.5 307.838   76.9595   76.9595
  1.0 467.580  193.8545  270.8140
  2.0 439.886  453.7330  724.5470
  4.0 489.621  929.5070 1654.0540
  8.0 306.011 1591.2640 3245.3180
 24.0  54.424 2883.4800 6128.7980
 48.0   3.628  696.6240 6825.4220
 72.0      NA        NA        NA

Confirmed that PKNCA is using the linear-up/log-down and you the linear trapezoidal in bloody Excel. If you really want to do that (I hope, you don’t), see there.
Note: If you give 0 at 72, you will get the same result with the lin-up/log-down, but 6868.958 with the linear rule. That’s Pharsight’s crappy AUCall.