From 220dbd4cfbbe6ca3241aeecab1ed47fec4977de5 Mon Sep 17 00:00:00 2001 From: Martin Elff Date: Sat, 14 Dec 2024 21:23:22 +0100 Subject: [PATCH] Use block inverse in random effects model fitting --- pkg/R/mmclogit-fitPQLMQL.R | 15 +++++++++++++-- 1 file changed, 13 insertions(+), 2 deletions(-) diff --git a/pkg/R/mmclogit-fitPQLMQL.R b/pkg/R/mmclogit-fitPQLMQL.R index b3eedad..bf0f2c6 100644 --- a/pkg/R/mmclogit-fitPQLMQL.R +++ b/pkg/R/mmclogit-fitPQLMQL.R @@ -454,7 +454,13 @@ PQLMQL_innerFit <- function(parms,aux,model.struct,method,estimator,control){ Phi <- lapply(Psi,safeInverse) ZWZiSigma <- ZWZ + iSigma - K <- solve(ZWZiSigma) + if(getOption("mclogit.use_blkinv", TRUE)) { + K <- blk_inv.squareBlockMatrix(ZWZiSigma) + } + else { + K <- solve(ZWZiSigma) + } + log.det.iSigma <- Lambda2log.det.iSigma(Lambda,m) @@ -558,7 +564,12 @@ PQLMQL_pseudoLogLik <- function(lambda, iSigma <- Psi2iSigma(Psi,m) H <- ZWZ + iSigma - K <- solve(H) + if(getOption("mclogit.use_blkinv", TRUE)) { + K <- blk_inv.squareBlockMatrix(H) + } + else { + K <- solve(H) + } XiVX <- XWX - fuseMat(bMatCrsProd(ZWX,bMatProd(K,ZWX))) XiVy <- XWy - fuseMat(bMatCrsProd(ZWX,bMatProd(K,ZWy)))