Hi everyone,
I am new to pystan (on MacOS X).
I am trying to fit a bayesian model developed in this paper:
https://static.googleusercontent.com/media/research.google.com/en//pubs/archive/46001.pdf
My model is the following:
   functions {
  real Hill(real t, real ec, real slope) {
    return 1 / (1 + (t / ec)^(-slope));
  }
  real Adstock(row_vector t, row_vector weights) {
    return dot_product(t, weights) / sum(weights);
  }
}
data {
  int<lower=1> N;
  real<lower=0> Y[N];
  int<lower=1> max_lag;
  row_vector[max_lag] lag_vec;
	row_vector[max_lag] x_tv[N];
	row_vector[max_lag] x_radio[N];
	row_vector[max_lag] x_ooh[N];
	row_vector[max_lag] x_incentive[N];
	real<lower=0> x_nb_users[N];
	real<lower=0> x_is_covid[N];
}
parameters {
	real<lower=0> beta_tv;
	real<lower=0,upper=1> ec_tv;
	real<lower=0> slope_tv;
	real<lower=0,upper=1> retain_rate_tv;
	real<lower=0,upper=max_lag-1> delay_tv;
	real<lower=0> beta_radio;
	real<lower=0,upper=1> ec_radio;
	real<lower=0> slope_radio;
	real<lower=0,upper=1> retain_rate_radio;
	real<lower=0,upper=max_lag-1> delay_radio;
	real<lower=0> beta_ooh;
	real<lower=0,upper=1> ec_ooh;
	real<lower=0> slope_ooh;
	real<lower=0,upper=1> retain_rate_ooh;
	real<lower=0,upper=max_lag-1> delay_ooh;
	real<lower=0> beta_incentive;
	real<lower=0,upper=1> ec_incentive;
	real<lower=0> slope_incentive;
	real<lower=0> gamma_ctrl_nb_users;
	real<lower=0> gamma_ctrl_is_covid;
  real<lower=0> tau;
  real<lower=0> noise_var;
}
transformed parameters {
  real mu[N];
  real cum_effect;
  row_vector[max_lag] lag_weights;
	row_vector[4] cum_effects_hill[N];
	row_vector[4] beta_medias;
	row_vector[2] cum_effects_linear[N];
	row_vector[2] gamma_ctrl;
	beta_medias[1] = beta_tv;
	beta_medias[2] = beta_radio;
	beta_medias[3] = beta_ooh;
	beta_medias[4] = beta_incentive;
	gamma_ctrl[1] = gamma_ctrl_nb_users;
	gamma_ctrl[2] = gamma_ctrl_is_covid;
  for (nn in 1:N) {
		for (lag in 1 : max_lag) {
			lag_weights[lag] <- pow(retain_rate_tv, (lag - 1 - delay_tv) ^ 2);
		}
		cum_effect <- Adstock(x_tv[nn], lag_weights);
		cum_effects_hill[nn, 1] <- Hill(cum_effect, ec_tv, slope_tv);
		for (lag in 1 : max_lag) {
			lag_weights[lag] <- pow(retain_rate_radio, (lag - 1 - delay_radio) ^ 2);
		}
		cum_effect <- Adstock(x_radio[nn], lag_weights);
		cum_effects_hill[nn, 2] <- Hill(cum_effect, ec_radio, slope_radio);
		for (lag in 1 : max_lag) {
			lag_weights[lag] <- pow(retain_rate_ooh, (lag - 1 - delay_ooh) ^ 2);
		}
		cum_effect <- Adstock(x_ooh[nn], lag_weights);
		cum_effects_hill[nn, 3] <- Hill(cum_effect, ec_ooh, slope_ooh);
		cum_effect <- x_incentive[nn, 1];
		cum_effects_hill[nn, 4] <- Hill(cum_effect, ec_incentive, slope_incentive);
		cum_effects_linear[nn, 1] <- x_nb_users[nn];
		cum_effects_linear[nn, 2] <- x_is_covid[nn];
    mu[nn] <- tau + dot_product(cum_effects_hill[nn], beta_medias) + dot_product(cum_effects_linear[nn], gamma_ctrl);
  }
}
model {
	slope_tv ~ gamma(3, 1);
	ec_tv ~ beta(2,2);
	beta_tv ~ normal(0.2, 1);
	retain_rate_tv ~ beta(3,3);
	delay_tv ~ uniform(0, max_lag - 1);
	slope_radio ~ gamma(3, 1);
	ec_radio ~ beta(2,2);
	beta_radio ~ normal(0.17, 1);
	retain_rate_radio ~ beta(3,3);
	delay_radio ~ uniform(0, max_lag - 1);
	slope_ooh ~ gamma(3, 1);
	ec_ooh ~ beta(2,2);
	beta_ooh ~ normal(0.05, 1);
	retain_rate_ooh ~ beta(3,3);
	delay_ooh ~ uniform(0, max_lag - 1);
	slope_incentive ~ gamma(3, 1);
	ec_incentive ~ beta(2,2);
	beta_incentive ~ normal(0.5, 1);
	gamma_ctrl_nb_users ~ normal(0, 0.5);
	gamma_ctrl_is_covid ~ normal(0, 0.25);
  tau ~ normal(0, 1);
  noise_var ~ inv_gamma(0.05, 0.05 * 0.01);
  Y ~ normal(mu, sqrt(noise_var));
}
I am getting the following error and my mu is getting nan value.
Rejecting initial value:
Error evaluating the log probability at the initial value.
Exception: normal_lpdf: Location parameter[1] is nan, but must be finite!
(in 'unknown file name' at line 141)
Initialization between (-2, 2) failed after 100 attempts.
Try specifying initial values, reducing ranges of constrained values, or reparameterizing the model.
My print statements in the loop are having the following outputs:
cum_effects_hill: [0,0.0265629,0,0.681911]
cum_effects_linear: [1,1]
dot_product(cum_effects_hill[nn], beta_medias): 1.17828
dot_product(cum_effects_linear[nn], gamma_ctrl): 4.99926
tau: 0.227982
mu (calculated): 6.40552
mu: nan
I am not really sure how to edit my code so that it doesn’t yield a nan for mu.
I would appreciate the community helps on this.
Thanks a lot!