Still not exactly sure [RSABE / ABEL]

posted by ElMaestro  – Denmark, 2013-04-09 12:12 (4418 d 18:00 ago) – Posting: # 10380
Views: 5,849

Thank you, Helmut and d_labes, for your input.

I am not planning to via the in-built confidence interval function because I can port it directly to C when I know exactly what's going on.

Helmut, I downloaded your spreadsheet, saved the dataset as csv with tab delimiter and adjusted the code to fit your format.

Here's some code:
A=read.table("EMAset1c.csv", header=T)

Model1=lm(log(A$Data)~0+factor(A$Formulation)+
factor(A$Period)+factor(A$Sequence)+factor(A$Subject))

lnPE=(coef(Model1)[2]-coef(Model1)[1])


B=subset(A, Formulation=='R')
Model2=lm(log(B$Data)~factor(B$Period)+factor(B$Sequence)+factor(B$Subject))
s=summary(Model2)$sigma
CV=sqrt(exp(s*s)-1.0)
DF=Model2$df.residual
CritT=qt(1-0.05, DF)

nTRTR=0
nRTRT=0

for (i in min(B$Subject):max(B$Subject))
if (i %in% B$Subject) ##if the subject exists
{
  C=subset(B, Subject==i)
  r=subset(C, Formulation=='R')
  if(length(r$Formulation)==2) ##if the subject ate pocus two times.
   {
    if (C$Sequence[1]=='TRTR') nTRTR=nTRTR+1
    if (C$Sequence[1]=='RTRT') nRTRT=nRTRT+1
   }
} ## there's prolly a more efficient way of doing it.

lnU=lnPE+s*CritT*sqrt(0.25*(1/nTRTR+1/nRTRT))
lnL=lnPE-s*CritT*sqrt(0.25*(1/nTRTR+1/nRTRT))
cat("PE=", exp(lnPE),"\n")
cat("Upper Limit=", exp(lnU),"\n")
cat("Lower Limit=", exp(lnL),"\n")
cat("CV=", CV, "\n")
cat("nTRTR , nRTRT:", nTRTR, ",", nRTRT, "\n")
 


I am still not really there :-D

Pass or fail!
ElMaestro

Complete thread:

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

Don’t undertake a project
unless it’s manifestly important
and nearly impossible.    Edwin H. Land

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