Generate new quantity using stanvars

Interfaces brms
1

Hi there,

I am confused on how to set up the names of the parameters for calculating new quantities. As you can see below, the parameter name I used ‘b_reactgm_c_group‘ is exactly the one stan used, then I still got this error. Please help!

> library(brms)

Loading 'brms' package (version 2.23.0). Useful instructions
can be found by typing help('brms'). A more detailed introduction
to the package is available through vignette('brms_overview').


The following object is masked from ‘package:stats’:

    ar
> f1=bf(nagatie_body_image_zjf~ c_judge + c_react+ c_time+react_gm+judge_gm+c_group+c_gender+c_age+(1|id))
> f2=bf(react_gm~c_group+c_qp+int)
> f3=bf(judge_gm~c_group+c_qp+int)
> 
> p=get_prior(f1+f2+f3+set_rescor(F),data=dat)
> # flat prior
> p
                  prior     class      coef group                resp dpar nlpar lb ub tag       source
                 (flat)         b                             judgegm                           default
                 (flat)         b   c_group                   judgegm                      (vectorized)
                 (flat)         b      c_qp                   judgegm                      (vectorized)
                 (flat)         b       int                   judgegm                      (vectorized)
   student_t(3, 0, 2.5) Intercept                             judgegm                           default
   student_t(3, 0, 2.5)     sigma                             judgegm             0             default
                 (flat)         b                 nagatiebodyimagezjf                           default
                 (flat)         b     c_age       nagatiebodyimagezjf                      (vectorized)
                 (flat)         b  c_gender       nagatiebodyimagezjf                      (vectorized)
                 (flat)         b   c_group       nagatiebodyimagezjf                      (vectorized)
                 (flat)         b   c_judge       nagatiebodyimagezjf                      (vectorized)
                 (flat)         b   c_react       nagatiebodyimagezjf                      (vectorized)
                 (flat)         b    c_time       nagatiebodyimagezjf                      (vectorized)
                 (flat)         b  judge_gm       nagatiebodyimagezjf                      (vectorized)
                 (flat)         b  react_gm       nagatiebodyimagezjf                      (vectorized)
 student_t(3, 3.2, 2.5) Intercept                 nagatiebodyimagezjf                           default
   student_t(3, 0, 2.5)        sd                 nagatiebodyimagezjf             0             default
   student_t(3, 0, 2.5)        sd              id nagatiebodyimagezjf             0        (vectorized)
   student_t(3, 0, 2.5)        sd Intercept    id nagatiebodyimagezjf             0        (vectorized)
   student_t(3, 0, 2.5)     sigma                 nagatiebodyimagezjf             0             default
                 (flat)         b                             reactgm                           default
                 (flat)         b   c_group                   reactgm                      (vectorized)
                 (flat)         b      c_qp                   reactgm                      (vectorized)
                 (flat)         b       int                   reactgm                      (vectorized)
   student_t(3, 0, 2.5) Intercept                             reactgm                           default
   student_t(3, 0, 2.5)     sigma                             reactgm             0             default

> fit_flat_prior=brm(f1+f2+f3+set_rescor(F),data=dat,iter=10000,cores = 12,backend='cmdstanr',silent=2) 

> fit_flat_prior
 Family: MV(gaussian, gaussian, gaussian) 
  Links: mu = identity
         mu = identity
         mu = identity 
Formula: nagatie_body_image_zjf ~ c_judge + c_react + c_time + react_gm + judge_gm + c_group + c_gender + c_age + (1 | id) 
         react_gm ~ c_group + c_qp + int 
         judge_gm ~ c_group + c_qp + int 
   Data: dat (Number of observations: 86) 
  Draws: 4 chains, each with iter = 10000; warmup = 5000; thin = 1;
         total post-warmup draws = 20000

Multilevel Hyperparameters:
~id (Number of levels: 43) 
                                  Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(nagatiebodyimagezjf_Intercept)     0.48      0.08     0.34     0.65 1.00     5826    10265

Regression Coefficients:
                              Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
