ModelFitting#

ModelFitting is a module for fitting flavor models of Leptons to experimental data. By inputting the mass matrices, the module can determine the masses and mixing angles of a model. Further specifying a parameter space and the experimental data, it can fit the model and test its compatibility with the experimental data.

As a quick start guide, take a look at this Simple example.

For the more advanced features, consider this Detailed example.

The Model classes#

exception modelfitting.model.FlavorPyError[source]#

Bases: Exception

class modelfitting.model.LeptonModel(mass_matrix_e=None, mass_matrix_n=None, parameterspace=None, ordering='NO', experimental_data=NuFit v5.2 NO with SK chisquare profiles, fitted_observables='all', name='Lepton Model', comments='', fit_results=None)[source]#

Bases: object

A model of leptons with its mass matrices, parameters, ordering, and experimental data.

Parameters:
  • mass_matrix_e (function) –

    The charged lepton mass matrix M, for Phi_left M Phi_right. This should be a function with one argument, namely a dictionary of parameters and the corresponding value. Example:

    def m_e(params):
        return params['c'] * numpy.identity(3)
    M = LeptonModel(mass_matrix_e=m_e, ...)
    

    Default is a function that returns ‘numpy.identity(3)’. The program in its current state only gives the dimensionless mass ratios m_e/m_mu and m_mu/m_tau. For fitting, it is advisable to only use dimensionless parameters and simply ignore the overall mass scale, as it would only confuse the fit with an additional unnecessary parameter. The convention Phi_left M Phi_right is used, where left and right indicates left- and right-handed chiral fields, respectively. I.e. L_i^c M_ij e_j, where L refers to the left-handed lepton doublet and e is the right-handed charged lepton field, and i,j=1,2,3. If you use the other convention, i.e. left-handed fields on the right-hand-side, simply transpose your mass matrix.

  • mass_matrix_n (function) – The light neutrino mass matrix M, for Phi_left M Phi_right. Should be a function with one argument, namely a dictionary of parameters. Default is a function that returns ‘numpy.identity(3)’. It is !strongly! recommended to only use dimensionless parameters and one overall mass scale with the exact name ‘n_scale’, i.e. mass_matrix_n = params[‘n_scale’]*dimensionless_mass_matrix(params). Otherwise, the program will add a parameter ‘n_scale’ itself to the ParameterSpace and mass_matrix_n.

  • parameterspace (ParameterSpace()) – The parameterspace of the model. See documentation of ‘ParameterSpace’. Default is an empty parameter space.

  • ordering (str) – Specify whether the neutrino spectrum is normal or inverted ordered. Has to be either ‘NO’ or ‘IO’. Default is ‘NO’.

  • experimental_data (ExperimentalData()) – The experimental data used when fitting the model. For more information on the structure please confer the documentation of the ‘ExperimentalData’ class. The default is ‘NuFit52_NO’, i.e. the data of NuFit v5.2 for Normal Ordering taking into account results of SK. Please consider citing NuFit (www.nu-fit.org) when using this data.

  • fitted_observables (list, optional) – A list of the observables that is considered when calculating chi-square and, hence, when making the fit. Default is ‘[‘me/mu’, ‘mu/mt’, ‘s12^2’, ‘s13^2’, ‘s23^2’, ‘d/pi’, ‘m21^2’, ‘m3l^2’]’.

  • name (str, optional) – If you want, you can give the model a name. This name will be used as __repr__.

  • comments (str, optional) – If you want, you can write some comments here.

  • fit_results (list, optional) – It is possible to store the results of make_fits() in this list. It is of course also possible to load the results from an external calculation into this list.

copy()[source]#
get_dimless_obs(params: dict) dict[source]#
get_obs(params: dict) dict[source]#

Get a dictionary of all lepton observables for a point in parameterspace.

Parameters:

params (dict) – The point in parameter space.

Returns:

All lepton observables, i.e. {‘me/mu’:0.0048, …}

Return type:

dict

residual(params: dict) list[source]#

The residual used to make the dimensionless fit.

Parameters:

params (dict) – The point in parameterspace.

Returns:

A list of values of individual chis (not chi-squares!). Only dimensionless observables are being considered.

Return type:

list

get_chisq(params=None, observables=None) float[source]#

Returns the value of chi-square for a given point in parameter space.

