Abbreviations

Scope

We try to reproduce the study presented in (Rosner,2015, page 713-721) in the data tennis2.dat: A Cross-Over study of N1=N2=88 of a Motrin treatment,an Ibuprofen-like against a Placebo.Usually Xover design are used for bioequivalence assesment. No direct references was found and nothing left in the read me file.Because a lack a data reference protocole annexed to the datafile,we have to take time to understand the how the variables were measured with different exercises settings(see summary and Data).

Data

The data were found on the web but were published with Fundamentals of Biostatistics’, 7th Edition, by Bernard Rosner, ISBN-10 0-538-73349-7; ISBN-13: 978-0-538-73349-6. We based our theory analysis on the 8th.edition. Changes has occured in the dataset as Rosner where it is stated that two drop-out occurs during the 8th edition.

The raw data are wide coded as treatment appeared as 2 columns for each ID(1 row|Id).

We choose the long format coding for LMM modelling as grouping factor is id and Mor and placebo columns merge as one,that is the grouping factor as two Yi,j=1,2,k.

load("TENNIS2.DAT.rdata")
summary(tennis2)
##        Id             age             sex           drg_ord       painmx_2    
##  Min.   :701.0   Min.   :24.00   Min.   :1.000   Min.   :1.0   Min.   :1.000  
##  1st Qu.:740.8   1st Qu.:39.00   1st Qu.:1.000   1st Qu.:1.0   1st Qu.:2.000  
##  Median :775.5   Median :43.50   Median :1.000   Median :1.5   Median :3.000  
##  Mean   :774.8   Mean   :43.56   Mean   :1.568   Mean   :1.5   Mean   :3.125  
##  3rd Qu.:810.2   3rd Qu.:48.25   3rd Qu.:2.000   3rd Qu.:2.0   3rd Qu.:4.000  
##  Max.   :848.0   Max.   :64.00   Max.   :9.000   Max.   :2.0   Max.   :9.000  
##     pain12_2        painav_2        painov_2        painmx_3    
##  Min.   :1.000   Min.   :1.000   Min.   :1.000   Min.   :1.000  
##  1st Qu.:2.000   1st Qu.:2.000   1st Qu.:2.000   1st Qu.:2.000  
##  Median :3.000   Median :3.000   Median :3.000   Median :2.000  
##  Mean   :3.261   Mean   :3.386   Mean   :3.409   Mean   :2.932  
##  3rd Qu.:4.000   3rd Qu.:4.250   3rd Qu.:5.000   3rd Qu.:4.000  
##  Max.   :9.000   Max.   :9.000   Max.   :9.000   Max.   :9.000  
##     pain12_3        painav_3        painov_3       painmx_4        pain12_4    
##  Min.   :1.000   Min.   :1.000   Min.   :1.00   Min.   :0.000   Min.   :1.000  
##  1st Qu.:1.000   1st Qu.:2.000   1st Qu.:1.00   1st Qu.:3.000   1st Qu.:2.750  
##  Median :2.000   Median :2.000   Median :2.00   Median :4.000   Median :4.000  
##  Mean   :2.148   Mean   :2.841   Mean   :2.33   Mean   :4.114   Mean   :4.057  
##  3rd Qu.:2.000   3rd Qu.:4.000   3rd Qu.:2.00   3rd Qu.:5.000   3rd Qu.:5.000  
##  Max.   :9.000   Max.   :9.000   Max.   :9.00   Max.   :9.000   Max.   :9.000  
##     painav_4        painov_4    
##  Min.   :1.000   Min.   :1.000  
##  1st Qu.:3.000   1st Qu.:3.000  
##  Median :4.000   Median :4.500  
##  Mean   :4.136   Mean   :4.216  
##  3rd Qu.:5.000   3rd Qu.:5.000  
##  Max.   :9.000   Max.   :9.000
tennis2$Id=factor(tennis2$Id)
tennis2$drg_ord=factor(tennis2$drg_ord)#other factor will be ignored yet
str(tennis2)
## 'data.frame':    88 obs. of  16 variables:
##  $ Id      : Factor w/ 88 levels "701","702","704",..: 1 14 16 17 18 19 21 22 23 24 ...
##  $ age     : num  42 45 43 48 56 44 31 49 44 38 ...
##  $ sex     : num  2 1 1 1 2 2 2 2 2 2 ...
##  $ drg_ord : Factor w/ 2 levels "1","2": 1 2 1 2 1 1 2 2 2 1 ...
##  $ painmx_2: num  5 2 4 1 5 5 3 2 2 2 ...
##  $ pain12_2: num  5 2 3 1 1 5 3 2 4 2 ...
##  $ painav_2: num  5 2 4 1 5 5 3 3 4 3 ...
##  $ painov_2: num  5 2 3 2 5 5 3 2 4 2 ...
##  $ painmx_3: num  2 4 3 2 6 5 3 3 3 2 ...
##  $ pain12_3: num  1 4 2 2 5 5 2 3 2 2 ...
##  $ painav_3: num  2 4 3 2 5 5 3 2 3 2 ...
##  $ painov_3: num  1 4 2 2 5 5 2 2 2 2 ...
##  $ painmx_4: num  3 6 2 3 2 5 4 5 5 3 ...
##  $ pain12_4: num  3 5 2 3 1 5 4 5 5 3 ...
##  $ painav_4: num  3 5 2 3 1 6 4 5 5 3 ...
##  $ painov_4: num  3 6 2 3 2 6 4 5 5 3 ...

