madminer.likelihood package
Submodules
madminer.likelihood.base module
- class madminer.likelihood.base.BaseLikelihood(filename, disable_morphing=False, include_nuisance_parameters=True)[source]
Bases:
DataAnalyzer
Methods
event_loader
([start, end, batch_size, ...])Yields batches of events in the MadMiner file.
weighted_events
([theta, nu, start_event, ...])Returns all events together with the benchmark weights (if theta is None) or weights for a given theta.
xsec_gradients
(thetas[, nus, partition, ...])Returns the gradient of total cross sections with respect to parameters.
xsecs
([thetas, nus, partition, test_split, ...])Returns the total cross sections for benchmarks or parameter points.
create_expected_negative_log_likelihood
create_negative_log_likelihood
madminer.likelihood.histo module
- class madminer.likelihood.histo.HistoLikelihood(filename, disable_morphing=False, include_nuisance_parameters=True)[source]
Bases:
BaseLikelihood
Methods
Returns a function which calculates the expected negative log likelihood for a given parameter point, evaluated with test data sampled according to theta_true.
create_negative_log_likelihood
(x_observed[, ...])Returns a function which calculates the negative log likelihood for a given parameter point, evaluated with a dataset (x_observed,n_observed,x_observed_weights).
event_loader
([start, end, batch_size, ...])Yields batches of events in the MadMiner file.
weighted_events
([theta, nu, start_event, ...])Returns all events together with the benchmark weights (if theta is None) or weights for a given theta.
xsec_gradients
(thetas[, nus, partition, ...])Returns the gradient of total cross sections with respect to parameters.
xsecs
([thetas, nus, partition, test_split, ...])Returns the total cross sections for benchmarks or parameter points.
- create_expected_negative_log_likelihood(theta_true, nu_true, observables=None, score_components=None, include_xsec=True, luminosity=300000.0, n_asimov=None, mode='sampled', n_histo_toys=100000, model_file=None, hist_bins=None, thetas_binning=None, test_split=None)[source]
Returns a function which calculates the expected negative log likelihood for a given parameter point, evaluated with test data sampled according to theta_true.
- Parameters
- theta_truendarray
Specifies the physical parameters according to which the test data is sampled.
- nu_truendarray
Specifies the nuisance parameters according to which the test data is sampled.
- observableslist of str or None , optional
Kinematic variables used in the histograms. The names are the same as used for instance in DelphesReader.
- score_componentsNone or list of int, optional
Defines the score components used. Default value: None.
- include_xsecbool, optional
Whether the Poisson likelihood representing the total number of events is included in the analysis. Default value: True.
- luminosityfloat, optional
Integrated luminosity in pb^{-1} assumed in the analysis. Default value: 300000.
- n_asimovint or None, optional
Size of the Asimov sample. If None, all weighted events in the MadMiner file are used. Default value: None.
- mode{“weighted” , “sampled”} , optional
If “sampled”, for each evaluation of the likelihood function, a separate set of events are sampled and histogram is created to construct the likelihood function. If “weighted”, first a set of weighted events is sampled which is then used to create histograms. Default value: “sampled”
- n_histo_toysint or None, optional
Number of events drawn to construct the histograms used. If None and weighted_histo is True, all events in the training fraction of the MadMiner file are used. If None and weighted_histo is False, 100000 events are used. Default value: 100000.
- model_filestr or None, optional
Filename of a saved neural network estimating the score. Required if score_components is not None. Default value: None.
- hist_binsint or list of (int or ndarray) or None, optional
Defines the histogram binning. If int, gives the number of bins automatically chosen for each summary statistic. If list, each entry corresponds to one summary statistic (e.g. kinematic variable specified by hist_vars); an int entry corresponds to the number of automatically chosen bins, an ndarray specifies the bin edges along this dimension explicitly. If None, the bins are chosen according to the defaults: for one summary statistic the default is 25 bins, for 2 it’s 8 bins along each direction, for more it’s 5 per dimension. Default value: None.
- thetas_binninglist of ndarray or None
Specifies the parameter points used to determine the optimal binning. If none, theta_true will be used. Default value : None
- test_split
- Returns
- negative_log_likelihoodlikelihood
Function that evaluates the negative log likelihood for a given parameter point
- create_negative_log_likelihood(x_observed, observables=None, score_components=None, n_observed=None, x_observed_weights=None, include_xsec=True, luminosity=300000.0, mode='sampled', n_histo_toys=100000, model_file=None, hist_bins=None, thetas_binning=None, test_split=None)[source]
Returns a function which calculates the negative log likelihood for a given parameter point, evaluated with a dataset (x_observed,n_observed,x_observed_weights).
- Parameters
- x_observedlist of ndarray
Set of event observables with shape (n_events, n_observables).
- observableslist of str or None , optional
Kinematic variables used in the histograms. The names are the same as used for instance in DelphesReader.
- score_componentsNone or list of int, optional
Defines the score components used. Default value: None.
- n_observedint or None , optional
If int, number of observed events. If None, n_observed is defined by the length of x_observed. Default: None.
- x_observed_weightslist of float or None , optional
List of event weights with shape (n_events). If None, all events have equal weights. Default: None.
- include_xsecbool, optional
Whether the Poisson likelihood representing the total number of events is included in the analysis. Default value: True.
- luminosityfloat, optional
Integrated luminosity in pb^{-1} assumed in the analysis. Default value: 300000.
- mode{“weighted” , “sampled”, “histo”} , optional
If “sampled”, for each evaluation of the likelihood function, a separate set of events are sampled and histogram is created to construct the likelihood function. If “weighted”, first a set of weighted events is sampled which is then used to create histograms. Default value: “sampled”
- n_histo_toysint or None, optional
Number of events drawn to construct the histograms used. If None and weighted_histo is True, all events in the training fraction of the MadMiner file are used. If None and weighted_histo is False, 100000 events are used. Default value: 100000.
- model_filestr or None, optional
Filename of a saved neural network estimating the score. Required if score_components is not None. Default value: None.
- hist_binsint or list of (int or ndarray) or None, optional
Defines the histogram binning. If int, gives the number of bins automatically chosen for each summary statistic. If list, each entry corresponds to one summary statistic (e.g. kinematic variable specified by hist_vars); an int entry corresponds to the number of automatically chosen bins, an ndarray specifies the bin edges along this dimension explicitly. If None, the bins are chosen according to the defaults: for one summary statistic the default is 25 bins, for 2 it’s 8 bins along each direction, for more it’s 5 per dimension. Default value: None.
- thetas_binninglist of ndarray or None
Specifies the parameter points used to determine the optimal binning. This is requires if hist_bins doesn’t already fully specify the binning of the histogram. Default value : None
- test_split
- Returns
- negative_log_likelihoodlikelihood
Function that evaluates the negative log likelihood for a given parameter point
madminer.likelihood.manipulate module
- madminer.likelihood.manipulate.fix_params(negative_log_likelihood, theta, fixed_components=None)[source]
Function that reduces the dimensionality of a likelihood function by fixing some of the components.
- Parameters
- negative_log_likelihoodlikelihood
Function returned by Likelihood class (for example NeuralLikelihood.create_expected_negative_log_likelihood()`) which takes an n-dimensional input parameter.
- thetalist of float
m-dimensional vector of coordinate which will be fixed.
- fixed_componentslist of int or None, optional.
m-dimensional vector of coordinate indices provided in theta. fixed_components=[0,1] will fix the 1st and 2nd component of the parameter point. If None, uses [0, …, m-1].
- Returns
- constrained_nll_negative_log_likelihoodlikelihood
Constrained likelihood function which takes an n-m dimensional input parameter.
- madminer.likelihood.manipulate.profile_log_likelihood(negative_log_likelihood, remaining_components=None, grid_ranges=None, grid_resolutions=25, thetas_eval=None, theta_start=None, dof=None, method='TNC')[source]
Takes a likelihood function depending on N parameters, and evaluates for a set of M-dimensional parameter points (either grid or explicitly specified) while the remaining N-M parameters are profiled over.
- Parameters
- negative_log_likelihoodlikelihood
Function returned by Likelihood class (for example NeuralLikelihood.create_expected_negative_log_likelihood()`).
- remaining_componentslist of int or None , optional
List with M entries, each an int with 0 <= remaining_components[i] < N. Denotes which parameters are kept, and their new order. All other parameters or projected out (set to zero). If None, all components are kept. Default: None
- grid_rangeslist of (tuple of float) or None, optional
Specifies the boundaries of the parameter grid on which the p-values are evaluated. It should be [(min, max), (min, max), …, (min, max)], where the list goes over all parameters and min and max are float. If None, thetas_eval has to be given. Default: None.
- grid_resolutionsint or list of int, optional
Resolution of the parameter space grid on which the p-values are evaluated. If int, the resolution is the same along every dimension of the hypercube. If list of int, the individual entries specify the number of points along each parameter individually. Doesn’t have any effect if grid_ranges is None. Default value: 25.
- thetas_evalndarray or None , optional
Manually specifies the parameter point at which the likelihood and p-values are evaluated. If None, grid_ranges and resolution are used instead to construct a regular grid. Default value: None.
- theta_startndarray or None , optional
Manually specifies a parameter point which is the starting point for the minimization algorithm which find the maximum likelihood point. If None, theta_start = 0 is used. Default is None.
- dofint or None, optional
If not None, sets the number of parameters for the calculation of the p-values. If None, the overall number of parameters is used. Default value: None.
- method{“TNC”, ” L-BFGS-B”} , optional
Minimization method used. Default value: “TNC”
- Returns
- parameter_gridndarray
Parameter points at which the p-values are evaluated with shape (n_grid_points, n_parameters).
- p_valuesndarray
Observed p-values for each parameter point on the grid, with shape (n_grid_points,).
- mleint
Index of the parameter point with the best fit (largest p-value / smallest -2 log likelihood ratio).
- log_likelihood_rationdarray or None
log likelihood ratio based only on kinematics for each point of the grid, with shape (n_grid_points,).
- madminer.likelihood.manipulate.project_log_likelihood(negative_log_likelihood, remaining_components=None, grid_ranges=None, grid_resolutions=25, dof=None, thetas_eval=None)[source]
Takes a likelihood function depending on N parameters, and evaluates for a set of M-dimensional parameter points (either grid or explicitly specified) while the remaining N-M parameters are set to zero.
- Parameters
- negative_log_likelihoodlikelihood
Function returned by Likelihood class (for example NeuralLikelihood.create_expected_negative_log_likelihood()`).
- remaining_componentslist of int or None , optional
List with M entries, each an int with 0 <= remaining_components[i] < N. Denotes which parameters are kept, and their new order. All other parameters or projected out (set to zero). If None, all components are kept. Default: None
- grid_rangeslist of (tuple of float) or None, optional
Specifies the boundaries of the parameter grid on which the p-values are evaluated. It should be [(min, max), (min, max), …, (min, max)], where the list goes over all parameters and min and max are float. If None, thetas_eval has to be given. Default: None.
- grid_resolutionsint or list of int, optional
Resolution of the parameter space grid on which the p-values are evaluated. If int, the resolution is the same along every dimension of the hypercube. If list of int, the individual entries specify the number of points along each parameter individually. Doesn’t have any effect if grid_ranges is None. Default value: 25.
- dofint or None, optional
If not None, sets the number of parameters for the calculation of the p-values. If None, the overall number of parameters is used. Default value: None.
- thetas_evalndarray or None , optional
Manually specifies the parameter point at which the likelihood and p-values are evaluated. If None, grid_ranges and resolution are used instead to construct a regular grid. Default value: None.
- Returns
- parameter_gridndarray
Parameter points at which the p-values are evaluated with shape (n_grid_points, n_parameters).
- p_valuesndarray
Observed p-values for each parameter point on the grid, with shape (n_grid_points,).
- mleint
Index of the parameter point with the best fit (largest p-value / smallest -2 log likelihood ratio).
- log_likelihood_rationdarray or None
log likelihood ratio based only on kinematics for each point of the grid, with shape (n_grid_points,).
madminer.likelihood.neural module
- class madminer.likelihood.neural.NeuralLikelihood(filename, disable_morphing=False, include_nuisance_parameters=True)[source]
Bases:
BaseLikelihood
Methods
event_loader
([start, end, batch_size, ...])Yields batches of events in the MadMiner file.
weighted_events
([theta, nu, start_event, ...])Returns all events together with the benchmark weights (if theta is None) or weights for a given theta.
xsec_gradients
(thetas[, nus, partition, ...])Returns the gradient of total cross sections with respect to parameters.
xsecs
([thetas, nus, partition, test_split, ...])Returns the total cross sections for benchmarks or parameter points.
create_expected_negative_log_likelihood
create_negative_log_likelihood