How to perform MRP using stan

I have 4 demographic features: state, age , gender and urban. i want to perform a multi-level regression and poststratification. The codes is as follows,

mrp <- '
data {
  int<lower=0> N;
  array[N] int<lower=11, upper=65> state; 
  array[N] int<lower=1, upper=6> age;
  array[N] int<lower=1, upper=2> gender;
  array[N] int<lower=1, upper=2> urban;
  array[N] int<lower=0> y;
  array[31,6,2,2] int<lower=0> P;
}

parameters {
  real alpha;
  real<lower=0> sigma_beta;
  vector<multiplier=sigma_beta>[31] beta;
  real<lower=0> sigma_gamma;
  vector<multiplier=sigma_gamma>[6] gamma;
  real delta;
  real epsilon;
}

model {
  y ~ bernoulli_logit(alpha + beta[state] + gamma[age] + [delta,-delta][gender]' + [epsilon, -epsilon][urban]');
  alpha ~ normal(0, 2);
  beta ~ normal(0, sigma_beta);
  gamma ~ normal(0, sigma_gamma);
  {sigma_beta,sigma_gamma} ~ normal(0, 2);
  delta ~ normal(0, 2);
  epsilon ~ normal(0, 2);
}

generated quantities {
  real expect_pos = 0;
  int total = 0;
  for (b in 1:31) {
    for (c in 1:6) {
      for (d in 1:2) {
        for (e in 1:2) {
        total += P[b, c, d, e];
        expect_pos += P[b, c,d,e] * inv_logit(alpha + beta[b] + gamma[c] + delta[d] + epsilon[e] );
      }
    }

  real<lower=0, upper=1> phi = expect_pos / total;
}

' 

There is a problems. when adding ’ at the end of two binary groups,gender and urban, Rstudio give the message that there is unexpected [,] in the following line.

y ~ bernoulli_logit(alpha + beta[state] + gamma[age] + [delta,-delta][gender]' + [epsilon, -epsilon][urban]');

any help is appreciated!

i fixed my codes and set two coefficients for dichotomous variables,gender and urban. The Stan can run. the codes are as follows.

parameters {
  real alpha;
  real<lower=0> sigma_beta;
  real<lower=0> sigma_gamma;
  real<lower=0> sigma_delta;
  real<lower=0> sigma_epsilon;
  vector<multiplier=sigma_beta>[31] beta; 
  vector<multiplier=sigma_gamma>[6] gamma;
  vector<multiplier=sigma_delta>[2] delta;
  vector<multiplier=sigma_epsilon>[2] epsilon;
}

model {
  y ~ bernoulli_logit(alpha + beta[state] + gamma[age] + delta[gender] + epsilon[urban]);
  alpha ~ normal(0, 2);
  beta ~ normal(0, sigma_beta);
  gamma ~ normal(0, sigma_gamma);
  delta ~ normal(0, sigma_delta);
  epsilon ~ normal(0, sigma_epsilon);
{sigma_beta,sigma_gamma,sigma_delta,sigma_epsilon} ~ normal(0, 1);
}

however,i think the discussion in the stan user is right. there should be only one coefficient for two groups.
how to settle the issue?
any help is appreciated.

If there are two groups, you can give both groups a coefficient, but you either want to pin one value to zero or you want to enforce that they sum to zero—either way you get the desired one degree of freedom out.

This is just the “MR” part of “MRP”.

thank you for helping me out

i want to add group-level predictors in MRP, such as average income in a state.
in the data for stan, i specify income is a vector with size of 31,

income = gdp$log_gdp
 [1] 12.0 11.5 10.8 10.8 11.2 11.0 10.8
 [8] 10.7 12.0 11.7 11.5 11.1 11.6 10.9
[15] 11.2 10.9 11.2 11.0 11.4 10.7 10.9
[22] 11.3 11.0 10.7 10.9 10.9 11.1 10.5
[29] 10.8 10.9 10.9

in the data block, i add a line,

array[31] real<lower=0> income;

in the data block, i add a line,

real psi;

in the model block, i add a line,

y ~ bernoulli_logit(alpha + beta[age] + ... + delta[state]
+ income[state] * psi);
psi ~ normal(0, 2);

the result shows that the first line in the model block is wrong.
the message is that Ill-typed arguments supplied to infix operator *.
what’s wrong with my code?
thanks in advane.

with the help of AI,doubao, i had settled the issue.
AI use a for loop to calculate the linear predictor for each observation explicitly.

vector[N] linear_predictor;
  for (i in 1:N) {
    linear_predictor[i] = alpha + beta[state[i]] + gamma[age[i]] + 
                          delta[gender[i]] + epsilon[urban[i]] + 
                          income[state[i]] * psi + edu[state[i]] * eta ;
  }
  
  y ~ bernoulli_logit(linear_predictor);