Step size continues to change during the sampling period when stepsize!=1

When a sampler is run with an option of stepsize !=1, the stepsize seems to continue to change during the sampling period. I am not sure why this is the case or whether this is an intended behavior? This also seems to happen when engine=static as well, and/or when engagement option is turned off with adapt engaged=0 and num_warmup=0

Example 1: run bernoulli example in cmdstan with stepsize=1. It seems after the warmup period the (adapted) constant step size = 0.620072 is used.

./bernoulli sample num_samples=10 num_warmup=10 save_warmup=1 \
  adapt engaged=1 \
  algorithm=hmc engine=nuts stepsize=1 \
  data file=bernoulli.data.json \
  random seed=1357911 \
  output file=nuts_test1.csv refresh=1000  

Output:

# stan_version_major = 2# stan_version_minor = 33
# stan_version_patch = 0
# model = bernoulli_model
# start_datetime = 2023-09-12 04:29:51 UTC
# method = sample (Default)
#   sample
#     num_samples = 10
#     num_warmup = 10
#     save_warmup = 1
#     thin = 1 (Default)
#     adapt
#       engaged = 1 (Default)
#       gamma = 0.050000000000000003 (Default)
#       delta = 0.80000000000000004 (Default)
#       kappa = 0.75 (Default)
#       t0 = 10 (Default)
#       init_buffer = 75 (Default)
#       term_buffer = 50 (Default)
#       window = 25 (Default)
#     algorithm = hmc (Default)
#       hmc
#         engine = nuts (Default)
#           nuts
#             max_depth = 10 (Default)
#         metric = diag_e (Default)
#         metric_file =  (Default)
#         stepsize = 1 (Default)
#         stepsize_jitter = 0 (Default)
#     num_chains = 1 (Default)
# id = 1 (Default)
# data
#   file = bernoulli.data.json
# init = 2 (Default)
# random
#   seed = 1357911
# output
#   file = nuts_test1.csv
#   diagnostic_file =  (Default)
#   refresh = 1000
#   sig_figs = -1 (Default)
#   profile_file = profile.csv (Default)
# num_threads = 10 (Default)
# stanc_version = stanc3 v2.33.0-rc1-21-g934496d9
# stancflags = --name=bernoulli_model
lp__,accept_stat__,stepsize__,treedepth__,n_leapfrog__,divergent__,energy__,theta
-11.6842,1.28295e-07,2,1,3,0,11.8514,0.691228
-11.6842,1.6848e-11,2.33507,1,1,0,12.0242,0.691228
-11.3994,1,0.230236,1,1,0,11.7439,0.679474
-11.3438,1,0.239791,1,1,0,11.5614,0.677113
-10.9485,1,0.324333,3,15,0,11.3599,0.028347
-8.86777,1,0.507206,2,3,0,10.9974,0.0633098
-8.49079,1,0.863377,1,1,0,9.10726,0.0743976
-6.87192,1,1.54496,1,1,0,7.64316,0.191467
-6.87192,5.37683e-34,2.84484,1,1,0,7.05756,0.191467
-6.75068,0.990502,0.226075,2,7,0,7.6069,0.240958
# Adaptation terminated
# Step size = 0.620072
# Diagonal elements of inverse mass matrix:
# 1
-7.6013,0.807688,0.620072,1,3,0,8.05102,0.114168
-7.96722,0.95234,0.620072,2,3,0,8.54374,0.094668
-6.79277,0.973714,0.620072,1,3,0,8.19835,0.213904
-6.75457,0.999022,0.620072,2,3,0,6.80227,0.235883
-6.74846,0.824895,0.620072,2,3,0,8.27092,0.246328
-6.74915,0.999851,0.620072,1,1,0,6.74915,0.244089
-6.84431,0.987084,0.620072,2,3,0,6.85664,0.197994
-7.22132,0.94955,0.620072,1,3,0,7.24814,0.143217
-9.08743,0.804713,0.620072,1,3,0,9.49342,0.0578167
-8.14034,1,0.620072,1,1,0,9.01039,0.0871781
# 
#  Elapsed Time: 0 seconds (Warm-up)
#                0 seconds (Sampling)
#                0 seconds (Total)

Example 2: run same code except setting stepsize=0.5. During the sampling period, sample size varies.

./bernoulli sample num_samples=10 num_warmup=10 save_warmup=1 \
>   adapt engaged=1 \
>   algorithm=hmc engine=nuts stepsize=0.5 \
>   data file=bernoulli.data.json \
>   random seed=1357911 \
>   output file=nuts_test1.csv refresh=1000  

Output:

# stan_version_major = 2
# stan_version_minor = 33
# stan_version_patch = 0
# model = bernoulli_model
# start_datetime = 2023-09-12 04:39:29 UTC
# method = sample (Default)
#   sample
#     num_samples = 10
#     num_warmup = 10
#     save_warmup = 1
#     thin = 1 (Default)
#     adapt
#       engaged = 1 (Default)
#       gamma = 0.050000000000000003 (Default)
#       delta = 0.80000000000000004 (Default)
#       kappa = 0.75 (Default)
#       t0 = 10 (Default)
#       init_buffer = 75 (Default)
#       term_buffer = 50 (Default)
#       window = 25 (Default)
#     algorithm = hmc (Default)
#       hmc
#         engine = nuts (Default)
#           nuts
#             max_depth = 10 (Default)
#         metric = diag_e (Default)
#         metric_file =  (Default)
#         stepsize = 0.5
#         stepsize_jitter = 0 (Default)
#     num_chains = 1 (Default)
# id = 1 (Default)
# data
#   file = bernoulli.data.json
# init = 2 (Default)
# random
#   seed = 1357911
# output
#   file = nuts_test1.csv
#   diagnostic_file =  (Default)
#   refresh = 1000
#   sig_figs = -1 (Default)
#   profile_file = profile.csv (Default)
# num_threads = 10 (Default)
# stanc_version = stanc3 v2.33.0-rc1-21-g934496d9
# stancflags = --name=bernoulli_model
lp__,accept_stat__,stepsize__,treedepth__,n_leapfrog__,divergent__,energy__,theta
-11.6842,2.70997e-11,2.22077,1,1,0,11.7777,0.691228
-11.6842,0.669326,1.334,1,3,0,11.8694,0.691228
-8.41406,1,0.55738,1,1,0,10.9489,0.507899
-7.33319,0.965726,0.47472,2,3,0,8.88201,0.398243
-6.83317,1,1.03159,1,1,0,7.10452,0.303782
-6.83317,0.000129385,1.68004,2,3,0,9.802,0.303782
-7.04704,0.998552,0.170997,4,15,0,7.082,0.353827
-6.961,0.991063,0.323372,2,3,0,7.25201,0.33682
-6.87354,0.973395,0.247474,4,15,0,8.30392,0.315817
-6.78384,0.99654,0.807968,2,3,0,6.91221,0.21758
# Adaptation terminated
# Step size = 0.622228
# Diagonal elements of inverse mass matrix:
# 1
-7.1544,0.749035,0.718792,2,3,0,8.52368,0.149961
-6.98185,1,0.546031,1,1,0,7.13381,0.171638
-7.03157,0.988452,0.745944,1,1,0,7.08104,0.164566
-7.87966,0.895127,0.635039,1,3,0,7.9396,0.098828
-10.1041,0.808039,0.749397,1,3,0,10.1779,0.0388005
-10.3424,0.999612,0.324777,3,7,0,10.351,0.0354672
-10.0319,0.961604,0.552549,2,3,0,12.7201,0.0398808
-7.5994,0.913906,0.636393,1,3,0,10.2428,0.114286
-6.82905,1,0.60127,2,3,0,7.55226,0.202064
-6.82905,0.821467,0.834338,1,3,0,7.50771,0.202064
# 
#  Elapsed Time: 0 seconds (Warm-up)
#                0 seconds (Sampling)
#                0 seconds (Total)

Example3: even with num_warmup=0, adapt engaged=0, the step size continues to vary.

./bernoulli sample num_samples=10 num_warmup=0 save_warmup=1 \
>   adapt engaged=0 \
>   algorithm=hmc engine=nuts stepsize=0.5 \
>   data file=bernoulli.data.json \
>   random seed=1357911 \
>   output file=nuts_test1.csv refresh=1000  

output:

# stan_version_major = 2
# stan_version_minor = 33
# stan_version_patch = 0
# model = bernoulli_model
# start_datetime = 2023-09-12 04:41:34 UTC
# method = sample (Default)
#   sample
#     num_samples = 10
#     num_warmup = 0
#     save_warmup = 1
#     thin = 1 (Default)
#     adapt
#       engaged = 0
#       gamma = 0.050000000000000003 (Default)
#       delta = 0.80000000000000004 (Default)
#       kappa = 0.75 (Default)
#       t0 = 10 (Default)
#       init_buffer = 75 (Default)
#       term_buffer = 50 (Default)
#       window = 25 (Default)
#     algorithm = hmc (Default)
#       hmc
#         engine = nuts (Default)
#           nuts
#             max_depth = 10 (Default)
#         metric = diag_e (Default)
#         metric_file =  (Default)
#         stepsize = 0.5
#         stepsize_jitter = 0 (Default)
#     num_chains = 1 (Default)
# id = 1 (Default)
# data
#   file = bernoulli.data.json
# init = 2 (Default)
# random
#   seed = 1357911
# output
#   file = nuts_test1.csv
#   diagnostic_file =  (Default)
#   refresh = 1000
#   sig_figs = -1 (Default)
#   profile_file = profile.csv (Default)
# num_threads = 10 (Default)
# stanc_version = stanc3 v2.33.0-rc1-21-g934496d9
# stancflags = --name=bernoulli_model
lp__,accept_stat__,stepsize__,treedepth__,n_leapfrog__,divergent__,energy__,theta
# Adaptation terminated
# Step size = 0.5
# Diagonal elements of inverse mass matrix:
# 1
-7.18016,1,0.4926,2,3,0,11.6683,0.376189
-6.88342,1,0.634186,2,3,0,7.10301,0.318471
-6.74821,0.970537,0.743492,2,3,0,7.01738,0.2476
-6.8329,0.982611,0.457896,2,5,0,6.94226,0.303695
-6.75158,0.985678,0.383565,2,5,0,7.0913,0.239549
-7.42011,0.848019,0.712129,1,3,0,7.57903,0.126474
-7.05328,0.997039,0.518224,2,3,0,7.51217,0.161721
-7.1152,0.996135,0.395159,1,1,0,7.11882,0.154265
-7.39386,0.987633,0.605211,2,3,0,7.39784,0.12846
-6.77238,0.94836,0.584992,2,7,0,7.96497,0.278244
# 
#  Elapsed Time: 0 seconds (Warm-up)
#                0 seconds (Sampling)
#                0 seconds (Total)
2 Likes

Yep, this is a bug! The stepsize_jitter argument was accidentally being overwritten by the stepsize argument.

Thank you for reporting!

2 Likes

Are you using CmdStan from github or a released version?

Thank you. I installed from the github, running the following code:

git clone https://github.com/stan-dev/cmdstan.git --recursive

following the instruction in 1 CmdStan Installation | CmdStan User’s Guide

Should I install a different version?

Thank you,

The latest released version is available here: Release v2.33.0 (5 September 2023) · stan-dev/cmdstan · GitHub