glmmrBase
R package to support the specification of generalised linear mixed models using the R6 objectorientated class system.
Generalised linear mixed models
A generalised linear mixed model (GLMM) has a mean function for observation
where
where
Generating data
The package includes the function nelder()
, which we use to generate data for the examples below. Nelder (1965) suggested a simple notation that could express
a large variety of different blocked designs. The notation was proposed in the context of splitplot experiments for agricultural research, where researchers
often split areas of land into blocks, subblocks, and other smaller divisions, and apply different combinations of treatments. However, the notation is useful
for expressing a large variety of experimental designs with correlation and clustering including cluster trials, cohort studies, and spatial and temporal
prevalence surveys. We have included the function \code{nelder()} that generates a data frame of a design using the notation.
There are two operations:

>
(or$\to$ in Nelder's notation) indicates "clustered in". 
*
(or$\times$ in Nelder's notation) indicates a crossing that generates all combinations of two factors.
The function takes a formula input indicating the name of the variable and a number for the number of levels, such as abc(12)
.
So for example ~cl(4) > ind(5)
means in each of five levels of cl
there are five levels of ind
, and the individuals are different between clusters.
The formula ~cl(4) * t(3)
indicates that each of the four levels of cl
are observed for each of the three levels of t
. Brackets are used to indicate the
order of evaluation.
Specifying covariance
The specification of a covariance object requires three inputs: a formula, data, and parameters. A new instance of each class can be generated with the
$new()
function, for example Covariance\$new(...)
.
A covariance function is specified as an additive formula made up of components with structure (1f(j))
. The left side of the vertical bar specifies the
covariates in the model that have a random effects structure. The right side of the vertical bar specify the covariance function f
for that term using variable
named in the data j
. Covariance functions on the right side of the vertical bar are multiplied together, i.e., (1f(j)*g(t))
. The table below shows the currently
implemented covariance functions
Function  Code 


Identity/ Group membership  gr(x) 

Exponential  fexp(x) 

fexp0(x) 

Squared Exponential  sqexp(x) 

sqexp0(x) 

Autoregressive order 1  ar1(x) 

Bessel 

bessel(x) 

Matern  matern(x) 

Compactly supported*  
Wendland 0  wend0(x) 

Wendland 1  wend1(x) 

Wendland 2  wend1(x) 

WhittleMatern 
prodwm(x) 

Cauchy 
prodcb(x) 

Exponential 
prodek(x) 

x
.
**Permissible in one or two dimensions. ***Only permissible in one dimension. ****Permissible in up to three dimensions.
One combines functions to provide the desired covariance function. For example, for a steppedwedge cluster trial we could consider the
standard specification with an exchangeable random effect for the cluster level, and a separate exchangeable random effect for the clusterperiod,
which would be ~(1gr(j))+(1gr(j,t))
or ~(1gr(j))+(1gr(j)*gr(t))
. Alternatively, we could consider an autoregressive clusterlevel random effect
that decays exponentially over time so that, for persons ~(1gr(j)*ar1(t))
.
Parameters are provided to the covariance function as a vector. The covariance functions described in the Table have different parameters
A new covariance object can then be created as such
R> df < nelder(~ (j(10)* t(5)) > ind(10))
R> cov < Covariance$new(formula = ~(1gr(j)*ar1(t)),
R> parameters = c(0.05,0.8),
R> data= df)
where a compactly supported function is used, then the effective range parameters should be provided in the order the function appears in the formula.
R> cov < Covariance$new(formula = ~(1prodwm(x,y)),
R> parameters = c(0.25,0.5),
R> eff_range = 0.5,
R> data= df)
Mean function specification
Specification of the mean function follows standard model formulae in R. A vector of values of the mean function parameters is required to complete the model specification along with the distribution as a standard R family object. A complete specification is thus:
R> mf < MeanFunction$new(formula = ~ factor(t)+ int  1,
R> data = df,
R> parameters = rep(0,6),
R> family = gaussian())
Note that factor
in this function does not drop one level, unlike standard R formulae, so removing the intercept is required to prevent a collinearity problem.
Model specification
A model is simply a covariance object and a mean function object:
R> model < Model$new(covariance = cov,
R> mean.function = mf,
R> var_par = 1)
For Gaussian models, and other distributions requiring an additional scale parameter var_par
which is the
conditional variance covariance
and mean.function
instead of model objects.
Accessing computed elements
Each class holds associated matrices and has member functions to compute basic summaries and analyses. The Matrix
package is used
for matrix operations in R. For example, a Covariance
object holds the matrix
R> cov$D[1:10,1:10]
10 x 10 sparse Matrix of class "dsCMatrix"
[1,] 0.002500 0.00200 0.0016 0.00128 0.001024 . . . . .
[2,] 0.002000 0.00250 0.0020 0.00160 0.001280 . . . . .
[3,] 0.001600 0.00200 0.0025 0.00200 0.001600 . . . . .
[4,] 0.001280 0.00160 0.0020 0.00250 0.002000 . . . . .
[5,] 0.001024 0.00128 0.0016 0.00200 0.002500 . . . . .
[6,] . . . . . 0.002500 0.00200 0.0016 0.00128 0.001024
[7,] . . . . . 0.002000 0.00250 0.0020 0.00160 0.001280
[8,] . . . . . 0.001600 0.00200 0.0025 0.00200 0.001600
[9,] . . . . . 0.001280 0.00160 0.0020 0.00250 0.002000
[10,] . . . . . 0.001024 0.00128 0.0016 0.00200 0.002500
glmmrBase
in other packages
Use of This package provides a foundation for other packages providing analysis or estimation of generalised linear mixed models. For example, we have the glmmrMCML
package,
which provides Markov Chain Monte Carlo Maximum Likelihood estimations for these models, and glmmrOptim
, which finds approximate optimal designs based on these
models. New classes can be defined that inherit from the classes included in this package. glmmrMCML
defines the modelMCML
class that adds the member
function MCML
. Then the new functions can access the model elements, such as covariance matrices, and benefit from automatic updating of these elements when specifications or parameters change.
As an example we can define a new class that has a member function that returns the determinant of the matrix
R> CovDet < R6::R6Class("CovDet",
R> inherit = Covariance,
R> public = list(
R> det = function(){
R> return(Matrix::determinant(self$D))
R> }))
R> cov < CovDet$new(formula = ~(1gr(j)*ar1(t)),
R> parameters = c(0.05,0.8),
R> data= df)
R> cov$det()
$modulus
[1] 340.4393
attr(,"logarithm")
[1] TRUE
$sign
[1] 1
attr(,"class")
[1] "det"