From 4c87e49630dd5c88411771edb8abb6c7834ea3f5 Mon Sep 17 00:00:00 2001 From: Daniel Date: Thu, 13 Jun 2024 14:48:18 +0200 Subject: [PATCH] update --- R/compute_variances.R | 4 ++-- WIP/nakagawa_suppl.R | 45 +++++++++++++++++++++++++++++++++++++++++++ 2 files changed, 47 insertions(+), 2 deletions(-) diff --git a/R/compute_variances.R b/R/compute_variances.R index a675a5abe..309529151 100644 --- a/R/compute_variances.R +++ b/R/compute_variances.R @@ -767,7 +767,7 @@ # transform expected mean of the null model if (is.null(faminfo$family)) { - mu <- exp(mu) + mu <- link_inverse(x)(mu) } else { # transform mu mu <- switch(faminfo$family, @@ -868,7 +868,7 @@ switch(approx_method, delta = cvsquared, - trimamma = trigamma(cvsquared), + trigamma = trigamma(cvsquared), log1p(cvsquared) ) } diff --git a/WIP/nakagawa_suppl.R b/WIP/nakagawa_suppl.R index eb6bbaeb8..642df9d67 100644 --- a/WIP/nakagawa_suppl.R +++ b/WIP/nakagawa_suppl.R @@ -384,3 +384,48 @@ c(R2glmmM = R2glmmM, R2glmmC = R2glmmC) performance::r2_nakagawa(parmodGLMERf) MuMIn::r.squaredGLMM(parmodGLMERf) + + + + + + +# Fit null model without fixed effects (but including all random effects) +glmmTMBr <- glmmTMB::glmmTMB( + count ~ (1 | site), + family = glmmTMB::nbinom1(), + data = Salamanders, REML = TRUE +) +glmmTMBf <- glmmTMB::glmmTMB( + count ~ mined + spp + (1 | site), + family = glmmTMB::nbinom1(), + data = Salamanders, REML = TRUE +) +# Calculation of the variance in fitted values +VarF <- var(as.vector(get_modelmatrix(glmmTMBf) %*% fixef(glmmTMBf)$cond)) +# getting the observation-level variance Null model +thetaN <- sigma(glmmTMBr) +lambda <- as.numeric(exp(fixef(glmmTMBr)$cond + 0.5 * (as.numeric(VarCorr(glmmTMBr)$cond[1])))) +# lambda2 <- mean(DataAll$Parasite) +VarOdN <- 1 / lambda + 1 / thetaN # the delta method +VarOlN <- log(1 + (1 / lambda) + (1 / thetaN)) # log-normal approximation +VarOtN <- trigamma((1 / lambda + 1 / thetaN)^(-1)) # trigamma function +# comparing the three +c(VarOdN = VarOdN, VarOlN = VarOlN, VarOtN = VarOtN) + +thetaF <- sigma(glmmTMBf) # note that theta is called alpha in glmmadmb +VarOdF <- 1 / lambda + 1 / thetaF # the delta method +VarOlF <- log(1 + (1 / lambda) + (1 / thetaF)) # log-normal approximation +VarOtF <- trigamma((1 / lambda + 1 / thetaF)^(-1)) # trigamma function +# comparing the three +c(VarOdF = VarOdF, VarOlF = VarOlF, VarOtF = VarOtF) + +# R2[GLMM(m)] - marginal R2[GLMM] +R2glmmM <- VarF / (VarF + sum(as.numeric(VarCorr(glmmTMBf)$cond)) + VarOlF) +# R2[GLMM(c)] - conditional R2[GLMM] for full model +R2glmmC <- (VarF + sum(as.numeric(VarCorr(glmmTMBf)$cond))) / (VarF + sum(as.numeric(VarCorr(glmmTMBf)$cond)) + VarOlF) +c(R2glmmM = R2glmmM, R2glmmC = R2glmmC) + +performance::r2_nakagawa(glmmTMBf, null_model = glmmTMBr) +MuMIn::r.squaredGLMM(glmmTMBf) +get_variance(glmmTMBf, null_model = glmmTMBr)