Option to thin

[luiarthur]

Is there currently an option to thin samples by a certain factor? (I've looked in sample.jl and there doesn't seem to be such an option.)

I'm thinking of something like:

number_of_samples = 100
chain = sample(some_model, some_sampler, number_of_iterations; thin=3, discard_initial=50)

would run a Markov chain for 50 iterations (discard_initial), then another 300 iterations (number_of_samples * thin) but only every third (thin) sample is kept. Thus, 100 samples are returned.

This would be useful for certain samplers, like MH, where samples tend to be more correlated (than say those from NUTS) and the user wants to keep only every thin-th sample so as to obtain less correlated samples while reducing memory consumption.

I'd be happy to implement this if there isn't currently an option to do so.

[devmotion]

I guess it could be useful (I think we had this discussion at some point at least) but it hasn't been implemented yet mainly for two reasons: Turing does not use AbstractMCMC 2 but the old interface (this will change soon though) and we noticed that one has to be careful with changing things in AbstractMCMC since it affects many also Turing-independent downstream packages and hence it is not a good place to play around with ideas.

That being said, I think it could useful, not sure if it should be called thin or thinning or something else. Also at some point the future (not the current) code in Turing for resuming sampling should be moved to AbstractMCMC I assume.

By the way, if you are concerned about memory, you can use the iteration and transducer interfaces.

[luiarthur]

Yes, I recall having that conversation but couldn't find the thread...

one has to be careful with changing things in AbstractMCMC since it affects many also Turing-independent downstream packages and hence it is not a good place to play around with ideas.

Makes sense. Is this then perhaps something that's more suited to be implemented in Turing? Though, my naive view is that this is mainly an accounting problem that involves something like the following change in mcmcsample (and perhaps other relevant places):

    Ntotal = N * thin + discard_initial

    # stuff ...

        # Step through the sampler.
        for i in 2:(N*thin)
            # Obtain the next sample and state.
            sample, state = step(rng, model, sampler, state; kwargs...)

            # Run callback.
            callback === nothing || callback(rng, model, sampler, sample, i)

            # Save the sample.
            (i % thin == 0) && (samples = save!!(samples, sample, i, model, sampler, N; kwargs...))

            # Update the progress bar.
            progress && ProgressLogging.@logprogress (i + discard_initial) / Ntotal
        end

Besides the benefit of the added functionality, it's not apparent to me what might inadvertently be affected downstream. Am I missing something?

By the way, if you are concerned about memory, you can use the iteration and transducer interfaces.

Thanks for the pointer. Could you elaborate?

[devmotion]

Am I missing something?

Yeah, that's it basically - but there are some caveats and it requires some design choices. In your sketch the initial sample would always be saved (i = 1) but then the next one would be i = thin. Rather it should be i = thin + 1. So the check should be changed or one should not save the first sample. The latter is problematic, however, since then samples is not defined. Alternatively, everything could be handled in save!! (that's also how you could implement it in Turing) but that's probably less stable and would easily break if someone defines a custom save!! implementation.

Could you elaborate?

You can use steps (https://github.com/TuringLang/AbstractMCMC.jl/blob/faf1d73433734f53bb24452ae273401e0969f799/src/stepper.jl#L17-L27) or Sample (https://github.com/TuringLang/AbstractMCMC.jl/blob/faf1d73433734f53bb24452ae273401e0969f799/src/transducer.jl#L8-L10) if you want to sample while iterating.