Minor inconsistency in the vignette #9

DanChaltiel · 2021-12-13T10:31:10Z

Hi,

Thanks for this great package.

I'm following the vignette found in https://cran.r-project.org/web/packages/trialr/vignettes/A200-CRM.html.

I ran its code with no problem (sorry, I couldn't find how to manage the verbosity):

library(trialr)
#> Warning: le package 'trialr' a été compilé avec la version R 4.1.2
#> Le chargement a nécessité le package : Rcpp
target <- 0.25
skeleton <- c(0.05, 0.15, 0.25, 0.4, 0.6)
fit1 <- stan_crm("2NN 3NN 4TT", skeleton=skeleton, target=target,
                model='empiric', beta_sd=sqrt(1.34), seed=123, verbose=FALSE)
#> 
#> SAMPLING FOR MODEL 'CrmEmpiricNormalPrior' NOW (CHAIN 1).
#> Chain 1: 
#> Chain 1: Gradient evaluation took 0 seconds
#> Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
#> Chain 1: Adjust your expectations accordingly!
#> Chain 1: 
#> Chain 1: 
#> Chain 1: Iteration:    1 / 2000 [  0%]  (Warmup)
#> Chain 1: Iteration:  200 / 2000 [ 10%]  (Warmup)
#> Chain 1: Iteration:  400 / 2000 [ 20%]  (Warmup)
#> Chain 1: Iteration:  600 / 2000 [ 30%]  (Warmup)
#> Chain 1: Iteration:  800 / 2000 [ 40%]  (Warmup)
#> Chain 1: Iteration: 1000 / 2000 [ 50%]  (Warmup)
#> Chain 1: Iteration: 1001 / 2000 [ 50%]  (Sampling)
#> Chain 1: Iteration: 1200 / 2000 [ 60%]  (Sampling)
#> Chain 1: Iteration: 1400 / 2000 [ 70%]  (Sampling)
#> Chain 1: Iteration: 1600 / 2000 [ 80%]  (Sampling)
#> Chain 1: Iteration: 1800 / 2000 [ 90%]  (Sampling)
#> Chain 1: Iteration: 2000 / 2000 [100%]  (Sampling)
#> Chain 1: 
#> Chain 1:  Elapsed Time: 0.023 seconds (Warm-up)
#> Chain 1:                0.021 seconds (Sampling)
#> Chain 1:                0.044 seconds (Total)
#> Chain 1: 
#> 
#> SAMPLING FOR MODEL 'CrmEmpiricNormalPrior' NOW (CHAIN 2).
#> Chain 2: 
#> Chain 2: Gradient evaluation took 0 seconds
#> Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
#> Chain 2: Adjust your expectations accordingly!
#> Chain 2: 
#> Chain 2: 
#> Chain 2: Iteration:    1 / 2000 [  0%]  (Warmup)
#> Chain 2: Iteration:  200 / 2000 [ 10%]  (Warmup)
#> Chain 2: Iteration:  400 / 2000 [ 20%]  (Warmup)
#> Chain 2: Iteration:  600 / 2000 [ 30%]  (Warmup)
#> Chain 2: Iteration:  800 / 2000 [ 40%]  (Warmup)
#> Chain 2: Iteration: 1000 / 2000 [ 50%]  (Warmup)
#> Chain 2: Iteration: 1001 / 2000 [ 50%]  (Sampling)
#> Chain 2: Iteration: 1200 / 2000 [ 60%]  (Sampling)
#> Chain 2: Iteration: 1400 / 2000 [ 70%]  (Sampling)
#> Chain 2: Iteration: 1600 / 2000 [ 80%]  (Sampling)
#> Chain 2: Iteration: 1800 / 2000 [ 90%]  (Sampling)
#> Chain 2: Iteration: 2000 / 2000 [100%]  (Sampling)
#> Chain 2: 
#> Chain 2:  Elapsed Time: 0.021 seconds (Warm-up)
#> Chain 2:                0.02 seconds (Sampling)
#> Chain 2:                0.041 seconds (Total)
#> Chain 2: 
#> 
#> SAMPLING FOR MODEL 'CrmEmpiricNormalPrior' NOW (CHAIN 3).
#> Chain 3: 
#> Chain 3: Gradient evaluation took 0 seconds
#> Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
#> Chain 3: Adjust your expectations accordingly!
#> Chain 3: 
#> Chain 3: 
#> Chain 3: Iteration:    1 / 2000 [  0%]  (Warmup)
#> Chain 3: Iteration:  200 / 2000 [ 10%]  (Warmup)
#> Chain 3: Iteration:  400 / 2000 [ 20%]  (Warmup)
#> Chain 3: Iteration:  600 / 2000 [ 30%]  (Warmup)
#> Chain 3: Iteration:  800 / 2000 [ 40%]  (Warmup)
#> Chain 3: Iteration: 1000 / 2000 [ 50%]  (Warmup)
#> Chain 3: Iteration: 1001 / 2000 [ 50%]  (Sampling)
#> Chain 3: Iteration: 1200 / 2000 [ 60%]  (Sampling)
#> Chain 3: Iteration: 1400 / 2000 [ 70%]  (Sampling)
#> Chain 3: Iteration: 1600 / 2000 [ 80%]  (Sampling)
#> Chain 3: Iteration: 1800 / 2000 [ 90%]  (Sampling)
#> Chain 3: Iteration: 2000 / 2000 [100%]  (Sampling)
#> Chain 3: 
#> Chain 3:  Elapsed Time: 0.025 seconds (Warm-up)
#> Chain 3:                0.02 seconds (Sampling)
#> Chain 3:                0.045 seconds (Total)
#> Chain 3: 
#> 
#> SAMPLING FOR MODEL 'CrmEmpiricNormalPrior' NOW (CHAIN 4).
#> Chain 4: 
#> Chain 4: Gradient evaluation took 0 seconds
#> Chain 4: 1000 transitions using 10 leapfrog steps per transition would take 0 seconds.
#> Chain 4: Adjust your expectations accordingly!
#> Chain 4: 
#> Chain 4: 
#> Chain 4: Iteration:    1 / 2000 [  0%]  (Warmup)
#> Chain 4: Iteration:  200 / 2000 [ 10%]  (Warmup)
#> Chain 4: Iteration:  400 / 2000 [ 20%]  (Warmup)
#> Chain 4: Iteration:  600 / 2000 [ 30%]  (Warmup)
#> Chain 4: Iteration:  800 / 2000 [ 40%]  (Warmup)
#> Chain 4: Iteration: 1000 / 2000 [ 50%]  (Warmup)
#> Chain 4: Iteration: 1001 / 2000 [ 50%]  (Sampling)
#> Chain 4: Iteration: 1200 / 2000 [ 60%]  (Sampling)
#> Chain 4: Iteration: 1400 / 2000 [ 70%]  (Sampling)
#> Chain 4: Iteration: 1600 / 2000 [ 80%]  (Sampling)
#> Chain 4: Iteration: 1800 / 2000 [ 90%]  (Sampling)
#> Chain 4: Iteration: 2000 / 2000 [100%]  (Sampling)
#> Chain 4: 
#> Chain 4:  Elapsed Time: 0.024 seconds (Warm-up)
#> Chain 4:                0.022 seconds (Sampling)
#> Chain 4:                0.046 seconds (Total)
#> Chain 4:

prob_tox_samp1 = as.data.frame(fit1, 'prob_tox') 
colMeans(prob_tox_samp1)
#> prob_tox[1] prob_tox[2] prob_tox[3] prob_tox[4] prob_tox[5] 
#>   0.1082452   0.2160202   0.3098107   0.4442895   0.6232141
fit1$prob_tox
#> [1] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141
fit1$median_prob_tox
#> [1] 0.07277617 0.19025108 0.29742550 0.44866472 0.63965960
unname(colMeans(prob_tox_samp1 > target))
#> [1] 0.11825 0.35625 0.60250 0.86175 0.99175

^{Created on 2021-12-13 by the reprex package (v2.0.1)}

However, the output is slightly different from the vignette, which reads:

#> prob_tox[1] prob_tox[2] prob_tox[3] prob_tox[4] prob_tox[5] 
#>   0.1081169   0.2159618   0.3098591   0.4444842   0.6235105
prob_too_toxic1 <- unname(colMeans(prob_tox_samp1 > target))
prob_too_toxic1
#> [1] 0.11700 0.35725 0.60075 0.86500 0.99150

Could there be a minor regression?

I'm using trialr v0.1.5 and rstan v2.21.2

The text was updated successfully, but these errors were encountered:

brockk · 2021-12-13T10:51:42Z

Hi Dan, MCMC fitting uses random numbers so without setting the RNG seed, you will get slightly different results. In fact, even with setting the seed, you could still get different results when comparing different versions of R or results from different OSs. As long as the results are not substantively different, you are good to go. Is that helpful? Thanks, Kristian

…

On Mon, 13 Dec 2021 at 10:31, Dan Chaltiel ***@***.***> wrote: Hi, Thanks for this great package. I'm following the vignette found in https://cran.r-project.org/web/packages/trialr/vignettes/A200-CRM.html. I ran its code with no problem (sorry, I couldn't find how to manage the verbosity): library(trialr) #> Warning: le package 'trialr' a été compilé avec la version R 4.1.2 #> Le chargement a nécessité le package : Rcpp target <- 0.25 skeleton <- c(0.05, 0.15, 0.25, 0.4, 0.6) fit1 <- stan_crm("2NN 3NN 4TT", skeleton=skeleton, target=target, model='empiric', beta_sd=sqrt(1.34), seed=123, verbose=FALSE) #> #> SAMPLING FOR MODEL 'CrmEmpiricNormalPrior' NOW (CHAIN 1). #> Chain 1: #> Chain 1: Gradient evaluation took 0 seconds #> Chain 1: 1000 transitions using 10 leapfrog steps per transition would take 0 seconds. #> Chain 1: Adjust your expectations accordingly! #> Chain 1: #> Chain 1: #> Chain 1: Iteration: 1 / 2000 [ 0%] (Warmup) #> Chain 1: Iteration: 200 / 2000 [ 10%] (Warmup) #> Chain 1: Iteration: 400 / 2000 [ 20%] (Warmup) #> Chain 1: Iteration: 600 / 2000 [ 30%] (Warmup) #> Chain 1: Iteration: 800 / 2000 [ 40%] (Warmup) #> Chain 1: Iteration: 1000 / 2000 [ 50%] (Warmup) #> Chain 1: Iteration: 1001 / 2000 [ 50%] (Sampling) #> Chain 1: Iteration: 1200 / 2000 [ 60%] (Sampling) #> Chain 1: Iteration: 1400 / 2000 [ 70%] (Sampling) #> Chain 1: Iteration: 1600 / 2000 [ 80%] (Sampling) #> Chain 1: Iteration: 1800 / 2000 [ 90%] (Sampling) #> Chain 1: Iteration: 2000 / 2000 [100%] (Sampling) #> Chain 1: #> Chain 1: Elapsed Time: 0.023 seconds (Warm-up) #> Chain 1: 0.021 seconds (Sampling) #> Chain 1: 0.044 seconds (Total) #> Chain 1: #> #> SAMPLING FOR MODEL 'CrmEmpiricNormalPrior' NOW (CHAIN 2). #> Chain 2: #> Chain 2: Gradient evaluation took 0 seconds #> Chain 2: 1000 transitions using 10 leapfrog steps per transition would take 0 seconds. #> Chain 2: Adjust your expectations accordingly! #> Chain 2: #> Chain 2: #> Chain 2: Iteration: 1 / 2000 [ 0%] (Warmup) #> Chain 2: Iteration: 200 / 2000 [ 10%] (Warmup) #> Chain 2: Iteration: 400 / 2000 [ 20%] (Warmup) #> Chain 2: Iteration: 600 / 2000 [ 30%] (Warmup) #> Chain 2: Iteration: 800 / 2000 [ 40%] (Warmup) #> Chain 2: Iteration: 1000 / 2000 [ 50%] (Warmup) #> Chain 2: Iteration: 1001 / 2000 [ 50%] (Sampling) #> Chain 2: Iteration: 1200 / 2000 [ 60%] (Sampling) #> Chain 2: Iteration: 1400 / 2000 [ 70%] (Sampling) #> Chain 2: Iteration: 1600 / 2000 [ 80%] (Sampling) #> Chain 2: Iteration: 1800 / 2000 [ 90%] (Sampling) #> Chain 2: Iteration: 2000 / 2000 [100%] (Sampling) #> Chain 2: #> Chain 2: Elapsed Time: 0.021 seconds (Warm-up) #> Chain 2: 0.02 seconds (Sampling) #> Chain 2: 0.041 seconds (Total) #> Chain 2: #> #> SAMPLING FOR MODEL 'CrmEmpiricNormalPrior' NOW (CHAIN 3). #> Chain 3: #> Chain 3: Gradient evaluation took 0 seconds #> Chain 3: 1000 transitions using 10 leapfrog steps per transition would take 0 seconds. #> Chain 3: Adjust your expectations accordingly! #> Chain 3: #> Chain 3: #> Chain 3: Iteration: 1 / 2000 [ 0%] (Warmup) #> Chain 3: Iteration: 200 / 2000 [ 10%] (Warmup) #> Chain 3: Iteration: 400 / 2000 [ 20%] (Warmup) #> Chain 3: Iteration: 600 / 2000 [ 30%] (Warmup) #> Chain 3: Iteration: 800 / 2000 [ 40%] (Warmup) #> Chain 3: Iteration: 1000 / 2000 [ 50%] (Warmup) #> Chain 3: Iteration: 1001 / 2000 [ 50%] (Sampling) #> Chain 3: Iteration: 1200 / 2000 [ 60%] (Sampling) #> Chain 3: Iteration: 1400 / 2000 [ 70%] (Sampling) #> Chain 3: Iteration: 1600 / 2000 [ 80%] (Sampling) #> Chain 3: Iteration: 1800 / 2000 [ 90%] (Sampling) #> Chain 3: Iteration: 2000 / 2000 [100%] (Sampling) #> Chain 3: #> Chain 3: Elapsed Time: 0.025 seconds (Warm-up) #> Chain 3: 0.02 seconds (Sampling) #> Chain 3: 0.045 seconds (Total) #> Chain 3: #> #> SAMPLING FOR MODEL 'CrmEmpiricNormalPrior' NOW (CHAIN 4). #> Chain 4: #> Chain 4: Gradient evaluation took 0 seconds #> Chain 4: 1000 transitions using 10 leapfrog steps per transition would take 0 seconds. #> Chain 4: Adjust your expectations accordingly! #> Chain 4: #> Chain 4: #> Chain 4: Iteration: 1 / 2000 [ 0%] (Warmup) #> Chain 4: Iteration: 200 / 2000 [ 10%] (Warmup) #> Chain 4: Iteration: 400 / 2000 [ 20%] (Warmup) #> Chain 4: Iteration: 600 / 2000 [ 30%] (Warmup) #> Chain 4: Iteration: 800 / 2000 [ 40%] (Warmup) #> Chain 4: Iteration: 1000 / 2000 [ 50%] (Warmup) #> Chain 4: Iteration: 1001 / 2000 [ 50%] (Sampling) #> Chain 4: Iteration: 1200 / 2000 [ 60%] (Sampling) #> Chain 4: Iteration: 1400 / 2000 [ 70%] (Sampling) #> Chain 4: Iteration: 1600 / 2000 [ 80%] (Sampling) #> Chain 4: Iteration: 1800 / 2000 [ 90%] (Sampling) #> Chain 4: Iteration: 2000 / 2000 [100%] (Sampling) #> Chain 4: #> Chain 4: Elapsed Time: 0.024 seconds (Warm-up) #> Chain 4: 0.022 seconds (Sampling) #> Chain 4: 0.046 seconds (Total) #> Chain 4: prob_tox_samp1 = as.data.frame(fit1, 'prob_tox') colMeans(prob_tox_samp1) #> prob_tox[1] prob_tox[2] prob_tox[3] prob_tox[4] prob_tox[5] #> 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141 fit1$prob_tox #> [1] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141 fit1$median_prob_tox #> [1] 0.07277617 0.19025108 0.29742550 0.44866472 0.63965960 unname(colMeans(prob_tox_samp1 > target)) #> [1] 0.11825 0.35625 0.60250 0.86175 0.99175 Created on 2021-12-13 by the reprex package <https://reprex.tidyverse.org> (v2.0.1) However, the output is slightly different from the vignette, which reads: #> prob_tox[1] prob_tox[2] prob_tox[3] prob_tox[4] prob_tox[5] #> 0.1081169 0.2159618 0.3098591 0.4444842 0.6235105 prob_too_toxic1 <- unname(colMeans(prob_tox_samp1 > target)) prob_too_toxic1 #> [1] 0.11700 0.35725 0.60075 0.86500 0.99150 Could there be a minor regression? I'm using trialr v0.1.5 and rstan v2.21.2 — You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub <#9>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ACNU3D5ZLVJIWPXBKZOKVPLUQXDPVANCNFSM5J5YKV2A> .

