Quickstart ========== A short walk-through of fitting a gamma-ray spectrum with *BaySpec*. It covers the three stages of a typical session: - **Data** — joint spectra from Fermi/GBM's NaI and BGO detectors. - **Model** — a cutoff power-law (``cpl``). - **Inference** — Bayesian posterior sampling via ``multinest`` (``emcee`` works as a drop-in replacement). .. code:: ipython3 import numpy as np from bayspec.model.local import * from bayspec import DataUnit, Data, BayesInfer, Plot, MaxLikeFit .. code:: ipython3 savepath = 'quickstart' 1. Load spectra and response data. .. code:: ipython3 nai = DataUnit( src='./ME/me.src', bkg='./ME/me.bkg', rsp='./ME/me.rsp', notc=[8, 900], stat='pgstat', rebn={'min_sigma': 2, 'max_bin': 10}) bgo = DataUnit( src='./HE/he.src', bkg='./HE/he.bkg', rsp='./HE/he.rsp', notc=[300, 38000], stat='pgstat', rebn={'min_sigma': 2, 'max_bin': 10}) data = Data([('nai', nai), ('bgo', bgo)]) data.save(savepath) data .. raw:: html
Data
Name Noticing Statistic Grouping Time
nai [[8, 900]] pgstat None None
bgo [[300, 38000]] pgstat None None
Data Parameters
par# Component Parameter Value Prior Frozen
1 nai sf 1 None True
2 nai bf 1 None True
3 nai rf 1 None True
4 nai ra 0 None True
5 nai dec 0 None True
6 bgo sf 1 None True
7 bgo bf 1 None True
8 bgo rf 1 None True
9 bgo ra 0 None True
10 bgo dec 0 None True
2. Define the spectral model. .. code:: ipython3 model = cpl() model.save(savepath) model .. raw:: html
Model

cpl [add]

power-law model with high-energy cutoff

Model Configurations
cfg# Component Parameter Value
1 cpl redshift 0.000
2 cpl pivot_energy 1.000
3 cpl vfv_peak True
Model Parameters
par# Component Parameter Value Prior
1 cpl $\alpha$ -1 unif(-2, 2)
2 cpl log$E_p$ 2 unif(0, 4)
3 cpl log$A$ 0 unif(-10, 10)
3. Run Bayesian inference. .. code:: ipython3 infer = BayesInfer([(data, model)]) infer.save(savepath) infer .. raw:: html
Bayesian Inference
Configurations
cfg# Class Expression Component Parameter Value
1 model cpl cpl redshift 0.000
2 model cpl cpl pivot_energy 1.000
3 model cpl cpl vfv_peak True
Parameters
par# Class Expression Component Parameter Value Prior
1* model cpl cpl $\alpha$ -1 unif(-2, 2)
2* model cpl cpl log$E_p$ 2 unif(0, 4)
3* model cpl cpl log$A$ 0 unif(-10, 10)
.. code:: ipython3 post = infer.multinest(nlive=400, resume=True, verbose=False, savepath=savepath) post.save(savepath) post .. raw:: html
Posterior Results
Parameters
par# Class Expression Component Parameter Mean Median Best 1sigma Best 1sigma CI
1 model cpl cpl $\alpha$ -1.562 -1.562 -1.561 -1.561 [-1.572, -1.551]
2 model cpl cpl log$E_p$ 2.331 2.331 2.330 2.330 [2.321, 2.341]
3 model cpl cpl log$A$ 2.353 2.354 2.352 2.352 [2.338, 2.368]
Statistics
Data Model Statistic Value Bins
nai cpl pgstat 405.171 119
bgo cpl pgstat 149.560 117
Total Total stat/dof 554.731/233 236
Information Criterias
AIC AICc BIC lnZ
560.731560.834571.123-295.844
.. code:: ipython3 fig = Plot.infer(post, style='CE') fig.save(f'{savepath}/ctsspec') .. raw:: html .. code:: ipython3 fig = Plot.post_corner(post, ploter='plotly') fig.save(f'{savepath}/corner') .. raw:: html .. code:: ipython3 earr = np.logspace(1, 3, 100) modelplot = Plot.model(ploter='plotly', style='vFv', post=True) modelplot.add_model(model, E=earr) fig = modelplot.get_fig() fig.save(f'{savepath}/model') .. raw:: html .. code:: ipython3 ergflux = model.best_ergflux(emin=10, emax=1000, ngrid=1000) ergflux_sample = model.ergflux_sample(emin=10, emax=1000, ngrid=1000) print(ergflux) print(ergflux_sample) .. parsed-literal:: 8.367943905909297e-06 {'mean': 8.369501130327865e-06, 'median': 8.369420051840691e-06, 'Isigma': array([8.31271099e-06, 8.42362544e-06]), 'IIsigma': array([8.25891689e-06, 8.48377819e-06]), 'IIIsigma': array([8.20528996e-06, 8.54385474e-06])}