This switch means almost a 20-fold speedup (for c++ routines) when using analytical first derivatives (i.e., need to calculate "trace" terms)
Speedup comes from only calculating elements of inverse matrix (C-inverse) following the non-zero pattern determined for the Cholesky decomposition of the C matrix.
follows the SuiteSparse Matlab_Tools sparseinv
by Tim Davis, but sparseinv works on LDL' factorization of C whereas I changed this to work on LL' factorization of C.
gremlinR()
now uses far less RAM per iteration of the model (previously was forming the entire C-inverse)
add finite difference algorithm to obtain first derivatives of likelihood function
h
) inside gremlinControl()
to set the "difference" or amount to alter parameters to calculate change in log-likelihood.created a REML function inside the c++ code (reml
) to calculate log-likelihood
error()
in c++
now issues with matrix singularities etc. do not stop code without returning model so far
should now be possible to use update()
to get "through" trouble spots
also allows for user interruptions to c++ code from terminal
Changed default parameterization so lambda transformation is no longer the default
Changed default convergence check criteria (cctol
)
deltaSE()
function to calculate approximate standard errors for functions of (co)variance parameters (e.g., h2, standard deviations of variances, or correlations)
this can take a formula for the function or a character expression
also allows for a list of formulas or character expressions e.g., calculate all variance components as proportions of total variance
Introduce Gcon
and Rcon
arguments to gremlin()
for constraining parameters
Gstart
and Rstart
argumentssire
model, we could restrain the sire
variance =0.38
.grSf <- gremlin(WWG11 ~ sex,
random= ~ sire,
data = Mrode11,
Gstart = list(matrix(0.38)),
Gcon = list("F"),
control = gremlinControl(lambda = FALSE))
Similar to above change (Gcon
/Rcon
), introduced steps to deal with parameters outside of the boundaries of their parameter space (e.g., variance < 0).
change version numbering to just 3 numbers (instead of 4)
update()
function
can now continue a model where it left off or change the structure (e.g., drop a single variance component for likelihood ratio test)
Implement "step-halving" algorithm for AI updates
gremlinControl()
using the step
argumentgremlinControl()
function for advanced changes to the way gremlin runsem
, ai
, and elsewhere (where relevant) in next versionmade a "modular" series of functions for setting up the model and optimizing the REML likelihood
new grMod
and gremlinR
classes.
grMod
is the model structure for which a log-likelihood can be calculatedgremlinR
class distinguishes from gremlin
class in that gremlinR
objects will only use R
code written by the package in order to run the model. Class gremlin
will execute underlying c++ code written in the package.Average Information algorithm has been vastly improved
ai()
efficiently calculates the AI matrix without directly computing several matrix inverses (as previously coded)lambda
and alternative parameterizations now possible and executed by the same code
lambda
parameterization is the REML likelihood of the variance ratios after factoring out a residual variance from the Mixed Model Equations.lambda
models).M
) matrix from which the Cholesky factorization (and logDetC
and tyPy
calculations are made)
C
) and obtain tyPy
and logDetC
using thisM
and C
, now do a solve
with Cholesky of C
(sLc
/Lc
in R
/c++
code) to calculate tyPy
based off Boldman and Van Vleckgremlin
objects
AIC
, residuals
, anova
, and nobs
summary
, print
, and logLik
methods as wellImproved algorithm that reduces computational resources and time! Also implemented c++ code in gremlin()
, while keeping gremlinR()
purely the R implementation (at least from the package writing standpoint).
Documentation has switched from filling out the .Rd
files manually to providing
documentation next to the function code in the .R
files using roxygen2
Congratulations, its a gremlin!