I fit an ordinal regression model with brms and ran loo using add_criterion. This returned 132 observations with pareto k values greater than 0.7 (out of 217,883 observations). I had fit the model with save_pars(all=TRUE), so I re-ran loo with moment matching, but got the following error: Error in if (quantities_i$ki < ki) { : missing value where TRUE/FALSE needed.
I don’t know what this means and haven’t found a discussion of this error. I’ve included my code and the warnings and errors below, along with session information. I can’t provide a reproducible example, as the data set is not shareable and too large in any case, but please let me know if there is additional information that would be helpful to have.
I don’t think multiple refits with reloo will be a realistic option, as this model took about 2.7 hours to fit, it is much simpler than the models I actually want to fit, and the data set I actually want to use is larger. Is there a way I can get moment matching to work?
library(brms)
library(cmdstanr)
# Model formula
form = course.grade ~ instruction.mode + incoming.gpa +
(1|instructor.id) + (1|student.id)
# Set priors
n.cat = 13
threshold.priors = qlogis(cumsum(rep(1/n.cat, n.cat)))[1:(n.cat-1)] %>%
round(., 3)
threshold.priors = map_df(seq_along(threshold.priors),
~prior_string(paste0("normal(",threshold.priors[.x],", 1)"),
class="Intercept", coef=.x))
priors = c(prior(normal(0,0.5), class="b"),
prior(normal(0,1), class="b", coef="incoming.gpa"),
threshold.priors)
# Fit model
m = brm(data = mdat, file="models/learning-mode/two-terms/tt1",
family = cumulative(link="logit", threshold="flexible"),
formula = form, prior = priors,
chains = 3, cores = 3, threads=threading(3),
backend="cmdstanr", sample_prior=TRUE, save_pars=save_pars(all=TRUE))
m = add_criterion(m, criterion="loo", moment_match=FALSE)
Automatically saving the model object in ‘models/learning-mode/two-terms/tt1.rds’
Warning messages:
1: The global prior ‘student_t(3, 0, 2.5)’ of class ‘Intercept’ will not be used in the model as all related coefficients have individual priors already. If you did not set those priors yourself, then maybe brms has assigned default priors. See ?set_prior and ?get_prior for more details.
2: Found 132 observations with a pareto_k > 0.7 in model ‘model’. It is recommended to set ‘moment_match = TRUE’ in order to perform moment matching for problematic observations.
loo(m)
Computed from 3000 by 217883 log-likelihood matrix Estimate SE elpd_loo -373018.3 497.2 p_loo 27621.2 72.0 looic 746036.6 994.3 ------ Monte Carlo SE of elpd_loo is NA. Pareto k diagnostic values: Count Pct. Min. n_eff (-Inf, 0.5] (good) 215836 99.1% 186 (0.5, 0.7] (ok) 1915 0.9% 68 (0.7, 1] (bad) 130 0.1% 88 (1, Inf) (very bad) 2 0.0% 15 See help('pareto-k-diagnostic') for details.
# Given the warnings above, redo loo with moment_match=TRUE
m = add_criterion(m, criterion="loo", moment_match=TRUE, overwrite=TRUE)
Recompiling the model with ‘rstan’
Recompilation done
Error in if (quantities_i$ki < ki) { :
missing value where TRUE/FALSE needed
In addition: Warning message:
The global prior ‘student_t(3, 0, 2.5)’ of class ‘Intercept’ will not be used in the model as all related coefficients have individual priors already. If you did not set those priors yourself, then maybe brms has assigned default priors. See ?set_prior and ?get_prior for more details.
Error: Moment matching failed. Perhaps you did not set ‘save_pars = save_pars(all = TRUE)’ when fitting your model? If you are running moment matching on another machine than the one used to fit the model, you may need to set recompile = TRUE.
Session Info
Macbook Pro M1 Max, Ventura 13.2.1
R 4.2.2
loo 2.6.0
m$version
$brms
[1] ‘2.19.0’
$rstan
[1] ‘2.26.13’
$stanHeaders
[1] ‘2.26.13’
$cmdstanr
[1] ‘0.5.2’
$cmdstan
[1] ‘2.31.0’