nagatiebodyimagezjf_Intercept     3.13      0.08     2.97     3.30 1.00     9062    12030
reactgm_Intercept                 0.02      0.05    -0.07     0.11 1.00    37909    14542
judgegm_Intercept                 0.01      0.05    -0.09     0.11 1.00    32985    12754
nagatiebodyimagezjf_c_judge       0.32      0.12     0.09     0.55 1.00    32814    15552
nagatiebodyimagezjf_c_react       0.27      0.10     0.07     0.47 1.00    25462    16128
nagatiebodyimagezjf_c_time        0.26      0.10     0.05     0.46 1.00    24076    15933
nagatiebodyimagezjf_react_gm      0.74      0.19     0.36     1.11 1.00     9343    11736
nagatiebodyimagezjf_judge_gm      0.28      0.19    -0.09     0.65 1.00    10472    12396
nagatiebodyimagezjf_c_group      -0.02      0.18    -0.38     0.34 1.00     9890    12395
nagatiebodyimagezjf_c_gender     -0.31      0.20    -0.70     0.07 1.00     9665    12039
nagatiebodyimagezjf_c_age        -0.01      0.04    -0.08     0.06 1.00     9950    11906
reactgm_c_group                   0.35      0.09     0.16     0.53 1.00    35374    14338
reactgm_c_qp                     -0.17      0.04    -0.25    -0.08 1.00    33561    14843
reactgm_int                      -0.29      0.09    -0.47    -0.12 1.00    34503    13811
judgegm_c_group                   0.34      0.10     0.14     0.53 1.00    32713    14003
judgegm_c_qp                     -0.12      0.05    -0.22    -0.03 1.00    33935    14375
judgegm_int                      -0.19      0.09    -0.38    -0.01 1.00    34639    13808

Further Distributional Parameters:
                          Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sigma_nagatiebodyimagezjf     0.36      0.04     0.29     0.46 1.00     7180     9966
sigma_reactgm                 0.42      0.03     0.36     0.50 1.00    31294    14195
sigma_judgegm                 0.46      0.04     0.39     0.54 1.00    31911    14828

Draws were sampled using sample(hmc). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).

> med_stanvars=stanvar(scode='real med=b_reactgm_c_group * b_nagatiebodyimagezjf_react_gm;',block='genquant')

> fit_flat_prior=brm(f1+f2+f3+set_rescor(F),stanvars=med_stanvars,data=dat,iter=10000,cores = 12,backend='cmdstanr',silent=2) 
Error : An error occured during compilation! See the message above for more information.
Compiling Stan program...
Semantic error in 'C:/Users/HUAWEI/AppData/Local/Temp/RtmpC2svtE/model-c80778b34d1.stan', line 104, column 11 to column 28:
   -------------------------------------------------
   102:    // actual population-level intercept
   103:    real b_judgegm_Intercept = Intercept_judgegm - dot_product(means_X_judgegm, b_judgegm);
   104:    real med=b_reactgm_c_group * b_nagatiebodyimagezjf_react_gm;
                    ^
   105:  }
   106:  
   -------------------------------------------------

Identifier 'b_reactgm_c_group' not in scope.
make: *** [make/program:66: C:/Users/HUAWEI/AppData/Local/Temp/RtmpC2svtE/model-c80778b34d1.hpp] Error 1
Error: An error occured during compilation! See the message above for more information.

> names(as.data.frame(fit_flat_prior))
 [1] "b_nagatiebodyimagezjf_Intercept"         "b_reactgm_Intercept"                     "b_judgegm_Intercept"                    
 [4] "b_nagatiebodyimagezjf_c_judge"           "b_nagatiebodyimagezjf_c_react"           "b_nagatiebodyimagezjf_c_time"           
 [7] "b_nagatiebodyimagezjf_react_gm"          "b_nagatiebodyimagezjf_judge_gm"          "b_nagatiebodyimagezjf_c_group"          
[10] "b_nagatiebodyimagezjf_c_gender"          "b_nagatiebodyimagezjf_c_age"             "b_reactgm_c_group"                      
[13] "b_reactgm_c_qp"                          "b_reactgm_int"                           "b_judgegm_c_group"                      
[16] "b_judgegm_c_qp"                          "b_judgegm_int"                           "sd_id__nagatiebodyimagezjf_Intercept"   
[19] "sigma_nagatiebodyimagezjf"               "sigma_reactgm"                           "sigma_judgegm"                          
[22] "Intercept_nagatiebodyimagezjf"           "Intercept_reactgm"                       "Intercept_judgegm"     

