# Problem with application of integrate_1d

Dear all,
In my model, I need to compute an integral (using integrate_1d) with the same function with several different upper limits. But, I encountered with following error in rstan:

## ’ SYNTAX ERROR, MESSAGE(S) FROM PARSER: Expressions for elements of array must have the same or promotable types; found first element type=real; type at position 2=vector error in ‘model390c40e3ef1_Onicesco_interval_independentmodel_Log_logistic34744’ at line 70, column 74

``````70:     log_integ[i]=log(1-integrate_1d(integrand, lb, ub[i],{ alpha, lambdaa[i]}, {0.0} , {0}));
71:     }
``````

'.

This is my codes in the RSTAN . Is there any idea how to solve this problem?

``````functions{
real integrand(real u, real xc, real[] theta, real[] x_r, int[] x_i)
{     // data (integer)

real df_dx = ((theta[2] *theta[1]*u^(theta[1]-1))/((1+(theta[2] *u^theta[1]))^2));
return(df_dx);
}
}
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;
real lb;
}
transformed data {
real x_r[0];
int x_i[0];
vector <lower=lb>[N] ub;
for (i in 1:N){
ub[i]=respon[i,4];
}
}
parameters{
vector [num] beta;
real <lower=0> sigma;
//real <lower=0> sigma1;
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
vector[N] log_integ;
time_ratio  = exp(beta);
exp_x = exp(x);
lambdaa = exp ((-1 * respon[,7:18] * beta) / sigma);
for (i in 1:N) {
log_integ[i]=log(1-integrate_1d(integrand, lb, ub[i],{ alpha, lambdaa[i]}, {0.0} , {0}));
}
{ // `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;
target += normal_lpdf(x|0,5);
target += normal_lpdf(beta | 0, 1000);
target += gamma_lpdf(alpha | 0.001, 0.001);
for (i in 1:N) {
f[i]= a[county[i]]*((lambdaa [i]*alpha*respon[i,4]^(alpha-1))/((1+(lambdaa [i]*respon[i,4]^alpha))^2));
//s[i]= a[county[i]]*(1- ((lambdaa [i]*respon[i,4]^alpha)/ (1+(lambdaa [i]*respon[i,4]^alpha)))); s[i]= a[county[i]]  * (log_integ[i]);
target += respon[i,6] *  (log(f[i])-log(s[i])) +  log(s[i]);
}
}
.`````````

Your model code compiled successfully for me after I explicitly specified the default relative tolerance:

``````log_integ[i]=log(1-integrate_1d(integrand, lb, ub[i],{ alpha, lambdaa[i]}, {0.0} , {0}, 1e-8));
``````
1 Like

Thanks a ton for your reply and help. Actually, I’m new to Stan and I had trouble running my hierarchical survival model (about ~ 8 days for iter=10000, warm up=5000, N=35000, province=31, num=12).

How can I run my model using within-chain parallelization method(25.1 Reduce-sum | Stan User’s Guide)? Can you help me to change my codes in RSTAN?

``````functions{
real icar_normal_lpdf(vector phi, int province, int[] node1, int[] node2) {
return -0.5 * dot_self(phi[node1] - phi[node2]);
}
// defines the log survival
vector log_S (vector t, real shape,vector scale){
vector[num_elements(t)] log_S;
for (i in 1:num_elements(t)){
log_S[i] =  log(1- ((scale [i]*(t[i])^shape)/ (1+(scale [i]*(t[i])^shape)))); // county is a vector that shows the location of patients with dim = num_elements(t)
}
return log_S;
}
//defines the log pdf
vector log_f (vector t,real shape,vector scale){
vector[num_elements(t)] log_f;
for (i in 1:num_elements(t)){
log_f[i] =log((scale[i]*shape*t[i]^(shape-1))/((1+(scale[i]*t[i]^shape))^2));
}
return log_f;
}
//defines the log-likelihood for right censored data
real surv_loglogistic_lpdf(vector t,vector d, real shape,vector scale){
vector[num_elements(t)] log_lik;
real prob;
log_lik = d .* log_f(t, shape, scale)+ (1 - d) .* log_S(t, shape, scale);
prob = sum(log_lik);
return prob;
}
real partial_sum(array[] int y_slice,
int start, int end,vector d,
real shape, vector scale) {
return surv_loglogistic_lpdf(y_slice | d, shape, scale);
}
}
//data block
data{
int <lower=0> province;
int <lower=0> N;
int <lower=0> num;
int <lower=0> number;
int <lower =0 , upper = province> county[N];
matrix [N,number] respon;
int<lower=0> N_edges;
int<lower=1, upper=province> node1[N_edges];  // node1[i] adjacent to node2[i]
int<lower=1, upper=province> node2[N_edges];  // and node1[i] < node2[i]
}
//parameters block
parameters{
vector [num] beta;
real <lower=0> sigma;//scale parameter sigma=1/shape
vector <lower=0> [province] sigma1;
vector[province] w1;  // random effect (CAR)
vector[province] theta;       // heterogeneous effects
}
// transformed parameters block
transformed parameters{
vector [num] time_ratio;
//vector[N] linpred;
vector[N] lambdaa;
real <lower=0> alpha = 1 / sigma;
time_ratio  = exp(beta);
lambdaa = exp ((-respon[,7:18] * beta - w1[county] - theta[county]) / sigma);
}
// model block
model{
target += normal_lpdf(beta | 0, 1000);
//target += student_t_lpdf(sigma| 3,0,1);
target += gamma_lpdf(alpha | 0.001, 0.001);
w1 ~ icar_normal_lpdf(province, node1, node2);
sum(w1) ~ normal(0, 0.001 * province);
target += normal_lpdf(theta| 0, sigma1);
target += student_t_lpdf(sigma1| 3,0,1);
respon[,4] ~ surv_loglogistic//model for the data
}
// generated quantities block
generated quantities{
vector[N] log_lik;//log-likelihood

{ for(n in 1:N){
log_lik[n] = respon[n,6] * (log((exp ((-respon[n,7:18] * beta - w1[county[n]] - theta[county[n]])*alpha)*alpha*respon[n,4]^(alpha-1)) / ((1+(exp ((-respon[n,7:18] * beta - w1[county[n]] - theta[county[n]])*alpha)*(respon[n,4])^(alpha)))^2))) + (1-respon[n,6]) * ( log(1- ((exp ((-respon[n,7:18] * beta - w1[county[n]] - theta[county[n]])*alpha)*(respon[n,4])^alpha) / (1+(exp ((-respon[n,7:18] * beta - w1[county[n]] - theta[county[n]])*alpha)*(respon[n,4])^alpha)))));}
}
}
.`````````