Validation

Validation provided with 3 sections:

• REML validation for public datasets with Metida & SPSS
• Parameters validation for public datasets Metida & SPSS & MixedModels
• Validation with bioequivalence datasets with Metida & SPSS

To run validation:

using Metida; include(joinpath(dirname(pathof(Metida)), "..", "test", "validation.jl"))

Section 1: REML validation for public datasets Metida & SPSS

REML result table

sleepstudy.csv

Model 1

lmm = LMM(@formula(Reaction~Days), df;
  random = VarEffect(@covstr(1|Subject), SI),
  )
  fit!(lmm)

SPSS:

MIXED Reaction BY Days
  /CRITERIA=DFMETHOD(SATTERTHWAITE) CIN(95) MXITER(100) MXSTEP(10) SCORING(1)
    SINGULAR(0.000000000001) HCONVERGE(0, ABSOLUTE) LCONVERGE(0, ABSOLUTE) PCONVERGE(0.000001, ABSOLUTE)
  /FIXED=Days | SSTYPE(3)
  /METHOD=REML
  /RANDOM=INTERCEPT | SUBJECT(Subject) COVTYPE(ID).

Model 2

lmm = LMM(@formula(Reaction~1), df;
  random = VarEffect(Metida.@covstr(Days|Subject), CS),
  )
  fit!(lmm)

SPSS:

MIXED Reaction BY Days
  /CRITERIA=DFMETHOD(SATTERTHWAITE) CIN(95) MXITER(100) MXSTEP(10) SCORING(1)
    SINGULAR(0.000000000001) HCONVERGE(0, ABSOLUTE) LCONVERGE(0, ABSOLUTE) PCONVERGE(0.000001, ABSOLUTE)
  /FIXED=INTERCEPT | SSTYPE(3)
  /METHOD=REML
  /RANDOM=Days | SUBJECT(Subject) COVTYPE(CS).

Model 3

lmm = LMM(@formula(Reaction~1), df;
  random = VarEffect(@covstr(Days|Subject), CSH)
  )
  fit!(lmm)

SPSS:

MIXED Reaction BY Days
  /CRITERIA=DFMETHOD(SATTERTHWAITE) CIN(95) MXITER(100) MXSTEP(10) SCORING(1)
    SINGULAR(0.000000000001) HCONVERGE(0, ABSOLUTE) LCONVERGE(0, ABSOLUTE) PCONVERGE(0.000001, ABSOLUTE)
  /FIXED=INTERCEPT | SSTYPE(3)
  /METHOD=REML
  /RANDOM=Days | SUBJECT(Subject) COVTYPE(CSH).

Model 4

lmm = LMM(@formula(Reaction~1), df;
  random = VarEffect(@covstr(Days|Subject), ARH)
  )
   fit!(lmm)

SPSS:

MIXED Reaction BY Days
  /CRITERIA=DFMETHOD(SATTERTHWAITE) CIN(95) MXITER(100) MXSTEP(10) SCORING(1)
    SINGULAR(0.000000000001) HCONVERGE(0, ABSOLUTE) LCONVERGE(0, ABSOLUTE) PCONVERGE(0.000001, ABSOLUTE)
  /FIXED=INTERCEPT | SSTYPE(3)
  /METHOD=REML
  /RANDOM=Days | SUBJECT(Subject) COVTYPE(ARH1).

pastes.csv

Model 5

lmm =  LMM(@formula(strength~1), df;
random = [VarEffect(@covstr(1|batch), SI),  VarEffect(@covstr(1|batch & cask), SI)]
)
fit!(lmm)

SPSS:

MIXED strength
  /CRITERIA=DFMETHOD(SATTERTHWAITE) CIN(95) MXITER(100) MXSTEP(10) SCORING(1)
    SINGULAR(0.000000000001) HCONVERGE(0, ABSOLUTE) LCONVERGE(0, ABSOLUTE) PCONVERGE(0.000001, ABSOLUTE)
  /FIXED=| SSTYPE(3)
  /METHOD=REML
  /RANDOM=INTERCEPT | SUBJECT(batch) COVTYPE(ID)
  /RANDOM=INTERCEPT | SUBJECT(cask * batch) COVTYPE(ID).

