Homecooked solution? [🇷 for BE/BA]

posted by ElMaestro  – Denmark, 2017-08-30 01:18 (2851 d 15:00 ago) – Posting: # 17753
Views: 8,847

Hi all,

since no standard function exists which helps me along, I wrote the stuff below.
The structure of the input file is like in the appendices for Schütz, Labes, Fuglsang's paper in the AAPS Journal.
I am sure the code police will tell me this isn't efficient :-D because the same could be achieved with fewer code lines or some apply family calls. I am sure that is correct. Does the code achieve its goal properly, though? That is a much more interesting question. Not at all throughly tested and certainly not validated.
The purpose is to allow some CI info where we don't ask for "T/R" but perhaps for "R/T" or for "B/D" or "Foo/Bar", whatever.




Extract.Hot.Info=function(TreatmentRatio, Z)
{
  i=regexpr("/", TreatmentRatio)
  T1=substring(TreatmentRatio, 1, i-1)
  T2=substring(TreatmentRatio, i+1, nchar(TreatmentRatio)) ##heck!
  for (i in 1:nrow(Z))
   {
    if (as.character(Z[i,1]) == T1) E1=as.numeric(Z[i,2])
    if (as.character(Z[i,1]) == T2) E2=as.numeric(Z[i,2])
   }
  return(E1-E2)
}

Get.CI.BE=function(FName,  alpha, TreatmentRatio, Stage.Term=F)
{
  A=read.table(FName, header=T) ##Note Var is untransformed
  if (Stage.Term)
   M=lm(log(Var)~factor(Trt)+factor(Seq)+factor(Per)+factor(Subj)+factor(Stg), data=A)
  else
   M=lm(log(Var)~factor(Trt)+factor(Seq)+factor(Per)+factor(Subj), data=A)
 
  df=M$df.residual
  Sc=summary(M)$coefficients
  Trt.Lvls=levels(A$Trt)
  Eff.Lvls=rep(0, length(Trt.Lvls)) ##pseudo effect levels. Usefyl for for extr. of differences
 
  for (i in 1:nrow(Sc))
  for (j in 1:length(Trt.Lvls))
   {
    x=paste("factor(Trt)", as.character(Trt.Lvls[j]), sep="")
    if (rownames(Sc)[i]==x) { Eff.Lvls[j]=Sc[i,1]; SEdiffs=Sc[i,2]}
   }
  Z=data.frame(Trt.Lvls, Eff.Lvls)
  Diff.LogScale=Extract.Hot.Info(TreatmentRatio,Z)
  logL=Diff.LogScale-SEdiffs*qt(p=1-alpha, df=df)
  logU=Diff.LogScale+SEdiffs*qt(p=1-alpha, df=df)
  L=list(PE=exp(Diff.LogScale), Lower=exp(logL), Upper=exp(logU))
  return(L)
}
## call it like Get.CI.BE("Cojones.txt",  0.05, "T/R", F)


Pass or fail!
ElMaestro

Complete thread:

UA Flag
Activity
 Admin contact
23,424 posts in 4,927 threads, 1,677 registered users;
33 visitors (0 registered, 33 guests [including 5 identified bots]).
Forum time: 16:18 CEST (Europe/Vienna)

Complex, statistically improbable things are by their nature
more difficult to explain than
simple, statistically probable things.    Richard Dawkins

The Bioequivalence and Bioavailability Forum is hosted by
BEBAC Ing. Helmut Schütz
HTML5