The problem we are facing is between ed 7th & 8th and without DOE description (readme variables),is as follows; we don’t know which variables are effectively take for the analysis ( diff-baseline, max, min ). We identified them according to their mean value published in book 8 th ed page 716:table 13:32 This identified the possible (?!) variables painov-pain_av4:

#identifying trt and placbo form analysis in the book
#makin a new df and names variable
tennis2$placebo = tennis2$painov_2##identified placebo
tennis2$Mortin = tennis2$painav_4#trt Mortin
ten= tennis2[,c("Id", "sex", "drg_ord", "placebo", "Mortin")]
ten$difftrt=ten$Mortin-ten$placebo#make the difference needed for testing trt
summary(ten)
##        Id          sex        drg_ord    placebo          Mortin     
##  701    : 1   Min.   :1.000   1:44    Min.   :1.000   Min.   :1.000  
##  702    : 1   1st Qu.:1.000   2:44    1st Qu.:2.000   1st Qu.:3.000  
##  704    : 1   Median :1.000           Median :3.000   Median :4.000  
##  705    : 1   Mean   :1.568           Mean   :3.409   Mean   :4.136  
##  706    : 1   3rd Qu.:2.000           3rd Qu.:5.000   3rd Qu.:5.000  
##  707    : 1   Max.   :9.000           Max.   :9.000   Max.   :9.000  
##  (Other):82                                                          
##     difftrt       
##  Min.   :-4.0000  
##  1st Qu.: 0.0000  
##  Median : 1.0000  
##  Mean   : 0.7273  
##  3rd Qu.: 2.0000  
##  Max.   : 7.0000  
## 
sum(is.na(ten))
## [1] 0
#no NA that is 88 Same ID for/in each period=BALANCED no D.OUT
xyplot(Mortin+placebo~Id|drg_ord,data=ten,type=c("p"),auto.key = TRUE,xlab = list(label="ID",cex=0.6))

Note on the lattice plot that score dont enable clear cut between treatment, this visually. Note that drop out in Xover trials are catastrophic for the study: Loosing 1 Id is like twice a // Drop out should not be linked with trt (or placebo): In general,if Drop-out is linked with improvement so the analysis of the effect can be adversely and reversly biased.Fortunately drop-out are more likely in second order Xover a 2X3,2x4.

Note also that trt is embeded in within patients data and the between patient variability mask the true effect (This can be assessed by a random effect variance [l]) In Xover each Id act as its own control hence limiting effect of possible cofoudners, and some researcher advocated “it removes the between-patient” variability", which we found partially agree: It just reveal a part of it: As it still exist and to asses it one might use LMM models presented at end and compare: In that instance the covariances of two measures is the variance at the ID levels, that in a regression settings.

