Bayesian Sampler Examples

Examples of running each sampler avaiable in 3ML.

Before, that, let’s discuss setting up configuration default sampler with default parameters. We can set in our configuration a default algorithm and default setup parameters for the samplers. This can ease fitting when we are doing exploratory data analysis.

With any of the samplers, you can pass keywords to access their setups. Read each pacakges documentation for more details.

[1]:

from threeML import *
from threeML.plugins.XYLike import XYLike

from packaging.version import Version
import numpy as np
import dynesty
from jupyterthemes import jtplot

%matplotlib inline
jtplot.style(context="talk", fscale=1, ticks=True, grid=False)
silence_warnings()
set_threeML_style()

[2]:

threeML_config.bayesian.default_sampler

[2]:

<Sampler.emcee: 'emcee'>

[3]:

threeML_config.bayesian.emcee_setup

[3]:

{'n_burnin': None, 'n_iterations': 500, 'n_walkers': 50, 'seed': 5123}

If you simply run bayes_analysis.sample() the default sampler and its default parameters will be used.

Let’s make some data to fit.

[4]:

sin = Sin(K=1, f=0.1)
sin.phi.fix = True
sin.K.prior = Log_uniform_prior(lower_bound=0.5, upper_bound=1.5)
sin.f.prior = Uniform_prior(lower_bound=0, upper_bound=0.5)

model = Model(PointSource("demo", 0, 0, spectral_shape=sin))

x = np.linspace(-2 * np.pi, 4 * np.pi, 20)
yerr = np.random.uniform(0.01, 0.2, 20)


xyl = XYLike.from_function("demo", sin, x, yerr)
xyl.plot()

bayes_analysis = BayesianAnalysis(model, DataList(xyl))

../_images/notebooks_sampler_docs_5_0.png

emcee

[5]:

bayes_analysis.set_sampler("emcee")
bayes_analysis.sampler.setup(n_walkers=20, n_iterations=500)
bayes_analysis.sample()

xyl.plot()
bayes_analysis.results.corner_plot()

Maximum a posteriori probability (MAP) point:

	result	unit
parameter
demo.spectrum.main.Sin.K	(9.95 -0.17 +0.19) x 10^-1	1 / (keV s cm2)
demo.spectrum.main.Sin.f	(1.006 -0.007 +0.005) x 10^-1	rad / keV

Values of -log(posterior) at the minimum:

	-log(posterior)
demo	-5.745958
total	-5.745958

Values of statistical measures:

	statistical measures
AIC	16.197799
BIC	17.483382
DIC	15.531233
PDIC	2.015471

[5]:

../_images/notebooks_sampler_docs_7_8.png

../_images/notebooks_sampler_docs_7_9.png

../_images/notebooks_sampler_docs_7_10.png

multinest

[6]:

bayes_analysis.set_sampler("multinest")
bayes_analysis.sampler.setup(n_live_points=400, resume=False, auto_clean=True)
bayes_analysis.sample()

xyl.plot()
bayes_analysis.results.corner_plot()

 *****************************************************
 MultiNest v3.10
 Copyright Farhan Feroz & Mike Hobson
 Release Jul 2015

 no. of live points =  400
 dimensionality =    2
 *****************************************************
  analysing data from chains/fit-.txt ln(ev)=  -14.485838790708810      +/-  0.13918331226363556
 Total Likelihood Evaluations:         5660
 Sampling finished. Exiting MultiNest

Maximum a posteriori probability (MAP) point:

	result	unit
parameter
demo.spectrum.main.Sin.K	(9.95 -0.17 +0.16) x 10^-1	1 / (keV s cm2)
demo.spectrum.main.Sin.f	(1.006 +/- 0.006) x 10^-1	rad / keV

Values of -log(posterior) at the minimum:

	-log(posterior)
demo	-5.751789
total	-5.751789

Values of statistical measures:

	statistical measures
AIC	16.209460
BIC	17.495042
DIC	15.557033
PDIC	2.032081
log(Z)	-6.291120

WARNING:root:Too few points to create valid contours

[6]:

../_images/notebooks_sampler_docs_9_8.png

../_images/notebooks_sampler_docs_9_9.png

../_images/notebooks_sampler_docs_9_10.png

dynesty

[7]:

bayes_analysis.set_sampler("dynesty_nested")
bayes_analysis.sampler.setup(nlive=400)
bayes_analysis.sample()

xyl.plot()
bayes_analysis.results.corner_plot()

4055it [00:04, 977.82it/s, +400 | bound: 11 | nc: 1 | ncall: 19770 | eff(%): 22.999 | loglstar:   -inf < -5.751 <    inf | logz: -14.839 +/-  0.142 | dlogz:  0.001 >  0.409]

Maximum a posteriori probability (MAP) point:

	result	unit
parameter
demo.spectrum.main.Sin.K	(9.95 -0.17 +0.18) x 10^-1	1 / (keV s cm2)
demo.spectrum.main.Sin.f	(1.006 +/- 0.006) x 10^-1	rad / keV

Values of -log(posterior) at the minimum:

	-log(posterior)
demo	-5.746066
total	-5.746066

Values of statistical measures:

	statistical measures
