# Ragged array of simplexes

**URL:** <https://discourse.mc-stan.org/t/ragged-array-of-simplexes/1382>\
**Category:** Modeling\
**Created:** [July 24, 2017, 11:41pm UTC](https://discourse.mc-stan.org/t/ragged-array-of-simplexes/1382 "2017-07-24T23:41:54Z")\
**Posts on this page:** 1\
**Showing post:** 4

<div class="post-metadata">

**Author:** ![bgoodri](https://yyz2.discourse-cdn.com/flex030/user_avatar/discourse.mc-stan.org/bgoodri/32/4451_2.png) [@bgoodri](https://discourse.mc-stan.org/u/bgoodri)\
**Post date:** [July 25, 2017, 3:49pm UTC](https://discourse.mc-stan.org/t/ragged-array-of-simplexes/1382/4 "2017-07-25T15:49:28Z")

</div>

That is quite possibly correct but more complicated than necessary. Dealing with simplexes of different sizes is annoying, but it it actually the least annoying of the ragged arrays. The trick to avoiding having to deal with Jacobians is to remember that a Dirichlet distribution over simplexes can be constructed by normalizing a set of unit-scale but independent Gamma random variables:

> **[Dirichlet distribution | Gamma distribution](https://en.wikipedia.org/wiki/Dirichlet_distribution#Gamma_distribution)**
>
> With a source of Gamma-distributed random variates, one can easily sample a random vector 
>   
>     
>       
> x
> =
> (
>         
> x
>           
> 1
>           
>         
> ,
> …
> ,
>         
> x
>           
> K
>           
>         
> )
>       
>     
> {\\displaystyle x=(x\_{1},\\ldots ,x\_{K})}
>   
> from the K-dimensional Dirichlet distribution with parameters 
>   
>     
>       
> (
>         
> α
>           
> 1...

  
Thus, in Stan you declare a looooong vector of non-negative parameters

```
parameters {
  vector<lower=0>[sum(K)] gamma; // K is an int[] declared in data
  ...
}

```

and in the model block you break off however many of those you need at the moment, divide them by their sum, and use the resulting simplex in your likelihood, as in

```
model {
  int pos = 1;
  for (j in 1:J) {
    vector[K[j]] pi = segment(gamma, pos, K[j]);
    pi = pi / sum(pi); // simplex
    target += // some function of pi
    pos = pos + K[j];
  }
}

```

Then you put a Gamma prior on gamma

```
  target += gamma_lpdf(gamma | alpha, 1);
   // implies pi ~ dirichlet(segment(alpha, pos, K[j]))

```

where `alpha` is a looooong vector of known shape hyperparameters for the Dirichlet distributions that are concatenated together. If `alpha = rep_vector(1, sum(K));` then each `pi` is uniform _a priori_ over simplexes of size `K[j]` and you can get away with

```
  target += exponential_lpdf(gamma | 1);
  // implies pi ~ uniform over simplexes

```

---

_[View the full topic](https://discourse.mc-stan.org/t/ragged-array-of-simplexes/1382)._