# Calculate R2_lik Calculate partial and total R2s for LMM, GLMM, PGLS, and PGLMM using R2_lik, an R2 based on the likelihood of the fitted model. ```r R2_lik(mod = NULL, mod.r = NULL) ``` ## Arguments - `mod`: A regression model with one of the following classes: 'lm', 'glm', 'lmerMod', 'glmerMod', 'phylolm', 'phyloglm', 'gls', 'pglmm', pglmm_compare' or 'communityPGLMM'. - `mod.r`: A reduced model; if not provided, the total R2 will be given by setting 'mod.r' to the model corresponding to 'mod' with the intercept as the only predictor. ## Returns R2_lik value. ## Details `R2_lik()` is implemented as $$ partial R2 = 1 - exp(-2/n * (logLik(mod.f) - logLik(mod.r))) $$ where 'mod.f' and 'mod.r' are the full and reduced models, respectively. The total R2 is given when 'mod.r' is the model corresponding to mod.f that contains only the intercept. For GLMMs and PGLMMs, `R2_lik()` is standardized to have a maximum of one following Nagelkerke (1991). Although you can use R2_lik with models fit with REML, you really shouldn't, because this makes it impossible to compare values when reduced models differ in independent variables (fixed effects). `R2_lik()` is also computed for LMMs and GLMMs in the MuMIn package. ## Examples ```r library(ape) library(phylolm) library(lme4) library(nlme) library(phyr) set.seed(12345) p1 <- 10 nsample <- 10 n <- p1 * nsample d <- data.frame(x1 = 0, x2 = 0, u1 = rep(1:p1, each = nsample), u2 = rep(1:p1, times = nsample)) d$u1 <- as.factor(d$u1) d$u2 <- as.factor(d$u2) b1 <- 1 b2 <- -1 sd1 <- 1.5 d$x1 <- rnorm(n = n) d$x2 <- rnorm(n = n) d$y.lmm <- b1 * d$x1 + b2 * d$x2 + rep(rnorm(n = p1, sd = sd1), each = nsample) + rep(rnorm(n = p1, sd = sd1), times = nsample) + rnorm(n = n) prob <- inv.logit(b1 * d$x1 + rep(rnorm(n = p1, sd = sd1), each = nsample)) d$y.glmm <- rbinom(n = n, size = 1, prob = prob) # LMM with two fixed and two random effects ---- z.f <- lmer(y.lmm ~ x1 + x2 + (1 | u1) + (1 | u2), data = d, REML = FALSE) z.x <- lmer(y.lmm ~ x1 + (1 | u1) + (1 | u2), data = d, REML = FALSE) z.v <- lmer(y.lmm ~ 1 + (1 | u2), data = d, REML = FALSE) z.0 <- lm(y.lmm ~ 1, data = d) R2_lik(z.f, z.x) R2_lik(z.f, z.v) R2_lik(z.f) # GLMM with one fixed and one random effect ---- z.f <- glmer(y.glmm ~ x1 + (1 | u1), data = d, family = 'binomial') z.x <- glmer(y.glmm ~ 1 + (1 | u1), data = d, family = 'binomial') z.v <- glm(y.glmm ~ x1, data = d, family = 'binomial') R2_lik(z.f, z.x) R2_lik(z.f, z.v) R2_lik(z.f) # PGLS with a single fixed effect ---- n <- 100 d <- data.frame(x = array(0, dim = n), y = 0) b1 <- 1.5 signal <- 0.5 phy <- compute.brlen(rtree(n = n), method = 'Grafen', power = 1) phy.x <- compute.brlen(phy, method = 'Grafen', power = .0001) # Generate random data x <- rTraitCont(phy.x, model = 'BM', sigma = 1) e <- signal^0.5 * rTraitCont(phy, model = 'BM', sigma = 1) + (1-signal)^0.5 * rnorm(n = n) d$x <- x[match(names(e), names(x))] d$y <- b1 * x + e rownames(d) <- phy$tip.label d$sp <- phy$tip.label z.f <- pglmm_compare(y ~ x, data = d, phy = phy, REML=FALSE) z.x <- pglmm_compare(y ~ 1, data = d, phy = phy, REML=FALSE) z.v <- glm(y ~ x, data = d) R2_lik(z.f, z.x) R2_lik(z.f, z.v) R2_lik(z.f) # These data can also be fit with phylolm() in {phylolm} z.f <- phylolm(y ~ x, phy = phy, data = d, model = 'lambda') z.x <- phylolm(y ~ 1, phy = phy, data = d, model = 'lambda') z.v <- lm(y ~ x, data = d) R2_lik(z.f, z.x) R2_lik(z.f, z.v) R2_lik(z.f) # This also works for models fit with gls() in {nlme} z.f <- gls(y ~ x, data = d, correlation = corPagel(1, phy, form = ~sp), method = "ML") z.x <- gls(y ~ 1, data = d, correlation = corPagel(1, phy, form = ~sp), method = "ML") z.v <- lm(y ~ x, data = d) R2_lik(z.f, z.x) R2_lik(z.f, z.v) R2_lik(z.f) # PGLMM with one fixed effect ---- n <- 100 b1 <- 1.5 signal <- 2 phy <- compute.brlen(rtree(n = n), method = 'Grafen', power = 1) phy.x <- compute.brlen(phy, method = 'Grafen', power = .0001) # Generate random data x <- rnorm(n) d <- data.frame(x = x, y = 0) e <- signal * rTraitCont(phy, model = 'BM', sigma = 1) e <- e[match(phy$tip.label, names(e))] d$y <- rbinom(n = n, size = 1, prob = inv.logit(b1 * d$x + e)) rownames(d) <- phy$tip.label z.f <- phyloglm(y ~ x, data = d, start.alpha = 1, phy = phy) z.x <- phyloglm(y ~ 1, data = d, phy = phy, start.alpha = min(20, z.f$alpha)) z.v <- glm(y ~ x, data = d, family = 'binomial') R2_lik(z.f, z.x) R2_lik(z.f, z.v) R2_lik(z.f) # These data can also be fit with pglmm_compare(), although note that # this is a different model from phyloglm() z.f <- pglmm_compare(y ~ x, data = d, family = "binomial", phy = phy, REML=FALSE) z.x <- pglmm_compare(y ~ 1, data = d, family = "binomial", phy = phy, REML=FALSE) z.v <- glm(y ~ x, data = d, family = "binomial") R2_lik(z.f, z.x) R2_lik(z.f, z.v) R2_lik(z.f) ``` ## References Ives A.R. and Li D. 2018. rr2: An R package to calculate R2s for regression models. Journal of Open Source Software. DOI:10.21105/joss.01028 Ives A.R. 2018. R2s for Correlated Data: Phylogenetic Models, LMMs, and GLMMs. Systematic Biology, Volume 68, Issue 2, March 2019, Pages 234-251. DOI:10.1093/sysbio/syy060 Nagelkerke 1991. A note on a general definition of the coefficient of determination. Biometrika 78:691–692. ## See Also MuMIn, lme4, ape, phylolm, phyr, pez ## Author(s) Anthony R. Ives

R2_lik function