RStan: may I submit some code?

[sp_r...@yahoo.it]

Prequel 1: sometimes it might be useful to try non-default sampler parameters (see: https://groups.google.com/forum/#!topic/stan-users/Fl3RBPscUXo)
Prequel 2: one could merge chains coming from different algorithms and/or iter/warmup/thin settings (as long as n_save and warmup2 are equal) and/or "controls", as they are named now (see: https://groups.google.com/forum/#!topic/stan-users/sCLQYj02tB8)

So I'm writing/rewriting some rstan functions. A few examples follow.

Let sm be a compiled stan model.

> iter <- c(10, 100, 500)
> warmup <- c(5, 50, 250)
> thin <- c(1, 10, 50)
> algorithm <- c("HMC", "NUTS", "NUTS")
> control <- list(list(adapt_delta = 0.61, int_time = 3 * pi,
+                      metric = "unit_e", stepsize = 0.6),
+                 list(adapt_delta = 0.62, max_treedepth = 12,
+                      metric = "unit_e", stepsize = 0.7),
+                 list(adapt_delta = 0.63, max_treedepth = 15,
+                      metric = "diag_e", stepsize = 0.5))
> sflist <- mclapply(1:3, mc.cores = detectCores(),
+                    function(i) sampling(sm, chains = 1, chain_id = i, seed = 12345,
+                        iter = iter[i], warmup = warmup[i], thin = thin[i],
+                        algorithm = algorithm[i], control = control[[i]]))
> sf <- sflist2stanfit(sflist)

I'd love to add a fourth case, but setting metric = dense_e doesn't work now (https://groups.google.com/forum/#!topic/stan-users/AJXmiruMiSA).

I'd like to look at something like:

> get_sampler_config(sf)
chain_id                1      2      3
algorithm             HMC   NUTS   NUTS
- int_time           9.42     NA     NA
- max_treedepth        NA     12     15
- metric               NA unit_e diag_e
- stepsize             NA   0.70   0.50
- stepsize_jitter      NA   0.00   0.00
iter                   10    100    500
warmup                  5     50    250
thin                    1     10     50
Adaptation:
- engaged            TRUE   TRUE   TRUE
- gamma              0.05   0.05   0.05
- delta              0.61   0.62   0.63
- kappa              0.75   0.75   0.75
- t0                10.00  10.00  10.00

I can now notice that metric and stepsize are not changed when algorithm is "HMC", but it'll be soon fixed.
Otherwise:

> sampler_config <- get_sampler_config(sf, print = FALSE)
> sampler_config
  chain_id sampler int_time max_treedepth metric stepsize stepsize_jitter iter
1        1     HMC 9.424778            NA   <NA>       NA              NA   10
2        2    NUTS       NA            12 unit_e      0.7               0  100
3        3    NUTS       NA            15 diag_e      0.5               0  500
  warmup thin engaged gamma delta kappa t0
1      5    1    TRUE  0.05  0.61  0.75 10
2     50   10    TRUE  0.05  0.62  0.75 10
3    250   50    TRUE  0.05  0.63  0.75 10

A data frame. As to get_sampler_params():

> get_sampler_params(sf)

Chain 1 (HMC):
              Min. 1st Qu. Median   Mean 3rd Qu.  Max.
accept_stat 0.0000  0.2173 0.7094 0.5975   1.000 1.000
stepsize    0.3383  1.5630 1.5630 1.7300   1.671 4.800
int_time    9.4250  9.4250 9.4250 9.4250   9.425 9.425

Chain 2 (NUTS):
                 Min. 1st Qu. Median   Mean 3rd Qu. Max.
accept_stat 2.232e-32  0.2365 0.4477 0.5134   0.927  1.0
stepsize    4.814e-01  1.5500 1.5500 1.6700   1.572  2.8
treedepth   0.000e+00  0.0000 0.0000 0.2000   0.000  1.0

Chain 3 (NUTS):
              Min. 1st Qu. Median   Mean 3rd Qu.  Max.
accept_stat  0.000 0.02386 0.2003 0.3332  0.5172  1.00
stepsize     1.896 2.02900 2.1150 3.7140  2.1150 11.78
treedepth   -1.000 0.00000 0.0000 0.0000  0.0000  1.00

Or:

> sampler_params <- get_sampler_params(sf, print = FALSE)
> sampler_params
   chain algorithm   accept_stat int_time treedepth
1      1       HMC  0.000000e+00 9.424778        NA
2      1       HMC  1.805636e-42 9.424778        NA
...
29     3      NUTS  5.472281e-01       NA         0
30     3      NUTS  2.905030e-02       NA         0

Another data frame. Furthermore,

> print(sf)
Inference for Stan model: normal.
3 chains, each with iter=10; warmup=5; thin=1;
...

is wrong, because n_save and warmup2 are equal, but iter/warmup/thin are not.

If someone thinks that such pedantic "improvements" could be useful, I'll polish up and submit my code.

Sergio

[Bob Carpenter]

If the RStan developers are interested, the thing to do with it is to
make a pull request through GitHub.

Then we (and by that I mean Jiqiang or Ben or someone else who
knows R better than me) can review the code and decide whether to
merge or to require fixes or other further work.

We want tests and documentation along with pull requests whenever
possible. Nothing should get merged into the main development branch
until that's done.

- Bob

> --
> You received this message because you are subscribed to the Google Groups "stan users mailing list" group.
> To unsubscribe from this group and stop receiving emails from it, send an email to stan-users+...@googlegroups.com.
> For more options, visit https://groups.google.com/groups/opt_out.