Parameters:
  • params (dict) – The point in parameterspace. Alternatively, you can also insert the observables. Default is None.

  • observables (dict) – Alternatively, you can directly insert a dictionary with the observables. Note that the value of ‘params’ is ignored one you define ‘observables’! Also you need to define either ‘params’ or ‘observables’ Default is None.

Returns:

The value of chi-square.

Return type:

float

print_chisq(params: dict)[source]#

Prints the value of all observables and the associated contribution to chi-square. Also prints total chi-square.

Parameters:

params (dict) – The point in parameterspace

make_fit(points: int, **fitting_kwargs) DataFrame[source]#

Does the fit for a specific number of random points in parameterspace.

Parameters:
  • points (int) – The number of random points in parameter space you want to fit. If you want to fit a specific starting point in parameter space, adjust the ‘sampling_fct’ in your ParameterSpace.

  • fitting_kwargs – properties of the Fit class. You can add keyword arguments that will be passed down to the Fit object used to make the fit. Please see the documentation of the Fit class for the specific keyword arguments. Of course, the keywords ‘model’ and ‘params’ can not be passed down to Fit.

Returns:

The result of the fit is returned in form of a pandas.DataFrame. Note that several (default:4) minimization algorithms are applied consecutively to one random point. Since the results of the intermediate steps are also written into the resulting DataFrame, it has more rows than the number entered as ‘points’.

Return type:

pandas.DataFrame

dimless_fit(points: int, **fitting_kwargs) DataFrame[source]#

Calling this function fits the dimensionless parameters of the model. The procedure of Model.make_fit() can be split into the fitting of dimensionless parameters (with Model.dimless_fit) and the fitting of the dimensionful ones (with Model.complete_fit), where Model.complete_fit() also adds the observables to the resulting pandas.DataFrame. This function comes in handy when running fits on an external machine, e.g. a cluster, since the result of dimless_fit() uses a smaller amount of memory when stored into a file compared to the result of complete_fit(). You can transfer the files from dimless_fit() to your local machine and easily run complete_fit() since complete_fit() does not take a lot of time compared to dimless_fit().

Parameters:
  • points (int) – The number of random points in parameter space you want to fit. If you want to fit a specific starting point in parameter space, adjust the ‘sampling_fct’ in your ParameterSpace.

  • fitting_kwargs – properties of the Fit class. You can add keyword arguments that will be passed down to the Fit object used to make the fit. Please see the documentation of the Fit class for the specific keyword arguments. Of course, the keywords ‘model’ and ‘params’ can not be passed down to Fit.

Returns:

The result of the fit is returned in form of a pandas.DataFrame. Note that several (default:4) minimization algorithms are applied consecutively to one random point. Since the results of the intermediate steps are also written into the resulting DataFrame, it has more rows than the number entered as ‘points’.

Return type:

pandas.DataFrame

complete_fit(df: DataFrame) DataFrame[source]#

Use this function to add the observables (while simultaneously fitting the dimensionful ‘n_scale’ parameter) to a pandas.DataFrame containing points in parameterspace.

Parameters:

df (pandas.DataFrame) – Points in parameterspace. E.g. the result of Model.dimless_fit().

Returns:

The as ‘df’ entered points in parameterspace plus their corresponding observables and chi-square value.

Return type:

pandas.DataFrame

class modelfitting.model.QuarkModel[source]#

Bases: object

The Parameterspace class#

class modelfitting.parameterspace.ParameterDimension(name, sample_fct=None, vary=True, min=-inf, max=inf, expr=None, brute_step=None)[source]#

Bases: object

A parameter dimension, i.e. one direction in a multidimensional parameter space

Parameters:
  • name (str) – The name of the dimension.

  • sample_fct (function) – A function that when evaluated as function() returns a float.

  • vary (bool, optional) – Specifies whether the parameter is varied during a fit.

  • min (float, optional) – Lower bound. Default is -numpy.inf

  • max (float, optional) – Upper bound. Default is numpy.inf

  • expr (str, optional) – A mathematical expression used to constrain the value during the fit. Default is ‘None’

  • brute_step (float, optional) – Step size for grid points, if you use the ‘brute’ method when fitting.

class modelfitting.parameterspace.ParameterSpace(name='Parameter space')[source]#

Bases: dict

