select_nested_models#
- gammapy.modeling.select_nested_models(datasets, parameters, null_values, n_sigma=2, n_free_parameters=None, fit=None)[source]#
Compute the test statistic (TS) between two nested hypotheses.
This function evaluates whether adding one or more parameters (alternative hypothesis) provides a statistically significant improvement (larger than the given threshold) over a simpler model (null hypothesis), where those parameters are fixed to given values.
The model is updated to the alternative hypothesis if there is a significant improvement (larger than the given threshold,
n_sigma
).- Parameters:
- datasets
Datasets
Datasets.
- parameters
Parameters
or list ofParameter
List of parameters frozen for the null hypothesis but free for the test hypothesis.
- null_valueslist of float or
Parameters
Values of the parameters frozen under the null hypothesis. If a
Parameters
object or a list ofParameters
is provided, the null hypothesis assumes these values. This allows testing of linked versus unlinked parameters.- n_sigmafloat, optional
Threshold in number of sigma for rejecting the null hypothesis in favour of the alternative. Default is 2. The TS is converted to sigma assuming that the Wilks’ theorem is verified.
- n_free_parametersint, optional
Number of free parameters to consider between the two hypotheses when estimating the
ts_threshold
from then_sigma
threshold. Default islen(parameters)
.- fit
Fit
, optional Fit instance specifying the backend and fit options. Default is None, which utilises the “minuit” backend with tol=0.1 and strategy=1.
- datasets
- Returns:
- resultdict
Dictionary with the TS of the best fit value compared to the null hypothesis and fit results for the two hypotheses. Entries are:
“ts” : fit statistic difference with null hypothesis
“fit_results” : results for the best fit
“fit_results_null” : fit results for the null hypothesis
Notes
This is useful for testing the significance of adding model components, such as a spectral cutoff or a source detection. It assumes the two models are nested — the null hypothesis is a special case of the alternative with fixed parameters.
See below for two simple examples and Compute upper limits for a source for a more in-depth example.
Examples
from gammapy.modeling.selection import select_nested_models from gammapy.datasets import Datasets, SpectrumDatasetOnOff from gammapy.modeling.models import SkyModel # Test if cutoff is significant dataset = SpectrumDatasetOnOff.read("$GAMMAPY_DATA/joint-crab/spectra/hess/pha_obs23523.fits") datasets = Datasets(dataset) model = SkyModel.create(spectral_model="ecpl", spatial_model="point", name='hess') datasets.models = model # If there is no cutoff (null hypothesis) then the spectral cutoff (`lambda`) is zero result = select_nested_models(datasets, parameters=[model.spectral_model.lambda_], null_values=[0], ) # Test if source is significant filename = "$GAMMAPY_DATA/fermi-3fhl-crab/Fermi-LAT-3FHL_datasets.yaml" filename_models = "$GAMMAPY_DATA/fermi-3fhl-crab/Fermi-LAT-3FHL_models.yaml" fermi_datasets = Datasets.read(filename=filename, filename_models=filename_models) model = fermi_datasets.models["Crab Nebula"] # Number of parameters previously fit for the source of interest n_free_parameters = len(model.parameters.free_parameters) # Freeze spatial parameters to ensure another weaker source does not move from its position # to replace the source of interest during the null hypothesis test. # (with all parameters free you test N vs. N+1 models and not the detection of a specific source.) fermi_datasets.models.freeze(model_type='spatial') # Null hypothesis value for the amplitude should be 0, i.e. the flux of the source is 0 results = select_nested_models(fermi_datasets, parameters=[model.spectral_model.amplitude], null_values=[0], n_free_parameters=n_free_parameters, n_sigma=4, )