Some estimators are presented in first instance:

MM=ten %>% 
    group_by(drg_ord) %>% 
  summarise(Mort=mean(Mortin),Plac=mean(placebo),Difference=mean(difftrt))
kable(MM,digits=3)
drg_ord Mort Plac Difference
1 3.864 3.977 -0.114
2 4.409 2.841 1.568 
SD= ten %>% 
   group_by(drg_ord) %>% 
  summarise(Mortsd=sd(Mortin),Placsd=sd(placebo),Differencesd=sd(difftrt))
kable(SD,digits =3)
drg_ord Mortsd Placsd Differencesd
1 2.030 1.406 2.003
2 1.716 1.140 1.717

Data analysis

Now we need to analize the data (treatment ,carryover, period).Recall that in DOE treatment is a mean of an orthogonal design (i.e High and Low).In biostatistics we’ d better analyse the difference between two measurements ,either in the period, in the group or in the treatment: this is the difficult and tricky part of the analysis. Recall that,

Lets check the Normality assumptions:

par(mfrow=c(1,3))
hist(ten$Mortin,main="MORTIN",col=8)
hist(ten$placebo,main="Placebo",col=8)
hist(as.numeric(ten$Mortin-ten$placebo),col=6,main="diffrence M-p")

#Although talking about the choosen scale one note the diffenece trt-placebo is N Distributed.
#Recall that T test is resistant to midl N departure!

It is not clear in this data how the author decipher between Periods/Group: Might be misleading in such study. As data is code as wide format one should found a Grand mean difference (of treatment) value of 0.723