A parameter space. This object is a dictionary that contains ParameterDimension() objects.

Parameters:

name (str, optional) – You can give your parameter space a name.

add_dim(name, sample_fct=None, vary=True, min=-inf, max=inf, expr=None, brute_step=None)[source]#

Adds a dimension to your parameter space. Can also be used to update or overwrite an existing dimension.

Parameters:
  • name (str) – The name of the dimension.

  • sample_fct (function) – A function that when evaluated as function() returns a float. Default is numpy.random.uniform.

  • vary (bool, optional) – Specifies whether the parameter is varied during a fit.

  • min (float, optional) – Lower bound. Default is -numpy.inf

  • max (float, optional) – Upper bound. Default is numpy.inf

  • expr (str, optional) – A mathematical expression used to constrain the value during the fit. Default is ‘None’

  • brute_step (float, optional) – Step size for grid points, if you use the ‘brute’ method when fitting.

random_pt() Parameters[source]#

Draws a sample in your parameter space.

Returns:

A lmfit.Parameters object.

The Experimentaldata class#

class modelfitting.experimentaldata.ExperimentalData(name='ExpData', data_table=None, data=None)[source]#

Bases: object

An experimental data set.

Parameters:
  • name (str) – The name of the set

  • data_table (dict, optional) –

    A dictionary that contains the experimental best fit value and 1 sigma errors. The dict needs to be in a very specific form! It is very important that you use the keys ‘best’, ‘1sig_min’, and ‘1sig_max’! An example is:

    data_table = {'me/mu': {'best':0.0048, '1sig_min':0.0046, '1sig_max':0.0050}, 'mu/mt': ...}
    

    The data entered in ‘data_table’ will be used to calculate a chi-square value assuming a gaussian error, i.e.:

    chisq = ( (model_value - best) / (0.5*(1sig_max - 1sig_min)) )^2
    

    Alternatively, especially for non-gaussian errors, one can enter the chisquare profiles using the ‘data’ parameter. Mixed entries are also possible.

  • data (dict, optional) –

    A dictionary, whose keys are the names of the experimental observables and the corresponding values are functions, typically chisquare profiles. So for example:

    data = {"me/mu": memuchisq, "mu/mt": mumtchisq}
    

    where memuchisq(value) is defined to return ((value - exp_best)/(exp_1sig_error))^2 or something equivalent.

get_chisq(values: dict, considered_obs='auto') float[source]#

Returns chi-square for a set of values. The values are added quadratically, so chisq = sum_i chisq_i.

Parameters:
  • values (dict) – The values that you want to compare to the experimental data. ‘values’ should be in the form of a dictionary like for example ‘values = {‘me/mu’:0.0045, ‘mu/mt’:0.546, …}’, or any object that returns a callable via a key, i.e. ‘values[‘me/mu’]’ returns 0.0045, for example a pandas dataframe.

  • considered_obs (list, optional) – A list that contains the keys of all observables that should be considered to calculate the chisquare. So if ‘considered_obs=[‘me/mu’, ‘s12^2’]’, then the returned chisquare only contains contributions from ‘me/mu’ and ‘s12^2’, even though values also has an entry ‘mu/mt’. The default value is ‘auto’, which will consider all keys present in ‘values’.

Returns:

List of chi-square values associated with a specific observable, i.e. [chisq_me/mu, chisq_s12^2, …].

Return type:

list

get_chisq_list(values: dict, considered_obs='auto') list[source]#

Returns a list that contains the contributions to chi-square.

Parameters:
  • values (dict) – The values that you want to compare to the experimental data. ‘values’ should be in the form of a dictionary like for example ‘values = {‘me/mu’:0.0045, ‘mu/mt’:0.546, …}’, or any object that returns a callable via a key, i.e. ‘values[‘me/mu’]’ returns 0.0045, for example a pandas dataframe.

  • considered_obs (list, optional) – A list that contains the keys of all observables that should be considered to calculate the chisquare. So if ‘considered_obs=[‘me/mu’, ‘s12^2’]’, then the returned chisquare only contains contributions from ‘me/mu’ and ‘s12^2’, even though values also has an entry ‘mu/mt’. The default value is ‘auto’, which will consider all keys present in ‘values’.

Returns:

