Num_warmup / discard_initial / thinning error (?!)

[penelopeysm]

Same setup as #195

using AbstractMCMC

struct M <: AbstractMCMC.AbstractModel end
struct S <: AbstractMCMC.AbstractSampler end
# The samples are a tuple of "warmup"/"not warmup" and the number of times
# step_warmup / step have been called
function AbstractMCMC.step_warmup(rng, ::M, ::S; kwargs...)
    sample = ("Warmup", 1)
    @info sample
    return sample, 1
end
function AbstractMCMC.step_warmup(rng, ::M, ::S, i::Int; kwargs...)
    sample = ("Warmup", i + 1)
    @info sample
    return sample, i + 1
end
function AbstractMCMC.step(rng, ::M, ::S; kwargs...)
    sample = ("Not warmup", 1)
    @info sample
    return sample, 1
end
function AbstractMCMC.step(rng, ::M, ::S, i::Int; kwargs...)
    sample = ("Not warmup", i + 1)
    @info sample
    return sample, i + 1
end

sample(M(), S(), 5; num_warmup=10, discard_initial=5, thinning=5)

Here, I would expect:

• 5 warmup steps to be discarded
• we then run 5 more warmup steps (for a total of 10), then 16 normal steps, for a total of 21 actual steps
• we thin those 21 steps down to 5 by taking samples[1:5:21]

Out of the final thinned 5 steps, only the first should be a warmup step.

But, voila, they are all warmup steps:

julia> sample(M(), S(), 5; num_warmup=10, discard_initial=5, thinning=5)
[ Info: ("Warmup", 1)
[ Info: ("Warmup", 2)
[ Info: ("Warmup", 3)
[ Info: ("Warmup", 4)
[ Info: ("Warmup", 5)
[ Info: ("Warmup", 6)
[ Info: ("Warmup", 7)
[ Info: ("Warmup", 8)
[ Info: ("Warmup", 9)
[ Info: ("Warmup", 10)
[ Info: ("Warmup", 11)
[ Info: ("Warmup", 12)
[ Info: ("Warmup", 13)
[ Info: ("Warmup", 14)
[ Info: ("Warmup", 15)
[ Info: ("Warmup", 16)
[ Info: ("Warmup", 17)
[ Info: ("Warmup", 18)
[ Info: ("Warmup", 19)
[ Info: ("Warmup", 20)
[ Info: ("Warmup", 21)
[ Info: ("Warmup", 22)
[ Info: ("Warmup", 23)
[ Info: ("Warmup", 24)
[ Info: ("Warmup", 25)
[ Info: ("Warmup", 26)
Sampling 100%|████████████████████████████████████████████| Time: 0:00:00
5-element Vector{Tuple{String, Int64}}:
 ("Warmup", 6)
 ("Warmup", 11)
 ("Warmup", 16)
 ("Warmup", 21)
 ("Warmup", 26)

[nsiccha]

That doesn't look good IMO.

[nsiccha]

I feel like "all" of the logic here is questionable and unnecessarily convoluted.

[nsiccha]

This (IMO weird) behaviour is due to warm-up samples counting as a single step to be discarded, but (like post-warm-up samples) as thinning steps to be made after the discarded steps, compare e.g.

AbstractMCMC.jl/src/sample.jl

Line 147 in a346c51

Ntotal = thinning * (N - 1) + discard_initial + 1

In the first loop the first num_warmup+1 steps (see #195 (comment)) are done as a warm-up step, while in the second loop the next (keep_from_warmup-1)*thinning == (num_warmup-min(num_warmup, discard_initial)-1)*thinning steps (compare

AbstractMCMC.jl/src/sample.jl

Lines 161 to 162 in a346c51

	discard_from_warmup = min(num_warmup, discard_initial)
	keep_from_warmup = num_warmup - discard_from_warmup

) are done as warm-up steps.

[penelopeysm]

With #194:

julia> sample(M(), S(), 5; num_warmup=10, discard_initial=5, thinning=5)
[ Info: ("Warmup", 1)
[ Info: ("Warmup", 2)
[ Info: ("Warmup", 3)
[ Info: ("Warmup", 4)
[ Info: ("Warmup", 5)
[ Info: ("Warmup", 6)
[ Info: ("Warmup", 7)
[ Info: ("Warmup", 8)
[ Info: ("Warmup", 9)
[ Info: ("Warmup", 10)
[ Info: ("Warmup", 11)
[ Info: ("Not warmup", 12)
[ Info: ("Not warmup", 13)
[ Info: ("Not warmup", 14)
[ Info: ("Not warmup", 15)
[ Info: ("Not warmup", 16)
[ Info: ("Not warmup", 17)
[ Info: ("Not warmup", 18)
[ Info: ("Not warmup", 19)
[ Info: ("Not warmup", 20)
[ Info: ("Not warmup", 21)
[ Info: ("Not warmup", 22)
[ Info: ("Not warmup", 23)
[ Info: ("Not warmup", 24)
[ Info: ("Not warmup", 25)
[ Info: ("Not warmup", 26)
Sampling 100%|█████████████████████████████████████████████████████████| Time: 0:00:01
5-element Vector{Tuple{String, Int64}}:
 ("Warmup", 6)
 ("Warmup", 11)
 ("Not warmup", 16)
 ("Not warmup", 21)
 ("Not warmup", 26)