FitResult#
In this tutorial, we will explore the FitResult
of zfit. Specifically, we will examine the error methods hesse and errors as well as attributes like info, valid etc. We will also provide an example with weighted data to demonstrate how FitResult works with weighted datasets.
We will start out by creating a simple gaussian model and sampling some data from it. We will then fit the data with the same model and explore the FitResult
.
import numpy as np
import zfit
import zfit.z.numpy as znp
obs = zfit.Space('x', 0, 10)
mu = zfit.Parameter('mu', 5, 0, 10)
sigma = zfit.Parameter('sigma', 1, 0, 10)
nsig = zfit.Parameter('nsig', 1000, 0, 10000)
gauss = zfit.pdf.Gauss(obs=obs, mu=mu, sigma=sigma,extended=nsig)
data = gauss.sample()
print(f"The sampled data (poisson fluctuated) has {data.nevents} events.")
The sampled data (poisson fluctuated) has 964 events.
We use an extended likelihood to fit the data.
nll = zfit.loss.ExtendedUnbinnedNLL(model=gauss, data=data)
minimizer = zfit.minimize.Minuit()
result = minimizer.minimize(nll)
Simply printing the result will give you a beautified overview of the fit result.
print(result)
FitResult
of
<ExtendedUnbinnedNLL model=[<zfit.<class 'zfit.models.dist_tfp.Gauss'> params=[mu, sigma]] data=[<zfit.Data: Data obs=('x',) shape=(964, 1)>] constraints=[]>
with
<Minuit Minuit tol=0.001>
╒═════════╤═════════════╤══════════════════╤════════╤══════════════════════════════╕
│ valid │ converged │ param at limit │ edm │ approx. fmin (full | opt.) │
╞═════════╪═════════════╪══════════════════╪════════╪══════════════════════════════╡
│
True
│ True
│ False
│ 0.0002 │ -4306.73 | 10013.85 │
╘═════════╧═════════════╧══════════════════╧════════╧══════════════════════════════╛
Parameters
name value (rounded) at limit
------ ------------------ ----------
nsig 963.567 False
mu 5.01469 False
sigma 0.984388 False
What happened#
First and foremost, the FitResult contains all the information about what happened with the minimization, most notably the loss
that was minimized, the minimizer
that was used and the params
that were fitted (the latter has a beautified presentation).
print(f"""
loss: {result.loss}
minimizer: {result.minimizer}
params: {result.params}
""")
loss: <ExtendedUnbinnedNLL model=[<zfit.<class 'zfit.models.dist_tfp.Gauss'> params=[mu, sigma]] data=[<zfit.Data: Data obs=('x',) shape=(964, 1)>] constraints=[]>
minimizer: <Minuit Minuit tol=0.001>
params: name value (rounded) at limit
------ ------------------ ----------
nsig 963.567 False
mu 5.01469 False
sigma 0.984388 False
params#
params
contains all the information of the parameters that was ever added to them. This includes the output of uncertainty methods, limits and much more.
The actual content looks like this:
print(f"params raw: {repr(result.params)}")
params raw: {<zfit.Parameter 'nsig' floating=True value=963.6>: {'value': 963.566634312413}, <zfit.Parameter 'mu' floating=True value=5.015>: {'value': 5.014692431429814}, <zfit.Parameter 'sigma' floating=True value=0.9844>: {'value': 0.9843877678152884}}
The FitResult
has a lot of attributes and methods. We will now explore some of them.
All the displayed information can be accessed via the attributes of the FitResult
object, namely
valid: whether the fit converged and is in general valid
converged: whether the fit converged
param at limit: whether any parameter is at its limit (approximate, hard to estimate)
edm: estimated distance to minimum
fmin: the minimum of the function, i.e. the negative log likelihood
values: the parameter values at the minimum in an array-like object
print(f"""
valid: {result.valid}
converged: {result.converged}
param at limit: {result.params_at_limit}
edm: {result.edm}
fmin: {result.fmin}
optimal values: {result.values}
""")
valid: True
converged: True
param at limit: False
edm: 0.00020095214326310348
fmin: -4306.733397818645
optimal values: [963.566634312413, 5.014692431429814, 0.9843877678152884]
Error methods#
There are two main ways to estimate the uncertainties: Either using a profiling method that varies the parameters one by one and finds the point of 1 sigma (or the specified n sigma), resulting in asymmetric errors, or using a matrix inversion method that calculates an approximation of the former by using a second derivative matrix.
The first method is called errors
and the second hesse
. Both methods are available in the FitResult
object.
Pitfall weights#
For weighted likelihoods, the errors
method will not report the correct uncertainties. Instead, hesse
should be used
as it will, by default, calculate the asymptotic correct approximations for weights as we will see a few lines below.
Arguments#
Both methods take some common arguments:
params
: the parameters to calculate the errors for. IfNone
, all parameters will be used. (this can be expensive!)name
: the name of the new result. IfNone
, the name will be chosen automatically.cl
: the confidence level for the errors. The default is 0.68, which corresponds to 1 sigma.method
: the method to use. The default isNone
which will use the default method of the uncertainty estimator.
errors, new_result = result.errors(name="errors")
print(f"New result: {new_result}")
print(result)
New result: None
FitResult
of
<ExtendedUnbinnedNLL model=[<zfit.<class 'zfit.models.dist_tfp.Gauss'> params=[mu, sigma]] data=[<zfit.Data: Data obs=('x',) shape=(964, 1)>] constraints=[]>
with
<Minuit Minuit tol=0.001>
╒═════════╤═════════════╤══════════════════╤════════╤══════════════════════════════╕
│ valid │ converged │ param at limit │ edm │ approx. fmin (full | opt.) │
╞═════════╪═════════════╪══════════════════╪════════╪══════════════════════════════╡
│
True
│ True
│ False
│ 0.0002 │ -4306.73 | 10013.85 │
╘═════════╧═════════════╧══════════════════╧════════╧══════════════════════════════╛
Parameters
name value (rounded) errors at limit
------ ------------------ ------------------- ----------
nsig 963.567 - 30 + 32 False
mu 5.01469 - 0.032 + 0.032 False
sigma 0.984388 - 0.022 + 0.023 False
The uncertainties are added to the fit result. The new_result
is usually None
but in case a new minimum was found, it will be returned
as the new result. In this case, the old result will be rendered invalid.
There are currently two implementations, the minos method from iminuit
(as minuit_minos
) and a completely independent implementation
(zfit_errors
).
More information#
To find more information about the uncertainty estimation, the return value can be inspected. This is though also automatically added to the result
to each parameter. Looking again at the raw params
attribute, we find that all the information is there:
Note: this part is still under WIP and future plans are to standardize these attributes as well. Any ideas or inputs are very welcome!
print(f"params raw: {repr(result.params)}")
params raw: {<zfit.Parameter 'nsig' floating=True value=963.6>: {'value': 963.566634312413, 'errors': {'lower': -30.28867825766317, 'upper': 31.822274981786542, 'is_valid': True, 'upper_valid': True, 'lower_valid': True, 'at_lower_limit': False, 'at_upper_limit': False, 'nfcn': 49, 'original': <MError number=0 name='nsig' lower=-30.28867825766317 upper=31.822274981786542 is_valid=True lower_valid=True upper_valid=True at_lower_limit=False at_upper_limit=False at_lower_max_fcn=False at_upper_max_fcn=False lower_new_min=False upper_new_min=False nfcn=49 min=963.566634312413>, 'cl': 0.68268949}}, <zfit.Parameter 'mu' floating=True value=5.015>: {'value': 5.014692431429814, 'errors': {'lower': -0.03173019813255402, 'upper': 0.031730107868945555, 'is_valid': True, 'upper_valid': True, 'lower_valid': True, 'at_lower_limit': False, 'at_upper_limit': False, 'nfcn': 20, 'original': <MError number=1 name='mu' lower=-0.03173019813255402 upper=0.031730107868945555 is_valid=True lower_valid=True upper_valid=True at_lower_limit=False at_upper_limit=False at_lower_max_fcn=False at_upper_max_fcn=False lower_new_min=False upper_new_min=False nfcn=20 min=5.014692431429814>, 'cl': 0.68268949}}, <zfit.Parameter 'sigma' floating=True value=0.9844>: {'value': 0.9843877678152884, 'errors': {'lower': -0.021690467331561568, 'upper': 0.02318825957509713, 'is_valid': True, 'upper_valid': True, 'lower_valid': True, 'at_lower_limit': False, 'at_upper_limit': False, 'nfcn': 45, 'original': <MError number=2 name='sigma' lower=-0.021690467331561568 upper=0.02318825957509713 is_valid=True lower_valid=True upper_valid=True at_lower_limit=False at_upper_limit=False at_lower_max_fcn=False at_upper_max_fcn=False lower_new_min=False upper_new_min=False nfcn=45 min=0.9843877678152884>, 'cl': 0.68268949}}}
errors2, _ = result.errors(name="zfit_unc", method="zfit_errors")
print(result)
FitResult
of
<ExtendedUnbinnedNLL model=[<zfit.<class 'zfit.models.dist_tfp.Gauss'> params=[mu, sigma]] data=[<zfit.Data: Data obs=('x',) shape=(964, 1)>] constraints=[]>
with
<Minuit Minuit tol=0.001>
╒═════════╤═════════════╤══════════════════╤════════╤══════════════════════════════╕
│ valid │ converged │ param at limit │ edm │ approx. fmin (full | opt.) │
╞═════════╪═════════════╪══════════════════╪════════╪══════════════════════════════╡
│
True
│ True
│ False
│ 0.0002 │ -4306.73 | 10013.85 │
╘═════════╧═════════════╧══════════════════╧════════╧══════════════════════════════╛
Parameters
name value (rounded) errors zfit_unc at limit
------ ------------------ ------------------- ------------------- ----------
nsig 963.567 - 30 + 32 - 30 + 32 False
mu 5.01469 - 0.032 + 0.032 - 0.032 + 0.032 False
sigma 0.984388 - 0.022 + 0.023 - 0.022 + 0.023 False
As we see, they both agree well. We can also change the confidence level to 0.95, which corresponds to 2 sigma and recalculate the errors.
errors3, _ = result.errors(name="zfit_2sigma", method="zfit_errors", cl=0.95)
print(result)
FitResult
of
<ExtendedUnbinnedNLL model=[<zfit.<class 'zfit.models.dist_tfp.Gauss'> params=[mu, sigma]] data=[<zfit.Data: Data obs=('x',) shape=(964, 1)>] constraints=[]>
with
<Minuit Minuit tol=0.001>
╒═════════╤═════════════╤══════════════════╤════════╤══════════════════════════════╕
│ valid │ converged │ param at limit │ edm │ approx. fmin (full | opt.) │
╞═════════╪═════════════╪══════════════════╪════════╪══════════════════════════════╡
│
True
│ True
│ False
│ 0.0002 │ -4306.73 | 10013.85 │
╘═════════╧═════════════╧══════════════════╧════════╧══════════════════════════════╛
Parameters
name value (rounded) errors zfit_unc zfit_2sigma at limit
------ ------------------ ------------------- ------------------- ------------------- ----------
nsig 963.567 - 30 + 32 - 30 + 32 - 59 + 63 False
mu 5.01469 - 0.032 + 0.032 - 0.032 + 0.032 - 0.062 + 0.062 False
sigma 0.984388 - 0.022 + 0.023 - 0.022 + 0.023 - 0.042 + 0.046 False
Hesse#
The hesse method approximates the errors by calculating the second derivative matrix of the function and inverting it.
As for errors
there are two implementations, one from iminuit
(minuit_hesse
) and one from zfit
(hesse_np
).
Additionally, the hesse
has a third option, approx
: this is the approximation of the hessian estimated by the minimizer
during the minimization procedure. This however can be None
! Also, the accuracy can be low, especially if the
fit converged rapidly.
hesse = result.hesse(name="h minuit", method="minuit_hesse", cl=0.95) # can also take the cl argument
hesse2 = result.hesse(name="h zfit", method="hesse_np")
hesse3 = result.hesse(name="h approx", method="approx")
print(result)
FitResult
of
<ExtendedUnbinnedNLL model=[<zfit.<class 'zfit.models.dist_tfp.Gauss'> params=[mu, sigma]] data=[<zfit.Data: Data obs=('x',) shape=(964, 1)>] constraints=[]>
with
<Minuit Minuit tol=0.001>
╒═════════╤═════════════╤══════════════════╤════════╤══════════════════════════════╕
│ valid │ converged │ param at limit │ edm │ approx. fmin (full | opt.) │
╞═════════╪═════════════╪══════════════════╪════════╪══════════════════════════════╡
│
True
│ True
│ False
│ 0.0002 │ -4306.73 | 10013.85 │
╘═════════╧═════════════╧══════════════════╧════════╧══════════════════════════════╛
Parameters
name value (rounded) errors zfit_unc zfit_2sigma h minuit h zfit h approx at limit
------ ------------------ ------------------- ------------------- ------------------- ----------- ----------- ----------- ----------
nsig 963.567 - 30 + 32 - 30 + 32 - 59 + 63 +/- 61 +/- 31 +/- 31 False
mu 5.01469 - 0.032 + 0.032 - 0.032 + 0.032 - 0.062 + 0.062 +/- 0.062 +/- 0.032 +/- 0.032 False
sigma 0.984388 - 0.022 + 0.023 - 0.022 + 0.023 - 0.042 + 0.046 +/- 0.044 +/- 0.022 +/- 0.022 False
Internally, zfit uses by default a numerical approximation of the hessian, which is usually sufficient and good for one-time use.
However, if you want to use the hessian for multiple fits, it is recommended to force it to use the exact gradient provided by the
backend. To make sure one or the other is used, you can set zfit.run.set_autograd_mode(False)
or zfit.run.set_autograd_mode(True)
.
with zfit.run.set_autograd_mode(True):
hesse4 = result.hesse(name="h autograd", method="hesse_np")
print(result)
FitResult
of
<ExtendedUnbinnedNLL model=[<zfit.<class 'zfit.models.dist_tfp.Gauss'> params=[mu, sigma]] data=[<zfit.Data: Data obs=('x',) shape=(964, 1)>] constraints=[]>
with
<Minuit Minuit tol=0.001>
╒═════════╤═════════════╤══════════════════╤════════╤══════════════════════════════╕
│ valid │ converged │ param at limit │ edm │ approx. fmin (full | opt.) │
╞═════════╪═════════════╪══════════════════╪════════╪══════════════════════════════╡
│
True
│ True
│ False
│ 0.0002 │ -4306.73 | 10013.85 │
╘═════════╧═════════════╧══════════════════╧════════╧══════════════════════════════╛
Parameters
name value (rounded) errors zfit_unc zfit_2sigma h minuit h zfit h approx h autograd at limit
------ ------------------ ------------------- ------------------- ------------------- ----------- ----------- ----------- ------------ ----------
nsig 963.567 - 30 + 32 - 30 + 32 - 59 + 63 +/- 61 +/- 31 +/- 31 +/- 31 False
mu 5.01469 - 0.032 + 0.032 - 0.032 + 0.032 - 0.062 + 0.062 +/- 0.062 +/- 0.032 +/- 0.032 +/- 0.032 False
sigma 0.984388 - 0.022 + 0.023 - 0.022 + 0.023 - 0.042 + 0.046 +/- 0.044 +/- 0.022 +/- 0.022 +/- 0.022 False
Weighted uncertainties#
A weighted likelihood is technically not a likelihood anymore and the errors are not calculated correctly. However, the hesse method can be corrected for weights, which is done automatically as soon as the dataset is weighted.
The method used is the asymptotically correct
yet computationally expensive method described in Parameter uncertainties in weighted unbinned maximum likelihood fits.
The computation involves the jacobian of each event that can be expensive to compute. Again, zfit offers the possibility to use the autograd or the numerical jacobian.
weighted_data = zfit.Data.from_tensor(obs=obs, tensor=data.value(), weights=znp.random.uniform(0.1, 5, size=(data.nevents,)))
weighted_nll = zfit.loss.UnbinnedNLL(model=gauss, data=weighted_data)
weighted_result = minimizer.minimize(weighted_nll)
weighted_result.errors(name="errors")
({<zfit.Parameter 'mu' floating=True value=5.031>: {'lower': -0.019707933988665726,
'upper': 0.019708444538968355,
'is_valid': True,
'upper_valid': True,
'lower_valid': True,
'at_lower_limit': False,
'at_upper_limit': False,
'nfcn': 12,
'original': ┌──────────┬───────────────────────┐
│ │ mu │
├──────────┼───────────┬───────────┤
│ Error │ -0.02 │ 0.02 │
│ Valid │ True │ True │
│ At Limit │ False │ False │
│ Max FCN │ False │ False │
│ New Min │ False │ False │
└──────────┴───────────┴───────────┘,
'cl': 0.68268949},
<zfit.Parameter 'sigma' floating=True value=0.9827>: {'lower': -0.013903836442438032,
'upper': 0.013968707147453756,
'is_valid': True,
'upper_valid': True,
'lower_valid': True,
'at_lower_limit': False,
'at_upper_limit': False,
'nfcn': 12,
'original': ┌──────────┬───────────────────────┐
│ │ sigma │
├──────────┼───────────┬───────────┤
│ Error │ -0.014 │ 0.014 │
│ Valid │ True │ True │
│ At Limit │ False │ False │
│ Max FCN │ False │ False │
│ New Min │ False │ False │
└──────────┴───────────┴───────────┘,
'cl': 0.68268949}},
None)
with zfit.run.set_autograd_mode(True):
weighted_result.hesse(name="hesse autograd")
weighted_result.hesse(name="hesse autograd np", method="hesse_np")
with zfit.run.set_autograd_mode(False):
weighted_result.hesse(name="hesse numeric")
weighted_result.hesse(name="hesse numeric np", method="hesse_np")
print(weighted_result) # FIXME: the errors are not correct for the nsig
FitResult
of
<UnbinnedNLL model=[<zfit.<class 'zfit.models.dist_tfp.Gauss'> params=[mu, sigma]] data=[<zfit.Data: Data obs=('x',) shape=(964, 1)>] constraints=[]>
with
<Minuit Minuit tol=0.001>
╒═════════╤═════════════╤══════════════════╤═════════╤══════════════════════════════╕
│ valid │ converged │ param at limit │ edm │ approx. fmin (full | opt.) │
╞═════════╪═════════════╪══════════════════╪═════════╪══════════════════════════════╡
│
True
│ True
│ False
│ 4.5e-05 │ 3484.47 | 9999.646 │
╘═════════╧═════════════╧══════════════════╧═════════╧══════════════════════════════╛
Parameters
name value (rounded) errors hesse autograd hesse autograd np hesse numeric hesse numeric np at limit
------ ------------------ ------------------- ---------------- ------------------- --------------- ------------------ ----------
mu 5.03112 - 0.02 + 0.02 +/- 0.036 +/- 0.036 +/- 0.036 +/- 0.036 False
sigma 0.982724 - 0.014 + 0.014 +/- 0.025 +/- 0.025 +/- 0.025 +/- 0.025 False
As we can see, the errors are underestimated for the nuisance parameters using the minos method while the hesse method is correct.
Standardized minimizer information#
Some of the minimizers collect information about the loss during the minimization process, such as an approximation of the hessian, inverse hessian, gradient etc. They can be retrieved via approx
, note however that they can be None
.
hessian
and inv_hessian
have an invert
argument: if True and only one of the two is available, the other one will be inverted to obtain the request.
print(f"Approx gradient: {result.approx.gradient()}") # gradient approx not available in iminuit
print(f"Approx hessian (no invert): {result.approx.hessian(invert=False)}") # hessian approximation is also not available
print(f"Approx inverse hessian: {result.approx.inv_hessian(invert=False)}") # inv_hessian is available
print(f"Approx hessian (can invert): {result.approx.hessian(invert=True)}") # allowing the invert now inverts the inv_hessian
Approx gradient: None
Approx hessian (no invert): None
Approx inverse hessian: [[ 9.63326575e+02 -3.63930352e-10 2.53892648e-10]
[-3.63930352e-10 1.00521956e-03 1.13533110e-06]
[ 2.53892648e-10 1.13533110e-06 5.02250933e-04]]
Approx hessian (can invert): [[ 1.03806957e-03 3.76417028e-10 -5.25604975e-10]
[ 3.76417028e-10 9.94810084e+02 -2.24875407e+00]
[-5.25604975e-10 -2.24875407e+00 1.99104170e+03]]
info#
The information returned by the minimizer. CAREFUL! This is a dictionary and can be different for different minimizers. The standardized keys can always be accessed in other ways, such as the approximations of the hessian, the covariance matrix etc.
result.info.keys()
dict_keys(['n_eval', 'minuit', 'original', 'inv_hessian'])
This can be helpful if underlying information from a specific minimizer should be retrieved. For example, the original
key contains the original result from the minimizer while “minuit” is the actual iminuit
minimizer that was used.
result.info.get("original", f"Not available for the minimizer: {result.minimizer}")
Migrad | |
---|---|
FCN = 1.001e+04 | Nfcn = 150 |
EDM = 0.000201 (Goal: 0.001) | |
Valid Minimum | Below EDM threshold (goal x 10) |
No parameters at limit | Below call limit |
Hesse ok | Covariance accurate |
result.info.get("minuit", "Not available, not iminuit used in minimization?")
Migrad | |
---|---|
FCN = 1.001e+04 | Nfcn = 166 |
EDM = 0.000201 (Goal: 0.001) | |
Valid Minimum | Below EDM threshold (goal x 10) |
No parameters at limit | Below call limit |
Hesse ok | Covariance accurate |
Name | Value | Hesse Error | Minos Error- | Minos Error+ | Limit- | Limit+ | Fixed | |
---|---|---|---|---|---|---|---|---|
0 | nsig | 964 | 31 | -30 | 32 | 0 | 1E+04 | |
1 | mu | 5.015 | 0.032 | -0.032 | 0.032 | 0 | 10 | |
2 | sigma | 0.984 | 0.022 | -0.022 | 0.023 | 0 | 10 |
nsig | mu | sigma | ||||
---|---|---|---|---|---|---|
Error | -30 | 32 | -0.032 | 0.032 | -0.022 | 0.023 |
Valid | True | True | True | True | True | True |
At Limit | False | False | False | False | False | False |
Max FCN | False | False | False | False | False | False |
New Min | False | False | False | False | False | False |
nsig | mu | sigma | |
---|---|---|---|
nsig | 963 | -0.0000 | -0 |
mu | -0.0000 | 0.00101 | 0 |
sigma | -0 | 0 | 0.000502 |
Finding problems#
If the fit failed for some reason, valid
may be False. To find the actual reason, message
should be human-readable information about what went wrong. If everything went well, the message will be empty
result.message
''