Bioequivalence and Bioavailability Forum • handling of missing data

gracehung
☆

Taiwan,
2013-05-20 09:35
(4355 d 06:26 ago)

Posting: # 10603
Views: 22,403

handling of missing data [NCA / SHAM]

I have one question about the handling of the missing data, and need your advice.
My question is:
As in BE study, One sample’s actual sampling time was miss recorded. However as according to protocol, this sample was analyzed by our analytical laboratory.

Considering from PK personnel’s aspect, as it was not specified in protocol, how would you handling this missed time point at statistical analysis?
Will this data be included in statistical analysis? If yes, what is the handling procedure?
If the data is excluded, what kind of method would you use in calculate your AUC?
1. Would you adapt interpolation method?
2. If yes, would you adapt logarithmic interpolation or linear interpolation?
Is there any guideline or reference document we could follow?

Helmut
★★★

Vienna, Austria,
2013-05-20 16:22
(4354 d 23:40 ago)

@ gracehung
Posting: # 10604
Views: 20,044

Lin-up/log-down trapezoidal avoids most trouble

Post reply

Hi Gracehung,

please follow the Forum’s Policy in the future. THX.

❝ As in BE study, One sample’s actual sampling time was miss recorded. However as according to protocol, this sample was analyzed by our analytical laboratory.

Perform an audit. Why was an existing sample recorded as missing?

❝ 1. it was not specified in protocol…

Don’t worry, it’s too late. :-D

❝ … how would you handling this missed time point at statistical analysis?

In the future describe how you would handle missing samples in an SOP or – much better IMHO – in the protocol. In the latter case any procedure is more transparent to an assessor (i.e., approved by the IEC and the local authority).

❝ 2. Will this data be included in statistical analysis? If yes, what is the handling procedure?