The overall mean (no design) is higher for Mortin resulting in a positive difference of 0.723 (Score under mortin is improved in OVERALL : That is it is score about 3/4 higher under Mortin (less pain).

Recall that taking a difference remove partial of 2*COV(XY) [Paired t test vs t test!!] . Lets review some estimators before running the test of trt effect:

# DIfference is N Distributed
#some estimators before inference
KG=ten %>% 
  group_by(drg_ord) %>% 
  summarize(Mgroupdiff=mean(difftrt),sdgroupdif=sd(difftrt),vargroup=var(difftrt)) 
kable(KG)
drg_ord Mgroupdiff sdgroupdif vargroup
1 -0.1136364 2.002509 4.010042
2 1.5681818 1.717187 2.948732 
mean(ten$difftrt)#check ok
## [1] 0.7272727
#Cross over particular:As the DIFF RV is also a Grand mean over the 2Group/Period
#we score slighty highe sd with this data of 7ed.
#The question of equal variance is left to the reader

t.test(ten$difftrt)
## 
##  One Sample t-test
## 
## data:  ten$difftrt
## t = 3.3471, df = 87, p-value = 0.001208
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  0.2953991 1.1591464
## sample estimates:
## mean of x 
## 0.7272727
# DF seems corract although we prefer 86 (n1+n2-2)
#because we use two mean in the difference debated in Stat.!

#NOTE:now t test with paired = TRUE should be done not on the difference but on treatment itself:
t.test(ten$Mortin,ten$placebo,paired=TRUE,var.equal = F)
## 
##  Paired t-test
## 
## data:  ten$Mortin and ten$placebo
## t = 3.3471, df = 87, p-value = 0.001208
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  0.2953991 1.1591464
## sample estimates:
## mean of the differences 
##               0.7272727
#same result

We replicate (Z,p-val,diff) about the same result as in text book (Confidence int. slightly different but precision not affected).

Same T test calculation but “By hand” Like in the formula 13:51:

sd(ten$difftrt)##overall difference sd
## [1] 2.038296
#dont make a mean of sd in 1 and 2 this gives you wrong variance estimates Use the raw datapoint yij
 
#VAR of the two 1/2 become 1/4 as coeff. X
#polling the sd
Sdpool=(43*4.01+43*2.95)/86
SEDiff=sqrt(1/4*Sdpool*(1/44+1/44))
T_Zstat=0.7272/SEDiff
T_Zstat#about the same of t test (reader should choose variance qual of not pooling of not...)
## [1] 3.656837
qt(0.975,86)#5% level with N-2(parameters)
## [1] 1.987934
CONFINT_INTERVAL=c(LowerB=0.7272-1.988*SEDiff,diff=0.7272,HighB=0.7272+1.988*SEDiff)
#in accordance with t test
CONFINT_INTERVAL
##    LowerB      diff     HighB 
## 0.3318655 0.7272000 1.1225345

Some un-significant variation occurs in automated T test and by hand probably from continuity/and or pooling /df estimates., but significance isn’t affected, that is R reproduce with precision the result of the author.

Testing the carry over effect

COE are persistents (Long Last) effect(s) occuring on period 2 after Wash-out as biophysical effect of inherited treatment.The problem in RCT Xover is that Wash-out period are estimated form Bio-pharmco-kinetics protocols: They might not reflect the entire description of treatment persistence: Luckily enough if it persistent, the better it is for our patient,isnt’it?. Further discussion see (Grizzle,1965;Senn,2002)

  • COE effects:

Lets see if summation of trt+placebo in P1 is different form summation mean trt+placebo in P2 (Note: that is not one way to consider COE!)

ten$Mpmean=(ten$Mortin+ten$placebo)
hist(ten$Mpmean,main="mean of MOR+PLAC",col=6)#Ok N

shapiro.test(ten$Mpmean)#ok for t test visually
## 
##  Shapiro-Wilk normality test
## 
## data:  ten$Mpmean
## W = 0.94258, p-value = 0.0007015
#The carry over is the mean effect in P1 compared to mean effect P2 (defined M1+p1)/2 in P1 (M2+p2/2 in P2)
t.test(ten$Mpmean~ten$drg_ord) # testing the two mean summation in periods
## 
##  Welch Two Sample t-test
## 
## data:  ten$Mpmean by ten$drg_ord
## t = 1.058, df = 82.925, p-value = 0.2931
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  -0.5199204  1.7017386
## sample estimates:
## mean in group 1 mean in group 2 
##        7.840909        7.250000
#No carry effect is noticeable at first sight (but Placebo effect might not be discarded)
#testing the mean of summation might not reflect the true COE because we neglect the PLACEBO effect in P1

The carry-over effect seem not prevalent but caution in the placebo effect in P1 sometimes forgotten in the analysis.

Note you can also test the sum of M1p1 vs M2p2 with a Non parametric test.

wilcox.test(ten$Mpmean~ten$drg_ord)
## 
##  Wilcoxon rank sum test with continuity correction
## 
## data:  ten$Mpmean by ten$drg_ord
## W = 1066, p-value = 0.4114
## alternative hypothesis: true location shift is not equal to 0

Note that paring in the carryover might not be fully tested or partially: we dont consider any justifying idea. Note that adding 2 RV (sum) as one if data are correlated that we add **2COV(M-P)** see below-,hence correlation is partly accounted in testing, that is COE:

sd(ten$Mpmean)##Sum
## [1] 2.62135
sd(ten$difftrt)#diff
## [1] 2.038296

More on carry over estimate

Now we have assesed the trt effect over Placebo we might question about the magnitude in Period1 (P1) and in Period2 (P2): The following scheme resume the table with P as periods,p as placebo,M as trt: M1_P1-p_P1 / P2_P2-M2_P2 This can be seen as a two way table like. - If there is a carry over-effect (COE) it is manifested for M in period2 in p2 (because ID received in 1 M1): Taht is the mean of M1+p1 should be lower that M2p2. - Conversely if there is a Placebo effect It will occurs mainly in p_P1, reducing they effects with Mortin in P1 (this effect is overlooked).

Lets test first between the Periods the effect: These two T-tests should give you insight (direction , magnitude) of assimetric trt effect.Lets test it:

#mean in P1 vs P2 by test
MM[1,c(2,3)]
## # A tibble: 1 x 2
##    Mort  Plac
##   <dbl> <dbl>
## 1  3.86  3.98
meanP1=(3.86+3.98)/2
meanP1
## [1] 3.92
MM[2,c(2,3)]
## # A tibble: 1 x 2
##    Mort  Plac
##   <dbl> <dbl>
## 1  4.41  2.84
meanp2=(4.41+2.81)/2
meanp2
## [1] 3.61
#the mean is higher in 1 which is counter-intuitive of due to Placebo or Mortin effect of both
#Howver M1 and M2 are less different than p_p1-p_p2 similar
#Testing the effect of Mortin in P1-P2 and Placebo P1-P2
ten %>% 
  group_by(drg_ord) %>% 
  summarise(mean(Mortin))
## # A tibble: 2 x 2
##   drg_ord `mean(Mortin)`
##   <fct>            <dbl>
## 1 1                 3.86
## 2 2                 4.41
t.test(ten$Mortin~ten$drg_ord)#t test in two groups MOR
## 
##  Welch Two Sample t-test
## 
## data:  ten$Mortin by ten$drg_ord
## t = -1.3612, df = 83.684, p-value = 0.1771
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  -1.3423839  0.2514748
## sample estimates:
## mean in group 1 mean in group 2 
##        3.863636        4.409091
##MORTIN Mean is not significance between P1 P2
t.test(ten$placebo~ten$drg_ord)
## 
##  Welch Two Sample t-test
## 
## data:  ten$placebo by ten$drg_ord
## t = 4.1646, df = 82.482, p-value = 7.628e-05
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  0.5936069 1.6791204
## sample estimates:
## mean in group 1 mean in group 2 
##        3.977273        2.840909
#placebo is!

Some insight might also be take form p1-p2 M1-M2 and look for the slope difference, as indication of COE.

library(gplots)
par(mfrow=c(1,2))
plotmeans(ten$Mortin ~ ten$drg_ord,ylim=c(2.5,6))
segments(2.1, 3.8, 2.1, 4.4, col = 2)
text(2.2, 3.7,"MORT sequence effect", col = 2, cex = 0.8)
plotmeans(ten$placebo ~ ten$drg_ord,ylim=c(2.5,6))

Parallel plot of ID might also help to capture the phenomena direction and Magnitude (not shown).

According to SENN,2002 there might be another way to test the period effect or to “control for it” in the t.test and it is as follows : Testing the mean trt difference in each period respectively:Therefore the trt is nested within the period:

t.test(ten$difftrt~ten$drg_ord)
## 
##  Welch Two Sample t-test
## 
## data:  ten$difftrt by ten$drg_ord
## t = -4.229, df = 84.045, p-value = 5.945e-05
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
##  -2.4726537 -0.8909827
## sample estimates:
## mean in group 1 mean in group 2 
##      -0.1136364       1.5681818
#compare with no controlling for Period:
t.test(ten$difftrt)
## 
##  One Sample t-test
## 
## data:  ten$difftrt
## t = 3.3471, df = 87, p-value = 0.001208
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  0.2953991 1.1591464
## sample estimates:
## mean of x 
## 0.7272727

With this nesting t.test technich, Period account for about 1 Z_value score with 3df:We cannot assume significance but 1 Z difference is far enough to consider trend in the Period as influential too.

The placebo effect:

A placebo effect occur in reverse order as usually expected. It take place prevalently in period 1 when a group received placebo (but believing in trt) they feel any “artificial” effect: This is a real effect in the measurements not an artefact of the data! We also draw attention of the Nocebo effect particularly effective in psychological drugs.

The problem in this study that PLACEBO effect and CARRY OVER might confounded in the analysis of data and it is terribly difficult to noticed it (see PlotMeans plot).

However there is a catch: make a baseline before P1 & P2 any departure form these differences might be in magnitude the COE estimate & Pacebo effect.

Regression to the Mean (RTM)

Not accounted in this study, RTM occurs when several measurements are made on the same ID given a selected strata of population (i.e. DBP measurement according BDP> 100):I.e when having too high score score regress to the overall mean,similarly when too low value are measurements in first place it might be possible that next one settle close to the mean.It is like scoring a 6/6 again in maths is less likely that another student who score a 2/6: Each one are more likely to score close to the mean on the 2nd exams. Some Statisticians believes it is an artefact of Bivariate Normal DIstribution; We argue that physical and biological process under action respond to the Normal law which tell that probability of scoring mean is higher that an extreme value. RTM effect can be visually observed when asking parallel plot. Note that RTM is also confounded with healing if healing is going closer to the population Mean. For further discussion See:(Greenland et al. 2016) and https://twitter.com/mudryjm/status/1394255580798365698?s=20

Power

It has been established that // study required 4 times more sample size:Let’s see:

#POWER T TEST 

#PAIRING Power 0.8,crit=0.05,balanced, Effect size at mean diff.
library(samplesize)
n.ttest(
  0.8,
  0.05,
  0.722,
  sd1 = 1.7,
  design = "paired",
  fraction = "balanced",
  variance = "equal"
)
## $`Total sample size`
## [1] 46
#REGULAR POWER T TEST
n.ttest(
  0.8,
  0.05,
  0.722,
  sd1 = 1.7,
  design = "unpaired",
  fraction = "balanced",
  variance = "equal"
)
## $`Total sample size`
## [1] 178
## 
## $`Sample size group 1`
## [1] 89
## 
## $`Sample size group 2`
## [1] 89

LMM (See Galecki,2013)

#making format as long with trt variables (1,2)

names(ten)
## [1] "Id"      "sex"     "drg_ord" "placebo" "Mortin"  "difftrt" "Mpmean"
str(ten)
## 'data.frame':    88 obs. of  7 variables:
##  $ Id     : Factor w/ 88 levels "701","702","704",..: 1 14 16 17 18 19 21 22 23 24 ...
##  $ sex    : num  2 1 1 1 2 2 2 2 2 2 ...
##  $ drg_ord: Factor w/ 2 levels "1","2": 1 2 1 2 1 1 2 2 2 1 ...
##  $ placebo: num  5 2 3 2 5 5 3 2 4 2 ...
##  $ Mortin : num  3 5 2 3 1 6 4 5 5 3 ...
##  $ difftrt: num  -2 3 -1 1 -4 1 1 3 1 1 ...
##  $ Mpmean : num  8 7 5 5 6 11 7 7 9 5 ...
attach(ten)
idvar=c("Id","drg_ord")
mes.var=c("Mortin","placebo")
longten=melt(ten,id.var=idvar,measure.var=mes.var)
summary(longten)#check doubled rows as 2MEasures/id grouping facto
##        Id      drg_ord    variable      value      
##  701    :  2   1:88    Mortin :88   Min.   :1.000  
##  702    :  2   2:88    placebo:88   1st Qu.:2.000  
##  704    :  2                        Median :4.000  
##  705    :  2                        Mean   :3.773  
##  706    :  2                        3rd Qu.:5.000  
##  707    :  2                        Max.   :9.000  
##  (Other):164
modlme0=lme(value~variable,random=~1|Id,data=longten)
modlme0
## Linear mixed-effects model fit by REML
##   Data: longten 
##   Log-restricted-likelihood: -336.8637
##   Fixed: value ~ variable 
##     (Intercept) variableplacebo 
##       4.1363636      -0.7272727 
## 
## Random effects:
##  Formula: ~1 | Id
##         (Intercept) Residual
## StdDev:   0.8241395 1.441293
## 
## Number of Observations: 176
## Number of Groups: 88
modlme1=lme(value~variable+drg_ord,random=~1|Id,data=longten)
modlme2=lme(value~variable,random=~drg_ord|Id,data=longten)
modlme3=lme(value~variable+variable*drg_ord,random=~1|Id,data=longten)
modlme4=lme(value~variable+variable*drg_ord,random=~drg_ord|Id,data=longten)

summary(modlme1)#recoreing the effect diff trt =contrast trt in regression R
## Linear mixed-effects model fit by REML
##  Data: longten 
##        AIC     BIC    logLik
##   683.3215 699.088 -336.6608
## 
## Random effects:
##  Formula: ~1 | Id
##         (Intercept) Residual
## StdDev:    0.822709 1.441293
## 
## Fixed effects: value ~ variable + drg_ord 
##                     Value Std.Error DF   t-value p-value
## (Intercept)      4.284091 0.2253706 87 19.009094  0.0000
## variableplacebo -0.727273 0.2172831 87 -3.347121  0.0012
## drg_ord2        -0.295455 0.2792451 86 -1.058047  0.2930
##  Correlation: 
##                 (Intr) vrblpl
## variableplacebo -0.482       
## drg_ord2        -0.620  0.000
## 
## Standardized Within-Group Residuals:
##         Min          Q1         Med          Q3         Max 
## -2.02660306 -0.69450125  0.04210248  0.63920976  2.96371992 
## 
## Number of Observations: 176
## Number of Groups: 88
modlme1
## Linear mixed-effects model fit by REML
##   Data: longten 
##   Log-restricted-likelihood: -336.6608
##   Fixed: value ~ variable + drg_ord 
##     (Intercept) variableplacebo        drg_ord2 
##       4.2840909      -0.7272727      -0.2954545 
## 
## Random effects:
##  Formula: ~1 | Id
##         (Intercept) Residual
## StdDev:    0.822709 1.441293
## 
## Number of Observations: 176
## Number of Groups: 88
head(coef(modlme1),5)#note that if you want ref Trt Mortin set levels of teh factor
##     (Intercept) variableplacebo   drg_ord2
## 701    4.315475      -0.7272727 -0.2954545
## 702    3.526382      -0.7272727 -0.2954545
## 704    5.499116      -0.7272727 -0.2954545
## 705    4.907296      -0.7272727 -0.2954545
## 706    3.642952      -0.7272727 -0.2954545
getVarCov(modlme1)
## Random effects variance covariance matrix
##             (Intercept)
## (Intercept)     0.67685
##   Standard Deviations: 0.82271
modlme1$sigma
## [1] 1.441293
ICC=0.68*0.68/(0.68*0.68+0.83*0.83)
ICC#part of the between patient 40% of total variability
## [1] 0.4016329
#something occuring within the patieant that is within the treatment

var(longten$value)#variance recovered by model
## [1] 2.873766
#the slope is trt difference

##TESTING THE INTERCATIOSN TRT-PERIOD
#

#considering mod1 simple with 3 random effect subject plus intercations FIXEF
modlme3
## Linear mixed-effects model fit by REML
##   Data: longten 
##   Log-restricted-likelihood: -328.4482
##   Fixed: value ~ variable + variable * drg_ord 
##              (Intercept)          variableplacebo                 drg_ord2 
##                3.8636364                0.1136364                0.5454545 
## variableplacebo:drg_ord2 
##               -1.6818182 
## 
## Random effects:
##  Formula: ~1 | Id
##         (Intercept) Residual
## StdDev:   0.9196021 1.318974
## 
## Number of Observations: 176
## Number of Groups: 88
#the coefficient and mean has change P1 Placebo has positive coefficient in P1 compare to Mortin which we failed to give an explanation
#however with the strong NEG SIGN of intercations gives a dreacrase of Placebo mean in P2 which is conform to data
tapply(longten$value,longten$drg_ord,mean)
##        1        2 
## 3.920455 3.625000
anova(modlme1)#fixef same as tteat due to orthogonality design
##             numDF denDF  F-value p-value
## (Intercept)     1    87 730.1289  <.0001
## variable        1    87  11.2032  0.0012
## drg_ord         1    86   1.1195  0.2930
library(car)
Anova(modlme3,type="II")
## Analysis of Deviance Table (Type II tests)
## 
## Response: value
##                    Chisq Df Pr(>Chisq)    
## variable         13.3775  1  0.0002547 ***
## drg_ord           1.1195  1  0.2900341    
## variable:drg_ord 17.8846  1  2.347e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova(modlme3,type="sequential")
##                  numDF denDF  F-value p-value
## (Intercept)          1    86 730.1283  <.0001
## variable             1    86  13.3775  0.0004
## drg_ord              1    86   1.1195  0.2930
## variable:drg_ord     1    86  17.8846  0.0001
anova(modlme1,modlme3,type="Chisq")
##         Model df      AIC      BIC    logLik   Test  L.Ratio p-value
## modlme1     1  5 683.3215 699.0880 -336.6608                        
## modlme3     2  6 668.8965 687.7814 -328.4482 1 vs 2 16.42506   1e-04
1-pchisq(8,1)#pval of the 2 nested models retained the intercation (test lower tail=F
## [1] 0.004677735

Testing a random effect (RI_RIAS) model via LRT/Anova is not quite easy.Especially with Interactions (knowing the F/Chis approx distribution and or Type of ANOVA with Interactions) and two points measures in time! Our guidelines is : - Interaction is tested with type II and III (marginality and hierarchical principal): Sequential usually not appropriate with interactions but might be questionable with COE as it is confounded with effect?

  • LMM have an approximate (KEN/ROG/SAITH. DF approx) Chi Sq Distribution
  • LMM have nesting properties test with LRT.

All models proposed with that interactions gives evidence of COE.(even as random coefficient too) Now, How to explain the non significant T test?

Well, possibly that T test test only the mean of PERIOD while regression find the E(Y|X) with included “slope”= the interaction in such instance(?).In the T test, the Placebo effect cannot be tested as it masking the effect difference in P1 (?).

We claim that the interaction:Period must always been tested in order to locate any trend difference/or COE.

More over, the within variance of the model (1-ICC) is far two high for an orthogonal design giving evidence of too high variance left for treatment part: Another effect might be involved.

Note:Contrast SUM will be tested too in orthogonal design to see if any insight/diff possible.

REGRESSION HIGLIGHT

We highlight contradictory analysis of DOE with regression and especially with the interactions will looks an “awsome modeling” but we failed to expain the consistency of the coefficients: Really the interpretation of interactions-slope-regressino-lmm is now not fully understood and need more robust study data.

CONCLUSION

Testing the validity of either:no carry effect , no period effects, no trt interactions effects are not easy tasks in Xover Trials. All these are confounded / masking / with the treatment and result as bias The< all might have also difference/similar causes in common.

We recommend the technic of mean of Summation of P1M1/P2M2 of (SENN2002) to test the COE effects

We also recommend the use of plot to estimate the direction & magnitude of COE and Placebo effect.

Testing the Period effect is like “controlling” for period the main effect (Like a regression covariate,see (*SENN,2002*).

The LMM analysis report a unsually low ICC not capturing all between patients variability and let part of the variance in the within measurments: this might be attributed to the scale of scoring which is too greedy. We definitely do not recommend to use score or discrete ordered factors in such analysis at the cost of adversely affect the precision and the true estimate of effects.

BIBLIOGRAPHY

ROSNER, Bernard. Fundamentals of biostatistics. Cengage learning, 2015.

GRIZZLE, James E. “The two-period change-over design and its use in clinical trials.” Biometrics (1965): 467-480.

BONATE, Peter L. Analysis of pretest-posttest designs. CRC Press, 2000.

GAŁECKI, Andrzej, BURZYKOWSKI, Tomasz, GAŁECKI, Andrzej, et al. Linear mixed-effects model. Springer New York, 2013.

SENN, Stephen S. Cross-over trials in clinical research. John Wiley & Sons, 2002.

Greenland, Sander, Stephen J. Senn, Kenneth J. Rothman, John B. Carlin, Charles Poole, Steven N. Goodman, and Douglas G. Altman. 2016. “Statistical Tests, P Values, Confidence Intervals, and Power: A Guide to Misinterpretations.” European Journal of Epidemiology 31 (4): 337–50. https://doi.org/10.1007/s10654-016-0149-3.