I don't know how to do this in R, because I'm more of a WinBUGS person for this sort of thing. There you'd just monitor the precision (or variance) and whatever ratio you're interested in, and it would just another parameter that was reported, with credible interval. So if you don't get any better suggestions, WinBUGS links with R very nicely. Somebody else might ask whether you really need the R-squared, and whether something else might be better.

Thank you so much, Shinchi! That's incredibly helpful!

Hi Rafter and Shinichi Nakagawa

Thanks for this question and response above here.

I'm trying to do something similar to this at the moment for a multiple response GLMM in mcmcGLMM in R, with the residual error in the model accounting for the covariance of the response variables. It has a form roughly like this is R.

mod_fin_1

You may be interested in this

http://www.stat.columbia.edu/~gelman/research/unpublished/bayes_R2.pdf

R-squared for Bayesian regression models∗

Andrew Gelman† Ben Goodrich‡ Jonah Gabry‡ Imad Ali

8 Nov 2017

Abstract The usual definition of R2 (variance of the predicted values divided by the variance of the data) has a problem for Bayesian fits, as the numerator can be larger than the denominator. We propose an alternative definition similar to one that has appeared in the survival analysis literature: the variance of the predicted values divided by the variance of predicted values plus the variance of the errors. This summary is computed automatically for linear and generalized linear regression models fit using rstanarm, our R package for fitting Bayesian applied regression models with Stan.

Thanks Kelvyn,

Sorry for the delayed reply, this paper and function are extremely useful to me. I like the approach very much and can write a similar function to work with univariate response models fitting with mcmcGLMM.

However, as my model is a multi-response model, by which I mean that there are multiple response variables being fitted simultaneously. I have a further issue with applying this. In my model, I have multiple-responses which are fitted simultaneously, which means the residual error structure contains a matrix of covariances which represent the covariances between the response traits. This approach is useful where there are known reasons why responses would covary (in my case the response are life-history metrics, so we expect trade-offs/constraint between reproduction and lifespan of species for example. Both of which are responses in the model, and modelling them separately in a univariate model, would miss out on this important piece of the picture) (see https://www.jstatsoft.org/article/view/v033i02/v33i02.pdf for details of the maths and fitting).

The reason that dealing with this residual error matrix needs extra consideration is that it is unclear how to treat these covariances when partitioning variance in the model. Note here that the covariances can frequently be negative, leading to odd results when one sums these, or adds them to other variances in the model. One could take the magnitude of the covariances as the important value, or the square of the covariance estimates (as is used in the calculation of 'trace'), or perhaps ignore the covariances and just partition the explained variance by the fitted variables relative to the residual error in terms of the unexplained variance (i.e. the residual variance within each response trait, i.e. the diagonal of the residual covariance matrix, ignoring the covariance in the off-diagonals).

I can calculate all of these relatively easily, but what I am need is further guidance on what is the most meaningful way to deal with these covariances, and what other people have done.

The article you passed on was really useful, in that it is the closest I've seen to what I'm attempting, and really helped with understanding why traditional frequentist r2 approaches won't work here. But it stops one step short of where I want to be. Any thoughts?

Many thanks again for your help.

Best regards,

Ruth

Ruth Kelly

I am a fan of multivariate models (see my long answer to https://www.researchgate.net/post/Can_you_share_any_examples_of_the_use_of_multivariate_regression)

but I know of no developments of Bayesian r-squared for these, just to say that if if you have two variables X and Y, the total variance is

Var(X*X) + 2*covar(X*Y) + Var(Y*Y)

It might be useful to post your follow up question on the Gelman website

https://statmodeling.stat.columbia.edu/2017/12/21/r-squared-bayesian-regression-models/