Depends on the intellectual horsepower you are willing/able to invest. Missings are generally not problematic if unambiguously located in the absorption¹ or distribution/elimination² phase. If the missing might be C_max³ data imputation is problematic. Regulators don’t like PK modeling, so I have tried regressive splines with not very convincing results (references there). :-(

❝ 3. If the data is excluded, what kind of method would you use in calculate your AUC?

❝ a. Would you adapt interpolation method?

❝ b. If yes, would you adapt logarithmic interpolation or linear interpolation?

I always suggest to use the lin-up/log-down trapezoidal method. I think that Pharsight will make it even the default in the next version of Phoenix/WinNonlin. No need to impute data. Fine unless the missing might be C_max… See also this presentation about different integration algos (especially slides 16–17 on what happens with missing values).

❝ 4. Is there any guideline or reference document we could follow?

Not that I recall.

Let’s denote denote the missing as C_i. Then if

C_i+1 > C_i-1 ~ absorption
C_i+1 < C_i-1 ~ distribution/elimination
C_i+1 = C_i-1 ~ C_max or pure chance? Any algo will come up with the linear interpolation C_i = C_i-1 + |(t_i – t_i-1)/(t_i+1 – t_i-1)|(C_i+1 – C_i-1) – which likely underestimates the true value.

—
Dif-tor heh smusma 🖖🏼 Довге життя Україна!
Helmut Schütz

The quality of responses received is directly proportional to the quality of the question asked. 🚮
Science Quotes

yjlee168
★★★

Kaohsiung, Taiwan,
2013-05-22 10:25
(4353 d 05:36 ago)

@ Helmut
Posting: # 10614
Views: 19,647

Lin-up/log-down trapezoidal avoids most trouble

Post reply

Dear Helmut,

Great. However is there any special consideration when using lin-up/log-down trapezoidal rule to calculate AUCs if there is more than one absorption peak? Thanks in advance.

❝ ...

❝ I always suggest to use the lin-up/log-down trapezoidal method. I think that Pharsight will make it even the default in the next version of Phoenix/WinNonlin. No need to impute data. Fine unless the missing might be C_max… See also this presentation about different integration algos (especially slides 16–17 on what happens with missing values).

—
All the best,
-- Yung-jin Lee
bear v2.9.2:- created by Hsin-ya Lee & Yung-jin Lee
Kaohsiung, Taiwan https://www.pkpd168.com/bear
Download link (updated) -> here

Helmut
★★★

Vienna, Austria,
2013-05-22 22:25
(4352 d 17:36 ago)

@ yjlee168
Posting: # 10620
Views: 20,013

Multiple peaks: fallback to linear trapezoidal

Post reply

Hi Yung-jin!

❝ However is there any special consideration when using lin-up/log-down trapezoidal rule to calculate AUCs if there is more than one absorption peak?

Not that I am aware of. If values increase again the algo falls back to the linear trapezoidal. Try this clumsy code (a bimodal release formulation with two lag-times, two absorption rate constants, and a common elimination):

#set.seed(29081957) # keep this line to compare results only

require(truncnorm)

n       <- 48        # no of profiles

t       <- c(0,0.25,0.5,0.75,1,1.25,1.5,2,2.5,3,3.5,

             4,4.25,4.5,4.75,5,5.25,5.5,6,6.5,7,7.75,9,12,16,24)

D       <- 40

D.f     <- c(0.5, 0.5) # dose fractions of components

# Index ".d" denotes default values, ".c" CV

F.d     <- 0.5      # fraction absorbed (BA)

F.c     <- 0.1      # ±0.1, uniform

V.d     <- 1        # volume of distribution

V.c     <- 0.1      # CV 10%, lognormal

k01.1.d <- 1.4      # absorption rate constant (1st component)

k01.2.d <- 0.7      # absorption rate constant (2nd component)

k01.c   <- 0.2      # CV 20%, lognormal

k10.d   <- log(2)/3 # elimination rate constant (common to both components)

k10.c   <- 0.2      # CV 20%, lognormal

tlag1.d <- 0.5      # 1st lag-time

tlag1.l <- 0        # lower truncation 1

tlag1.u <- 1.5      # upper truncation 1

tlag2.d <- 4        # 2nd lag-time

tlag2.l <- 3        # lower truncation 2

tlag2.u <- 6        # upper truncation 2

tlag.c  <- 0.5      # CV 50%, truncated normal

tmax.th <- max((log(k01.1.d/k10.d)/(k01.1.d-k10.d)+tlag1.d),

               (log(k01.2.d/k10.d)/(k01.2.d-k10.d)+tlag2.d-tlag1.d))

Cmax.th <- F.d*D*D.f[1]/V.d*k01.1.d/(k01.1.d-k10.d)*

             (exp(-k10.d*(tmax.th-tlag1.d))-exp(-k01.1.d*(tmax.th-tlag1.d)))+

           F.d*D*D.f[2]/V.d*k01.2.d/(k01.2.d-k10.d)*

             (exp(-k10.d*(tmax.th-tlag2.d))-exp(-k01.2.d*(tmax.th-tlag2.d)))

x       <- seq(min(t), max(t), length.out=((max(t)-min(t))*4+1))

C1      <- F.d*D*D.f[1]/V.d*k01.1.d/(k01.1.d-k10.d)*

             (exp(-k10.d*(x-tlag1.d))-exp(-k01.1.d*(x-tlag1.d)))

C1[C1<0]<- 0        # set negative values of 1st component to zero

C2      <- F.d*D*D.f[2]/V.d*k01.2.d/(k01.2.d-k10.d)*

            (exp(-k10.d*(x-tlag2.d))-exp(-k01.2.d*(x-tlag2.d)))

C2[C2<0]<- 0        # set negative values of 2nd component to zero

C       <- C1+C2    # add components

AErr    <- 0.10     # analytical error CV 10%, lognormal

LLOQ    <- 0.01*Cmax.th # 1% of theoretical Cmax

split.screen(c(2,1))

for (s in 1:2) {

  screen(s)

  if (s == 1){

    plot(x, C, type="l", lwd=3, col="red", ylim=c(0, Cmax.th*1.75),

      main="linear trapezoidal", cex.main=1.1,

      xlab="time", ylab="concentration")

  } else {

    plot(x, C, type="l", lwd=3, col="red", ylim=c(0, Cmax.th*1.75),

      main="lin-up/log-down trapezoidal", cex.main=1.1,

      xlab="time", ylab="concentration")

  }

  rug(t, side=3)

  abline(v=0, lty=3)

  abline(h=0, lty=3)

}

AUCt.th <- NULL; AUCt.lin <- NULL; AUCt.linlog <- NULL

for (i in 1:n) {

  F      <- runif(n=1, min=F.d-F.c, max=F.d+F.c)

  V      <- rlnorm(n=1, meanlog=log(V.d)-0.5*log(V.c^2+1),

              sdlog=sqrt(log(V.c^2+1)))

  k01.1  <- rlnorm(n=1, meanlog=log(k01.1.d)-0.5*log(k01.c^2+1),

              sdlog=sqrt(log(k01.c^2+1)))

  k01.2  <- rlnorm(n=1, meanlog=log(k01.2.d)-0.5*log(k01.c^2+1),

              sdlog=sqrt(log(k01.c^2+1)))

  k10    <- rlnorm(n=1, meanlog=log(k10.d)-0.5*log(k10.c^2+1),

              sdlog=sqrt(log(k10.c^2+1)))

  tlag1  <- rtruncnorm(n=1, a=tlag1.l, b=tlag1.u, mean=tlag1.d, sd=tlag.c)

  tlag2  <- rtruncnorm(n=1, a=tlag2.l, b=tlag2.u, mean=tlag2.d, sd=tlag.c)

  C0.1   <- F*D*D.f[1]/V  # fictive iv concentration at t=0

  B1.1   <- C0.1*k01.1/(k01.1-k10)*exp(k10*tlag1)

  B2.1   <- C0.1*k01.1/(k01.1-k10)*exp(k01.1*tlag1) 

  C1     <- B1.1*exp(-k10*t)-B2.1*exp(-k01.1*t)

  C1[C1<0] <- 0

  C0.2   <- F*D*D.f[2]/V  # fictive iv concentration at t=0

  B1.2   <- C0.2*k01.2/(k01.2-k10)*exp(k10*tlag2)

  B2.2   <- C0.2*k01.2/(k01.2-k10)*exp(k01.2*tlag2) 

  C2     <- B1.2*exp(-k10*t)-B2.2*exp(-k01.2*t)

  C2[C2<0] <- 0

  C <- C1+C2

  for (j in 1:length(t)) { # analytical error (multiplicative 10% C)

    if (j==1)

      AErr1 <- rnorm(n=1, mean=0, sd=abs(C[j]*AErr))

    else

      AErr1 <- c(AErr1, rnorm(n=1, mean=0, sd=abs(C[j]*AErr)))

  }

  C <- C+AErr1 # add analytical error

  C[C<LLOQ] <- NA  # assign NAs to Cs below LLOQ

  tfirst <- (t[!is.na(C)])[1]              # 1st C >= LLOQ

  t0     <- t[length(head((t[is.na(C)])))] # last LLOQ before tfirst

  tlast  <- rev(t[!is.na(C)])[1]           # last C >= LLOQ

  screen(1, FALSE)

  plot(t, C, type="l", axes=FALSE,

    ylim=c(0, Cmax.th*1.75), xlab="", ylab="") # linear spaghetti plots

  screen(2, FALSE)

  for (j in 1:(length(t)-1)) {

    if (is.na(C[j]) == FALSE & is.na(C[j+1]) == FALSE) {

      if (C[j+1] >= C[j]) {

        lines(x=c(t[j],t[j+1]), y=c(C[j],C[j+1])) # increasing segments

      } else {

        tk <- NULL; Ck <- NULL # interpolated segments

        for (k in seq(t[j], t[j+1], length.out=10)) {

          tk <- c(tk, k)

          Ck <- c(Ck, exp(log(C[j])+abs((k-t[j])/(t[j+1]-t[j]))*(log(C[j+1])-log(C[j]))))

        }

        lines(tk, Ck) # decreasing segments

      }

    }

  }

  Cwo.err1 <- function(x) {B1.1*exp(-k10*x)-B2.1*exp(-k01.1*x)}

  Cwo.err2 <- function(x) {B1.2*exp(-k10*x)-B2.2*exp(-k01.2*x)}

  AUCt.1 <- integrate(Cwo.err1, lower=tlag1, upper=max(t))

  AUCt.2 <- integrate(Cwo.err2, lower=tlag2, upper=max(t))

  AUCt <- AUCt.1$value+AUCt.2$value # AUCt without analytical error

  AUCt.th <- c(AUCt.th, AUCt)

  pAUC.lin <- 0 # linear trapezidal rule

  for (j in seq(which(t == t0), which(t == tlast)-1)){

    if (t[j] == t0) { # 1st triangle

      pAUC.lin <- pAUC.lin + 0.5*(t[j+1]-t[j])*C[j+1]


    } else {

      pAUC.lin <- pAUC.lin + 0.5*(t[j+1]-t[j])*(C[j+1]+C[j])

    }

  }

  AUCt.lin <- c(AUCt.lin, pAUC.lin)

  pAUC.linlog <- 0 # lin-up/log-down trapezidal rule

  for (j in seq(which(t == t0), which(t == tlast)-1)){

    if (t[j] == t0) { # 1st triangle

      pAUC.linlog <- pAUC.linlog + 0.5*(t[j+1]-t[j])*C[j+1]


      } else {

      if (C[j+1] >= C[j]) { # linear trap.

        pAUC.linlog <- pAUC.linlog + 0.5*(t[j+1]-t[j])*(C[j+1]+C[j])

      } else { # loglinear trap.

        pAUC.linlog <- pAUC.linlog + (t[j+1]-t[j])*(C[j+1]-C[j])/(log(C[j+1]/C[j]))

      }

    }

  }

  AUCt.linlog <- c(AUCt.linlog, pAUC.linlog)

}

close.screen(all = TRUE)

par(ask=T)

RE.AUCt.lin <- 100*(AUCt.lin-AUCt.th)/AUCt.th

RE.AUCt.linlog <- 100*(AUCt.linlog-AUCt.th)/AUCt.th

boxplot(RE.AUCt.lin, RE.AUCt.linlog, ylab="% RE",

  ylim=c(-max(RE.AUCt.lin, RE.AUCt.linlog)*1.1, max(RE.AUCt.lin, RE.AUCt.linlog)*1.1),

  main="Bias of algos", cex.main=1.1,

  names=c("linear", "lin-up/log-down"))

abline(h=0, lty=3)

cat("\nCVtotal; linear trap.:", round(sd(AUCt.lin)/mean(AUCt.lin),4),

  "– lin-up/log-down trap.:", round(sd(AUCt.linlog)/mean(AUCt.linlog),4), "\n\n")

wilcox.test(RE.AUCt.lin, mu=0, conf.int=T)

wilcox.test(RE.AUCt.linlog, mu=0, conf.int=T)

wilcox.test(RE.AUCt.lin, RE.AUCt.linlog, conf.int=T)

par(ask=F)

No big deal, but the lin-up/log-down should perform better if there are missings.

—
Dif-tor heh smusma 🖖🏼 Довге життя Україна!
Helmut Schütz

The quality of responses received is directly proportional to the quality of the question asked. 🚮
Science Quotes

yjlee168
★★★

Kaohsiung, Taiwan,
2013-05-22 23:57
(4352 d 16:04 ago)

@ Helmut
Posting: # 10622
Views: 19,705

Multiple peaks: fallback to linear trapezoidal

Post reply

Dear Helmut,

Thank you for your messages and codes. I tried to run your codes without luck, but I will try it later. I will change the way of AUC calculation with bear later. I remembered that I read something about method comparisons for trapezoidal rule in one textbook, but just cannot remember which book.

❝ Not that I am aware of. If values increase again the algo falls back to the linear trapezoidal. Try this clumsy code (a bimodal release formulation with two lag-times, two absorption rate constants and a common elimination):

❝

❝ ...

❝ B1.1 <- C0.1*k01.1/(k01.1-K)*exp(k10*tlag1)

❝ B2.1 <- C0.1*k01.1/(k01.1-K)*exp(k01.1*tlag1)

Here I got Error: object 'K' not found and some error messages after this.

One more question, in your above presentation (slide#17[edited]; missing data occurred in R at 12h), it was a lin-up/log-down plot. As you said that no need to do data imputation if using lin-up/log-down. Do you do any line smoothing with your data first? Otherwise, the line should not be connected in that way. It should be connected as the line that the arrow points. Do I miss anything here?

[image]

—
All the best,
-- Yung-jin Lee
bear v2.9.2:- created by Hsin-ya Lee & Yung-jin Lee
Kaohsiung, Taiwan https://www.pkpd168.com/bear
Download link (updated) -> here

Helmut
★★★

Vienna, Austria,
2013-05-23 01:25
(4352 d 14:36 ago)

@ yjlee168
Posting: # 10623
Views: 19,836

No interpolation

Post reply

Hi Yung-jin,

❝ I will change the way of AUC calculation with bear later.

Relax. Would make only sense if the user has a choice. Some people still get crazy if one tries anything else than the linear trapezoidal.
Don’t use my code without improving it! In its current form it will not work with NAs and/or isolated BQLs within the profile.

❝ I remembered that I read something about method comparisons for trapezoidal rule in one textbook, but just cannot remember which book.

Yes. I think I have posted even references somewhere. To lazy to search now.

❝ ❝ ...

❝ ❝ B1.1 <- C0.1*k01.1/(k01.1-K)*exp(k10*tlag1)

❝ ❝ B2.1 <- C0.1*k01.1/(k01.1-K)*exp(k01.1*tlag1)

❝ Here I got Error: object 'K' not found and some error messages after this.

Sorry, I made a copy/paste error. Should work now (use k10 instead of K).

❝ One more question, in your above presentation (slide#17[edited]; missing data occurred in R at 12h), it was a lin-up/log-down plot.

Nope. 16/17 are linear plots. In slide 16 I connected all points by straight lines (what the linear trap. would do), and in slide 17 by straight lines if ascending and by an exponential if descending.

❝ As you said that no need to do data imputation if using lin-up/log-down. Do you do any line smoothing with your data first?

No.

❝ Otherwise, the line should not be connected in that way. It should be connected as the line that the arrow points.

Why?

❝ Do I miss anything here?

Yes, you do. ;-)

If the 12 h sample is missing the linear trapezoidal works as if we would impute a value by linear interpolation. The lin-up/log-down works as if we would impute a value by logarithmic interpolation. But the trick is that nothing has to be imputed in the lin-up/log-down at all!

  t    C

  8  33.17

 12  miss.

 16   7.86

Interpolated values for imputation:

linear: C₁₂^*=33.17+|(12-8)/(16-8)|(7.86-33.17)=20.52

linlog: C₁₂^*=exp(log(33.17)+|(12-8)/(16-8)|log(7.86/33.17))=16.15

Calculation based on imputed values:

linear/lin imp.: pAUC_8-16=0.5(12- 8)(33.17+20.52)+

                          0.5(16-12)(20.52+ 7.86)=164.1

linear/log imp.: pAUC_8-16=0.5(12- 8)(33.17+20.52)+

                          0.5(16-12)(20.52+ 7.86)=146.7

linlog/log imp.: pAUC_8-16=(12- 8)(16.15-33.17)/log(16.15/33.17)+

                          (16-12)( 7.86-16.15)/log( 7.86/16.15)=140.6

Direct calculation without imputation:

linear: pAUC_8-16=0.5(16-8)(33.17+7.86)=164.1

linlog: pAUC_8-16=(16-8)(7.86-33.17)/log(7.86/33.17)=140.6

BTW, the theoretical pAUC_8-16 based on the model is 140.62…

To make a long story short: Only if one has to use the linear trapezoidal (for whatever wacky reasons) I would impute an estimate based on the log/linear interpolation (second variant above). Otherwise the bias is substantial. Q.E.D.

—
Dif-tor heh smusma 🖖🏼 Довге життя Україна!
Helmut Schütz

The quality of responses received is directly proportional to the quality of the question asked. 🚮
Science Quotes

yjlee168
★★★

Kaohsiung, Taiwan,
2013-05-23 13:43
(4352 d 02:18 ago)

@ Helmut
Posting: # 10626
Views: 19,318

No interpolation

Post reply

Dear Helmut,

Thanks for your messages, and the codes works well now.

❝ ...

❝ Relax. Would make only sense if the user has a choice. Some people still get crazy if one tries anything else than the linear trapezoidal.

Yes, that's my concern too. I should provide users options for this step in NCA. It should be easy to implement such a function with bear at the moment.

❝ Don’t use my code without improving it! In its current form it will not work with NAs and/or isolated BQLs within the profile.

Sure and thanks for your reminding.

❝ ...

❝ ❝ Otherwise, the line should not be connected in that way. It should be connected as the line that the arrow points.

❝

❝ Why?

Sorry, I thought it was using lin-up/log-down algo and plots in slide#17.

—
All the best,
-- Yung-jin Lee
bear v2.9.2:- created by Hsin-ya Lee & Yung-jin Lee
Kaohsiung, Taiwan https://www.pkpd168.com/bear
Download link (updated) -> here

Helmut
★★★

Vienna, Austria,
2013-05-23 17:15
(4351 d 22:46 ago)

@ yjlee168
Posting: # 10630
Views: 19,693

Spaghetti & other pasta

Post reply

Dear Yung-jin,

❝ ❝ ❝ Otherwise, the line should not be connected in that way. It should be connected as the line that the arrow points.

❝ ❝ Why?

❝ Sorry, I thought it was using lin-up/log-down algo and plots in slide#17.

We are all guilty in not representing in spaghettis plots the method we use for calculating the AUC. For simplicity we connect point with straight lines, regardless whether we use a linear or semilogarithmic scale (all left panels below).
If we use the linear trapezoidal, the semilog plot is wrong. What we are actually calculating would be reflected in the lower right.

[image]

Weird, isn’t it? But then hopefully everybody would realize that the algo overestimates AUCs in the distribution / elimination phase.

On the other hand if we use the lin-up/log-down trapezoidal the linear plot is wrong; should use the upper right one instead.

[image]

Have you ever seen something like this before? I didn’t. Would be a tough job to come up with such a setup in commercial software. :-D

PS: Very rarely you find publications where in mean plots the points are not connected by lines at all.
PPS: What do you think about adding rug(t, side=3) to plots in bear?

—
Dif-tor heh smusma 🖖🏼 Довге життя Україна!
Helmut Schütz

The quality of responses received is directly proportional to the quality of the question asked. 🚮
Science Quotes

yjlee168
★★★

Kaohsiung, Taiwan,
2013-05-23 22:47
(4351 d 17:14 ago)

@ Helmut
Posting: # 10632
Views: 19,411

Spaghetti & other pasta

Post reply

Dear Helmut,

Wow! That really quite surprises me. I know that all errors may result from "the approximation of trapezoid" (or assumption of trapezoid area between two points) when calculating AUC. Yes, right now we just use line connection between points to plot linear/semilog graphs. From your plots, it looks like lin-up/log-down method is more close to actual AUC. In this case, it does no matter what there is any missing data or not. Lin-up/log-down is more accurate than linear, isn't?

❝ If we use the linear trapezoidal, the semilog plot is wrong. What we are actually calculating would be reflected in the lower right.

❝

❝ [image]

❝

❝ Weird, isn’t it? But then hopefully everybody would realize that the algo overestimates AUCs in the distribution / elimination phase.

Weird indeed but very convincible. Great demo. How do you plot the low-right graph? Is it possible to overlap both the lower-left and the lower-right graph as one plot? That make differences more clear.

❝ Have you ever seen something like this before? I didn’t. Would be a tough job to come up with such a setup in commercial software. :-D

Me either. Why do you think that it would be a tough job for commercial software?

❝ PPS: What do you think about adding rug(t, side=3) to plots in bear?

Should be do-able in bear. But I don't quite understand why you want to add rug() plot? We don't have many data points in each subject (usually less than 50's in each period) for a BE/BA study.

—
All the best,
-- Yung-jin Lee
bear v2.9.2:- created by Hsin-ya Lee & Yung-jin Lee
Kaohsiung, Taiwan https://www.pkpd168.com/bear
Download link (updated) -> here

Helmut
★★★

Vienna, Austria,
2013-05-23 23:29
(4351 d 16:32 ago)

@ yjlee168
Posting: # 10633
Views: 19,685

Spaghetti & other pasta

Post reply

Dear Yung-jin!

❝ […] From your plots, it looks like lin-up/log-down method is more close to actual AUC. In this case, it does no matter what there is any missing data or not. Lin-up/log-down is more accurate than linear, isn't?

I think so. Whatever PK a drug follows, distribution/elimination likely is first order (well, leaving :party:

and Michaelis-Menten aside).

❝ Weird indeed but very convincible. […] How do you plot the low-right graph? Is it possible to overlap both the lower-left and the lower-right graph as one plot?

Sure. Just interpolate linear within each decreasing segment. Alternatively to lines 12–17 you could set up a function and use curve(). Was too lazy.

t <- c(0,0.5,0.75,1,1.5,2,2.5,3,4,6,8,16,24)

C <- c(0,69.59,84.64,92.61,96.73,93.45,87.51,

       80.85,68.02,47.54,33.17,7.86,1.86)

plot(t, C, type="p", ylim=c(1, max(C)*1.1), log="y", las=1,

  main="semilogarithmic", cex.main=1.1,

  xlab="time", ylab="concentration")

lines(t, C, col="blue")

for (j in 1:(length(t)-1)) {

  if (C[j+1] >= C[j]) {

    lines(x=c(t[j], t[j+1]), y=c(C[j], C[j+1]), col="blue")

  } else {

    tk <- NULL; Ck <- NULL

    for (k in seq(t[j], t[j+1], length.out=10)) {

      tk <- c(tk, k)

      Ck <- c(Ck, C[j]+abs((k-t[j])/(t[j+1]-t[j]))*(C[j+1]-C[j]))

    }

    lines(tk, Ck, col="red")

  }

}

rug(t, side=1)

legend("topright",

  c("trapezoidal rules", "linear", "lin-up/log-down"),

  lty=c(NA, 1, 1), col=c(NA, "red", "blue"))

If you want to use the code in a linear plot, interpolate logarithmically, i.e., change line 15 to
Ck <- c(Ck, exp(log(C[j])+abs((k-t[j])/(t[j+1]-t[j]))*(log(C[j+1])-log(C[j]))))

❝ Me either. Why do you think that it would be a tough job for commercial software?

I have no experiences with SAS, but in Phoenix/WinNonlin it is tricky indeed.

❝ […] But I don't quite understand why you want to add rug() plot? We don't have many data points in each subject (usually less than 50's in each period) for a BE/BA study.

It was just an idea to display the sampling schedule. Forget it. ;-)

—
Dif-tor heh smusma 🖖🏼 Довге життя Україна!
Helmut Schütz

The quality of responses received is directly proportional to the quality of the question asked. 🚮
Science Quotes

yjlee168
★★★

Kaohsiung, Taiwan,
2013-05-24 00:46
(4351 d 15:16 ago)

@ Helmut
Posting: # 10634
Views: 19,289

Spaghetti & other pasta

Post reply

Dear Helmut,

Thanks again. Now I think that I have fully understood the differences between the linear and the lin-up/log-down method for AUC calculations. Actually, lin-up/log-down method is still a kind of data smoothing approach. It is applied to calculate more accurate AUC than the linear one. Some people should not get crazy when using lin-up/log-down if they can see the differences here. Just because the line is too smooth to be true? That may scare people a little bit. That's why I asked if you used some kinds of data smoothing methods in previous message. The key point is that involved calculation step still uses the raw data to calculate the interpolated data. As we can see, the errors between linear and lin-up/log-down become significantly apparent, especially when sampling time between two data points is getting greater. Therefore, when there is missing data occurring, the lin-up/log-down method shows more accurate than the linear.

❝ [...]

❝ I have no experiences with SAS, but in Phoenix/WinNonlin it is tricky indeed.

should be no problem at all. They are professional after all.

❝ [...]

❝ It was just an idea to display the sampling schedule. Forget it. ;-)

I see. However, it could be misunderstood or confused by people to interpret it (the rug) as a 'log scale' with x axis.

—
All the best,
-- Yung-jin Lee
bear v2.9.2:- created by Hsin-ya Lee & Yung-jin Lee
Kaohsiung, Taiwan https://www.pkpd168.com/bear
Download link (updated) -> here

Helmut
★★★

Vienna, Austria,
2013-05-24 03:24
(4351 d 12:37 ago)

@ yjlee168
Posting: # 10635
Views: 19,306

Spaghetti & other pasta

Post reply

Dear Yung-jin!

❝ Now I think that I have fully understood the differences between the linear and the lin-up/log-down method for AUC calculations.

Hhm.

❝ Actually, lin-up/log-down method is still a kind of data smoothing approach.

No, it isn’t.

❝ The key point is that involved calculation step still uses the raw data to calculate the interpolated data.

There is no interpolation done at all. Have a look at my example in the previous post. I gave the interpolation/imputation formulas only as a hint what might happen if we need an intermediate data point. If we want to calculate the area in the interval [t₂, t₂] no interpolation, smoothing, whatsoever is done. Both methods are straightforward:

linear: AUC_t1→t2=0.5(t₂–t₁)(C₂+C₁)

linlog: AUC_t1→t2=(t₂–t₁)(C₂–C₁)/log(C₂/C₁)

❝ As we can see, the errors between linear and lin-up/log-down become significantly apparent, especially when sampling time between two data points is getting greater.

Right.

❝ Therefore, when there is missing data occurring, the lin-up/log-down method shows more accurate than the linear.

Xactly. See also Nathan’s PK/PD Blog.

—
Dif-tor heh smusma 🖖🏼 Довге життя Україна!
Helmut Schütz

The quality of responses received is directly proportional to the quality of the question asked. 🚮
Science Quotes

yjlee168
★★★

Kaohsiung, Taiwan,
2013-05-24 11:17
(4351 d 04:45 ago)

@ Helmut
Posting: # 10638
Views: 19,371

Spaghetti & other pasta

Post reply

Dear Helmut,

❝ ❝ Actually, lin-up/log-down method is still a kind of data smoothing approach.

❝

❝ No, it isn’t.

You are right. It isn't. Sorry about this.

❝ [...] I gave the interpolation/imputation formulas only as a hint what might happen if we need an intermediate data point. If we want to calculate the area in the interval [t₂, t₂] no interpolation, smoothing, whatsoever is done. Both methods are straightforward:

I see. Thank you for your detailed explanation once again.

❝ linear: AUC_t1→t2=0.5(t₂–t₁)(C₂+C₁)

❝ linlog: AUC_t1→t2=(t₂–t₁)(C₂–C₁)/log(C₂/C₁)

—
All the best,
-- Yung-jin Lee
bear v2.9.2:- created by Hsin-ya Lee & Yung-jin Lee
Kaohsiung, Taiwan https://www.pkpd168.com/bear
Download link (updated) -> here

Ken Peh
★

Malaysia,
2013-05-30 21:55
(4344 d 18:06 ago)

@ Helmut
Posting: # 10688
Views: 19,182

No interpolation

Post reply

Dear Helmut,

❝ Interpolated values for imputation:

❝ linear: C₁₂^*=33.17+|(12-8)/(16-8)|(7.86-33.17)=20.52

ok.

❝ linlog: C₁₂^*=exp(log(33.17)+|(12-8)/(16-8)|log(7.86/33.17))=16.15

How did you get the readings of 16.15? I tried to follow your equation but I could not get the answer. :confused:

❝ Calculation based on imputed values:

❝ linear/lin imp.: pAUC_8-16=0.5(12- 8)(33.17+20.52)+

❝ 0.5(16-12)(20.52+ 7.86)=164.1

❝ linear/log imp.: pAUC_8-16=0.5(12- 8)(33.17+20.52)+

❝ 0.5(16-12)(20.52+ 7.86)=146.7

linear/lin imp and linear/log imp have the same data set. Why the output is different ? :confused:

❝ linlog/log imp.: pAUC_8-16=(12- 8)(16.15-33.17)/log(16.15/33.17)+

❝ (16-12)( 7.86-16.15)/log( 7.86/16.15)=140.6

❝

❝ Direct calculation without imputation:

❝ linear: pAUC_8-16=0.5(16-8)(33.17+7.86)=164.1

ok.

❝ linlog: pAUC_8-16=(16-8)(7.86-33.17)/log(7.86/33.17)=140.6

Similarly, I could not get the answer of 140.61. How did you calculate ? :confused:

Please kindly elaborate.

Thank you.

Ken

Helmut
★★★

Vienna, Austria,
2013-05-31 00:38
(4344 d 15:23 ago)

@ Ken Peh
Posting: # 10690
Views: 19,343

Different algos!

Post reply

Hi Ken!

❝ ❝ linlog: C₁₂^*=exp(log(33.17)+|(12-8)/(16-8)|log(7.86/33.17))=16.15

❝ How did you get the readings of 16.15? I tried to follow your equation but I could not get the answer.

round(exp(log(33.17)+abs((12-8)/(16-8))*log(7.86/33.17)), 2)

[1] 16.15

SAS

x=round(exp(log(33.17)+abs((12-8)/(16-8))*log(7.86/33.17)), 2);

16.15

Excel / OO Calc / Gnumeric

=ROUND(EXP(LN(33.17)+ABS((12-8)/(16-8))*LN(7.86/33.17)), 2)

16.15

❝ linear/lin imp and linear/log imp have the same data set. Why the output is different ?

Since I’m using different algorithms. Once more: Forget linear interpolation in the descending part of the profile.

[image]

Have you ever asked yourself why the arithmetic and geometric means of the same data are different?
Compare (7.86+33.17)*0.5
with (7.86*33.17)^0.5
Are these values vaguely reminiscent of something you have seen before? :lookaround:

❝ ❝ linlog: pAUC_8-16=(16-8)(7.86-33.17)/log(7.86/33.17)=140.6

❝ Similarly, I could not get the answer of 140.61.

round((16-8)*(7.86-33.17)/log(7.86/33.17), 1)

[1] 140.6

=ROUND((16-8)*(7.86-33.17)/LN(7.86/33.17), 1)

140.6

❝ How did you calculate ?

See above. And you?

—
Dif-tor heh smusma 🖖🏼 Довге життя Україна!
Helmut Schütz

The quality of responses received is directly proportional to the quality of the question asked. 🚮
Science Quotes

Ken Peh ★ Malaysia, 2013-06-03 19:52 (4340 d 20:09 ago) @ Helmut Posting: # 10717 Views: 18,961	Different algos! Post reply
	Dear Helmut, Thank you very much. ❝ ❝ How did you calculate ? I know where my mistake is. It is the positive and negative signs. Brain is jammed up. Need more practice. May I know which software in R you refer to that can be used to calculate the data ? Regards, Ken

Helmut
★★★

Vienna, Austria,
2013-06-03 20:11
(4340 d 19:51 ago)

@ Ken Peh
Posting: # 10718
Views: 19,197

Calculate what?

Post reply

Hi Ken!

❝ May I know which software in R you refer to that can be used to calculate the data ?

… to calculate what? For imputation use my code above (simply give the timepoints in the t-vector and concentrations in C). For the lin-up/log-down trapezoidal rule in the calculation of AUC see the end of my code in this post.
I guess it’s easier to wait for the next version of bear.* ;-)

there

—
Dif-tor heh smusma 🖖🏼 Довге життя Україна!
Helmut Schütz

The quality of responses received is directly proportional to the quality of the question asked. 🚮
Science Quotes

Helmut
★★★

Vienna, Austria,
2013-05-25 17:09
(4349 d 22:52 ago)

@ yjlee168
Posting: # 10656
Views: 19,595

Lin-up/log-down trapezoidal example

Post reply

Dear Yung-jin,

❝ […] is there any special consideration when using lin-up/log-down trapezoidal rule to calculate AUCs if there is more than one absorption peak?

I checked one study of bimodal release formulations. Imbalanced, n=15. I evaluated AUC_t based on the linear trapezoidal and the lin-up/log-down method. Then I introduced a failure of the centrifuge resulting in a complete loss of 12 h samples in the second period. Let’s see:

PE 90% CI %CV_intra %CV_inter complete, linear 99.89 95.23–104.78 7.38 33.75 complete, lin-up/log-down 99.78 95.10–104.69 7.42 33.40 incomplete, linear 100.08 95.64–104.72 7.00 33.84 incomplete, lin-up/log-down 99.93 95.37–104.72 7.23 33.39

With the linear trapezoidal results change whereas with the lin-up/log-down results are affected to a lesser degree. I didn’t check but expect to get identical results in a balanced case.

—
Dif-tor heh smusma 🖖🏼 Довге життя Україна!
Helmut Schütz

The quality of responses received is directly proportional to the quality of the question asked. 🚮
Science Quotes

yjlee168
★★★

Kaohsiung, Taiwan,
2013-05-25 21:45
(4349 d 18:17 ago)

@ Helmut
Posting: # 10658
Views: 19,255

Lin-up/log-down trapezoidal example

Post reply

Dear Helmut,

I think the differences at elimination phase between linear and log-down may results from the degree of concavity ([C1-C2]) within a know time interval (t2-t1). I also did a simple test with two different scenarios:

** if t2-t1 = 12 h; C1 = 4.1, C2 =1.1 (the curve is more concave between t1 and t2)



  AUC_(t1-t2) linear = 31.2; log-down = 27.362 (log-down should be more accurate in this situation)

  

** if t2-t1 = 12 h; C1 = 4.1, C2 =4.05 (the curve is less concave between t1 and t2; AUC_(t1-t2) should be more close to real trapezoidal shape)



  AUC_(t1-t2) linear = 48.9; log-down = 48.9 (before truncated was 48.89939) (the linear is supposed to be better than log-down; however log-down still keeps the accuracy as good as linear and almost as same as the linear.)

Therefore, you're right about log-down method. The log-down can be set as default method for AUC calculation safely, accurately and reliably.

—
All the best,
-- Yung-jin Lee
bear v2.9.2:- created by Hsin-ya Lee & Yung-jin Lee
Kaohsiung, Taiwan https://www.pkpd168.com/bear
Download link (updated) -> here

Ohlbe
★★★

France,
2013-05-20 23:45
(4354 d 16:16 ago)

@ gracehung
Posting: # 10609
Views: 19,872

handling of missing data

Post reply

Dear Gracehung,

❝ As in BE study, One sample’s actual sampling time was miss recorded. However as according to protocol, this sample was analyzed by our analytical laboratory.

I understand you question differently from Helmut. My understanding is that the sample was collected, but the phlebotomist failed to record the actual blood sampling time. In such a case I would just file a protocol deviation and re-train the staff, full stop. Even if there were a sampling time deviation, the influence on the AUC 90 % CI is so limited that it would only change the outcome of the trial in a very borderline situation.

—
Regards
Ohlbe

Helmut ★★★ Vienna, Austria, 2013-05-21 15:53 (4354 d 00:09 ago) @ Ohlbe Posting: # 10610 Views: 19,480	Sorry Post reply
	Hi Ohlbe, ❝ I understand you question differently from Helmut. Though I understood it correctly (see the start of my post), I am notorious for answering unasked questions as well. — Dif-tor heh smusma 🖖🏼 Довге життя Україна! Helmut Schütz The quality of responses received is directly proportional to the quality of the question asked. 🚮 Science Quotes

gracehung
☆

Taiwan,
2013-05-23 03:06
(4352 d 12:55 ago)

@ Helmut
Posting: # 10624
Views: 19,541

Sorry

Post reply

Dear all:

Thanks for your reply!
The dilemma is we have the concentration data but don't have the actual sampling time. If we treat it as missing data, don't use it in bioequivalence evaluation, would assessor think we are hiding data? Should we include the data in PK analysis and use schedule time to be the "best guess"? One literature suggest to use log linear interpolation (at elimination phase) to minimize the bias (also as your reply). But someone said the assessor will think you are manipulating data. :confused:

which one would be the best choice?

Grace

Helmut
★★★

Vienna, Austria,
2013-05-23 03:35
(4352 d 12:26 ago)

@ gracehung
Posting: # 10625
Views: 19,521

Uncertain time point

Post reply

Hi Grace!

❝ If we treat it as missing data, don't use it in bioequivalence evaluation, would assessor think we are hiding data?

Duno. But if you treat the value as missing – depending on the algo – you might get a substantially biased AUC (see my conversation with Yung-jin above).

❝ Should we include the data in PK analysis and use schedule time to be the "best guess"?

Definitely (as Ohlbe suggested). I once had to assess deficiencies where in two studies the scheduled time points (instead of the actual ones) were used. Impact on the AUC comparison was in the second decimal place. And you have only one (!) value in an entire study…
BTW, scheduled time points were my standard method in hundreds (!) of studies. So what?

❝ One literature suggest to use log linear interpolation (at elimination phase) to minimize the bias (also as your reply).

Yes, but only if the value is missing. You know the value – only its time point is uncertain.

❝ But someone said the assessor will think you are manipulating data.

I would not use an estimate. Use the measured one at the scheduled time point. Report the deviation from the protocol, don’t panic, & sleep well.

Another idea – only if you really have nothing better to do: Most protocols require a justification of a time deviation to be recorded in CRF only if a certain time interval is exceeded (f.i. ±2’ in the very early phase, ±5’ later on, and ±15’ for t ≥24 h). Let’s assume only the actual sampling was not recorded, but still was within the allowed interval (otherwise we have two errors instead of one). Set up two data sets: One assuming the sample was taken at the maximum allowed negative deviation (too early) and the other one with the max. positive deviation (too late). Go ahead with the BE and assess the impact on the outcome.

—
Dif-tor heh smusma 🖖🏼 Довге життя Україна!
Helmut Schütz

The quality of responses received is directly proportional to the quality of the question asked. 🚮
Science Quotes