DanChaltiel · 2021-12-13T11:41:39Z

Hi Kristian,

Yes, that's very helpful, thanks. I see that the result is very close so I'm not really worried about it.

However (sorry for splitting hairs), the vignette states that:

In this demonstration, we set the seed for reproducibility when sampling [...].

Indeed, the output is very reproducible if the seed is set, which makes your answer is a bit confusing.

For instance, this code always yearns the same results, even after restarting R and even in reprex():

library(trialr)
rslt = replicate(10, {
  sink(tempfile())#I found how to make rstan quiet
  x=stan_crm("2NN 3NN 4TT", skeleton=c(0.05, 0.15, 0.25, 0.4, 0.6), target=0.25, seed=123, 
             model='empiric', beta_sd=sqrt(1.34))$prob_tox
  sink()
  x
})
t(rslt)
#>            [,1]      [,2]      [,3]      [,4]      [,5]
#>  [1,] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141
#>  [2,] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141
#>  [3,] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141
#>  [4,] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141
#>  [5,] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141
#>  [6,] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141
#>  [7,] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141
#>  [8,] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141
#>  [9,] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141
#> [10,] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141

^{Created on 2021-12-13 by the reprex package (v2.0.1)}

brockk · 2021-12-13T12:06:27Z

Yeah, fair enough. I will add more nuanced prose at some point. Thanks for raising it. Kristian