List of chi-square values associated with a specific observable, i.e. [chisq_me/mu, chisq_s12^2, …].

Return type:

list

The Fit class#

class modelfitting.fit.Fit(model=None, params=None, methods=None, nr_methods=4, shuffle_methods=True, max_time=45, retry_time=5, dig_deeper=False, dig_deeper_threshold=1000)[source]#

Bases: object

This Class is supposed to represent the fitting of a single random point.

Parameters:
  • model (py:meth:~modelfitting.model.LeptonModel) – The Model whose parameters you want to fit

  • params (A lmfit.Parameters object.) – The parameters you want to fit

  • methods (list, optional) – A list of all methods that should be used for fitting. Note that only ‘nr_methods’ of these are actually used. Either the first ones or random ones, depending on ‘shuffle_methods’. Default is a mixture of ‘least_square’, ‘nelder’, and further more advanced algorithms.

  • nr_methods (int, optional) – Default is 4

  • shuffle_methods (bool, optional) – Default is True

  • max_time (int, optional) – The amount of time in seconds, that the minimizer-algorithm is allowed to run. After this time is elapsed the algorithm is aborted.

  • retry_time (int, optional) – If the minimization-algorithm is aborted for any reason, as a replacement ‘least-squares’ is used for ‘retry_time’ second.

  • dig_deeper (bool, optional) – If ‘dig_deeper’ is True, then a second round of methods is applied, if the usual calculation has lead to a chi-square less than ‘dig_deeper_threshold’. Default is False.

  • dig_deeper_threshold (float, optional) –

make_fit() list[source]#

Call this function to execute the fit.

Returns:

A list that contains the results of the fit in form of lmfit.MinimizerResult objects.

Return type:

list

fit_results_into_dataframe(fit_results: list) DataFrame[source]#

Convertes the result of Fit.make_fit() into a pandas.DataFrame.

Parameters:

fit_results (list) – A list that contains elements of the lmfit.MinimizerResult.

Returns:

A pandas.DataFrame object that contains the best-fit parameters as well as the value of chi-square of the elements of fit_results.

Return type:

pandas.DataFrame

class modelfitting.fit.MinimizeStopper(max_sec=60)[source]#

Bases: object

This object is able to stop a minimization procedure after max_sec seconds is elapsed.

Some mixingcalculations utils#

modelfitting.mixingcalculations.calculate_quark_observables(mass_matrix_u=array([[1., 0., 0.], [0., 1., 0.], [0., 0., 1.]]), mass_matrix_d=array([[1., 0., 0.], [0., 1., 0.], [0., 0., 1.]]), parameterization='standard') dict[source]#

Calculates the quark observables of up and down-type mass matrices.

Parameters:
  • mass_matrix_u (3x3 matrix) – The up-type quark mass matrix M, for Phi_left M Phi_right.

  • mass_matrix_d (3x3 matrix) – The down-type quark mass matrix M, for Phi_left M Phi_right.

  • parameterization (str, default:'standard') – Specify whether you want the result in standard or wolfenstein parametrization. Has to be either ‘standard’ or ‘wolfenstein’.

Returns:

dict Contains the standard or wolfenstein parameters as well as the quark mass ratios.

modelfitting.mixingcalculations.calculate_lepton_dimensionless_observables(mass_matrix_e=array([[1., 0., 0.], [0., 1., 0.], [0., 0., 1.]]), mass_matrix_n=array([[1., 0., 0.], [0., 1., 0.], [0., 0., 1.]]), ordering='NO') dict[source]#

Calculates the dimensionless observables of a charged lepton and neutrino mass matrix.

Parameters:
  • mass_matrix_e (3x3 matrix) – The charged lepton mass matrix M, for Phi_left M Phi_right, where left and right indicates left- and right-handed chiral fields, respectively. I.e. L_i^c M_ij e_j, where L refers to the left-handed lepton doublet and e is the right-handed charged lepton field, and i,j=1,2,3. If you use the other convention, i.e. left-handed fields on the right-hand-side, simply transpose your mass matrix.

  • mass_matrix_n (3x3 matrix) – The neutrino mass matrix M, for Phi_left M Phi_right.

  • ordering (str) – Specify whether the neutrino spectrum is normal or inverted ordered. Has to be either ‘NO’ or ‘IO’. Default is ‘NO’.

Returns:

