 Predict in generated quantities for mixture models

I use stan to fit zero-inflated beta regression, which is mixture models of logistic regression and beta regression,

for (n in 1:N) {
if (y[n] == 0) {
target += bernoulli_logit_lpmf(1 | zi);
} else {
target += bernoulli_logit_lpmf(0 | zi)
+ beta_lpdf(y[n] | mu[n] * phi[n], (1 - mu[n]) * phi[n]);
}
}

The complete code is as follows

data {
int<lower=0> N;
int<lower=0> K;
matrix[N, K] X;
vector[N] y;
}
parameters {
vector[K] beta;
vector[K] beta_phi;
vector[K] beta_zi;
}
model {
vector[N] mu  = X * beta;
vector[N] phi = X * beta_phi;
vector[N] zi  = X * beta_zi;

for (n in 1:N) {
mu[n] = inv_logit(mu[n]);
phi[n] = exp(phi[n]);
}

for (n in 1:N) {
if (y[n] == 0) {
target += bernoulli_logit_lpmf(1 | zi);
} else {
target += bernoulli_logit_lpmf(0 | zi)
+ beta_lpdf(y[n] | mu[n] * phi[n], (1 - mu[n]) * phi[n]);
}
}

target += student_t_lpdf(beta | 3, 0, 2.5);
target += student_t_lpdf(beta_phi | 3, 0, 2.5);
target += logistic_lpdf(beta_zi | 0, 1);

}

The code above works fine. My question is how to write code in generated quantities { } to predict on new data, such as X_new，I do not know how to mix the bernoulli_logit_rng() and beta_rng()

generated quantities {
vector[M] y_predict;

vector[M] mu  = X_new * beta;
vector[M] phi = X_new * beta_phi;
vector[M] zi  = X_new * beta_zi;

for (i in 1:M) {
mu[i] = inv_logit(mu[i]);
phi[i] = exp(phi[i]);
}

// how to code ?
for(i in 1:M) {
//  bernoulli_logit_rng(zi);
//  beta_rng(mu[i] * phi[i], (1 - mu[i]) * phi[i]);
}

}

I don’t know if it’s correct

data {
int<lower=0> N;
int<lower=0> K;
matrix[N, K] X;
vector[N] y;
matrix[2, K] new_X;
}

parameters {
vector[K] beta;
vector[K] beta_phi;
vector[K] beta_zi;
}
model {
vector[N] mu  = X * beta;
vector[N] phi = X * beta_phi;
vector[N] zi  = X * beta_zi;

for (n in 1:N) {
mu[n] = inv_logit(mu[n]);
phi[n] = exp(phi[n]);
}

for (n in 1:N) {
if (y[n] == 0) {
target += bernoulli_logit_lpmf(1 | zi[n]);
} else {
target += bernoulli_logit_lpmf(0 | zi[n])
+ beta_lpdf(y[n] | mu[n] * phi[n], (1 - mu[n]) * phi[n]);
}
}

target += student_t_lpdf(beta | 3, 0, 2.5);
target += student_t_lpdf(beta_phi | 3, 0, 2.5);
target += logistic_lpdf(beta_zi | 0, 1);

}

generated quantities {
vector y_predict;

vector mu  = new_X * beta;
vector phi = new_X * beta_phi;
vector zi  = new_X * beta_zi;

for (i in 1:2) {
mu[i] = inv_logit(mu[i]);
phi[i] = exp(phi[i]);
}

for(i in 1:2) {
real tmp = bernoulli_logit_rng(zi[i]);
if (tmp == 1) {
y_predict[i] = 0;
} else {
y_predict[i] = beta_rng(mu[i] * phi[i], (1 - mu[i]) * phi[i]);
}
}

}