Correct way to code interactions

591 views
Skip to first unread message

Daniel Padfield

unread,
Jun 1, 2017, 5:23:59 AM6/1/17
to stan-...@googlegroups.com
Hi all

Been looking for some documentation on coding interactions in Stan models and haven't found anything explicitly explaining things. Any recommendations for where to look?

My current model looks at how flux of an entire community is influenced by the individuals within it. Consequently I have multiple masses (K) for each rate (N).

Finally came back to this after 2 months and managed to get it working but if I have a group factor for each sample how would I code the interaction on each parameter, especially that a is wrapped up in log_sum_exp(a * segment(mass, pos, mass_n[n]))

Here is the model as it is, not overly complex apart from the log_sum_exp() part.

stan_model <- '
data{
  int N;                              // number of samples
  int K;                              // number of mass samples
  int mass_n[N];                      // number of masses per sample
  vector[N] ln_rate;                  // rate
  vector[N] temp;                     // temperature
  vector[K] mass;                     // masses
}
parameters{  
  real<lower=0> sigma;      // standard deviation
  real ln_r0;               // community normalisation constant
  real bT;                  // effect of temperature
  real a;                   // size scaling exponent
}
model{
  int pos;
  vector[N] mu;
  vector[N] ln_biomass;
  ln_r0 ~ normal(0,10);
  bT ~ normal(0,10);
  a ~ normal(0,10);
  sigma ~ cauchy(0,2);


  pos = 1;
  for(n in 1:N){
    ln_biomass[n] = log_sum_exp(a * segment(mass, pos, mass_n[n]));
    pos = pos + mass_n[n];
  }  


  # likelihood regression for size dependence of rate
  mu = ln_r0 + ln_biomass + bT*temp;
  ln_rate ~  normal(mu, sigma);


}'

For ln_r0, I can see I can just add another vector for the group (0 or 1) and add a

 
ln_r0 + ln_r0_inter * group+ ...

into the model block

And for bT, I have read that I should included a transformed data block:

transformed data{
vector
[N] temp_inter;
temp_inter
= temp .* ancestral;
}

and then add that into the model block as

... bT*temp + bT_inter*temp_inter...

I have ran this and I am happy these are doing what I want them to. However when it comes to adding in an interaction to allow an interaction between "a" and "group" I do not really know where to start.

I have a reproducible gist of this model using rstan. It has the data created for "group" but it has not been incorporated into the model yet...

Cheers
Dan

Bob Carpenter

unread,
Jun 2, 2017, 1:14:08 PM6/2/17
to stan-...@googlegroups.com
I didn't understand what "group" was in your notation. When you
say "specify another vector", what exactly did you mean---a parameter,
data, or what?

Can you write out the model in mathematical terms, specifically what
you want mu to be in

mu = ln_r0 + ln_biomass + bT*temp;
ln_rate ~ normal(mu, sigma);

Translating to Stan shouldn't be a problem.

- Bob

P.S. We're shutting this list down in a week or two; the list
has moved to:

http://discourse.mc-stan.org

You'll find more people there now than here to answer your question.
Feel free (nay, encouraged) to cross-post.

> On Jun 1, 2017, at 5:23 AM, Daniel Padfield <d.padfi...@gmail.com> wrote:
>
> Hi all
>
> Been looking for some documentation on coding interactions in Stan models and haven't found anything explicitly explaining things. Any recommendations for where to look?
>
> My current model looks at how flux of an entire community is influenced by the individuals within it. Consequently I have multiple masses (K) for each rate (N).
>
> Finally came back to this after 2 months and managed to get it working but if I have a group factor for each sample how would I code the interaction on each parameter, especially that a is wrapped up in log_sum_exp(a * segment(mass, pos, mass_n[n]))



You just need to code up the linear predictor. So I'm not sure what
you're having problem with doing this. I also have problems with
words like "factor" because it's all tied up with stats terminology I never
understood and R conventions.
> --
> 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.

Reply all
Reply to author
Forward
0 new messages