# Copyright (c) 2020 zfit import zfit # create space obs = zfit.Space("x", limits=(-10, 10)) # parameters mu = zfit.Parameter("mu", 1., -4, 6) sigma = zfit.Parameter("sigma", 1., 0.1, 10) lambd = zfit.Parameter("lambda", -1., -5., 0) frac = zfit.Parameter("fraction", 0.5, 0., 1.) # pdf creation gauss = zfit.pdf.Gauss(mu=mu, sigma=sigma, obs=obs) exponential = zfit.pdf.Exponential(lambd, obs=obs) sum_pdf = zfit.pdf.SumPDF([gauss, exponential], fracs=frac)
# Copyright (c) 2019 zfit # Functors are PDfs and Functions that depend on other PDFs or Functions. They can be used to define # in a custom way combinations of PDFs or wrapping a single PDF. # An example would be to create the sum of two PDFs. Of course this is already implemented in zfit # as the functor SumPDF([pdf1, pdf2], fracs=...). For advanced uses, you # can define your own functor as demonstrated below. import zfit import numpy as np class CombinePolynomials(zfit.pdf.BaseFunctor): """Example of a functor pdf that adds three polynomials. DEMONSTRATION PURPOSE ONLY, DO **NOT** USE IN REAL CASE. """ def __init__(self, angle1, angle2, angle3, name="3dPolynomial"): pdfs = [angle1, angle2, angle3] super().__init__(pdfs=pdfs, name=name) _N_OBS = 3 def _unnormalized_pdf(self, x): cosangle1_data, cosangle2_data, angle3_data = x.unstack_x() cosangle1, cosangle2, angle3 = self.pdfs pdf = cosangle1.pdf(cosangle1_data) + cosangle2.pdf(cosangle2_data) + angle3.pdf(angle3_data) return pdf def create_angular(): # create parameters # c00 = zfit.Parameter(...) # ... and so on cosangle1_space = zfit.Space('cos angle 1', limits=(-1, 1)) cosangle2_space = zfit.Space('cos angle 2', limits=(-1, 1)) angle3_space = zfit.Space('angle 3', limits=(-np.pi, np.pi)) # this part could also be moved inside the __init__, but then the init would need to take all # the coefficients and the spaces as arguments cosangle1_pdf = zfit.pdf.Chebyshev(obs=cosangle1_space, coeffs=[c00, c01, c02]) cosangle2_pdf = zfit.pdf.Chebyshev(obs=cosangle2_space, coeffs=[c10, c11, c12]) angle3_pdf = zfit.pdf.Chebyshev(obs=angle3_space, coeffs=[c20, c21, c22]) return CombinePolynomials(angle1=cosangle1_pdf, angle2=cosangle2_pdf, angle3=angle3_pdf)
# Copyright (c) 2020 zfit import zfit from zfit import z class CustomPDF(zfit.pdf.ZPDF): """1-dimensional PDF implementing the exp(alpha * x) shape.""" _PARAMS = ['alpha'] # specify which parameters to take def _unnormalized_pdf(self, x): # implement function data = z.unstack_x(x) alpha = self.params['alpha'] return z.exp(alpha * data) obs = zfit.Space("obs1", limits=(-4, 4)) custom_pdf = CustomPDF(obs=obs, alpha=0.2) integral = custom_pdf.integrate(limits=(-1, 2)) sample = custom_pdf.sample(n=1000) prob = custom_pdf.pdf(sample)
# Copyright (c) 2020 zfit import tensorflow as tf import zfit from zfit import z class CustomPDF2D(zfit.pdf.BasePDF): """My custom, 2 dimensional pdf. The axes are: Energy, Momentum. """ def __init__(self, param1, param2, param3, obs, name="CustomPDF", ): # we can now do complicated stuff here if needed # only thing: we have to specify explicitly here what is which parameter params = {'super_param': param1, # we can change/compose etc parameters 'param2': param2, 'param3': param3} super().__init__(obs, params, name=name) def _unnormalized_pdf(self, x): energy, momentum = x.unstack_x() param1 = self.params['super_param'] param2 = self.params['param2'] param3 = self.params['param3'] # just a fantasy function probs = param1 * tf.cos(energy ** 2) + tf.math.log(param2 * momentum ** 2) + param3 return probs # add an analytic integral # define the integral function def integral_full(x, limits, norm_range, params, model): lower, upper = limits.limit1d param1 = params['super_param'] param2 = params['param2'] param3 = params['param3'] lower = z.convert_to_tensor(lower) upper = z.convert_to_tensor(upper) # calculate the integral here, dummy integral, wrong! integral = param1 * param2 * param3 + z.reduce_sum([lower, upper]) return integral # define the space over which it is defined. Here, we use the axes lower_full = (-10, zfit.Space.ANY_LOWER) upper_full = (10, zfit.Space.ANY_UPPER) integral_full_limits = zfit.Space(axes=(0, 1), limits=(lower_full, upper_full)) CustomPDF2D.register_analytic_integral(func=integral_full, limits=integral_full_limits) # define the partial integral function def integral_axis1(x, limits, norm_range, params, model): data_0 = x.unstack_x() # data from axis 0 param1 = params['super_param'] param2 = params['param2'] param3 = params['param3'] lower, upper = limits.limit1d lower = z.convert_to_tensor(lower) # the limits are now 1-D, for axis 1 upper = z.convert_to_tensor(upper) # calculate the integral here, dummy integral integral = data_0 ** 2 * param1 * param2 * param3 + z.reduce_sum([lower, upper]) return integral # define the space over which it is defined. Here, we use the axes lower_axis1 = ((-5,),) upper_axis1 = ((zfit.Space.ANY_UPPER,),) integral_axis1_limits = zfit.Space(axes=(1,), limits=(lower_axis1, upper_axis1)) CustomPDF2D.register_analytic_integral(func=integral_axis1, limits=integral_axis1_limits) if __name__ == '__main__': obs = zfit.Space('obs1', (-10, 10)) * zfit.Space('obs2', (-3, 5)) pdf = CustomPDF2D(1, 2, 3, obs=obs) sample = pdf.sample(n=1000)
# Copyright (c) 2019 zfit import zfit from zfit import z import tensorflow as tf # IMPORTANT! The communication of which axis corresponds to which data point happens here. So the user knows now that # he should create this pdf with a space in the obs (x, y, z). class CustomPDF(zfit.pdf.ZPDF): """3-dimensional PDF calculating doing something fancy. Takes x, y, z as data.""" _PARAMS = ['alpha', 'beta'] # specify which parameters to take _N_OBS = 3 def _unnormalized_pdf(self, x): # implement function x, y, z = x.unstack_x() alpha = self.params['alpha'] beta = self.params['beta'] x_new = tf.math.cos(alpha) * x y_new = tf.math.sinh(beta) * y z_new = z + 4.2 return x_new ** 2 + y_new ** 2 + z_new ** 2 xobs = zfit.Space('xobs', (-4, 4)) yobs = zfit.Space('yobs', (-3, 3)) zobs = zfit.Space('z', (-2, 2)) obs = xobs * yobs * zobs alpha = zfit.Parameter("alpha", 0.2) # floating beta = zfit.Parameter("beta", 0.4, floating=False) # non-floating custom_pdf = CustomPDF(obs=obs, alpha=alpha, beta=beta) integral = custom_pdf.integrate(limits=obs) # = 1 since normalized sample = custom_pdf.sample(n=1000) # DO NOT USE THIS FOR TOYS! prob = custom_pdf.pdf(sample) # DO NOT USE THIS FOR TOYS! integral_np, sample_np, prob_np = [integral.numpy(), sample.numpy(), prob.numpy()]
# Copyright (c) 2020 zfit import numpy as np import zfit # create space xobs = zfit.Space('xobs', (-4, 4)) yobs = zfit.Space('yobs', (-3, 5)) zobs = zfit.Space('z', (-2, 4)) obs = xobs * yobs * zobs # parameters mu1 = zfit.Parameter("mu1", 1., -4, 6) mu23 = zfit.Parameter("mu_shared", 1., -4, 6) sigma12 = zfit.Parameter("sigma_shared", 1., 0.1, 10) sigma3 = zfit.Parameter("sigma3", 1., 0.1, 10) # model building, pdf creation gauss_x = zfit.pdf.Gauss(mu=mu1, sigma=sigma12, obs=xobs) gauss_y = zfit.pdf.Gauss(mu=mu23, sigma=sigma12, obs=yobs) gauss_z = zfit.pdf.Gauss(mu=mu23, sigma=sigma3, obs=zobs) product_gauss = zfit.pdf.ProductPDF([gauss_x, gauss_y, gauss_z]) # OR create directly from your 3 dimensional pdf as # model = MyPDF(obs=obs, param1=..., param2,...) # data normal_np = np.random.normal(loc=[2., 2.5, 2.5], scale=[3., 3, 1.5], size=(10000, 3)) data_raw = zfit.Data.from_numpy(obs=obs, array=normal_np) # or from anywhere else, e.g. root df = data_raw.to_pandas() # preprocessing here, rename things. Match column names with the observable names "xobs", "yobs", "z" (they have to be # contained, more columns in the df is not a problem) data = zfit.Data.from_pandas(df, obs=obs) # create NLL nll = zfit.loss.UnbinnedNLL(model=product_gauss, data=data) # create a minimizer minimizer = zfit.minimize.Minuit() result = minimizer.minimize(nll) print(result.params) # do the error calculations, here with minos param_errors, _ = result.errors() print(param_errors)
# Copyright (c) 2020 zfit import numpy as np import zfit import matplotlib.pyplot as plt # create space obs = zfit.Space("x", limits=(-10, 10)) # parameters mu = zfit.Parameter("mu", 1., -4, 6) sigma = zfit.Parameter("sigma", 1., 0.1, 10) lambd = zfit.Parameter("lambda", -0.06, -1, -0.01) frac = zfit.Parameter("fraction", 0.3, 0, 1) # model building, pdf creation gauss = zfit.pdf.Gauss(mu=mu, sigma=sigma, obs=obs) exponential = zfit.pdf.Exponential(lambd, obs=obs) model = zfit.pdf.SumPDF([gauss, exponential], fracs=frac) # data n_sample = 10000 exp_data = exponential.sample(n=n_sample * (1 - frac)).numpy() gauss_data = gauss.sample(n=n_sample * frac).numpy() data = model.create_sampler(n_sample, limits=obs) data.resample() mu.set_value(0.5) sigma.set_value(1.2) lambd.set_value(-0.05) frac.set_value(0.07) # plot the data data_np = data[:, 0].numpy() color = 'black' n_bins = 50 linewidth = 2.5 plot_scaling = n_sample / n_bins * obs.area() x = np.linspace(-10, 10, 1000) # plot the pdf BEFORE fitting plt.figure() plt.title("Before fitting") # plot the data plt.hist(data_np, color=color, bins=n_bins, histtype="stepfilled", alpha=0.1) plt.hist(data_np, color=color, bins=n_bins, histtype="step") # plot the pdfs y = model.pdf(x).numpy() y_gauss = (gauss.pdf(x) * frac).numpy() # notice the frac! y_exp = (exponential.pdf(x) * (1 - frac)).numpy() # notice the frac! plt.plot(x, y * plot_scaling, label="Sum - Model", linewidth=linewidth * 2) plt.plot(x, y_gauss * plot_scaling, '--', label="Gauss - Signal", linewidth=linewidth) plt.plot(x, y_exp * plot_scaling, '--', label="Exponential - Background", linewidth=linewidth) plt.xlabel("Physical observable") plt.legend() # create NLL nll = zfit.loss.UnbinnedNLL(model=model, data=data) # create a minimizer minimizer = zfit.minimize.Minuit() result = minimizer.minimize(nll) # do the error calculations, here with minos param_errors, _ = result.errors() plt.figure() plt.title("After fitting") # plot the data plt.hist(data_np, color=color, bins=n_bins, histtype="stepfilled", alpha=0.1) plt.hist(data_np, color=color, bins=n_bins, histtype="step") y = model.pdf(x).numpy() # rerun now after the fitting y_gauss = (gauss.pdf(x) * frac).numpy() y_exp = (exponential.pdf(x) * (1 - frac)).numpy() plt.plot(x, y * plot_scaling, label="Sum - Model", linewidth=linewidth * 2) plt.plot(x, y_gauss * plot_scaling, '--', label="Gauss - Signal", linewidth=linewidth) plt.plot(x, y_exp * plot_scaling, '--', label="Exponential - Background", linewidth=linewidth) plt.xlabel("Physical observable") plt.legend() plt.show()
# Copyright (c) 2020 zfit import numpy as np import zfit # create space obs = zfit.Space("x", limits=(-2, 3)) # parameters mu = zfit.Parameter("mu", 1.2, -4, 6) sigma = zfit.Parameter("sigma", 1.3, 0.1, 10) # model building, pdf creation gauss = zfit.pdf.Gauss(mu=mu, sigma=sigma, obs=obs) # data normal_np = np.random.normal(loc=2., scale=3., size=10000) data = zfit.Data.from_numpy(obs=obs, array=normal_np) # create NLL nll = zfit.loss.UnbinnedNLL(model=gauss, data=data) # create a minimizer minimizer = zfit.minimize.Minuit() result = minimizer.minimize(nll) # do the error calculations, here with minos param_errors, _ = result.errors()
# Copyright (c) 2020 zfit import zfit # create space obs = zfit.Space("x", limits=(-10, 10)) # parameters mu = zfit.Parameter("mu", 1., -4, 6) sigma = zfit.Parameter("sigma", 1., 0.1, 10) lambd = zfit.Parameter("lambda", -0.06, -1, -0.01) frac = zfit.Parameter("fraction", 0.3, 0, 1) # model building, pdf creation gauss = zfit.pdf.Gauss(mu=mu, sigma=sigma, obs=obs) exponential = zfit.pdf.Exponential(lambd, obs=obs) model = zfit.pdf.SumPDF([gauss, exponential], fracs=frac) # data n_sample = 10000 data = model.create_sampler(n_sample, limits=obs) data.resample() # set the values to a start value for the fit mu.set_value(0.5) sigma.set_value(1.2) lambd.set_value(-0.05) frac.set_value(0.07) # create NLL nll = zfit.loss.UnbinnedNLL(model=model, data=data) # create a minimizer minimizer = zfit.minimize.Minuit() result = minimizer.minimize(nll) print(result) # do the error calculations, here with minos param_hesse = result.hesse() param_errors, new_result = result.errors() print(result.params)
# Copyright (c) 2020 zfit import numpy as np import zfit import matplotlib.pyplot as plt n_bins = 50 # create space obs = zfit.Space("x", limits=(-10, 10)) # parameters mu = zfit.Parameter("mu", 1., -4, 6) sigma = zfit.Parameter("sigma", 1., 0.1, 10) lambd = zfit.Parameter("lambda", -0.06, -1, -0.01) # model building, pdf creation gauss = zfit.pdf.Gauss(mu=mu, sigma=sigma, obs=obs) exponential = zfit.pdf.Exponential(lambd, obs=obs) n_bkg = zfit.Parameter('n_bkg', 20000) n_sig = zfit.Parameter('n_sig', 1000) gauss_extended = gauss.create_extended(n_sig) exp_extended = exponential.create_extended(n_bkg) model = zfit.pdf.SumPDF([gauss_extended, exp_extended]) # data # n_sample = 10000 data = model.create_sampler(n=21200) data.resample() # set the values to a start value for the fit mu.set_value(0.5) sigma.set_value(1.2) lambd.set_value(-0.05) # create NLL nll = zfit.loss.ExtendedUnbinnedNLL(model=model, data=data) # create a minimizer minimizer = zfit.minimize.Minuit(verbosity=7) result = minimizer.minimize(nll) print(result.params) # do the error calculations, here with hesse, than with minos param_hesse = result.hesse() param_errors, _ = result.errors() # this returns a new FitResult if a new minimum was found print(result.valid) # check if the result is still valid
# Copyright (c) 2019 zfit # create space import zfit import numpy as np obs = zfit.Space("x", limits=(-10, 10)) # parameters mu_shared = zfit.Parameter("mu_shared", 1., -4, 6) # mu is a shared parameter sigma1 = zfit.Parameter("sigma_one", 1., 0.1, 10) sigma2 = zfit.Parameter("sigma_two", 1., 0.1, 10) # model building, pdf creation gauss1 = zfit.pdf.Gauss(mu=mu_shared, sigma=sigma1, obs=obs) gauss2 = zfit.pdf.Gauss(mu=mu_shared, sigma=sigma2, obs=obs) # data normal_np = np.random.normal(loc=2., scale=3., size=10000) data1 = zfit.Data.from_numpy(obs=obs, array=normal_np) # data normal_np = np.random.normal(loc=2., scale=4., size=10000) data2 = zfit.Data.from_numpy(obs=obs, array=normal_np) # create simultaenous loss, two possibilities nll_simultaneous = zfit.loss.UnbinnedNLL(model=[gauss1, gauss2], data=[data1, data2]) # OR, equivalently nll1 = zfit.loss.UnbinnedNLL(model=gauss1, data=data1) nll2 = zfit.loss.UnbinnedNLL(model=gauss2, data=data2) nll_simultaneous2 = nll1 + nll2 minimizer = zfit.minimize.Minuit() minimizer.minimize(nll_simultaneous2)
# Copyright (c) 2020 zfit import zfit # Addition of limit with the same observables simple_limit1 = zfit.Space(obs='obs1', limits=(-5, 1)) simple_limit2 = zfit.Space(obs='obs1', limits=(3, 7.5)) added_limits = simple_limit1 + simple_limit2 # OR equivalently added_limits = simple_limit1.add(simple_limit2) # multiplication of limits with different observables first_limit_lower = (-5, 6) first_limit_upper = (5, 10) second_limit_lower = (7, 12) second_limit_upper = (9, 15) space1 = zfit.Space(obs=['obs1', 'obs2'], limits=(first_limit_lower, first_limit_upper)) space2 = zfit.Space(obs=['obs3', 'obs4'], limits=(second_limit_lower, second_limit_upper)) space4 = space1 * space2 assert space4.obs == ('obs1', 'obs2', 'obs3', 'obs4') assert space4.n_obs == 4