Model 6

lmm =  LMM(@formula(strength~1), df;
random = VarEffect(Metida.@covstr(cask|batch),  ARMA),
)
fit!(lmm)

SPSS:

MIXED strength by cask
  /CRITERIA=DFMETHOD(SATTERTHWAITE) CIN(95) MXITER(100) MXSTEP(10) SCORING(1)
    SINGULAR(0.000000000001) HCONVERGE(0, ABSOLUTE) LCONVERGE(0, ABSOLUTE) PCONVERGE(0.000001, ABSOLUTE)
  /FIXED=| SSTYPE(3)
  /METHOD=REML
  /RANDOM=cask | SUBJECT(batch) COVTYPE(ARMA11).

penicillin.csv

Model 7

lmm =  LMM(@formula(diameter ~ 1), df;
random = [VarEffect(@covstr(1|plate), SI), VarEffect(@covstr(1|sample), SI)]
)
fit!(lmm)

SPSS:

MIXED diameter
  /CRITERIA=DFMETHOD(SATTERTHWAITE) CIN(95) MXITER(100) MXSTEP(10) SCORING(1)
    SINGULAR(0.000000000001) HCONVERGE(0, ABSOLUTE) LCONVERGE(0, ABSOLUTE) PCONVERGE(0.000001, ABSOLUTE)
  /FIXED=| SSTYPE(3)
  /METHOD=REML
  /RANDOM=INTERCEPT | SUBJECT(plate) COVTYPE(ID)
  /RANDOM=INTERCEPT | SUBJECT(sample) COVTYPE(ID).

RepeatedPulse.csv

Model 8

sort!(df, :Day)
lmm =  LMM(@formula(Pulse~1), df;
random =  VarEffect(Metida.@covstr(Time|Time),  SI),
repeated =  VarEffect(Metida.@covstr(Day|Time),  AR),
)
fit!(lmm)

SPSS:

MIXED Pulse BY Day Time
  /CRITERIA=DFMETHOD(SATTERTHWAITE) CIN(95) MXITER(100) MXSTEP(10) SCORING(1)
    SINGULAR(0.000000000001) HCONVERGE(0, ABSOLUTE) LCONVERGE(0, ABSOLUTE) PCONVERGE(0.000001, ABSOLUTE)
  /FIXED=| SSTYPE(3)
  /METHOD=REML
  /RANDOM=Time | SUBJECT(Time) COVTYPE(ID)
  /REPEATED = Day | SUBJCET(Time) COVTYPE(AR1).

Model 9

sort!(df, :Day)
lmm =  LMM(@formula(Pulse~1), df;
repeated = VarEffect(Metida.@covstr(Day|Time),  AR),
)
 fit!(lmm)

SPSS:

MIXED Pulse BY Day Time
  /CRITERIA=DFMETHOD(SATTERTHWAITE) CIN(95) MXITER(100) MXSTEP(10) SCORING(1)
    SINGULAR(0.000000000001) HCONVERGE(0, ABSOLUTE) LCONVERGE(0, ABSOLUTE) PCONVERGE(0.000001, ABSOLUTE)
  /FIXED=| SSTYPE(3)
  /METHOD=REML
  /REPEATED = Day | SUBJCET(Time) COVTYPE(AR1).

Model 10

sort!(df, :Day)
lmm =  LMM(@formula(Pulse~1), df;
random =  VarEffect(Metida.@covstr(Day|Time),  AR),
)
 fit!(lmm)

SPSS:

MIXED Pulse BY Day Time
  /CRITERIA=DFMETHOD(SATTERTHWAITE) CIN(95) MXITER(100) MXSTEP(10) SCORING(1)
    SINGULAR(0.000000000001) HCONVERGE(0, ABSOLUTE) LCONVERGE(0, ABSOLUTE) PCONVERGE(0.000001, ABSOLUTE)
  /FIXED=| SSTYPE(3)
  /METHOD=REML
  /RANDOM=Day | SUBJECT(Time) COVTYPE(AR1).

