Fit¶
-
class
gammapy.modeling.
Fit
(datasets, store_trace=False)[source]¶ Bases:
object
Fit class.
The fit class provides a uniform interface to multiple fitting backends. Currently available: “minuit”, “sherpa” and “scipy”
- Parameters
- datasets
Datasets
Datasets
- datasets
Methods Summary
confidence
(parameter[, backend, sigma, …])Estimate confidence interval.
covariance
([backend])Estimate the covariance matrix.
minos_contour
(x, y[, numpoints, sigma])Compute MINOS contour.
optimize
([backend])Run the optimization.
run
([backend, optimize_opts, covariance_opts])Run all fitting steps.
stat_profile
(parameter[, values, bounds, …])Compute fit statistic profile.
stat_surface
(x, y, x_values, y_values[, …])Compute fit statistic surface.
Methods Documentation
-
confidence
(parameter, backend='minuit', sigma=1, reoptimize=True, **kwargs)[source]¶ Estimate confidence interval.
Extra
kwargs
are passed to the backend. E.g.iminuit.Minuit.minos
supports amaxcall
option.For the scipy backend
kwargs
are forwarded tobrentq
. If the confidence estimation fails, the bracketing interval can be adapted by modifying the the upper bound of the interval (b
) value.- Parameters
- backendstr
Which backend to use (see
gammapy.modeling.registry
)- parameter
Parameter
Parameter of interest
- sigmafloat
Number of standard deviations for the confidence level
- reoptimizebool
Re-optimize other parameters, when computing the confidence region.
- **kwargsdict
Keyword argument passed ot the confidence estimation method.
- Returns
- resultdict
Dictionary with keys “errp”, ‘errn”, “success” and “nfev”.
-
covariance
(backend='minuit', **kwargs)[source]¶ Estimate the covariance matrix.
Assumes that the model parameters are already optimised.
- Parameters
- backendstr
Which backend to use (see
gammapy.modeling.registry
)
- Returns
- result
CovarianceResult
Results
- result
-
minos_contour
(x, y, numpoints=10, sigma=1.0)[source]¶ Compute MINOS contour.
Calls
iminuit.Minuit.mncontour
.This is a contouring algorithm for a 2D function which is not simply the fit statistic function. That 2D function is given at each point
(par_1, par_2)
by re-optimising all other free parameters, and taking the fit statistic at that point.Very compute-intensive and slow.
- Parameters
- x, y
Parameter
Parameters of interest
- numpointsint
Number of contour points
- sigmafloat
Number of standard deviations for the confidence level
- x, y
- Returns
- resultdict
Dictionary containing the parameter values defining the contour, with the boolean flag “success” and the info objects from
mncontour
.
-
optimize
(backend='minuit', **kwargs)[source]¶ Run the optimization.
- Parameters
- backendstr
Which backend to use (see
gammapy.modeling.registry
)- **kwargsdict
Keyword arguments passed to the optimizer. For the
"minuit"
backend see https://iminuit.readthedocs.io/en/latest/api.html#iminuit.Minuit for a detailed description of the available options. If there is an entry ‘migrad_opts’, those options will be passed toiminuit.Minuit.migrad()
.For the
"sherpa"
backend you can from the optionsmethod = {"simplex", "levmar", "moncar", "gridsearch"}
Those methods are described and compared in detail on http://cxc.cfa.harvard.edu/sherpa/methods/index.html. The available options of the optimization methods are described on the following pages in detail:For the
"scipy"
backend the available options are desribed in detail here: https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.minimize.html
- Returns
- fit_result
FitResult
Results
- fit_result
-
run
(backend='minuit', optimize_opts=None, covariance_opts=None)[source]¶ Run all fitting steps.
- Parameters
- backendstr
Backend used for fitting, default : minuit
- optimize_optsdict
Options passed to
Fit.optimize
.- covariance_optsdict
Options passed to
Fit.covariance
.
- Returns
- fit_result
FitResult
Results
- fit_result
-
stat_profile
(parameter, values=None, bounds=2, nvalues=11, reoptimize=False, optimize_opts=None)[source]¶ Compute fit statistic profile.
The method used is to vary one parameter, keeping all others fixed. So this is taking a “slice” or “scan” of the fit statistic.
See also:
Fit.minos_profile
.- Parameters
- parameter
Parameter
Parameter of interest
- values
Quantity
(optional) Parameter values to evaluate the fit statistic for.
- boundsint or tuple of float
When an
int
is passed the bounds are computed frombounds * sigma
from the best fit value of the parameter, wheresigma
corresponds to the one sigma error on the parameter. If a tuple of floats is given those are taken as the min and max values andnvalues
are linearly spaced between those.- nvaluesint
Number of parameter grid points to use.
- reoptimizebool
Re-optimize other parameters, when computing the fit statistic profile.
- parameter
- Returns
- resultsdict
Dictionary with keys “values”, “stat” and “fit_results”. The latter contains an empty list, if
reoptimize
is set to False
-
stat_surface
(x, y, x_values, y_values, reoptimize=False, **optimize_opts)[source]¶ Compute fit statistic surface.
The method used is to vary two parameters, keeping all others fixed. So this is taking a “slice” or “scan” of the fit statistic.
Caveat: This method can be very computationally intensive and slow
See also:
Fit.minos_contour
- Parameters
- x, y
Parameter
Parameters of interest
- x_values, y_valueslist or
numpy.ndarray
Parameter values to evaluate the fit statistic for.
- reoptimizebool
Re-optimize other parameters, when computing the fit statistic profile.
- **optimize_optsdict
Keyword arguments passed to the optimizer. See
Fit.optimize
for further details.
- x, y
- Returns
- resultsdict
Dictionary with keys “x_values”, “y_values”, “stat” and “fit_results”. The latter contains an empty list, if
reoptimize
is set to False