Question about Log probability evaluates to log(0), i.e. negative infinity

Modeling
,
1

Hello all!
I’m trying to fit recurrent events with a Bayesian survival model with frailty term which specifies the baseline hazard as Weibull distribution. But currently I meet a problem as the title. My STAN code is shown as following.

data {
int<lower=0> N; // number of subject
int<lower=0> N_row; // number of rows in the dataset
int<lower=0> K; // number of observed events
int<lower=0> Cn; // number of censored observations
int ncov; // number of covariates
matrix[K,ncov] Xobs; // covariate matrix for observed event
matrix[Cn,ncov] Xcen; // covariate matrix for censored event
vector[K] tobs; // time for observed event
vector[Cn] tcen; // time for censored event
int<lower=1,upper=N> IDobs[K];
int<lower=1,upper=N> IDcen[Cn];

}

parameters {
real<lower = 0> nu; // nu = shape
real<lower = 0> theta; // lambda=scale
real<lower = 0> sigmaw; // variance parameter for frailty term
vector[N] w; // frailty term
real beta1;
real beta2;
}

model {
real lambda_Tmp1;
real lambda_obs;
real lambda_Tmp2;
real lambda_cen;
// prior
nu ~ gamma(0.01,0.01);
theta ~ normal(0,sqrt(1000));
beta1 ~ normal(0,sqrt(1000));
beta2 ~ normal(0,sqrt(1000));
w ~ gamma(1/sigmaw,1/sigmaw); // my guess: the part caused the error message
sigmaw ~ inv_gamma(0.01, 0.01);

// likelihood - observed times:
for (i in 1:K){
lambda_Tmp1 = w[IDobs[K]] + beta1Xobs[i,1] + beta2Xobs[i,2] + log(theta);
lambda_obs = exp(-lambda_Tmp1/nu);
target += weibull_lpdf(tobs[i] | nu, lambda_obs);
}

// likelihood - censored times:
for (i in 1:Cn){
lambda_Tmp2 = w[IDcen[Cn]] + beta1Xcen[i,1] + beta2Xcen[i,2] + log(theta);
lambda_cen = exp(-lambda_Tmp2/nu);
target += weibull_lccdf(tcen[i] | nu, lambda_cen);
}

}

I don’t understand why it gives error message as below.
Chain 1: Log probability evaluates to log(0), i.e. negative infinity.
Chain 1: Stan can’t start sampling from this initial value.
Chain 1:
Chain 1: Initialization between (-2, 2) failed after 100 attempts.
Chain 1: Try specifying initial values, reducing ranges of constrained values, or reparameterizing the model.
[1] “Error in sampler$call_sampler(args_list[[i]]) : Initialization failed.”
[1] “error occurred during calling the sampler; sampling not done”

Even though I tried giving different initial values, it still doesn’t work. I changed the distribution of theta, it doesn’t work as well.(since I think the log(0) may be caused by theta.) However, interestingly, if I replaced the w ~ gamma(1/sigmaw,1/sigmaw) by w ~ normal(0,sqrt(sigmaw)). The whole code could work, but it will take around 8 to 12 hours to finish running (and a series of w values do not converge).

Could anyone help me with this question? Why does the gamma distribution of frailty term w cause the problem?

Thank you very much in advanced!

1 Like
2

Presumably some of your w's are getting initialized to be negative (by default they’re getting uniformly initialized on the interval [-2,2]), which doesn’t play nicely with a gamma prior but is fine for a normal prior. If you think w is gamma distributed, you need to declare it with a lower bound of zero.

2 Likes
3

Hi, jsocolar

Thanks for your response! Since I’m still thinking using Gamma distribution for w. I’m wondering how to declare the w with lower bound as 0 since I have made w as vector, can I define the lower bound for the vector?

Thanks!

1 Like
4

Yup! Just use the same syntax that you’d use for a scalar, i.e. vector<lower=0>[N] w;.

1 Like
5

Hi, jsocolar

Thank you so much for your quick response! I will let you know if it works after the definition of the lower bound!

Thank you! :)

6

Hi, jsocolar

After I add lower bound as you provided, it works now!!

Thank you so much!

1 Like