I have a simple problem whereby I want to put a prior on the covariance matrix of multivariate normal, for that I will use a LKJ prior on the corresponding correlation matrix, and pre and post multiply by a vector of standard deviations, as discussed here and elsewhere. Specifically I will use the lkj_corr_cholesky
prior the Cholesky factor of the correlation matrix.
Whenever I run this, even in simple toy problems, I get the same warnings: The largest R-hat is NA ...
, and Bulk Effective Samples Size (ESS) is too low ...
, and Tail Effective Samples Size (ESS) is too low ...
.
Here is an example of a simple model (adopted from code in this blog).
data {
int<lower=1> N;
int<lower=1> J;
vector[J] y[N];
vector[J] Zero;
}
parameters {
cholesky_factor_corr[J] L_omega;
vector<lower=0>[J] sigma;
}
model {
y ~ multi_normal_cholesky(Zero, diag_pre_multiply(sigma, L_omega));
sigma ~ cauchy(0, 5);
L_omega ~ lkj_corr_cholesky(1.0);
}
Assuming that the above code is stored in a string in R named stancode
, I can then run it with the following R code.
library(tidyverse)
library(rstan)
set.seed(10101)
J <- 3
N <- 100
Sigma <- rbind(c(1.58, 0.30, -0.40),
c(0.30, 0.67, -1.04),
c(-0.40, -1.04, 1.86)
)
Zero <- rep(0, J)
y <- mvtnorm::rmvnorm(N, Zero, Sigma)
M <- stan(model_code = stancode,
data = list(J = J, N = N, y = y, Zero = Zero))
The warnings I get here, which are typical of whenever I use the lkj_corr_cholesky
are as follows:
Warning messages:
1: The largest R-hat is NA, indicating chains have not mixed.
Running the chains for more iterations may help. See
http://mc-stan.org/misc/warnings.html#r-hat
2: Bulk Effective Samples Size (ESS) is too low, indicating posterior means and medians may be unreliable.
Running the chains for more iterations may help. See
http://mc-stan.org/misc/warnings.html#bulk-ess
3: Tail Effective Samples Size (ESS) is too low, indicating posterior variances and tail quantiles may be unreliable.
Running the chains for more iterations may help. See
http://mc-stan.org/misc/warnings.html#tail-ess
Looking at the summary information for NA
s, we can see the following, again which is typical.
summary(M)$summary %>%
as_tibble(rownames = 'par') %>%
filter_all(any_vars(is.na(.)))
# A tibble: 4 x 11
par mean se_mean sd `2.5%` `25%` `50%` `75%` `97.5%` n_eff Rhat
<chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
1 L_omega[1,1] 1 NaN 0 1 1 1 1 1 NaN NaN
2 L_omega[1,2] 0 NaN 0 0 0 0 0 0 NaN NaN
3 L_omega[1,3] 0 NaN 0 0 0 0 0 0 NaN NaN
4 L_omega[2,3] 0 NaN 0 0 0 0 0 0 NaN NaN
I presume that the warnings above can and should be ignored and are simply as a consequence of the fact that the Cholesky factor is a lower triangle matrix whose row vectors have unit length, and so the upper triangle is fixed at zero and the first element of the first row is fixed at 1.0. In other words, these elements have fixed values at 0 or 1 hence a sd
of 0 etc.
However, the lkj_corr_cholesky
prior is used routinely in Stan based multilevel models, such as those created by brms
, and I never see these warnings.
More generally, assuming that these warnings are to be expected and are inconsequential, it feels like they are necessarily false alarms and so it seems like it would be better if they could be suppressed, though I don’t know of any way to do that using e.g. rstan::stan
.