# Multilevel Model for Regression in STAN [closed]

I have, for N_countries countries and for each day during a certain period, a list of purchases from my website and how much each purchase was. I think that money spent as a function of days since the start of this period follows a sigmoid:

$$y \sim \mathcal N\left(\frac{\alpha}{1 + e^{-(\beta\cdot day + \gamma)}} + \delta, \sigma\right)$$

I want to model this according to a hierarchical or multilevel model, in which $\alpha$ and $\delta$ are drawn from distributions whose hyperparametres are drawn from another distribution (my population level is the country). The way I did it is, however, incredibly computationally expensive:

data {
int<lower = 1> N_people_day;
int<lower = 1> N_countries;
int<lower = 1, upper = N_countries> country[N_people_day];
int<lower = 0> dsi[N_people_day];
vector<lower = 0>[N_people_day] net_revenue;
}

parameters {
// parametres
vector<lower = 0>[N_countries] scale;                       // the sigmoid's scale parametre
real<lower = 0> rate;                                       // rate parametre
real<lower = 0> x_location;                                 // location parametre
vector<lower = 0>[N_countries] y_location;                  // location parametre
real<lower = 0> sigma;                                      // dispersion parametre

// hyperparametres
real<lower = 0> y_location_hyper_shape;
real<lower = 0> y_location_hyper_rate;
real<lower = 0> scale_hyper_shape;
real<lower = 0> scale_hyper_rate;
}

transformed parameters {
vector[N_people_day] mu;

for(i in 1:N_people_day){
mu[i] = scale[country[i]] / (1 + exp(-(rate * dsi[i] + x_location))) + y_location[country[i]];
}
}

model {
net_revenue ~ normal(mu, sigma);
y_location ~ gamma(y_location_hyper_shape, y_location_hyper_rate);
x_location ~ gamma(10, 0.5);
scale ~ gamma(scale_hyper_shape, scale_hyper_rate);
}


I'm almost certain there's a better way to do it and that the biggest problem is that one for in the transformed parameters block. How should I do it instead, though?

mu = scale[country] ./ (1 + exp(-(rate * dsi + x_location))) + y_location[country];