Bioequivalence and Bioavailability Forum

Bong sure,

referring to the other thread where the new BE guideline from EMEA is discussed and where HS proposed a simulation.

As usual I am several steps behind as I am not sure I understand all details of the discussion. Nevertheless, the word anticonservative was right in my face.
Can that claim be generalised in any way?

Following the philosophy of HS' single simulations experiment I think we could simulate different 3-per, 3-trt datasets under various conditions, analyse them by a mixed model and record sigma = the residual error (this is the s used in calculation of the conficence interval). Next, prune the dataset so it's possible to analyse it as a 2,2,2-BE dataset with a normal linear model, and then obtain the residual and compare. The philosophy here is that if the pruning results in a sigma that is on average larger than the sigma obtained through the full dataset, then the approach is conservative; otherwise it is not.
Here's some elmaestrolophystic code (probably bugged and wrong):


load(nlme)
nConservative=0
nAntiConservative=0

seq1 = as.factor(c(1,1,1, 2,2,2, 3,3,3, 4,4,4, 5,5,5, 6,6,6, 1,1,1, 2,2,2, 3,3,3, 4,4,4, 5,5,5, 6,6,6))
sub1 = as.factor(c(1,1,1, 2,2,2, 3,3,3, 4,4,4, 5,5,5, 6,6,6, 7,7,7, 8,8,8, 9,9,9, 10,10,10, 11,11,11, 12,12,12))
per1 = as.factor(c(1,2,3, 1,2,3, 1,2,3, 1,2,3, 1,2,3, 1,2,3, 1,2,3, 1,2,3,1,2,3, 1,2,3, 1,2,3, 1,2,3))
trt1 = as.factor(c(1,2,3, 2,3,1, 3,1,2, 1,3,2, 2,1,3, 3,2,1 ,1,2,3 ,2,3,1 ,3,1,2 ,1,3,2 ,2,1,3 ,3,2,1))

##massacrating out everything related to R2
##means that when T was before R1 we now have a seq 1 otherwise a seq 2
sub2 = as.factor(c(1,1,2,2,3,3,4,4,5,5,6,6,7,7,8,8,9,9,10,10,11,11,12,12))
per2 = as.factor(c(1,2,1,2,1,2,1,2,1,2,1,2,1,2,1,2,1,2,1,2,1,2,1,2))
trt2 = as.factor(c(1,2, 2,1 ,1,2, 1,2, 2,1, 2,1, 1,2, 2,1, 1,2, 1,2, 2,1, 2,1))
seq2 = as.factor(c(1,1, 2,2, 1,1, 1,1, 2,2, 2,2, 1,1, 2,2, 1,1, 1,1, 2,2, 2,2))

for (iter in 1:100)
{

## T is 1, R1 is 2 and R2 is 3
lny1 = rnorm(n=36, mean=95, sd=0.15*95)
for (i in 1:36) if (trt1[i]==2) lny1[i]=rnorm(n=1, mean = 90, sd = 0.10*90)
for (i in 1:36) if (trt1[i]==3) lny1[i]=rnorm(n=1, mean = 110, sd = 0.15*110)

ElMaestrosMixedModel= lme(lny1~0+seq1+per1+trt1, random = ~1|sub1)
Smm=ElMaestrosMixedModel$sigma

##do the pruning!
lny2=lny1
length(lny2)=24
j=0
for (i in 1:36) if (trt1[i]!=3)
  {
    j=j+1;
    lny2[j]=lny1[i]   }
ElMaestrosLinearModel= lm(lny2~0+seq2+per2+trt2+sub2)
Slm=summary(ElMaestrosLinearModel)$sigma
if (Slm<Smm) nAntiConservative=nAntiConservative+1
if (Slm>Smm) nConservative=nConservative+1
}

cat("n AntiConservative = ", nAntiConservative, "\n")
cat("n Conservative = ", nConservative, "\n")

On my machine the code above is a win for anticonservativism. But that is dependent on the settings of means and CV in the treatment groups.
HS: I got the impression that you simulated your dataset on the ordinary scale, then log-transformed the data. At least, if I log-transform your data then I can reproduce your result (with the code above). Here I am simulating normal dist data on the log scale. Seems more right to me.

Bugs galore....

EM.

Update a few mins after posting: Seems escape characters do not go well on this forum?! The last two lines starting with "cat" was intended to finish with a letter n preceded by a down-right slash but I cannot paste that for some reason.

---

Edit: Yeah, that's a bug (see here). Corrected in the v.1.7.7beta scripts - hopefully available here soon. [HS]