Bioequivalence and Bioavailability Forum

lmer: Method B (PE catched for imbalanced dataset!!!) and Method C [🇷 for BE/BA]

posted by mittyri – Russia, 2016-11-05 18:38 (2718 d 22:43 ago) – Posting: # 16777
Views: 19,385

Dear Astea!

so many questions
I'll just put here my code as a solution. Yes, we cannot get SAS results, sorry...
Please let me know if after investigation you still have some questions

library(readxl)
library(lmerTest)
options(contrasts=c("contr.treatment","contr.poly"))
methodC <- function(Dataset){
  muddle.lmer <- (lmer(log(Data)~Sequence+Period+Formulation + (Formulation-1|Subject), data=Dataset))
  cat("\n", "lmer treats Sequence, Period and Formulation as factors, Random is /~Formulation -1|Subject/ (close to Method C), then", "\n")
  FormulationEffect.lmer <- summary(muddle.lmer)$coefficients["FormulationT","Estimate"]   print(FormulationEffect.lmer)
  cat("\n", "PE is", round(exp(as.numeric(FormulationEffect.lmer))*100, digits = 2),
      "\n", "SAS gives 115.66", "\n")
 
  cat("\n", "what about LSMeans?", "\n", "PHX gives for Method C  7.67042954073919 and 7.81589387050252", "\n", "and lmer gives", "\n")
  LSM.lmer <- lmerTest::lsmeans(muddle.lmer, "Formulation")
  print(LSM.lmer$lsmeans.table[c(1,2), c(3,4)])
 
  cat("\n", "\n", "what about CI? SAS gives 107.10-124.89", "\n")
  round(exp(confint(muddle.lmer, "FormulationT", level = 0.90))*100, 2)
}

methodB <- function(Dataset){
  muddle.lmer <- (lmer(log(Data)~Sequence+Period+Formulation + (1|Subject), data=Dataset))
  cat("\n", "lmer treats Sequence, Period and Formulation as factors, Random is /~1|Subject/ (close to Method B), then", "\n")
  FormulationEffect.lmer <- summary(muddle.lmer)$coefficients["FormulationT","Estimate"]   print(FormulationEffect.lmer)
  cat("\n", "PE is", round(exp(as.numeric(FormulationEffect.lmer))*100, digits = 2),
      "\n", "SAS gives 115.73", "\n")
 
  cat("\n", "what about LSMeans?", "\n", "PHX gives for Method B 7.67001367911898 and 7.81610190985527",  "\n", "and lmer gives", "\n")
  LSM.lmer <- lmerTest::lsmeans(muddle.lmer)
  print(LSM.lmer$lsmeans.table[c(7,8), c(3,4)])
 
  cat( "\n", "\n", "what about CI? SAS gives 107.17-124.97", "\n")
  round(exp(confint(muddle.lmer, "FormulationT", level = 0.90))*100, 2)
}

# downloading the reference file from BEBAC; that's the most convenient way to get the file
download.file('http://bebac.at/downloads/Validation Replicate Design EMA.xls', "Dataset.xls", cacheOK = FALSE, mode="wb")
# read it without Perl and Java; thanks, Hadley!
Dataset<-read_excel("Dataset.xls", sheet = 1)
# prepare dataset

Dataset$Formulation<-factor(Dataset$Formulation, levels = c("R", "T")) # for LSMeans in lmer
Dataset$Sequence<-factor(Dataset$Sequence, levels = c("TRTR", "RTRT")) # for LSMeans in lmer
Dataset$Period<-factor(Dataset$Period) # for LSMeans in lmer

methodC(Dataset)
methodB(Dataset)



Complete thread:

• Bear vs. Phoenix & SAS Helmut 2015-04-20 17:34 
  • R vs. Phoenix & SAS? yjlee168 2015-04-20 19:36
    • R vs. Phoenix & SAS? Helmut 2015-04-21 01:02
      • lme() does not work with all fixed effects yjlee168 2015-04-21 23:41
        • lme() does not work with all fixed effects Astea 2016-11-04 00:13
          • lmer: Method B (PE catched for imbalanced dataset!!!) and Method Cmittyri 2016-11-05 17:38
            • lmer: Method B (PE catched for imbalanced dataset!!!) and Method C Astea 2016-11-05 19:27
              • lmer: Method B is ready for scaling mittyri 2016-11-05 20:01
                • lmer: Method B is ready for scaling Astea 2016-11-06 11:50
                  • lmer: Method B is ready for scaling mittyri 2016-11-07 06:07
  • info for lsmeans yjlee168 2015-04-20 21:34
    • info for lsmeans Helmut 2015-04-21 01:15
      • once more about R and replicate designes Astea 2016-11-02 23:43
        • once more about R and replicate designes VStus 2016-11-06 11:34
          • Getting variance components from the lmer output StatR 2017-02-03 13:53
            • Getting variance components from the lmer output VStus 2017-02-03 15:47
              • Getting variance components from the lmer output StatR 2017-02-03 17:12
                • Getting variance components d_labes 2017-02-07 11:16
                  • Getting variance components StatR 2017-02-07 11:36
                  • Getting variance components StatR 2017-02-08 08:41
                    • Getting variance components d_labes 2017-02-08 10:13
                      • Getting variance components StatR 2017-02-08 10:19
                    • Data structure Helmut 2017-02-08 10:33
                      • Data structure StatR 2017-02-08 10:49