##
## model 6:
##     dixon and Simon
##
##     theta_g ~ N(theta, sigma^2)
##         theta = b0+b1*X1+..+bp*Xp
##          b0  ~ N(0, 10^3)
##          b1..bp ~ N(0, omega^2)
##          omega^2 ~ HalfN(0,1)

data {
	int<lower=0> SIZE;
	int<lower=0> NX;
	vector[SIZE] Y;
	vector[SIZE] VY;
	matrix[SIZE, NX] X;
	real<lower=0> VTAU;
	real<lower=0> VW;
	real VRANGE[2];
}

parameters {
    vector[NX] beta;
    real b0;
    real<lower=0> omega;
	  real<lower=0, upper=1> uvs[SIZE];
}

transformed parameters{
	real<lower=0> vs[SIZE];
	vector[2] range;

	for (i in 1:SIZE) {
		range[1] = log(VY[i]) + VRANGE[1];
		range[2] = log(VY[i]) + VRANGE[2];
		vs[i]    = exp(range[1] + uvs[i]*(range[2] - range[1]));
	}
}

model {
  b0    ~ normal(0, sqrt(VTAU));
  beta  ~ normal(0, sqrt(omega));
  uvs   ~ uniform(0,1);

  if (0 == VW) {
    target += -log(omega);
  } else {
    omega ~ normal(0, sqrt(VW));
  }

  for (i in 1:SIZE) {
 		 Y[i] ~ normal(b0+X[i]*beta, sqrt(vs[i]));
  }
}

generated quantities {
	real mu[SIZE];
	for (i in 1:SIZE) {
		mu[i] = b0+X[i]*beta;
	}
}
