# Compute (Bayesian) R2 for Turing models
**URL:** https://discourse.julialang.org/t/compute-bayesian-r2-for-turing-models/104537
**Category:** Probabilistic Programming
**Tags:** turing, bayesian-inference
**Created:** [October 3, 2023, 10:49am UTC](https://discourse.julialang.org/t/compute-bayesian-r2-for-turing-models/104537 "2023-10-03T10:49:01Z")
**Posts on this page:** 3
**Page:** 1
### Author:  [@DominiqueMakowski](https://discourse.julialang.org/u/DominiqueMakowski)
#### Post date: [October 3, 2023, 10:49am UTC](https://discourse.julialang.org/t/compute-bayesian-r2-for-turing-models/104537/1 "2023-10-03T10:49:01Z")
Somehow my google search didn’t yield many results for how to compute R2 for Turing models.
StatsBase.jl [defines](https://juliastats.org/StatsBase.jl/stable/statmodels/#StatsAPI.r2) a few generic methods, but I was in particular interested in the so-called Bayesian R2 ([Gelman et al., 2019](https://www.tandfonline.com/doi/full/10.1080/00031305.2018.1549100)), defined as:
Would that be easy to implement in Julia/Turing? Its [R implementation](https://github.com/jgabry/bayes_R2/) seems pretty straightforward:
```r
bayes_R2 <- function(fit) {
y <- rstanarm::get_y(fit)
ypred <- rstanarm::posterior_linpred(fit, transform = TRUE)
if (family(fit)$family == "binomial" && NCOL(y) == 2) {
trials <- rowSums(y)
y <- y[, 1]
ypred <- ypred %*% diag(trials)
}
e <- -1 * sweep(ypred, 2, y)
var_ypred <- apply(ypred, 1, var)
var_e <- apply(e, 1, var)
var_ypred / (var_ypred + var_e)
}
```
---
### Author:  [@p-gw](https://discourse.julialang.org/u/p-gw)
#### Post date: [October 3, 2023, 12:35pm UTC](https://discourse.julialang.org/t/compute-bayesian-r2-for-turing-models/104537/2 "2023-10-03T12:35:08Z")
Yes, this is easily possible in Turing.
The simplest approach is probably to return the R² value from the model directly,
```julia
@model function mymodel()
# ...model definition
r2 = compute_r2(args)
return r2
end
```
where `compute_r2` is a function that computes the R² values given some parameters of the model.
Then you can call `generated_quantities` on the sampled model to get the R² values for each iteration.
```julia
model = mymodel()
chain = sample(model, NUTS(), 1000)
generated_quantities(model, chain)
```
---
### Author:  [@sethaxen](https://discourse.julialang.org/u/sethaxen)
#### Post date: [October 3, 2023, 5:17pm UTC](https://discourse.julialang.org/t/compute-bayesian-r2-for-turing-models/104537/3 "2023-10-03T17:17:47Z")
PosteriorStats.jl implements Gelman’s Bayesian R^2. See [API · PosteriorStats.jl](https://julia.arviz.org/PosteriorStats/stable/api/#PosteriorStats.r2_score).