…

On Mon, 13 Dec 2021 at 11:41, Dan Chaltiel ***@***.***> wrote: Hi Kristian, Yes, that's very helpful, thanks. I see that the result is very close so I'm not really worried about it. However (sorry for splitting hairs), the vignette states that: In this demonstration, we set the seed for reproducibility when sampling [...]. Indeed, the output is very reproducible if the seed is set, which makes your answer is a bit confusing. For instance, this code always yearns the same results, even after restarting R and even in reprex(): library(trialr)rslt = replicate(10, { sink(tempfile())#I found how to make rstan quiet x=stan_crm("2NN 3NN 4TT", skeleton=c(0.05, 0.15, 0.25, 0.4, 0.6), target=0.25, seed=123, model='empiric', beta_sd=sqrt(1.34))$prob_tox sink() x }) t(rslt)#> [,1] [,2] [,3] [,4] [,5]#> [1,] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141#> [2,] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141#> [3,] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141#> [4,] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141#> [5,] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141#> [6,] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141#> [7,] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141#> [8,] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141#> [9,] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141#> [10,] 0.1082452 0.2160202 0.3098107 0.4442895 0.6232141 Created on 2021-12-13 by the reprex package <https://reprex.tidyverse.org> (v2.0.1) — You are receiving this because you commented. Reply to this email directly, view it on GitHub <#9 (comment)>, or unsubscribe <https://github.com/notifications/unsubscribe-auth/ACNU3D7T5XJ6UYPVC75LO2DUQXLX5ANCNFSM5J5YKV2A> .

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

Minor inconsistency in the vignette #9

Minor inconsistency in the vignette #9

Minor inconsistency in the vignette #9

Minor inconsistency in the vignette #9

Comments