gqs() still has some issues with the names of posterior draws

[skolenik]

Summary:

gqs() has issues with "structured" versions of matrix-valued parameters and/or transformed parameters, I am not sure which.

Relates to #695 reported and closed in 2019 by @jgabry.

Description:

corr_matrix parameter declarations results in "subscript out of bounds error"

Reproducible Steps:

mc1 <- " parameters { vector[2] tau_u;
corr_matrix[2] rho_matrix; vector[2] y; } transformed parameters { cov_matrix[2] Sigma;

Sigma = quad_form_diag(rho_matrix,tau_u); } model { rho_matrix ~ lkj_corr(1); tau_u ~ cauchy(0,1); y ~ multi_normal([0,0]', Sigma); } " m1 <- stan_model(model_code = mc1) f1 <- sampling(m1, iter = 300)

mc2 <- " parameters {matrix[2,2] Sigma;} generated quantities {vector[2] yy_new = multi_normal_rng( [0,0]', Sigma );} " m2 <- stan_model(model_code = mc2) f2 <- gqs(m2, draws = as.matrix(f1))

mc3 <- " parameters { // variance components vector[2] tau_u;
// correlation matrix corr_matrix[2] rho_matrix;

} transformed parameters { cov_matrix[2] Sigma;

// covariance matrix Sigma = quad_form_diag(rho_matrix,tau_u); } generated quantities {vector[2] yy_new = multi_normal_rng( [0,0]', Sigma );} " m3 <- stan_model(model_code = mc3) f3 <- gqs(m3, draws = as.matrix(f1))

f4 <- gqs(m3, draws = as.matrix(f1)[, c("tau_u[1]", "tau_u[2]", "rho_matrix[1,1]", "rho_matrix[1,2]", "rho_matrix[2,1]", "rho_matrix[2,2]")] )

Current Output:

f3 <- gqs(m3, draws = as.matrix(f1)) Wrong number of parameter values in draws from fitted model. Expecting 6 columns, found 10 columns.

f4 <- gqs(m3, draws = as.matrix(f1)[, c("tau_u[1]", "tau_u[2]", "rho_matrix[1,1]", "rho_matrix[1,2]", "rho_matrix[2,1]", "rho_matrix[2,2]")] ) Error in draws[, draws_colnames, drop = FALSE] : subscript out of bounds

Expected Output:

Stuff should run?..

RStan Version:

packageVersion("rstan") [1] ‘2.21.3’

R Version:

R.version.string [1] "R version 4.1.2 (2021-11-01)"

Operating System:

Windows Server 2016 Standard OS Build 14393.5006

[skolenik]

I can work around in my application with just the second version that declares Sigma and picks it up. But I originally coded the third version... by copy-paste of the original model code... and it broke.

[wds15]

HI!

I can confirm that constrained types are a problem as it seems. Also arrays of matrices cause problems as it looks.

Here is an example:

library(OncoBayes2)

example_model("combo2")

## now try to run the gqs block of the model using only the actually
## sampled parameters

model <- OncoBayes2:::stanmodels$blrm_exnex

## extract what is actually sampled in the parameters block, not more
mpost <- as.matrix(blrmfit$stanfit, pars=c("log_beta_raw", "eta_raw", "mu_log_beta", "tau_log_beta_raw", "L_corr_log_beta", "mu_eta", "tau_eta_raw", "L_corr_eta"))

res <- gqs(model_hacked, list2env(blrmfit$standata), mpost, seed=2545)

This is with rstan 2.21.2