Credibility intervals for one group

Modeling
, , , ,
1

Hello all,
I would like to have a good interpret for the credible interval for one treatment group in the clinical trial, what does that mean? and how I can compute that for different sample sizes with rstanarm package.

I simulated the event and treatment variables then I applied stan_glm with different prior. Does this true?

Thanks in advance

pi = 0.2 # event rate
n= c(30, 50, 100)
j=1

event = rep(1, n[j])
n0 = round(n[j]*pi,0)
event[sample(n[j],n0)] = 0 
mydata= data.frame(event = event, trt= rep(1, length(event
                                                     )))

glm_prior <- normal(location = 0, scale = 0.2, autoscale = FALSE )
model.or <- stan_glm(event~trt, data = mydata,
                     family =  "binomial",
                     prior = glm_prior,
                     prior_intercept = NULL,
                     #prior_aux = prior_aux,
                     QR=FALSE,
                     mean_PPD = FALSE,
                     seed = 123450,
                     chains=10,
                     #iter= 500,
                     #cores=ncores,
                     refresh = 0)
post = as.data.frame(model.or)

post = apply(post, 2, function(x) exp(x))

ci.beta.l = round(apply(post, 2, function(x) quantile(x, 0.025) ), 2) 
ci.beta.u = round(apply(post, 2, function(x) quantile(x, 0.975) ), 2)

CI = cat("(", ci.beta.l, ",",ci.beta.u, ")" )
2 Likes
2

Hi,
I don’t think I understand your question very well. Do you have a specific dataset you are trying to analyze or are you trying to use simulations to determine sample size for a study that has not yet been run (i.e. some sort of power analysis)?

3

I want to study the behavior of credible intervals as the sample size increases when I have just one treatment group. As you saw, I simulated one treatment group (trt) and the event var (event). Does the credible interval I got for n= 30, e.g ( 1.82, 11.46 ), is true or not? I meant the code construction is true or no? you can use any simulated data. Is it clear?

4

So the goal is just understanding how the model behaves with increasing sample size?

The way you call the model looks mostly OK. You are putting prior on model coefficients other than the intercept, but your model has only the intercept, so the prior is not used. The credible interval for the odds looks correct to me.

Best of luck with your model.

1 Like
5

yes, but just for one treatment group? does my simulation is ok?
Yes, I missed adding the prior to the intercept, it is also in my interest to see the model behaves with different types of prior.

library(rstanarm)

pi = 0.2 # event rate
n= c(30, 50, 100)
res = NULL
for (i in 1:length(n)) {
  
  
  
  event = rep(1, n[i])
  n0 = round(n[i]*pi,0)
  event[sample(n[i],n0)] = 0 
  mydata= data.frame(event = event, trt= rep(1, length(event
  )))
  
  glm_prior <- normal(location = 0.3, scale = 0.2, autoscale = FALSE )
  
  model.or <- stan_glm(event~trt, data = mydata,
                       family =  "binomial",
                       prior = glm_prior,
                       prior_intercept = glm_prior,
                       QR=FALSE,
                       mean_PPD = FALSE,
                       seed = 123450,
                       chains=10,
                       refresh = 0)
  post = as.data.frame(model.or)
  
  post = apply(post, 2, function(x) exp(x))
  
  ci.beta.l = round(apply(post, 2, function(x) quantile(x, 0.025) ), 2) 
  ci.beta.u = round(apply(post, 2, function(x) quantile(x, 0.975) ), 2)
  
  res = rbind(res,c( n[i],  LCI = ci.beta.l, UCI = ci.beta.u))
  
}

>> res

      n              LCI           UCI
[1,]  30            1.18            2.36
[2,]  50            1.32            2.58
[3,] 100            1.61            2.90

```

I noticed there is no significant difference as the sample size increase. Do you agree with me?
6

I am sorry, but I don’t think I am able to understand what is your goal. Could you be more precise and provide some background? Does “one treatment group” mean that you have just one group of patients or does it mean that you have one control group and one treatment group? The simulation you’ve shown is OK for some goals but not for others, so it is hard to say without more information.

1 Like
7

I meant one group of patients.