Analyzing results of a scan¶
Tip
See the Example 9: Analyzing the Results of an NPTFit Run for an interactive version of the analysis options described here.
While the chain samples of a non-Poissonian fit performed using MultiNest can be readily
accessed, a basic analysis module dnds_analysis.py
is provided which includes the
basic functions described below.
Initializing the analysis module¶
Having performed a scan using an instance of nptfit.NPTF
, the first thing to do is load the scan parameters in. This is done with
>>> nptf.load_scan()
where nptf
is an instance of nptfit.NPTF
.
Note
The analysis can be performed on an already existing scan. To do this an instance of nptfit.NPTF
should be created with the same template and model configuration as used to perform the scan. Then instead of running the scan, simply load it and proceed with the analysis as below.
An instance of the analysis module can then be created as follows:
>>> an = dnds_analysis.Analysis(nptf, mask=None, pixarea=0.)
where mask
specifies the ROI used for the analysis if this is different from that used for the scan, and pixarea
is the area of a pixel in sr
if using non-HEALPix
maps.
Making triangle plots¶
Triangle/corner plots can be used to visualize multidimensional samples using a scatterplot matrix. A triangle plot with the default options using the corner package can be made using:
>>> an.make_triangle()
To use your own custom plotting options, use corner
as follows
>>> corner.corner(an.nptf.samples, labels=an.nptf.params, range=[1 for i in range(an.nptf.n_params)])
with additional arguments as specified in http://corner.readthedocs.io/en/latest/.
Getting template intensities¶
Template intensities [counts/cm2/s/sr] can be calculated with
>>> an.return_intensity_arrays_poiss(comp)
>>> an.return_intensity_arrays_non_poiss(comp)
for Poissonian and non-Poissonian templates with the key comp
respectively. This returns a list of intensities corresponding to the posterior parameters
for the given template.
The non-Poissonian templates (NPT) intensity is calculated by integrating up \(\int_{S_{min}}^{S_{max}} dS~S~dN/dS\). This is approximated as a sum between \(S_{min}\) and \(S_{max}\). The options associated with the non-Poissonian template intensity are:
Argument | Default | Purpose |
---|---|---|
comp |
- | The NPT key |
smin |
0.01 | Minimum counts to sum up from |
smax |
10000 | Maximum counts to sum up to |
nsteps |
10000 | Number of bins in s while summing up |
Getting source count distributions¶
The posterior arrays for the source count distributions \(dN/dF\) [counts-1 cm2 s deg-2] associated with a given template comp
at a given flux
(in counts/cm2/s) can be obtained using
>>> an.return_dndf_arrays(comp,flux)
The source count distribution can be plotted with
>>> an.plot_source_count_median(comp, smin, smax, nsteps, spow, **kwargs)
>>> an.plot_source_count_band(comp, smin, smax, nsteps, spow, qs, **kwargs)
The options being the same as for obtaining the NPT intensity above. Additionally, spow is the power \(n\) in \(F^ndN/dF\) to return while plotting, and qs is an array of quantiles for which to return the dN/dF band.
Plotting intensity fractions¶
Intensity fractions (fraction of template intensity to total intensity) for Poissonian and non-Poissonian templates respectively can be plotting using
>>> an.plot_intensity_fraction_poiss(comp, bins, **kwargs)
>>> an.plot_intensity_fraction_non_poiss(comp, bins, **kwargs)
where comp
is the template key, bins
is the number of bins between 0 and 100 and **kwargs
specify plotting options.
Accessing posteriors¶
While the posteriors can be accessed with nptf.samples
(or an.nptf.samples
) as above, the following functions provide a useful interfact to access individual parameters:
>>> an.return_poiss_parameter_posteriors(comp)
>>> an.return_non_poiss_parameter_posteriors(comp)
where comp
is the (non-)Poissonian template key.
For Poissonian models, this returns a list of posterior normalizaion parameters for that model. For non-Poissonian models, this returns three arrays:
>>> A_non_poiss_post, n_non_poiss_post, Sb_non_poiss_post = an.return_non_poiss_parameter_posteriors(comp)
where
A_non_poiss_post
is an array of non-Poissonian normalization parameter posteriorsn_non_poiss_post
is a 2-D array, each sub-array containing posteriors for a given slope parameter, starting from the highest to the lowestSb_non_poiss_post
is a 2-D array, each sub-array containing posteriors for a given break parameter, starting from the highest to the lowest
Getting Bayesian log-evidences¶
The Bayesian log-evidence and associated error can be accessed as follows:
>>> l_be, l_be_err = an.get_log_evidence()