Compute Bayesian R-squared from matrix of posterior replications

The approach is not fully Bayesian and provides a global estimate rather than an estimate for each sample (this is because the predictive means and residual variance are estimated from replications than given by the model).

compute_rsquared(yrep)

Arguments

yrep

Matrix with rows representing samples and columns representing observations

Value

Bayesian R-squared (scalar, between 0 and 1)

References

A. Gelman, B. Goodrich, J. Gabry, and A. Vehtari, “R-squared for Bayesian Regression Models,” Am. Stat., vol. 73, no. 3, pp. 307–309, Jul. 2019

Examples

N <- 50
N_sample <- 1e2
y <- runif(N, 0, 10)
yrep <- do.call(cbind,
                lapply(1:N,
                       function(i) {
                         y[i] + rnorm(N_sample)
                       }))
compute_rsquared(yrep)
#> [1] 0.8904137