Dear all,

I’m new to Stan. The following paragraph is my first stan code for the hierarchical survival model. But I had trouble introducing the attached one-dimensional integral used to calculate the survival function part of likelihood.

Can anyone help me to compute the attached Integral in RSTAN?

Thank you in advance!

My stan code is as below (N=35000, n_grid=350, province=31, num=12).

fit ← stan(“final model.stan”,data = da,iter=iter, warmup=warmup, cores =4,chains = 1, thin = 10, init = “random”,control=list(max_treedepth=15))

Difinition of survival function.pdf (192.9 KB)

“”

data {

int <lower=0> province;

int <lower=0> n_grid;

int <lower=0> N;

int <lower=0> num;

int <lower=0> number;

int <lower =0 , upper = province> county[N];

vector <lower=0> [province] p;

row_vector[n_grid] kernel[province];

matrix [N,number] respon;

}

parameters{

vector [num] beta;

real <lower=0> sigma;

vector [n_grid] x;

}

transformed parameters{

vector [num] time_ratio;

vector [N] lambdaa;

real alpha = 1 / sigma;

vector [n_grid] exp_x;

vector[province] a; // `a`

saved to output because it’s a toplevel transformed parameter

time_ratio = exp(beta);

exp_x = exp(x);

lambdaa = exp (((-1 * respon[,7:17] * beta) / sigma));

// }

{ // `z`

not part of output

vector[province] landa;

for (k in 1 : province) {

landa[k] = kernel[k] * exp_x * p[k];

}

a = landa / sum(landa); // assign `a`

}

}

model{

vector [N] f;

vector [N] s;

real log_alpha;

log_alpha = log (alpha);

target += normal_lpdf(x|0,5);

target += normal_lpdf(beta | 0, 1000);

target += normal_lpdf(log_alpha | 0, 5);

target += -1 * log_alpha;

for (i in 1:N) {

f[i]= a[county[i]]*((lambdaa [i] alpharespon[i,4]^(alpha-1))/((1+(lambdaa [i]respon[i,4]^alpha))^2));*((lambdaa [i]

// How can we introduce the integral related to the calculation of the mentioned survival function in the RSTAN?

s[i]= a[county[i]](1- ???);

target += respon[i,6] * (log(f[i])-log(s[i])) + log(s[i]);

}

}

generated quantities{

vector [N] log_lik;

vector [N] s;

vector [N] f;

real log_alpha;

log_alpha = log (alpha);

for (i in 1:N){

f[i]= a[county[i]]

*alpha*respon[i,4]^(alpha-1))/((1+(lambdaa [i]

*respon[i,4]^alpha))^2));*

// How can we introduce the integral related to the calculation of the mentioned survival function in the RSTAN?

s[i]= a[county[i]](1- ???);

// How can we introduce the integral related to the calculation of the mentioned survival function in the RSTAN?

s[i]= a[county[i]]

log_lik [i]= respon[i,6] * (log(f[i])-log(s[i])) + log(s[i]) - log_alpha;

}

}

“”