Hi everyone, I wanna reproduce the Blooic estimate and se of rstan::loo(), using the pointwise data, but I can´t manage to reproduce it? Or anyone knows how to make a credible interval of the looic?
Hi,
what have you tried in order to reproduce? I’ve previously used small bits of code from loo (which I’ve found reasonably readable) to reproduce some steps in the computation, so reading the code of loo is IMHO a good place to start. Or do you want to start with the paper and going from there?
2 Likes
Do you mean you want to take the pointwise estimates in the loo object (loo_object$pointwise) and reproduce the point estimate and (approximate) standard error? Or do you want to reproduce the pointwise estimates themselves?
If you just want to take the pointwise estimates already calculated by loo and calculate the single point estimate and standard error for yourself (instead of relying on the package) then that’s just a few lines:
where in your case you’d replace x in the code above with:
x <- loo_object$pointwise[, "elpd_loo", drop=FALSE]
That will give you the total elpd and approximate standard error. That’s for a single model. It’s a bit different if you want the se of the the difference of elpd between two models. (Let me know if that’s what you were actually looking for.)
Thanks. It’s getting trickier now that the loo package functionality is continuing to expand (and there’s more coming soon!), but I’ve tried to keep the code as readable as possible! Always happy to answer questions about the code too.
2 Likes
I shall make a read on it thank u
Is just a sample mean of the point estimates? Ok I will check my math.
I am just making a package and have the idea to include the looic estimates (mean sd, and CI95%) on the default model summary. So It will be analogous to the bic and aic estimates that some classic package prints
The elpd point estimate reported by loo is the sum of the pointwise elpd contributions. You could report something like elpd_loo (total elpd) ± a few se. If you want it on the deviance scale you can use looic which is just multiplying the pointwise elpds by -2 before computing the se. That said none of these numbers (including aic and the others) have any interpretation except for relative to other models so if you’re including this in a package I recommend focusing on / encouraging the loo::loo_compare() results rather than the results for individual models (although you can show both of course).
Also you can get all these point estimates precomputed for you via loo_object$estimates.
Hope that helps!
2 Likes
Ok thanks with that feedback I think I can solve it, yes I know that buy themselves they dont say anything at all, I just want to create some kind of equivalent forecast package, so I just wanna present, Blooic just as the forecast present bic and aic.