Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Marginal R2 is a misnomer in brms #526

Open
vincentarelbundock opened this issue Dec 30, 2022 · 12 comments
Open

Marginal R2 is a misnomer in brms #526

vincentarelbundock opened this issue Dec 30, 2022 · 12 comments
Labels
consistency 🍏 🍎 Expected output across functions could be more similar docs 📚 Something to be adressed in docs and/or vignettes

Comments

@vincentarelbundock
Copy link
Contributor

According to the brms author, performance should change its label:

https://twitter.com/paulbuerkner/status/1608896620699226113?s=46&t=3uzWoIWG2nAVyzeyGrxTEg

@DominiqueMakowski
Copy link
Member

If we can expect to obtain a better marginal R2, we could explicit that in the documentation (rather than changing the name which would be super breaking probably) and try to find a way to compute a better score - assuming this is a goal we can reasonably have

@strengejacke
Copy link
Member

It's a possible misnomer for Bayesian models, not for frequentist.

@vincentarelbundock
Copy link
Contributor Author

TBC, I have no expertise here. I just wanted to bring your attention to this because it seemed like a credible argument from authority. Feel free to close once you decide what should be done.

@strengejacke
Copy link
Member

I think you can use r2_nakagawa() with Bayesian models, too, so we would have "valid" marginal and conditional R2, but that function is also not fully validated for more exotic model families.

@bwiernik
Copy link
Contributor

bwiernik commented Jan 1, 2023

I'm not sure what you mean @strengejacke about the distinction between Bayesian and frequentist frameworks?

@strengejacke
Copy link
Member

For frequentist models, we follow the approach from Nakagawa et. al, computing the variances for the different components (via insight::get_variance(), and then - similar to ICC - divide fixed effects vs. total variance or FE+RE vs. total, see

performance/R/r2_nakagawa.R

Lines 107 to 110 in 9b95409

} else {
r2_marginal <- vars$var.fixed / (vars$var.fixed + vars$var.random + vars$var.residual)
r2_conditional <- (vars$var.fixed + vars$var.random) / (vars$var.fixed + vars$var.random + vars$var.residual)
}

This is the distinction of "marginal" vs. "conditional" proposed by Nakagawa et al., and which is in line with many other packages that compute R2 for mixed models.

For Bayesian mixed models, we rely on bayes_R2() and once set re.form = NA and the other time re.form = NULL, where bayes_R2() internally calls posterior_epred() with the re.form argument. We call these two "marginal" and "conditional", which may be a misnomer, because since we rely on "predictions", "marginalizing" over random effects is not achieved by setting re.form either to NA or NULL. Thus, "marginal" is probably not the correct term when thinking in predictions.

But, since we're interested in the variance, re.form = NA and re.form = NULL still might be an equivalent to marginal and conditional R2 in the frequentist framework, because results are fairly similar, see:

library(lme4)
library(brms)
library(performance)

data(sleepstudy)
data(cbpp)
cbpp$out <- ifelse(cbpp$incidence / cbpp$size > 0.2, 1L, 0L)

m1 <- lmer(Reaction ~ Days + (1 + Days | Subject), data = sleepstudy)
m2 <- brm(Reaction ~ Days + (1 + Days | Subject), data = sleepstudy, refresh = 0)

gm1 <- glmer(
  out ~ period + (1 | herd),
  data = cbpp,
  family = binomial
)

gm2 <- brm(
  out ~ period + (1 | herd),
  data = cbpp,
  family = bernoulli,
  refresh = 0
)


r2(m1)
r2(m2)

r2(gm1)
r2(gm2)

Hence, in terms of "predictions", we don't have predictions marginalized over random effects, but since we're not interested in point estimates, but rather in the variances, we still might produce the "correct" marginal and conditional R2. Not sure about this, though, and that's why there might be a misnomer.

@strengejacke
Copy link
Member

r2() calls r2_bayes() for Bayesian (mixed) models, while r2() calls r2_nakagawa() for frequentist mixed models. But I think r2_nakagawa() can also be used for Bayesian models, just not for too exotic families.

@bwiernik
Copy link
Contributor

bwiernik commented Jan 2, 2023

Okay yeah I agree. I think May not really be a Bayes/freq thing but rather than Paul would probably just object to the Nakagawa use of marginal there.

@mattansb
Copy link
Member

mattansb commented Jan 2, 2023

But I think r2_nakagawa() can also be used for Bayesian models, just not for too exotic families.

It already can:

library(brms)
library(performance)

data(sleepstudy)
data(cbpp)
cbpp$out <- ifelse(cbpp$incidence / cbpp$size > 0.2, 1L, 0L)

m2 <- brm(Reaction ~ Days + (1 + Days | Subject), 
          data = sleepstudy, 
          refresh = 0, backend = "cmdstanr")

r2_bayes(m2)
#> # Bayesian R2 with Compatibility Interval
#> 
#>   Conditional R2: 0.793 (95% CI [0.758, 0.824])
#>      Marginal R2: 0.287 (95% CI [0.160, 0.406])

r2_nakagawa(m2)
#> # R2 for Mixed Models
#> 
#>   Conditional R2: 0.822
#>      Marginal R2: 0.240

@mattansb
Copy link
Member

mattansb commented Jan 2, 2023

Am I understanding correctly that there isn't any problem with the code per se, but the issue is with the labeling?

So it's not marginal in the sense that the random effect's aren't integrated out, but are "ignored" (as in https://twitter.com/tjmahr/status/1581563839459385344). Correct?
If this is the case, I think clearing this up in the docs is enough?

How does this explain why the original tweet had the marginal larger than the conditional?

@strengejacke
Copy link
Member

How does this explain why the original tweet had the marginal larger than the conditional?

Maybe related to
easystats/insight#664

@bwiernik
Copy link
Contributor

bwiernik commented Jan 2, 2023

@strengejacke strengejacke added docs 📚 Something to be adressed in docs and/or vignettes consistency 🍏 🍎 Expected output across functions could be more similar labels Jan 4, 2023
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
consistency 🍏 🍎 Expected output across functions could be more similar docs 📚 Something to be adressed in docs and/or vignettes
Projects
None yet
Development

No branches or pull requests

5 participants