data {

   int N;
   int y[N];
   real roach1[N];
   real treatment[N];
   real senior[N];
   real exposure2[N];

}

transformed data {

   real offset[N];

   for (i in 1:N) {
       offset[i] <- log(exposure2[i]);
   }

}

parameters {

   real b0;
   real b_roach1;
   real b_treatment;
   real b_senior;
}

transformed parameters {

  real loglambda[N];

  for (i in 1:N) {
       loglambda[i] <- b0 + b_roach1*roach1[i] + b_treatment*treatment[i] +
                    b_senior*senior[i] + offset[i];
  }

}

model {

  for (i in 1:N) {
       y[i] ~ poisson_log(loglambda[i]);
  }

  b0          ~ normal(0,1e6);
  b_roach1    ~ normal(0,1e6);
  b_treatment ~ normal(0,1e6);
  b_senior    ~ normal(0,1e6);

}

generated quantities {

  real y_pred[N];

  for (i in 1:N) {
      y_pred[i] <- poisson_rng(exp(loglambda[i]));
  }

}

CXXFLAGS = -O2 -Wall $(DEBUGFLAG) -mfpmath=sse -msse2 -mstackrealign $(LTO)