Hello,
I’m having trouble implementing the sum-to-zero constraint of a dichotomous predicted variable with multiple nominal predictors model. I’m attempting to create something of a fusion of the models Kruschke implements in chapters 20 and 21 of Doing Bayesian Data Analysis.
What’s specifically giving me trouble is trying to take the mean of an N dimensional array, and I can’t seem to find a function or implementation that can do it without awkward looping or summing up the entire N dimensional array and diving by the number of elements.
I’ve attached my model code below. Thank you for the assistance!
data {
int<lower=1> nrows;
int<lower=1> x1_ncats;
int<lower=1, upper=x1_ncats> x1[nrows];
int<lower=1> x2_ncats;
int<lower=1, upper=x2_ncats> x2[nrows];
int<lower=1> x3_ncats;
int<lower=1, upper=x3_ncats> x3[nrows];
int<lower=0> y[nrows];
int<lower=1> N[nrows];
}
parameters {
real a0;
vector[x1_ncats] a1;
vector[x2_ncats] a2;
vector[x3_ncats] a3;
real sigma;
}
model {
a0 ~ normal(0, 2);
a1 ~ normal(0, sigma);
a2 ~ normal(0, sigma);
a3 ~ normal(0, sigma);
y ~ binomial_logit(N, a0 + a1[x1] + a2[x2] + a3[x3]);
}
generated quantities {
real m[x1_ncats, x2_ncats, x3_ncats];
real b0;
vector[x1_ncats] b1;
vector[x2_ncats] b2;
vector[x3_ncats] b3;
for (k1 in 1:x1_ncats)
for (k2 in 1:x2_ncats)
for (k3 in 1:x3_ncats)
m[k1, k2, k3] = a0 + a1[k1] + a2[k2] + a3[k3];
b0 = mean(m);
for (k1 in 1:x1_ncats) b1[k1] = mean(m[k1, 1:x2_ncats, 1:x3_ncats]) - b0;
for (k2 in 1:x2_ncats) b2[k2] = mean(m[1:x1_ncats, k2, 1:x3_ncats]) - b0;
for (k3 in 1:x3_ncats) b3[k3] = mean(m[1:x1_ncats, 1:x2_ncats, k3]) - b0;
}