Hans Dembinski | TU Dortmund
PyPI https://pypi.org/project/iminuit
Source https://github.com/scikit-hep/iminuit
Documentation http://iminuit.readthedocs.org
Latest release v2.16.0
import numpy as np
from matplotlib import pyplot as plt
from argparse import Namespace
# make 10 data points with scatter on y-coordinate
rng = np.random.default_rng(1)
x = np.linspace(0, 1, 10)
ye = np.ones_like(x) * 0.2
y = rng.normal(2 * x + 1, ye)
data = Namespace(x=x, y=y, ye=ye)
# draw that
plt.errorbar(data.x, data.y, data.ye, fmt="o")
plt.xlabel("x")
plt.ylabel("y");
To get parameters of the line:
iminuit.cost
contains all common onesfrom iminuit import Minuit, cost
# line model
def model(x, a, b):
return a + b * x
# see docs for iminuit.cost to pick correct cost function
lsq = cost.LeastSquares(data.x, data.y, data.ye, model)
# - initialize Minuit object by passing cost function
# - Minuit uses a local minimizer, we need to set starting values
m = Minuit(lsq, a=0, b=0)
m
Name | Value | Hesse Error | Minos Error- | Minos Error+ | Limit- | Limit+ | Fixed | |
---|---|---|---|---|---|---|---|---|
0 | a | 0.0 | 0.1 | |||||
1 | b | 0.0 | 0.1 |
lsq
has parameters named a
and b
Minuit
: stateful fitting object# call Migrad minimizer (for other options, see docs)
m.migrad()
Migrad | ||||
---|---|---|---|---|
FCN = 3.959 (chi2/ndof = 0.5) | Nfcn = 34 | |||
EDM = 2.65e-23 (Goal: 0.0002) | time = 0.2 sec | |||
Valid Minimum | No Parameters at limit | |||
Below EDM threshold (goal x 10) | Below call limit | |||
Covariance | Hesse ok | Accurate | Pos. def. | Not forced |
Name | Value | Hesse Error | Minos Error- | Minos Error+ | Limit- | Limit+ | Fixed | |
---|---|---|---|---|---|---|---|---|
0 | a | 1.05 | 0.12 | |||||
1 | b | 1.99 | 0.20 |
a | b | |
---|---|---|
a | 0.0138 | -0.0196 (-0.843) |
b | -0.0196 (-0.843) | 0.0393 |
See docs for more details on Hesse, Minos, EDM.
Output by parts:
m.params
Name | Value | Hesse Error | Minos Error- | Minos Error+ | Limit- | Limit+ | Fixed | |
---|---|---|---|---|---|---|---|---|
0 | a | 1.05 | 0.12 | |||||
1 | b | 1.99 | 0.20 |
m.fmin
Migrad | ||||
---|---|---|---|---|
FCN = 3.959 (chi2/ndof = 0.5) | Nfcn = 34 | |||
EDM = 2.65e-23 (Goal: 0.0002) | time = 0.2 sec | |||
Valid Minimum | No Parameters at limit | |||
Below EDM threshold (goal x 10) | Below call limit | |||
Covariance | Hesse ok | Accurate | Pos. def. | Not forced |
m.covariance
a | b | |
---|---|---|
a | 0.0138 | -0.0196 (-0.843) |
b | -0.0196 (-0.843) | 0.0393 |
# Minos Errors are not automatically computed, you need to run *Minos* explicitly
m.minos()
Migrad | ||||
---|---|---|---|---|
FCN = 3.959 (chi2/ndof = 0.5) | Nfcn = 60 | |||
EDM = 2.65e-23 (Goal: 0.0002) | ||||
Valid Minimum | No Parameters at limit | |||
Below EDM threshold (goal x 10) | Below call limit | |||
Covariance | Hesse ok | Accurate | Pos. def. | Not forced |
Name | Value | Hesse Error | Minos Error- | Minos Error+ | Limit- | Limit+ | Fixed | |
---|---|---|---|---|---|---|---|---|
0 | a | 1.05 | 0.12 | -0.12 | 0.12 | |||
1 | b | 1.99 | 0.20 | -0.20 | 0.20 |
a | b | |||
---|---|---|---|---|
Error | -0.12 | 0.12 | -0.2 | 0.2 |
Valid | True | True | True | True |
At Limit | False | False | False | False |
Max FCN | False | False | False | False |
New Min | False | False | False | False |
a | b | |
---|---|---|
a | 0.0138 | -0.0196 (-0.843) |
b | -0.0196 (-0.843) | 0.0393 |
# Minos errors can be read-out like this
me = m.merrors["a"] # returns a data struct with many fields
me.name, me.is_valid, me.lower, me.upper
('a', True, -0.11755076272896135, 0.11755076272914197)
# draw data and fitted model
plt.errorbar(data.x, data.y, data.ye, fmt="o", label="data")
plt.plot(data.x, model(data.x, *m.values[:]), label="model")
# m.fval is only chi2 if cost function supports that
# telltale sign: m.ndof returns something finite and not None
chi2 = m.fval
ndof = m.ndof
title = [
f"χ²/ndof = {chi2:.1f}/{ndof} = {chi2/ndof:.1f}",
]
for par in m.parameters:
title.append(
f"{par} = {m.values[par]:.2f} +/- {m.errors[par]:.2f}"
)
plt.legend(title="\n".join(title))
plt.xlabel("x")
plt.ylabel("y");
# or use builtin visualization (only for 1D data)
lsq.visualize(m.values)
m.draw_profile("a");
m.draw_profile("b");
m.draw_contour("a", "b");
Minuit.draw_mncontour
draws a 68 % confidence regionMinuit.draw_contour
m.draw_mncontour("a", "b");
# tilt of ellipse indicates anti-correlation
m.covariance
a | b | |
---|---|---|
a | 0.0138 | -0.0196 (-0.843) |
b | -0.0196 (-0.843) | 0.0393 |
lsq.visualize(m.values)
ab = m.mncontour("a", "b")
for a, b in ab:
plt.plot(data.x, model(data.x, a, b), color="k", alpha=0.01)
import jacobi
lsq.visualize(m.values)
y, ycov = jacobi.propagate(lambda p: model(x, *p), m.values, m.covariance)
ye = np.diag(ycov) ** 0.5
plt.fill_between(x, y - ye, y + ye, alpha=0.5);
m.draw_mnmatrix();
# numba_stats has faster versions of statistical distributions in scipy.stats
from numba_stats import norm
rng = np.random.default_rng(1)
x = rng.normal(size=100)
nll = cost.UnbinnedNLL(x, norm.pdf)
# uh oh, starting value for loc far away
m = Minuit(nll, loc=50, scale=1)
m.migrad()
Migrad | ||||
---|---|---|---|---|
FCN = 1.487e+05 | Nfcn = 137 | |||
EDM = nan (Goal: 0.0002) | ||||
INVALID Minimum | No Parameters at limit | |||
ABOVE EDM threshold (goal x 10) | Below call limit | |||
Covariance | Hesse ok | Accurate | Pos. def. | Not forced |
Name | Value | Hesse Error | Minos Error- | Minos Error+ | Limit- | Limit+ | Fixed | |
---|---|---|---|---|---|---|---|---|
0 | loc | 5e1 | nan | |||||
1 | scale | 1 | nan |
loc | scale | |
---|---|---|
loc | nan | nan |
scale | nan | nan |
plt.hist(x)
xm = np.linspace(30, 60, 1000)
ym = norm.pdf(xm, *m.values) * 50
plt.plot(xm, ym);
norm.pdf([1, 10, 100, 1000], 0, 1)
array([2.41970725e-01, 7.69459863e-23, 0.00000000e+00, 0.00000000e+00])
loc
we fix the issuem = Minuit(nll, loc=0, scale=1)
m.migrad()
Migrad | ||||
---|---|---|---|---|
FCN = 251.6 | Nfcn = 31 | |||
EDM = 2.2e-05 (Goal: 0.0002) | ||||
Valid Minimum | No Parameters at limit | |||
Below EDM threshold (goal x 10) | Below call limit | |||
Covariance | Hesse ok | Accurate | Pos. def. | Not forced |
Name | Value | Hesse Error | Minos Error- | Minos Error+ | Limit- | Limit+ | Fixed | |
---|---|---|---|---|---|---|---|---|
0 | loc | -0.07 | 0.09 | |||||
1 | scale | 0.85 | 0.06 |
loc | scale | |
---|---|---|
loc | 0.00725 | 2.67e-06 |
scale | 2.67e-06 | 0.00362 |
# Generate signal and background mixture
rng = np.random.default_rng(1)
s = rng.normal(0.5, 0.1, size=1000)
b = rng.exponential(0.8, size=2000)
x = np.append(s, b)
x = x[x < 1]
data = Namespace(x=x, range=(0, 1))
# Bin and draw data
import boost_histogram as bh
# or use numpy.histogram ...
h = bh.Histogram(bh.axis.Regular(50, 0, 1))
h.fill(data.x)
data.edges = h.axes[0].edges
data.n = h.values()
plt.stairs(data.n, data.edges, fill=True);
cost.ExtendedBinnedNLL?
from numba_stats import truncnorm, truncexpon
# scaled cdf
def model(x, ns, mu, sigma, nb, slope):
s = ns * truncnorm.cdf(x, *data.range, mu, sigma)
b = nb * truncexpon.cdf(x, *data.range, 0, slope)
return s + b
nll = cost.ExtendedBinnedNLL(data.n, data.edges, model)
m = Minuit(nll, ns=1, mu=5, sigma=0.1, nb=1, slope=1)
m.migrad()
Migrad | ||||
---|---|---|---|---|
FCN = nan | Nfcn = 441 | |||
EDM = nan (Goal: 0.0002) | ||||
INVALID Minimum | No Parameters at limit | |||
ABOVE EDM threshold (goal x 10) | Below call limit | |||
Covariance | Hesse ok | Accurate | Pos. def. | Not forced |
Name | Value | Hesse Error | Minos Error- | Minos Error+ | Limit- | Limit+ | Fixed | |
---|---|---|---|---|---|---|---|---|
0 | ns | 1 | nan | |||||
1 | mu | 5 | nan | |||||
2 | sigma | 1e-1 | nan | |||||
3 | nb | 1 | nan | |||||
4 | slope | 1 | nan |
ns | mu | sigma | nb | slope | |
---|---|---|---|---|---|
ns | nan | nan | nan | nan | nan |
mu | nan | nan | nan | nan | nan |
sigma | nan | nan | nan | nan | nan |
nb | nan | nan | nan | nan | nan |
slope | nan | nan | nan | nan | nan |
nll(1, 5, 0.1, 1, 1)
nan
m = Minuit(nll, ns=1, mu=0.5, sigma=0.1, nb=1, slope=1)
m.limits["ns", "nb", "sigma", "slope"] = (0, None)
m.limits["mu"] = (0, 1) # not strictly necessary
m.migrad()
Migrad | ||||
---|---|---|---|---|
FCN = 203.9 (chi2/ndof = 4.5) | Nfcn = 411 | |||
EDM = 1.41e-06 (Goal: 0.0002) | ||||
Valid Minimum | No Parameters at limit | |||
Below EDM threshold (goal x 10) | Below call limit | |||
Covariance | Hesse ok | Accurate | Pos. def. | Not forced |
Name | Value | Hesse Error | Minos Error- | Minos Error+ | Limit- | Limit+ | Fixed | |
---|---|---|---|---|---|---|---|---|
0 | ns | 919.1 | 2.0 | 0 | ||||
1 | mu | 0.484 | 0.005 | 0 | 1 | |||
2 | sigma | 0.094 | 0.004 | 0 | ||||
3 | nb | 1.5179e3 | 0.0020e3 | 0 | ||||
4 | slope | 1.3171e3 | 0.0020e3 | 0 |
ns | mu | sigma | nb | slope | |
---|---|---|---|---|---|
ns | 3.99 | -5.53e-05 (-0.006) | 0.000141 (0.017) | -0.00482 (-0.001) | 1.86e-07 |
mu | -5.53e-05 (-0.006) | 2.47e-05 | -5.01e-06 (-0.246) | 3.36e-05 (0.003) | -3.24e-08 |
sigma | 0.000141 (0.017) | -5.01e-06 (-0.246) | 1.68e-05 | -8.58e-05 (-0.010) | 1.19e-08 |
nb | -0.00482 (-0.001) | 3.36e-05 (0.003) | -8.58e-05 (-0.010) | 3.99 | -1.13e-07 |
slope | 1.86e-07 | -3.24e-08 | 1.19e-08 | -1.13e-07 | 4 |
nll.visualize(m.values)
m.draw_profile("ns");
m.draw_profile("slope");
m.limits["slope"] = (0, 2)
m.migrad()
Migrad | ||||
---|---|---|---|---|
FCN = 40.6 (chi2/ndof = 0.9) | Nfcn = 744 | |||
EDM = 3.19e-05 (Goal: 0.0002) | ||||
Valid Minimum | No Parameters at limit | |||
Below EDM threshold (goal x 10) | Below call limit | |||
Covariance | Hesse ok | Accurate | Pos. def. | Not forced |
Name | Value | Hesse Error | Minos Error- | Minos Error+ | Limit- | Limit+ | Fixed | |
---|---|---|---|---|---|---|---|---|
0 | ns | 950 | 50 | 0 | ||||
1 | mu | 0.498 | 0.005 | 0 | 1 | |||
2 | sigma | 0.093 | 0.005 | 0 | ||||
3 | nb | 1.49e3 | 0.05e3 | 0 | ||||
4 | slope | 0.81 | 0.07 | 0 | 2 |
ns | mu | sigma | nb | slope | |
---|---|---|---|---|---|
ns | 2.45e+03 | -0.0228 (-0.097) | 0.108 (0.484) | -1.5e+03 (-0.554) | -0.716 (-0.216) |
mu | -0.0228 (-0.097) | 2.26e-05 | -2.86e-06 (-0.134) | 0.0228 (0.088) | -6.01e-05 (-0.188) |
sigma | 0.108 (0.484) | -2.86e-06 (-0.134) | 2.02e-05 | -0.108 (-0.438) | -4.48e-05 (-0.149) |
nb | -1.5e+03 (-0.554) | 0.0228 (0.088) | -0.108 (-0.438) | 2.99e+03 | 0.716 (0.195) |
slope | -0.716 (-0.216) | -6.01e-05 (-0.188) | -4.48e-05 (-0.149) | 0.716 (0.195) | 0.0045 |
nll.visualize(m.values)
m.draw_mnmatrix(figsize=(10, 7));
# add another background component to model to create ambiguity
def model(x, ns, mu, sigma, nb1, slope1, nb2, slope2):
s = ns * truncnorm.cdf(x, *data.range, mu, sigma)
b1 = nb1 * truncexpon.cdf(x, *data.range, 0, slope1)
b2 = nb2 * truncexpon.cdf(x, *data.range, 0, slope2)
return s + b1 + b2
nll = cost.ExtendedBinnedNLL(data.n, data.edges, model)
m = Minuit(nll, ns=1, mu=0.5, sigma=0.1, nb1=1, slope1=1, nb2=1, slope2=0.1)
m.limits["ns", "nb1", "nb2", "sigma"] = (0, None)
m.limits["mu"] = (0, 1)
m.limits["slope1", "slope2"] = (0, 2)
m.migrad()
Migrad | ||||
---|---|---|---|---|
FCN = 694.8 (chi2/ndof = 16.2) | Nfcn = 1118 | |||
EDM = 3.37e-10 (Goal: 0.0002) | ||||
INVALID Minimum | No Parameters at limit | |||
Below EDM threshold (goal x 10) | Below call limit | |||
Covariance | Hesse FAILED | APPROXIMATE | NOT pos. def. | Not forced |
Name | Value | Hesse Error | Minos Error- | Minos Error+ | Limit- | Limit+ | Fixed | |
---|---|---|---|---|---|---|---|---|
0 | ns | 4.8263e-9 | 0.0000e-9 | 0 | ||||
1 | mu | 8.5436e-1 | 0.0000e-1 | 0 | 1 | |||
2 | sigma | 1.8477e5 | 0.0000e5 | 0 | ||||
3 | nb1 | 7.7481e2 | 0.0000e2 | 0 | ||||
4 | slope1 | 1.3316 | 0.0000 | 0 | 2 | |||
5 | nb2 | 1.6622e3 | 0.0000e3 | 0 | ||||
6 | slope2 | 1.3317 | 0.0000 | 0 | 2 |
nll.visualize(m.values)
# resolve ambiguity by fixing nb2 to 0... but forget to fix slope2, too
m.values = [p.value for p in m.init_params]
m.values["nb2"] = 0
m.fixed["nb2"] = True
m.migrad()
Migrad | ||||
---|---|---|---|---|
FCN = 40.6 (chi2/ndof = 0.9) | Nfcn = 1808 | |||
EDM = 3.68e-07 (Goal: 0.0002) | ||||
INVALID Minimum | No Parameters at limit | |||
Below EDM threshold (goal x 10) | Below call limit | |||
Covariance | Hesse FAILED | APPROXIMATE | NOT pos. def. | Not forced |
Name | Value | Hesse Error | Minos Error- | Minos Error+ | Limit- | Limit+ | Fixed | |
---|---|---|---|---|---|---|---|---|
0 | ns | 9.5004e2 | 0.0000e2 | 0 | ||||
1 | mu | 4.9776e-1 | 0.0000e-1 | 0 | 1 | |||
2 | sigma | 9.3454e-2 | 0.0000e-2 | 0 | ||||
3 | nb1 | 1.487e3 | 0.000e3 | 0 | ||||
4 | slope1 | 8.1152e-1 | 0.0000e-1 | 0 | 2 | |||
5 | nb2 | 0 | 0 | 0 | yes | |||
6 | slope2 | 1e-1 | 0e-1 | 0 | 2 |
nll.visualize(m.values)
slope2
not constrained by data now, since nb2 = 0
slope2
resolve the issuem.values = [p.value for p in m.init_params]
m.fixed["slope2"] = True
m.migrad()
Migrad | ||||
---|---|---|---|---|
FCN = 40.56 (chi2/ndof = 0.9) | Nfcn = 2136 | |||
EDM = 2.02e-05 (Goal: 0.0002) | ||||
Valid Minimum | No Parameters at limit | |||
Below EDM threshold (goal x 10) | Below call limit | |||
Covariance | Hesse ok | Accurate | Pos. def. | Not forced |
Name | Value | Hesse Error | Minos Error- | Minos Error+ | Limit- | Limit+ | Fixed | |
---|---|---|---|---|---|---|---|---|
0 | ns | 950 | 50 | 0 | ||||
1 | mu | 0.498 | 0.005 | 0 | 1 | |||
2 | sigma | 0.093 | 0.005 | 0 | ||||
3 | nb1 | 1.49e3 | 0.05e3 | 0 | ||||
4 | slope1 | 0.81 | 0.07 | 0 | 2 | |||
5 | nb2 | 1 | 0 | 0 | yes | |||
6 | slope2 | 1e-1 | 0e-1 | 0 | 2 | yes |
ns | mu | sigma | nb1 | slope1 | nb2 | slope2 | |
---|---|---|---|---|---|---|---|
ns | 2.45e+03 | -0.0228 (-0.097) | 0.108 (0.484) | -1.5e+03 (-0.554) | -0.716 (-0.215) | 0 | 0 |
mu | -0.0228 (-0.097) | 2.26e-05 | -2.87e-06 (-0.134) | 0.0229 (0.088) | -6.05e-05 (-0.189) | 0 | 0 |
sigma | 0.108 (0.484) | -2.87e-06 (-0.134) | 2.02e-05 | -0.108 (-0.439) | -4.48e-05 (-0.148) | 0 | 0 |
nb1 | -1.5e+03 (-0.554) | 0.0229 (0.088) | -0.108 (-0.439) | 2.99e+03 | 0.715 (0.194) | 0 | 0 |
slope1 | -0.716 (-0.215) | -6.05e-05 (-0.189) | -4.48e-05 (-0.148) | 0.715 (0.194) | 0.00454 | 0 | 0 |
nb2 | 0 | 0 | 0 | 0 | 0 | 0 | 0 |
slope2 | 0 | 0 | 0 | 0 | 0 | 0 | 0 |
def model(x, ns, mu, sigma, nb, slope):
s = ns * truncnorm.cdf(x, *data.range, mu, sigma)
b = nb * truncexpon.cdf(x, *data.range, 0, slope)
return s + b
b = rng.exponential(1, size=2000)
b = b[b < 1]
data.nb = np.histogram(b, bins=data.edges)[0]
nll = cost.ExtendedBinnedNLL(data.nb, data.edges, model)
m = Minuit(nll, ns=10, mu=0.5, sigma=0.1, nb=200, slope=1)
m.limits["ns", "nb", "sigma"] = (0, None)
m.limits["mu"] = (0, 1)
m.limits["slope"] = (0, 2)
m.migrad()
Migrad | ||||
---|---|---|---|---|
FCN = 33.56 (chi2/ndof = 0.7) | Nfcn = 679 | |||
EDM = 9.4e-05 (Goal: 0.0002) | ||||
Valid Minimum | SOME Parameters at limit | |||
Below EDM threshold (goal x 10) | Below call limit | |||
Covariance | Hesse ok | Accurate | Pos. def. | Not forced |
Name | Value | Hesse Error | Minos Error- | Minos Error+ | Limit- | Limit+ | Fixed | |
---|---|---|---|---|---|---|---|---|
0 | ns | 0.1e3 | 0.5e3 | 0 | ||||
1 | mu | 0.75 | 0.18 | 0 | 1 | |||
2 | sigma | 0.3 | 0.6 | 0 | ||||
3 | nb | 1.1e3 | 0.6e3 | 0 | ||||
4 | slope | 0.8 | 0.7 | 0 | 2 |
ns | mu | sigma | nb | slope | |
---|---|---|---|---|---|
ns | 3.58e+05 | -49 (-0.441) | 410 (0.952) | -3.58e+05 (-0.998) | -472 (-0.985) |
mu | -49 (-0.441) | 0.0345 | -0.0517 (-0.387) | 49 (0.441) | 0.0553 (0.372) |
sigma | 410 (0.952) | -0.0517 (-0.387) | 0.518 | -410 (-0.951) | -0.521 (-0.904) |
nb | -3.58e+05 (-0.998) | 49 (0.441) | -410 (-0.951) | 3.59e+05 | 472 (0.984) |
slope | -472 (-0.985) | 0.0553 (0.372) | -0.521 (-0.904) | 472 (0.984) | 0.64 |
nll.visualize(m.values)
mu
and sigma
nll_mod = nll + cost.NormalConstraint(("mu", "sigma"), (0.5, 0.1), (0.1, 0.05))
m = Minuit(nll_mod, ns=1, mu=0.5, sigma=0.1, nb=1, slope=1)
m.limits["ns", "nb", "sigma"] = (0, None)
m.limits["mu"] = (0, 1)
m.limits["slope"] = (0, 2)
m.migrad()
Migrad | ||||
---|---|---|---|---|
FCN = 33.89 (chi2/ndof = 0.7) | Nfcn = 732 | |||
EDM = 6.68e-06 (Goal: 0.0002) | ||||
Valid Minimum | SOME Parameters at limit | |||
Below EDM threshold (goal x 10) | Below call limit | |||
Covariance | Hesse ok | Accurate | Pos. def. | Not forced |
Name | Value | Hesse Error | Minos Error- | Minos Error+ | Limit- | Limit+ | Fixed | |
---|---|---|---|---|---|---|---|---|
0 | ns | 0 | 15 | 0 | ||||
1 | mu | 0.5 | 0.1 | 0 | 1 | |||
2 | sigma | 0.10 | 0.05 | 0 | ||||
3 | nb | 1.255e3 | 0.035e3 | 0 | ||||
4 | slope | 1.06 | 0.11 | 0 | 2 |
ns | mu | sigma | nb | slope | |
---|---|---|---|---|---|
ns | 0.217 | -0.000413 (-0.009) | 0.000106 (0.005) | -0.396 (-0.024) | -0.000339 (-0.006) |
mu | -0.000413 (-0.009) | 0.01 | -1.83e-07 | 0.000752 | 5.96e-07 |
sigma | 0.000106 (0.005) | -1.83e-07 | 0.0025 | -0.000193 | -1.65e-07 |
nb | -0.396 (-0.024) | 0.000752 | -0.000193 | 1.26e+03 | 0.000617 |
slope | -0.000339 (-0.006) | 5.96e-07 | -1.65e-07 | 0.000617 | 0.0126 |
nll_mod.visualize(m.values)
ns
is 0, so fitted value may be close to lower limitdef model(x, ns, mu, sigma, nb, slope):
s = ns * truncnorm.cdf(x, *data.range, mu, sigma)
b = nb * truncexpon.cdf(x, *data.range, 0, slope)
return s + b
nll = cost.ExtendedBinnedNLL(data.n, data.edges, model)
m1 = Minuit(nll, ns=1000, mu=0.5, sigma=0.1, nb=1000, slope=1)
m1.limits["ns", "nb", "sigma"] = (0, None)
m1.limits["mu"] = (0, 1)
m1.limits["slope"] = (0, 2)
m1.migrad()
Migrad | ||||
---|---|---|---|---|
FCN = 40.6 (chi2/ndof = 0.9) | Nfcn = 144 | |||
EDM = 6.35e-05 (Goal: 0.0002) | ||||
Valid Minimum | No Parameters at limit | |||
Below EDM threshold (goal x 10) | Below call limit | |||
Covariance | Hesse ok | Accurate | Pos. def. | Not forced |
Name | Value | Hesse Error | Minos Error- | Minos Error+ | Limit- | Limit+ | Fixed | |
---|---|---|---|---|---|---|---|---|
0 | ns | 950 | 50 | 0 | ||||
1 | mu | 0.498 | 0.005 | 0 | 1 | |||
2 | sigma | 0.093 | 0.005 | 0 | ||||
3 | nb | 1.49e3 | 0.05e3 | 0 | ||||
4 | slope | 0.81 | 0.07 | 0 | 2 |
ns | mu | sigma | nb | slope | |
---|---|---|---|---|---|
ns | 2.45e+03 | -0.0228 (-0.097) | 0.108 (0.484) | -1.5e+03 (-0.554) | -0.716 (-0.216) |
mu | -0.0228 (-0.097) | 2.26e-05 | -2.86e-06 (-0.134) | 0.0228 (0.088) | -6.01e-05 (-0.188) |
sigma | 0.108 (0.484) | -2.86e-06 (-0.134) | 2.02e-05 | -0.108 (-0.438) | -4.48e-05 (-0.149) |
nb | -1.5e+03 (-0.554) | 0.0228 (0.088) | -0.108 (-0.438) | 2.98e+03 | 0.716 (0.195) |
slope | -0.716 (-0.216) | -6.01e-05 (-0.188) | -4.48e-05 (-0.149) | 0.716 (0.195) | 0.00449 |
m2 = Minuit(nll, ns=1000, mu=0.5, sigma=0.1, nb=1000, slope=1)
m2.limits["ns", "nb"] = (0, None)
# set tight limits around previously fitted value of sigma
m2.limits["sigma"] = m1.values["sigma"] - 1e-3, m1.values["sigma"] + 1e-3
m2.limits["mu"] = (0, 1)
m2.limits["slope"] = (0, 2)
m2.migrad()
Migrad | ||||
---|---|---|---|---|
FCN = 40.65 (chi2/ndof = 0.9) | Nfcn = 141 | |||
EDM = 2.11e-06 (Goal: 0.0002) | ||||
Valid Minimum | SOME Parameters at limit | |||
Below EDM threshold (goal x 10) | Below call limit | |||
Covariance | Hesse ok | Accurate | Pos. def. | Not forced |
Name | Value | Hesse Error | Minos Error- | Minos Error+ | Limit- | Limit+ | Fixed | |
---|---|---|---|---|---|---|---|---|
0 | ns | 955.5 | 1.8 | 0 | ||||
1 | mu | 0.498 | 0.005 | 0 | 1 | |||
2 | sigma | 0.0945 | 0.0017 | 0.0925 | 0.0945 | |||
3 | nb | 1.4815e3 | 0.0018e3 | 0 | ||||
4 | slope | 0.81 | 0.06 | 0 | 2 |
ns | mu | sigma | nb | slope | |
---|---|---|---|---|---|
ns | 3.33 | -1.02e-05 (-0.001) | 4.54e-11 | -0.00281 (-0.001) | -0.000644 (-0.005) |
mu | -1.02e-05 (-0.001) | 2.24e-05 | -7.37e-13 | 6.55e-06 | -6.95e-05 (-0.226) |
sigma | 4.54e-11 | -7.37e-13 | 1.97e-13 | -2.93e-11 | -3.28e-12 |
nb | -0.00281 (-0.001) | 6.55e-06 | -2.93e-11 | 3.33 | 0.000416 (0.004) |
slope | -0.000644 (-0.005) | -6.95e-05 (-0.226) | -3.28e-12 | 0.000416 (0.004) | 0.00422 |
sigma
is now too smallns
, are artificially reduced 🙁ipywidgets.interact
offers an easy way to visualize how a function reacts to parameter changesUseful
m1.interactive()
HBox(children=(Output(), VBox(children=(HBox(children=(Button(description='Fit', style=ButtonStyle()), ToggleB…
cost.BarlowBeestonLite?
rng = np.random.default_rng(1)
# generate templates from "simulation"
nb = np.histogram(rng.normal(0.5, 0.1, size=1000), bins=data.edges)[0]
ns = np.histogram(rng.exponential(1, size=1000), bins=data.edges)[0]
nll = cost.BarlowBeestonLite(data.n, data.edges, (nb, ns), name=("ns", "nb"))
m = Minuit(nll, ns=1, nb=1)
m.migrad()
Migrad | ||||
---|---|---|---|---|
FCN = 31.92 (chi2/ndof = 0.7) | Nfcn = 113 | |||
EDM = 1.8e-06 (Goal: 0.0002) | ||||
Valid Minimum | No Parameters at limit | |||
Below EDM threshold (goal x 10) | Below call limit | |||
Covariance | Hesse ok | Accurate | Pos. def. | Not forced |
Name | Value | Hesse Error | Minos Error- | Minos Error+ | Limit- | Limit+ | Fixed | |
---|---|---|---|---|---|---|---|---|
0 | ns | 990 | 70 | |||||
1 | nb | 1.45e3 | 0.09e3 |
ns | nb | |
---|---|---|
ns | 4.58e+03 | -2.62e+03 (-0.448) |
nb | -2.62e+03 (-0.448) | 7.47e+03 |
ns
larger compared to fitting parametric model (70 vs 50), since uncertainty of template not negligiblenll.visualize(m.values)
m.interactive()
HBox(children=(Output(), VBox(children=(HBox(children=(Button(description='Fit', style=ButtonStyle()), ToggleB…
BarlowBeestonLite
also supports higher dimensional data (2D, 3D, ...)