Can rejected samples be plotted in parameter space?

[lm1909]

Is there a way to access the rejected particles and plot their positions in parameter space (similar to the upper right-hand corner in plot_kde_matrix_highlevel)? This would be very helpful to debug situations where the acceptance rate drops off quickly.

(According to the docstrings in pyabc/storage/history, rejected particles do not get stored in the database. However, can they somehow be retrieved directly from the SMCABC object?)

[yannikschaelte]

Hi @lm1909 , no, unfortunately rejected particles are currently not stored anywhere. If needed, this could be done as the sampler already has the ability to also return rejected particles (https://github.com/ICB-DCM/pyABC/blob/main/pyabc/population/population.py#L269), but the subsequent processing (akin to here https://github.com/ICB-DCM/pyABC/blob/main/pyabc/inference/smc.py#L566 all particles could be extracted and stored somewhere) and storage would need to be rewritten.

[lm1909]

Thanks a lot!

[lm1909]

If someone faces the same problem (displaying the rejected samples in parameter space to debug acceptance rate issues), I solved it (in a hacky way) with code that looks something like this:

def resample_with_rejected(abc, t, n, max_eval):
    resample = abc.sampler.sample_until_n_accepted(
        n=n,
        simulate_one=abc._create_simulate_function(t),
        t=t,
        max_eval=max_eval,
        ana_vars=abc._vars(t=t)
    )
    if not resample.ok:
        logger.critical('Sample not ok!')
    resample.normalize_weights()
    return pd.DataFrame([dict(particle.parameter) | {'accepted': particle.accepted, 'distance': particle.distance} for particle in resample.all_particles])

seaborn.pairplot(
    resample_with_rejected(abc_of_the_run_one_wants_to_debug, generation_that_needs_debugging, 100, 1000),
    kind='scatter',
    diag_kind=None,
    hue='accepted',
    palette='viridis',
    corner=True,
)

In the ABC run enable recording rejected particles, e.g. like:

sampler = pyabc.sampler.multicore_evaluation_parallel.MulticoreEvalParallelSampler(n_procs=5)
sampler.sample_factory = pyabc.population.SampleFactory(record_rejected=True)
   
abc = pyabc.ABCSMC(
   ...,
   sampler=sampler,
   ...)

Note that this

• does actually not display the particles that were rejected during the inference procedure - instead this is a resampling from the same transition distribution. However this should be enough to identify issues with the transition KDE for example.
• requires the argument abc to be the returned object of a e.g. abc.run() as compared to something loaded from a file (because rejected particles are not stored as of this version of pyABC as in the answer above)

(this code could be made much cleaner, but I just wanted to briefly copy it here in case it helps someone with the same issue :) )