Simplex to use in a dirichlet distribution (Hierarchical model)

1,918 views
Skip to first unread message

David Vivas

unread,
Mar 14, 2017, 1:50:50 AM3/14/17
to Stan users mailing list
Hi, I am trying to implement a hierarchical model (Multinomial - dirichlet), the dirichlet distribution receives as parameter a vector whose entries add up to one, I get this vector(to dirichlet) from a gamma (4,2) distribution but this values don´t add one (obviously), I tried to create a simplex to get this vector (add one) but I don't have the desire result.

Does anyone have any suggestions to solve this ?
Sorry if this is a simple question but I am a new user of Stan

Part of my model is:

parameters{
  real <lower=0> alfa11;
  real <lower=0> alfa12;
  real <lower=0> alfa13;
  real <lower=0> alfa21;
  real <lower=0> alfa22;
  real <lower=0> alfa23;
  real <lower=0> alfa31;
  real <lower=0> alfa32;
  real <lower=0> alfa33;
  real <lower=0> alfa41;
  real <lower=0> alfa42;
  real <lower=0> alfa43;
}
transformed parameters{

  matrix[4,3] alfas[R];   //array of matrices 4x3
  simplex [3]  alfarow1;  //vector simplex in R3
  simplex [3]  alfarow2;
  simplex [3]  alfarow3;
  simplex [3]  alfarow4;

  alfarow1[1]=alfa11;
  alfarow1[2]=alfa12;
  alfarow1[3]=alfa13;
  alfarow2[1]=alfa21;
  alfarow2[2]=alfa22;
  alfarow2[3]=alfa23;
  alfarow3[1]=alfa31;
  alfarow3[2]=alfa32;
  alfarow3[3]=alfa33;
  alfarow4[1]=alfa41;
  alfarow4[2]=alfa42;
  alfarow4[3]=alfa43;  
  for (i in 1:R){
    alfas[i,1]=alfarow1';
    alfas[i,2]=alfarow2';
    alfas[i,3]=alfarow3';
    alfas[i,4]=alfarow4';
  }
}

model{
  for (r in 1:R){
    for (i in 1:4){
      for (j in 1:3){
        alfas[r,i,j]~gamma(4,2);  // where "alfas" is an array with R matrices of size 4x3
      }
    } 
  }
  for (r in 1:R){
    for(i in 1:4){
      betas[r,i]'~dirichlet(alfas[r,i]');  //alfas[r,i] : simplex vector that not add one 
    }
  }

Thanks
David

Daniel Lee

unread,
Mar 14, 2017, 8:02:02 AM3/14/17
to stan-...@googlegroups.com
Hi David,

The easiest thing to do is to change alfas to a parameter like:
simplex[3] alfas[4]
And remove the alfa?? parameters. 

You might want to revisit the chapter on parameters and transformed parameters. In transformed parameters, the Stan programmer needs to ensure the constraints are satisfied. 



Daniel
--
You received this message because you are subscribed to the Google Groups "Stan users mailing list" group.
To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+...@googlegroups.com.
To post to this group, send email to stan-...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

Bob Carpenter

unread,
Mar 14, 2017, 1:30:30 PM3/14/17
to stan-...@googlegroups.com
As Daniel pointed out, you can try to code up the Dirichlet
yourself using a *normalized* set of gammas. But that's
sub-optimal for clarity and for efficiency in Stan.
We have everything you need built-in. A basic Dirichlet-multinomial
with N observations over K categories looks like this:

data {
int K; // num dims
int N; // num observations
int y[N, K]; // multinomial observations
}
parameters {
real<lower=0> alpha;
simplex[K] theta[N];
}
model {
alpha ~ lognormal(0, 5); // or some other prior on alpha
for (n in 1:N) {
theta[n] ~ dirichlet(alpha);
y[n] ~ multinomial(theta[n]);
}
}

It'd be much more efficient to code up the compound Dirichlet-multinomial
by marginalizing out the theta. That's what you get in the beta-binomial,
but we don't have Dirichlet-multinomial built in.

- Bob

P.S. You can use matrices rather than all those hand-coded variables
and do all the linear algebra you do row- or column-wise with matrix ops,
but you don't need to if you choose the dimensions caefully.

> On Mar 14, 2017, at 1:50 AM, David Vivas <vivad...@gmail.com> wrote:
>

David Vivas

unread,
Mar 14, 2017, 6:27:02 PM3/14/17
to stan-...@googlegroups.com
Dear Daniel, 

Thanks a lot ! I have solved my entire problem.
Thank you again!

David

To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+unsubscribe@googlegroups.com.

To post to this group, send email to stan-...@googlegroups.com.
For more options, visit https://groups.google.com/d/optout.

--
You received this message because you are subscribed to a topic in the Google Groups "Stan users mailing list" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/stan-users/oXqN8z77_xE/unsubscribe.
To unsubscribe from this group and all its topics, send an email to stan-users+unsubscribe@googlegroups.com.

David Vivas

unread,
Mar 14, 2017, 6:27:41 PM3/14/17
to stan-...@googlegroups.com
Dear Bob, thank you very much, I have implemented my model, 
Thank you for your support

David 

> To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+unsubscribe@googlegroups.com.
> To post to this group, send email to stan-...@googlegroups.com.
> For more options, visit https://groups.google.com/d/optout.

--

You received this message because you are subscribed to a topic in the Google Groups "Stan users mailing list" group.
To unsubscribe from this topic, visit https://groups.google.com/d/topic/stan-users/oXqN8z77_xE/unsubscribe.
To unsubscribe from this group and all its topics, send an email to stan-users+unsubscribe@googlegroups.com.
Reply all
Reply to author
Forward
0 new messages