0
$\begingroup$

I am interested to calculate the pseudo-r-squuared for a zero-truncated poisson mixel model (using glmmTMB). The r.squaredglmm (package MuMin) gives a message that it can not calculate pseudo-r2 for zero-truncated poisson. Is there an easy work around or quick-and-dirty alternative?

$\endgroup$
2
  • $\begingroup$ Welcome to Cross Validated! What exactly do you mean by a pseudo $R^2?$ Once you get away from linear models with least squares fits, a generalization of $R^2$ is not straightforward, and there are multiple ways that can be defended. For instance, UCLA lists many reasonable generalizations of $R^2$ to logistic regression. $\endgroup$
    – Dave
    Commented Jun 27, 2022 at 14:07
  • $\begingroup$ Thanks for welcome. Actually, I thought that in case of mixed models, you speak of pseude-R2 and that you might calculate R2, and indeed there are several approaches for that. I am not an expert, but just want to have a measure of explanation of variation of the fixed (marginal) and totla (conditional) effects if possible at all. $\endgroup$
    – Wil Tamis
    Commented Jun 28, 2022 at 6:01

1 Answer 1

2
$\begingroup$

I hope this answer is not too late.

The performance library offers a function r2_zeroinflated() which is supposed to compute "R2 for Zero-Inflated and Hurdle Regression". However, I couldn't get the feature to work.

If you'd settle for the pseudo-R2 based on a likelihood-ratio, you can use MuMIn::r.squaredLR(). For the hurdle model you will have to manually specify the null argument (null model). But otherwise it works fine.

Third option would be to somehow implement equations (1.17) or (1.18) from this dissertation. On page 10 the author refers to an article Zero tolerance ecology: improving ecological inference by modelling the source of zero observations by Martin et al. (2005). I looked it up, but unfortunately I couldn't find any measures of pseudo-R2 in it...

As the last one I'll share a link to this page discussing calculation of conditional and marginal pseudo-R2 for glmmTMB objects but I'm doing it with a big big warning because they discuss a model without a zero-inflation part and with an observation-level random effect (OLRE). It isn't particularly what you have asked for. But if you can handle it, below I attach the script which is a bit clearer than in the link above and from which I removed the part where the variance of the OLRE term was calculated.

## fixed effects variance (manually)
model_linear.predictor = model.matrix(model) %*% fixef(model)$cond  ## "$cond" will extract coefficients from the conditional model; if you'd like to obtain coefficients from the zero-inflation part, use "$zi" instead
model_fixed.var = var(model_linear.predictor)

## random effects variance
model_randeff.var = unlist(VarCorr(model)$cond)

## originally the OLRE variance (model_OLRE.var) was calculated here

## the observation-level variance (it is distinct from the OLRE variance; the observation-level variance is the variance added by the distribution itself (a.k.a. distribution-specific variance))
model_beta0 = fixef(update(model, ~ 1))$cond  ## beta0 is the intercept from an intercept-only (y ~ 1) refit of the model
model_lambda = exp(model_beta0 + model_randeff.var/2)
model_ol.var = trigamma(model_lambda)

## total variance
model_total.var = model_fixed.var + model_randeff.var + model_ol.var  ## + model_OLRE.var

## conditional pseudo-R2 (takes both the fixed and random effects into account)
(model_fixed.var + model_randeff.var) / model_total.var

## marginal pseudo-R2 (considers only the variance of the fixed effects)
model_fixed.var / model_total.var

But again, although the script works technically well, this version of the script is definitely not the right solution for zero-inflated and hurdle models. It'd definitely need a modification.

Good luck with your analysis.

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.