> sessionInfo()$running
[1] "Windows 11 x64 (build 26100)"

> packageVersion("brms")
[1] ‘2.23.0’

> packageVersion("cmdstanr")
[1] ‘0.9.0’

> library(cmdstanr)
This is cmdstanr version 0.9.0
- CmdStanR documentation and vignettes: mc-stan.org/cmdstanr
- CmdStan path: C:/Users/HUAWEI/.cmdstan/cmdstan-2.37.0
- CmdStan version: 2.37.0
2

Can show us the generated stancode that brms produces to for your model?

3

sure, I paste it below.

> stancode(fit_flat_prior)
// generated with brms 2.23.0
functions {
}
data {
  int<lower=1> N;  // total number of observations
  int<lower=1> N_nagatiebodyimagezjf;  // number of observations
  vector[N_nagatiebodyimagezjf] Y_nagatiebodyimagezjf;  // response variable
  int<lower=1> K_nagatiebodyimagezjf;  // number of population-level effects
  matrix[N_nagatiebodyimagezjf, K_nagatiebodyimagezjf] X_nagatiebodyimagezjf;  // population-level design matrix
  int<lower=1> Kc_nagatiebodyimagezjf;  // number of population-level effects after centering
  int<lower=1> N_reactgm;  // number of observations
  vector[N_reactgm] Y_reactgm;  // response variable
  int<lower=1> K_reactgm;  // number of population-level effects
  matrix[N_reactgm, K_reactgm] X_reactgm;  // population-level design matrix
  int<lower=1> Kc_reactgm;  // number of population-level effects after centering
  int<lower=1> N_judgegm;  // number of observations
  vector[N_judgegm] Y_judgegm;  // response variable
  int<lower=1> K_judgegm;  // number of population-level effects
  matrix[N_judgegm, K_judgegm] X_judgegm;  // population-level design matrix
  int<lower=1> Kc_judgegm;  // number of population-level effects after centering
  // data for group-level effects of ID 1
  int<lower=1> N_1;  // number of grouping levels
  int<lower=1> M_1;  // number of coefficients per level
  array[N_nagatiebodyimagezjf] int<lower=1> J_1_nagatiebodyimagezjf;  // grouping indicator per observation
  // group-level predictor values
  vector[N_nagatiebodyimagezjf] Z_1_nagatiebodyimagezjf_1;
  int prior_only;  // should the likelihood be ignored?
}
transformed data {
  matrix[N_nagatiebodyimagezjf, Kc_nagatiebodyimagezjf] Xc_nagatiebodyimagezjf;  // centered version of X_nagatiebodyimagezjf without an intercept
  vector[Kc_nagatiebodyimagezjf] means_X_nagatiebodyimagezjf;  // column means of X_nagatiebodyimagezjf before centering
  matrix[N_reactgm, Kc_reactgm] Xc_reactgm;  // centered version of X_reactgm without an intercept
  vector[Kc_reactgm] means_X_reactgm;  // column means of X_reactgm before centering
  matrix[N_judgegm, Kc_judgegm] Xc_judgegm;  // centered version of X_judgegm without an intercept
  vector[Kc_judgegm] means_X_judgegm;  // column means of X_judgegm before centering
  for (i in 2:K_nagatiebodyimagezjf) {
    means_X_nagatiebodyimagezjf[i - 1] = mean(X_nagatiebodyimagezjf[, i]);
    Xc_nagatiebodyimagezjf[, i - 1] = X_nagatiebodyimagezjf[, i] - means_X_nagatiebodyimagezjf[i - 1];
  }
  for (i in 2:K_reactgm) {
    means_X_reactgm[i - 1] = mean(X_reactgm[, i]);
    Xc_reactgm[, i - 1] = X_reactgm[, i] - means_X_reactgm[i - 1];
  }
  for (i in 2:K_judgegm) {
    means_X_judgegm[i - 1] = mean(X_judgegm[, i]);
    Xc_judgegm[, i - 1] = X_judgegm[, i] - means_X_judgegm[i - 1];
  }
}
parameters {
  vector[Kc_nagatiebodyimagezjf] b_nagatiebodyimagezjf;  // regression coefficients
  real Intercept_nagatiebodyimagezjf;  // temporary intercept for centered predictors
  real<lower=0> sigma_nagatiebodyimagezjf;  // dispersion parameter
  vector[Kc_reactgm] b_reactgm;  // regression coefficients
  real Intercept_reactgm;  // temporary intercept for centered predictors
  real<lower=0> sigma_reactgm;  // dispersion parameter
  vector[Kc_judgegm] b_judgegm;  // regression coefficients
  real Intercept_judgegm;  // temporary intercept for centered predictors
  real<lower=0> sigma_judgegm;  // dispersion parameter
  vector<lower=0>[M_1] sd_1;  // group-level standard deviations
  array[M_1] vector[N_1] z_1;  // standardized group-level effects
}
transformed parameters {
  vector[N_1] r_1_nagatiebodyimagezjf_1;  // actual group-level effects
  // prior contributions to the log posterior
  real lprior = 0;
  r_1_nagatiebodyimagezjf_1 = (sd_1[1] * (z_1[1]));
  lprior += student_t_lpdf(Intercept_nagatiebodyimagezjf | 3, 3.2, 2.5);
  lprior += student_t_lpdf(sigma_nagatiebodyimagezjf | 3, 0, 2.5)
    - 1 * student_t_lccdf(0 | 3, 0, 2.5);
  lprior += student_t_lpdf(Intercept_reactgm | 3, 0, 2.5);
  lprior += student_t_lpdf(sigma_reactgm | 3, 0, 2.5)
    - 1 * student_t_lccdf(0 | 3, 0, 2.5);
  lprior += student_t_lpdf(Intercept_judgegm | 3, 0, 2.5);
  lprior += student_t_lpdf(sigma_judgegm | 3, 0, 2.5)
    - 1 * student_t_lccdf(0 | 3, 0, 2.5);
  lprior += student_t_lpdf(sd_1 | 3, 0, 2.5)
    - 1 * student_t_lccdf(0 | 3, 0, 2.5);
}
model {
  // likelihood including constants
  if (!prior_only) {
    // initialize linear predictor term
    vector[N_nagatiebodyimagezjf] mu_nagatiebodyimagezjf = rep_vector(0.0, N_nagatiebodyimagezjf);
    mu_nagatiebodyimagezjf += Intercept_nagatiebodyimagezjf;
    for (n in 1:N_nagatiebodyimagezjf) {
      // add more terms to the linear predictor
      mu_nagatiebodyimagezjf[n] += r_1_nagatiebodyimagezjf_1[J_1_nagatiebodyimagezjf[n]] * Z_1_nagatiebodyimagezjf_1[n];
    }
    target += normal_id_glm_lpdf(Y_nagatiebodyimagezjf | Xc_nagatiebodyimagezjf, mu_nagatiebodyimagezjf, b_nagatiebodyimagezjf, sigma_nagatiebodyimagezjf);
    target += normal_id_glm_lpdf(Y_reactgm | Xc_reactgm, Intercept_reactgm, b_reactgm, sigma_reactgm);
    target += normal_id_glm_lpdf(Y_judgegm | Xc_judgegm, Intercept_judgegm, b_judgegm, sigma_judgegm);
  }
  // priors including constants
  target += lprior;
  target += std_normal_lpdf(z_1[1]);
}
generated quantities {
  // actual population-level intercept
  real b_nagatiebodyimagezjf_Intercept = Intercept_nagatiebodyimagezjf - dot_product(means_X_nagatiebodyimagezjf, b_nagatiebodyimagezjf);
  // actual population-level intercept
  real b_reactgm_Intercept = Intercept_reactgm - dot_product(means_X_reactgm, b_reactgm);
  // actual population-level intercept
  real b_judgegm_Intercept = Intercept_judgegm - dot_product(means_X_judgegm, b_judgegm);
}
1 Like
4

As you can see in the generated Stan code, there is no parameter, transformed parameter, or generated quantity named b_reactgm_c_group. The naming conventions in the prior summary are not the same as the internal variables that brms uses in the stan code. Presumably what you want is one of the elements of the vector parameter b_reactgm.

1 Like