Section 2: Parameters validation for public datasets Metida & SPSS & MixedModels

Model 7

Metida:

lmm = LMM(@formula(diameter ~ 1), df;
random = [VarEffect(@covstr(1|plate), SI), VarEffect(@covstr(1|sample), SI)]
)
fit!(lmm)

MixedModels:

fm = @formula(diameter ~ 1 + (1|plate) + (1|sample))
mm  = fit(MixedModel, fm, df, REML=true)

SPSS:

MIXED diameter
  /CRITERIA=DFMETHOD(SATTERTHWAITE) CIN(95) MXITER(100) MXSTEP(10) SCORING(1)
    SINGULAR(0.000000000001) HCONVERGE(0, ABSOLUTE) LCONVERGE(0, ABSOLUTE) PCONVERGE(0.000001, ABSOLUTE)
  /FIXED=| SSTYPE(3)
  /METHOD=REML
  /EMMEANS=TABLES(OVERALL)
  /RANDOM=INTERCEPT | SUBJECT(plate) COVTYPE(ID)
  /RANDOM=INTERCEPT | SUBJECT(sample) COVTYPE(ID).

Section 3: Validation with bioequivalence datasets with Metida & SPSS

Model BE-B

lmm =  LMM(@formula(lnpk~sequence+period+treatment), dfrds;
random =  VarEffect(Metida.@covstr(1|subject),  SI),
)
fit!(lmm)

Model BE-C

lmm =  LMM(@formula(lnpk~sequence+period+treatment), dfrds;
    random =  VarEffect(Metida.@covstr(treatment|subject),  CSH),
    repeated =  VarEffect(Metida.@covstr(treatment|subject),  DIAG),
    )
     fit!(lmm)

Typical SPSS code

MIXED lnpk BY period sequence treatment subject
  /CRITERIA=DFMETHOD(SATTERTHWAITE) CIN(95) MXITER(100) MXSTEP(10) SCORING(1)
    SINGULAR(0.000000000001) HCONVERGE(0, ABSOLUTE) LCONVERGE(0, ABSOLUTE) PCONVERGE(0.000001, ABSOLUTE)
  /FIXED=period sequence treatment | SSTYPE(3)
  /METHOD=REML
  /RANDOM= subject(sequence) | COVTYPE(ID)
  /EMMEANS=TABLES(treatment) COMPARE REFCAT(FIRST) ADJ(LSD).

MIXED lnpk BY period treatment sequence subject
  /CRITERIA=CIN(90) MXITER(200) MXSTEP(20) SCORING(2) SINGULAR(0.000000000001) HCONVERGE(0,
    RELATIVE) LCONVERGE(0.0000000000001, RELATIVE) PCONVERGE(0, RELATIVE)
  /FIXED=period treatment sequence | SSTYPE(3)
  /METHOD=REML
  /RANDOM=treatment | SUBJECT(subject) COVTYPE(CSH)
  /REPEATED=treatment | SUBJECT(subject*period) COVTYPE(DIAG)
  /EMMEANS=TABLES(treatment) COMPARE REFCAT(FIRST) ADJ(LSD).

Results

using Metida; # hide
include(joinpath(dirname(pathof(Metida)), "..", "test", "validation_init.jl")); # hide
include(joinpath(dirname(pathof(Metida)), "..", "test", "validation_s3.jl"))

Full SPSS code provided in validation folder (here).

SPSS output in DOCX format.

Validation dataset available here, 122482020427MOESM2ESM.xls.

Validation report

Validation and report can be done on local machine with Weave.jl and Pandoc.

using Weave, Metida
weave(joinpath(dirname(pathof(Metida)), "..", "test", "validation_report.jmd");
doctype = "pandoc2pdf",
out_path = :pwd,
pandoc_options=["--toc", "-V colorlinks=true" , "-V linkcolor=blue", "-V urlcolor=red", "-V toccolor=gray"])

Report will be saved in Julia working directory. For your own purpose you can edit validation_report.jmd template.