# How to write generated quantities syntax from a joint density function

**URL:** <https://discourse.mc-stan.org/t/how-to-write-generated-quantities-syntax-from-a-joint-density-function/28860>\
**Category:** Modeling\
**Created:** [September 15, 2022, 11:01pm UTC](https://discourse.mc-stan.org/t/how-to-write-generated-quantities-syntax-from-a-joint-density-function/28860 "2022-09-15T23:01:26Z")\
**Posts on this page:** 2\
**Page:** 1

<div class="post-metadata">

**Author:** ![Cengiz\_Zopluoglu](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/cengiz_zopluoglu/32/2342_2.png) [@Cengiz\_Zopluoglu](https://discourse.mc-stan.org/u/Cengiz_Zopluoglu)\
**Post date:** [September 15, 2022, 11:01pm UTC](https://discourse.mc-stan.org/t/how-to-write-generated-quantities-syntax-from-a-joint-density-function/28860/1 "2022-09-15T23:01:26Z")

</div>

Hi,

I have a model syntax that successfully compiles and fits simulated data with no issues. I am trying to add the generated quantities block to compute model fit indices with the loo package.

Does anyone have any suggestions? Thank you.

The model syntax is below. Please ignore any redundancy in the code pointed out by @Bob_Carpenter in another post. I will deal with it later.

```nohighlight
data{
  int I; // number of items
  int P; // number of individuals
  int n_obs; // number of observed responses
  int item[n_obs]; // item id
  int id[n_obs]; // person id
  real Y[n_obs]; // vector outcome
}

parameters {
  real beta[I]; // vector of beta parameters for
  real ln_alpha[I]; // vector of alpha parameters for I items     
  real omega[I]; // vector of omega parameters for I items     
  real gamma0[I]; // vector of gamma0 parameters for I items     
  real gamma1[I]; // vector of gamma1 parameters for I items     
  real theta[P]; // vector of theta parameters for P individuals      
}

transformed parameters{
  real<lower=0> alpha[I];
  for(i in 1:I){
       alpha[i] = exp(ln_alpha[i]);
  } 
}
   
model{
    
  beta ~ normal(0,10);
  ln_alpha ~ normal(0,10);
  omega ~ normal(0,10);
  gamma0 ~ normal(0,10);
  theta ~ normal(0,1);      
  
  for(i in 1:I){
    gamma1[i] ~ normal(0,10)T[gamma0[i],];
  }
  
  for(i in 1:n_obs) {
    
    real k0 = inv_logit(gamma0[item[i]] - alpha[item[i]]*theta[id[i]]);
    real k1 = inv_logit(gamma1[item[i]] - alpha[item[i]]*theta[id[i]]);
  
    if (Y[i] == 0) {
      target += log(k0);
    } else if (Y[i] == 1) {
      target += log1m(k1);                    
    } else {
      
      real sh1 = exp((alpha[item[i]]*theta[id[i]] + beta[item[i]] + omega[item[i]])/2);
      real sh2 = exp((-(alpha[item[i]]*theta[id[i]] + beta[item[i]]) + omega[item[i]])/2);

      target += log(k1 - k0) + beta_lpdf(Y[i] | sh1,sh2);   
    }	

  }
}

```

---

<div class="post-metadata">

**Author:** ![caesoma](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/caesoma/32/19767_2.png) [@caesoma](https://discourse.mc-stan.org/u/caesoma)\
**Post date:** [September 19, 2022, 7:48pm UTC](https://discourse.mc-stan.org/t/how-to-write-generated-quantities-syntax-from-a-joint-density-function/28860/2 "2022-09-19T19:48:41Z")

</div>

> [@Cengiz\_Zopluoglu](#):
>
> to compute model fit indices with the loo package

Can you be more specific on what the issue is? I am not sure what the loo package requires (or what is the package, for that matter), but you can compute any quantity in the generated quantities block, and since it use samples from the posterior it will be from a joint density (as your question states). So at taken at face value your answer is simply to put whatever function you need to retrieve (that is not already in the transformed parameters) in that block and use it for whatever you need outside of the Stan program.
