# Generating random numbers from a conditional distribution

Hello,

Everyone! I am new to Stan and am trying to use it to fit a Weibull distribution. I got my model correct to get posterior estimates. Now I need to get some predicted values based on a conditional Weibull, i.e., conditional on the predicted value is greater than say ti (i=1, …, q_n) for q_n subjects. I tried the generated quantities block below, but it takes forever to run, how can I simplify it?

``````data {
int<lower=0> n;                                // number of subjects
real adt[n];                                        // survival time
int<lower=0, upper=1> evt[n];          // event indicator
int<lower=0> q_n;                            // number of subjects with administrative censor
}
parameters {
real<lower=0> alpha;
real<lower=0> sigma;
}
model {
for (i in 1:n){
}
}
generated quantities{
for (j in 1:q_n){
adt_pred[j] = weibull_rng(alpha, sigma);
}
}
}
``````
1 Like

hello again,

I am able to create a function to make generate the conditional weibull, with the code below:

``````functions {
real weibull_trunc_rng(real alpha, real sigma, real lb){
real lp = weibull_cdf(lb, alpha, sigma);
// real up = weibull_cdf(ub, alpha, sigma);
real u = uniform_rng(lp, 1);
real y = sigma * (-log1m(u))^inv(alpha);
return y;
}
}
...
generated quantities{
// real predt;
for (j in 1:q_n){
}
}
``````

But when I add the commented out codes to generate more, i got the following error message:
Error in sampler\$call_sampler(args_list[[i]]) :
Exception: Exception: uniform_rng: Upper bound parameter is 1, but must be greater than 1 (in ‘model3e5f28701a7_evt_pred_weibull’ at line 5)
(in ‘model3e5f28701a7_evt_pred_weibull’ at line 33)

Any help on what is causing this and how to solve it will be appreciated!!

1 Like

Hi @Lei_Hua, welcome to the Stan forums! Would it be possible to share the data you’re using (or some simulated data that also works) as well as a the full Stan program that results in the error? I think that will help me or someone else debug this. Thanks!

Sure, I can provide the data. What is the easiest way of providing the data?

This happens when I work on something similar. I think the issue is about the function generating u.

Thanks. When responding to this post, one of the editing-related buttons at the top will let you upload a csv file. Or you can just share some R or Python code for generating data that demonstrates the same problem. Whatever you prefer.

Also one thing I just noticed is that `adt_q` doesn’t seem to be defined anywhere in the code you shared (unless I missed it). I see `adt` and `adt_pred`, but there’s an `adt_q` variabled used that I can’t find the definition of.

Thanks, @jonah, I have the adt_q defined the data block. I find a misspecification of my model. Fixing the models seems solved the problem now.

1 Like

Ok great, glad you solved it!

1 Like

My issue comes back, i’ve uploaded the data as an R objectmydata.csv (24.4 KB)

``````functions {
real exponential_trunc_rng(real lambda, real lb){// function to generate truncated exponential;
real lp = exponential_cdf(lb, lambda);
real u = uniform_rng(lp, 1);
real y = -log1m(u)/lambda;
return y;
}
int fst_index(vector x, real y) {// function to find the first positio of matches;
real tmp = x[1];
int pos = 1;
while(tmp<y){// to find location
pos += 1;
tmp = x[pos];
}
return pos;
}
}
data {
int<lower=0> n;
int<lower=-1, upper=1> evt[n];
int<lower=0> q_n;
real<lower=0> enrdt[q_n];
}
parameters {
// for weibull event time;
real<lower=0> lambda_e;
// for exponential discontinuation time;
real<lower=0> lambda_d;
}
model {
for (i in 1:n){
}
}
generated quantities{
real evt_pred[q_n];
real disc_pred[q_n];
real evtdt_pred[q_n];
real evtdt_pred_order[q_n];
vector[q_n] ind;
vector[q_n] ind_order;
real predt;
for (j in 1:q_n){
evt_pred[j] = exponential_trunc_rng(lambda_e, adt_q[j]);
disc_pred[j] = exponential_trunc_rng(lambda_d, adt_q[j]);
adt_pred[j] = fmin(evt_pred[j], disc_pred[j]);
ind[j] = (evt_pred[j]<=disc_pred[j]);
evtdt_pred[j] = enrdt[j] + adt_pred[j];
}
evtdt_pred_order=evtdt_pred[sort_indices_asc(evtdt_pred)];
ind_order=ind[sort_indices_asc(evtdt_pred)];
// predt = evtdt_pred_order[min(which_equal(cumulative_sum(ind_order), addevt))];
predt = evtdt_pred_order[fst_index(cumulative_sum(ind_order), addevt)];
}
``````

below is my R code to get this work in R

mydata_mod <- list(
n = nrow(mydata),
evt = mydata\$evt,
# to simulate posterior dist and get estimates;
q_n = sum(mydata\$evt==0),
enrdt = mydata\$enrdt[mydata\$evt==0],
)
fit.exp <- stan(
file = “/usrfiles/shared-projects/clinicalInformatics/rstan/LeiHua_demo/Lei’s code/event prediction/bayesian prediction/oneArm/draft/evt_pred_exp.stan”, # Stan program
data = mydata_mod, # named list of data
# seed=1234,
chains = 4, # number of Markov chains
warmup = 1000, # number of warmup iterations per chain
iter = 5000, # total number of iterations per chain
cores = 2, # number of cores (could use one per chain)
refresh = 0 # no progress shown
)

When I run this with 4 chains, 3 chains do not have samples, while 1 has. It gave me the following errors:

[1] “Error in sampler\$call_sampler(args_list[[i]]) : "
[2] " Exception: Exception: uniform_rng: Upper bound parameter is 1, but must be greater than 1 (in ‘model68ff229467f6_evt_pred_exp’ at line 4)”
[3] " (in ‘model68ff229467f6_evt_pred_exp’ at line 50)"

It seems the function I defined using runif is causing such. How can I modify to make it work?

Sorry, I just saw this. I think this was solved in Error with uniform_rng, right? If so, that’s great. If not then let me know and I can take another look.

Yes, thanks!