Contains the PMNS-parameters and the charged lepton mass matrices. It also contains wrongly scaled neutrino masses, that are needed for the function ‘calculate_lepton_observables’.

Return type:

dict

modelfitting.mixingcalculations.calculate_lepton_observables(mass_matrix_e=array([[1., 0., 0.], [0., 1., 0.], [0., 0., 1.]]), mass_matrix_n=array([[1., 0., 0.], [0., 1., 0.], [0., 0., 1.]]), ordering='NO', m21sq_best=None, m3lsq_best=None) dict[source]#

Calculates the observables of the lepton sector from charged lepton and neutrino mass matrices.

Parameters:
  • mass_matrix_e (3x3 matrix) – The charged lepton mass matrix M, for Phi_left M Phi_right, where left and right indicates left- and right-handed chiral fields, respectively. I.e. L_i^c M_ij e_j, where L refers to the left-handed lepton doublet and e is the right-handed charged lepton field, and i,j=1,2,3. If you use the other convention, i.e. left-handed fields on the right-hand-side, simply transpose your mass matrix.

  • mass_matrix_n (3x3 matrix) – The neutrino mass matrix M, for Phi_left M Phi_right.

  • ordering (str) – Specify whether the neutrino spectrum is normal or inverted ordered. Has to be either ‘NO’ or ‘IO’. Default is ‘NO’.

  • m21sq_best (float) – The best fit value for the squared neutrino mass difference m_1^2 - m_2^2. Default is None, which will yield 7.41e-05.

  • m3lsq_best (float) – The best fit value for the squared neutrino mass difference m_3^2 - m_l^2, where l=1 for NO and l=2 for IO. Default is None, yielding 2.507e-03 for NO and -2.486e-03 for IO.

Returns:

Contains the PMNS parameters as well as the neutrino masses and charged lepton mass ratios.

Return type:

dict

modelfitting.mixingcalculations.get_wolfenstein_parameters(CKM) dict[source]#

Get values of the Wolfenstein-parameterization of CKM matrix, according to PDG.

Parameters:

CKM (3x3 matrix) – The CKM matrix in a 3x3-matrix-shape.

Returns:

Dictionary that contains the parameters

Return type:

dict

modelfitting.mixingcalculations.get_standard_parameters_ckm(CKM) dict[source]#

Get the values of the standard-parametrization of CKM matrix.

Parameters:

CKM (3x3 matrix) – The CKM matrix in a 3x3-matrix-shape.

Returns:

Dictionary that contains the parameters

Return type:

dict

modelfitting.mixingcalculations.get_standard_parameters_pmns(PMNS) dict[source]#

Get the values of the standard-parameterization of PMNS matrix.

Parameters:

PMNS (3x3 matrix) – The PMNS matrix in a 3x3-matrix-shape.

Returns:

Dictionary that contains the parameters.

Return type:

dict

modelfitting.mixingcalculations.calculate_ckm(Mu=array([[1., 0., 0.], [0., 1., 0.], [0., 0., 1.]]), Md=array([[1., 0., 0.], [0., 1., 0.], [0., 0., 1.]])) ndarray[source]#

Calculates the CKM matrix out of the up and down-type quark mass matrices.

Parameters:
  • Mu (3x3 matrix) – The up-type quark mass matrix M, for Phi_left M Phi_right.

  • Md (3x3 matrix) – The down-type quark mass matrix M, for Phi_left M Phi_right.

Returns:

The CKM matrix

Return type:

3x3 matrix

modelfitting.mixingcalculations.calculate_pmns(Me=array([[1., 0., 0.], [0., 1., 0.], [0., 0., 1.]]), Mn=array([[1., 0., 0.], [0., 1., 0.], [0., 0., 1.]]), ordering='NO') ndarray[source]#

Calculates the PMNS matrix out of the charged lepton and light neutrino mass matrices.

Parameters:
  • Me (3x3 matrix) – The charged lepton mass matrix M, for Phi_left M Phi_right.

  • Mn (3x3 matrix) – The light neutrino mass matrix M, for Phi_left M Phi_right.

  • ordering (str) – Specify whether the neutrino spectrum is normal or inverted ordered. Has to be either ‘NO’ or ‘IO’.

Returns:

The PMNS matrix.

Return type:

3x3 matrix