Installation

import Pkg; Pkg.add("Metida")

Simple example

Step 1: Load data

Load provided data with CSV and DataFrames:

using Metida, CSV, DataFrames, CategoricalArrays

df = CSV.File(joinpath(dirname(pathof(Metida)), "..", "test", "csv", "df0.csv")) |> DataFrame;
Note

Check that all categorical variables are categorical.

transform!(df, :subject => categorical, renamecols=false)
transform!(df, :period => categorical, renamecols=false)
transform!(df, :sequence => categorical, renamecols=false)
transform!(df, :formulation => categorical, renamecols=false)

Step 2: Make model

lmm = LMM(@formula(var~sequence+period+formulation), df;
random = VarEffect(@covstr(formulation|subject), CSH),
repeated = VarEffect(@covstr(formulation|subject), DIAG));
Linear Mixed Model: var ~ sequence + period + formulation
Random 1: 
    Model: formulation|subject
    Type: CSH (3), Subjects: 5
Repeated: 
    Model: formulation|subject
    Type: DIAG (2)
Blocks: 5, Maximum block size: 4
Not fitted.

Step 3: Fit

fit!(lmm)
Linear Mixed Model: var ~ sequence + period + formulation
Random 1: 
    Model: formulation|subject
    Type: CSH (3), Subjects: 5
Repeated: 
    Model: formulation|subject
    Type: DIAG (2)
Blocks: 5, Maximum block size: 4
Status: converged See warnings in log.
    -2 logREML: 10.0652    BIC: 23.9282

    Fixed-effects parameters:
───────────────────────────────────────────────────────
                     Coef.  Std. Error      z  Pr(>|z|)
───────────────────────────────────────────────────────
(Intercept)      1.57749     0.334543    4.72    <1e-05
sequence: 2     -0.170833    0.384381   -0.44    0.6567
period: 2        0.195984    0.117228    1.67    0.0946
period: 3        0.145014    0.109171    1.33    0.1841
period: 4        0.157363    0.117228    1.34    0.1795
formulation: 2  -0.0791667   0.0903709  -0.88    0.3810
───────────────────────────────────────────────────────
    Variance components:
    θ vector: [0.455584, 0.367656, 1.0, 0.143682, 0.205657]
  Random 1   σ² formulation: 1   var   0.207557
  Random 1   σ² formulation: 2   var   0.135171
  Random 1   ρ                   rho   1.0
  Residual   σ² formulation: 1   var   0.0206445
  Residual   σ² formulation: 2   var   0.0422948
Check warnings and errors in log.
lmm.log
4-element Vector{Metida.LMMLogMsg}:
   INFO  : Initial θ: [0.21085850593184408, 0.21085850593184408, 0.0001, 0.21085850593184408, 0.21085850593184408]

   INFO  : Resulting θ: [0.45558429062258543, 0.3676561978014235, 0.9999999999971845, 0.14368242865088793, 0.20565739028506563]; 26 iterations.

   INFO  : Model fitted.

   WARN  : det(G) of random effect 1 is less 1e-08.
Type III Tests of Fixed Effects
Warning

Experimental

typeiii(lmm)
 ------------- ---------- ----- --------- -----------
         Name          F   ndf       ddf        pval
 ------------- ---------- ----- --------- -----------
  (Intercept)    22.2346   1.0   3.16597   0.0159885
     sequence   0.197525   1.0   3.00055    0.686828
       period     1.0789   3.0   8.66853    0.407683
  formulation   0.767409   1.0   5.46311    0.417867
 ------------- ---------- ----- --------- -----------

Model construction

Metida.LMMType
LMM(model, data; contrasts=Dict{Symbol,Any}(),  random::Union{Nothing, VarEffect, Vector{VarEffect}} = nothing, repeated::Union{Nothing, VarEffect} = nothing)

Make Linear-Mixed Model object.

model: is a fixed-effect model (@formula)

data: tabular data

contrasts: contrasts for fixed factors

random: vector of random effects or single random effect

repeated: is a repeated effect (only one)

source
  • model - example: @formula(var ~ sequence + period + formulation)

  • random - effects can be specified like this: VarEffect(@covstr(formulation|subject), CSH). @covstr is a effect model: @covstr(formulation|subject). CSH is a CovarianceType structure. Premade constants: SI, DIAG, AR, ARH, CS, CSH, ARMA. If not specified only repeated used.

  • repeated - can be specified like random effect. If not specified VarEffect(@covstr(1|1), SI) used. If no repeated effects specified vector of ones used.

Random/repeated model

Metida.@covstrMacro
@covstr(ex)

Macros for random/repeated effect model.

Example

@covstr(factor|subject)
source

Random/repeated effect construction

Metida.VarEffectType
VarEffect(formula, covtype::T, coding) where T <: AbstractCovarianceType

VarEffect(formula, covtype::T; coding = nothing) where T <: AbstractCovarianceType

VarEffect(formula; coding = nothing)

Random/repeated effect.

  • formula from @covstr(ex) macros.

  • covtype - covariance type (SI, DIAG, CS, CSH, AR, ARH, ARMA, TOEP, TOEPH, TOEPP, TOEPHP)

Note

Categorical factors are coded with FullDummyCoding() by default, use coding for other contrast codeing.

Example

VarEffect(@covstr(1+factor|subject), CSH)

VarEffect(@covstr(1 + formulation|subject), CSH; coding = Dict(:formulation => StatsModels.DummyCoding()))
source

Fitting

StatsBase.fit!Function
fit!(lmm::LMM{T}; kwargs...
) where T

Fit LMM model.

Keywords:

  • solver - :default / :nlopt for using with MetidaNLopt.jl/ :cuda for using with MetidaCu.jl
  • verbose - :auto / 1 / 2 / 3
  • varlinkf - :exp / :sq / :identity ref
  • rholinkf - :sigm / :atan / :sqsigm / :psigm
  • aifirst - first iteration with AI-like method - :default / :ai / :score
  • g_tol - absolute tolerance in the gradient
  • x_tol - absolute tolerance of theta vector
  • f_tol - absolute tolerance in changes of the REML
  • hes - calculate REML Hessian
  • init - initial theta values
  • refitinit - true/false - if true - use last values for initial condition
  • optmethod - Optimization method. Look at Optim.jl documentation. (Newton by default)
  • io - output IO
source
  • solver - :default solving with Optim.jl, for :nlopt and :cuda MetidaNLopt.jl and MetidaCu.jl should be installed.

  • verbose - 1 - only log, 2 - log and print, 3 - print only errors, other log, 0 (or any other value) - no logging.