AIC	16.198013
BIC	17.483596
DIC	15.412783
PDIC	1.958908
log(Z)	-6.444439

[7]:

../_images/notebooks_sampler_docs_11_7.png

../_images/notebooks_sampler_docs_11_8.png

../_images/notebooks_sampler_docs_11_9.png

[8]:

bayes_analysis.set_sampler("dynesty_dynamic")
bayes_analysis.sampler.setup()

if Version(dynesty.__version__) >= Version("3.0.0"):
    bayes_analysis.sample(n_effective=None)
else:
    bayes_analysis.sample(
        stop_function=dynesty.utils.old_stopping_function, n_effective=None
    )

xyl.plot()
bayes_analysis.results.corner_plot()

16305it [00:16, 1012.00it/s, batch: 8 | bound: 5 | nc: 1 | ncall: 38255 | eff(%): 42.526 | loglstar: -10.724 < -5.751 < -6.077 | logz: -14.684 +/-  0.072 | stop:  0.884]

Maximum a posteriori probability (MAP) point:

	result	unit
parameter
demo.spectrum.main.Sin.K	(9.95 -0.18 +0.17) x 10^-1	1 / (keV s cm2)
demo.spectrum.main.Sin.f	(1.006 +/- 0.006) x 10^-1	rad / keV

Values of -log(posterior) at the minimum:

	-log(posterior)
demo	-5.746067
total	-5.746067

Values of statistical measures:

	statistical measures
AIC	16.198017
BIC	17.483599
DIC	15.549755
PDIC	2.028647
log(Z)	-6.374337

[8]:

../_images/notebooks_sampler_docs_12_7.png

../_images/notebooks_sampler_docs_12_8.png

../_images/notebooks_sampler_docs_12_9.png

zeus

[9]:

bayes_analysis.set_sampler("zeus")
bayes_analysis.sampler.setup(n_walkers=20, n_iterations=500)
bayes_analysis.sample()

xyl.plot()
bayes_analysis.results.corner_plot()

The run method has been deprecated and it will be removed. Please use the new run_mcmc method.
Initialising ensemble of 20 walkers...
Sampling progress : 100%|██████████| 625/625 [00:04<00:00, 149.15it/s]
fit restored to maximum of posterior
fit restored to maximum of posterior
Summary
-------
Number of Generations: 625
Number of Parameters: 2
Number of Walkers: 20
Number of Tuning Generations: 28
Scale Factor: 1.410082
Mean Integrated Autocorrelation Time: 2.75
Effective Sample Size: 4540.19
Number of Log Probability Evaluations: 64636
Effective Samples per Log Probability Evaluation: 0.070242

None

Maximum a posteriori probability (MAP) point:

	result	unit
parameter
demo.spectrum.main.Sin.K	(9.95 -0.17 +0.18) x 10^-1	1 / (keV s cm2)
demo.spectrum.main.Sin.f	(1.006 +/- 0.006) x 10^-1	rad / keV

Values of -log(posterior) at the minimum:

	-log(posterior)
demo	-5.746018
total	-5.746018

Values of statistical measures:

	statistical measures
AIC	16.197919
BIC	17.483501
DIC	15.356045
PDIC	1.931703

[9]:

../_images/notebooks_sampler_docs_14_8.png

../_images/notebooks_sampler_docs_14_9.png

../_images/notebooks_sampler_docs_14_10.png

ultranest

[10]:

bayes_analysis.set_sampler("ultranest")
bayes_analysis.sampler.setup(
    min_num_live_points=400, frac_remain=0.5, use_mlfriends=False
)
bayes_analysis.sample()

xyl.plot()
bayes_analysis.results.corner_plot()

sampler set to [blue]ultranest[/blue]

[ultranest] Sampling 400 live points from prior ...

[ultranest] Explored until L=-6
[ultranest] Likelihood function evaluations: 9583
[ultranest]   logZ = -14.76 +- 0.1131
[ultranest] Effective samples strategy satisfied (ESS = 970.0, need >400)
[ultranest] Posterior uncertainty strategy is satisfied (KL: 0.46+-0.07 nat, need <0.50 nat)
[ultranest] Evidency uncertainty strategy is satisfied (dlogz=0.42, need <0.5)
[ultranest]   logZ error budget: single: 0.14 bs:0.11 tail:0.41 total:0.42 required:<0.50
[ultranest] done iterating.

fit restored to maximum of posterior
fit restored to maximum of posterior

Maximum a posteriori probability (MAP) point:

	result	unit
parameter
demo.spectrum.main.Sin.K	(9.96 -0.16 +0.18) x 10^-1	1 / (keV s cm2)
demo.spectrum.main.Sin.f	(1.006 +/- 0.006) x 10^-1	rad / keV

Values of -log(posterior) at the minimum:

	-log(posterior)
demo	-5.746389
total	-5.746389

Values of statistical measures:

	statistical measures
AIC	16.198660
BIC	17.484242
DIC	15.472174
PDIC	1.987420
log(Z)	-6.399863

[10]:

../_images/notebooks_sampler_docs_16_11.png

../_images/notebooks_sampler_docs_16_12.png

../_images/notebooks_sampler_docs_16_13.png