#  HEPLike: an open source framework for experimental likelihood evaluation

Jihyun Bhom  Marcin Chrzaszcz  Henryk Niewodniczanski Institute of Nuclear
Physics Polish Academy of Sciences, Krak´ow, Poland

###### Abstract

We present a computer framework to store and evaluate likelihoods coming from
High Energy Physics experiments. Due to its flexibility it can be interfaced
with existing fitting codes and allows to uniform the interpretation of the
experimental results among users. The code is provided with large open
database, which contains the experimental measurements. The code is of use for
users who perform phenomenological studies, global fits or experimental
averages.

###### keywords:

experimental high energy physics , likelihoods

††journal: Computer Physics Communications

PROGRAM SUMMARY/NEW VERSION PROGRAM SUMMARY

Program Title: HEPLike
Licensing provisions(please choose one): GPLv3
Programming language: C++

Supplementary material:
Journal reference of previous version: FIRST VERSION OF PROGRAM

Nature of problem(approx. 50-250 words): Provide a uniform way of store, share
and evaluate experimental likelihoods in a proper statistical manner. The code
can be easily interfaced with existing global fitting codes. In addition a
large database with the measurements is published. The program targets users
who perform in their scientific work: phenomenological studies, global fits or
measurements averages. The HEPLike has been created for FlavBit project[1],
which was used to perform several analysis[2,3] and here we present an updated
version, which can be used in standalone mode.
Solution method(approx. 50-250 words): C++ code that evaluates the statistical
properties of the measurements without user intervention. The large open
database is provided as well. The measurements are stored in YAML files
allowing for easy readability and extensions.

## References

  * [1] arXiv: 1705.07933
  * [2] arXiv: 1705.07935
  * [3] arXiv: 1705.07917

##  1 Introduction

In the High Energy Physics (HEP) the experimental measurements are performed
by several collaborations, which measure various different observables. The
experimental results are presented in various ways; some being as simple as a
measurement with an Gaussian error, some more complicated as multiple
correlated measurements with asymmetric errors or in some places even a full
likelihood function is being published. To make things more complicated in
some cases multiple representations of the same measurement are being
published. All of this makes it hard to directly use and compare various
different results. It also leaves a room for misinterpreting the results by
theorists, which use these inputs to their studies. It happens that the
asymmetric errors are being symmetrized, instead of using the full likelihood
only central value with approximated asymmetric error is being used.

The High Energy Physics Likelihoods (HEPLike) is a computer program that
allows to store and share the likelihoods of various measured quantities. The
published code can be useful for users performing phenomenological studies
using experimental results, global fitting collaborations or experimental
averages. Thanks to its structure it is easy to be interface with existing
codes. It simplifies the work of people as instead of looking up the
appropriate measurement and coding up their own likelihood they can download
the database of measurements and choose the one they need. Furthermore, it
shifts the burden of constructing the proper likelihood functions back to the
experimentalists, which performed the measurement at the first place and are
clearly the most appropriate people to handle this task.

The computer code described in this paper is written in C++, making it useful
for majority of fitting programs available on the market [1, 2, 3, 4, 5, 6].
The library can be used in both the $\chi^{2}$ and likelihood fits. Moreover,
it contains a statistical module with useful functions that can be used in the
statistical analysis. Besides the computer code a database with the
likelihoods is being published. The measurements are stored in the YAML files
making them easy to read by both the machine and human. This database can be
easily extended by adding new YAML files if new measurement becomes available.
With the software we provide useful utilities, which allows to perform
searches inside the database, create BiBtex containing publications, which
have been in the fit, etc.

The paper is organized as follows: in Sec. 2 construction of the likelihood
functions is presented. Sec. 3 explains the detailed code implementations and
data storage, while Sec. 4 describes how to install and use HEPLike software.

##  2 Likelihood constructions

In this section we will present how likelihoods in HEPLike are stored and
constructed. Each measurement is stored in separate YAML file. There are
several ways collaborations published their results depending on the
measurements itself:

  * 1.

Upper limits,

  * 2.

Single measurement with symmetric uncertainty,

  * 3.

Single measurement with asymmetric uncertainty,

  * 4.

Multiple measurements with symmetric uncertainty,

  * 5.

Multiple measurements with asymmetric uncertainty,

  * 6.

One dimensional likelihood function,

  * 7.

n-dimensional likelihood function.

In addition, there is growing interest from the community that the
experimental collaborations instead of only the results of the analysis
publish also the dataset that has been used to obtain the result. For this
future occasion we have also implement a way that this data can be directly
used in the fits.

Each of these cases has a different module of HEPLike that is designed to
evaluate the likelihood functions. In this section we will present the
statistical treatment of the above cases and the modules that are responsible
for their evaluation. Each of the YAML files is required to have the following
information (here for example we use the
$R_{\mathup{{{K}}^{\scriptstyle{\ast}}}}$ measurement [7]):

⬇

1BibCite: Aaij:2017vbb

2BibEntry: ’@article{Aaij:2017vbb,

3 author = ”Aaij, R. and others”,

4 title = ”{Test of lepton universality

5 with $B^{0} \rightarrow

6 K^{*0}\ell^{+}\ell^{-}$ decays}”,

7 collaboration = ”LHCb”,

8 journal = ”JHEP”,

9 volume = ”08”,

10 year = ”2017”,

11 pages = ”055”,

12 doi = ”10.1007/JHEP08(2017)055”,

13 eprint = ”1705.05802”,

14 archivePrefix = ”arXiv”,

15 primaryClass = ”hep-ex”,

16 reportNumber = ”LHCB-PAPER-2017-013,

17 CERN-EP-2017-100”,

18 SLACcitation = ”%%CITATION = ARXIV:1705.05802;%%”

19 }

20 ’

21DOI: 10.1007/JHEP08(2017)055

22Process: R_{Kstar^{*}}

23FileName: RKstar.yaml

24Name: RKstar

25Source: HEPDATA

26SubmissionYear: 2017

27PublicationYear: 2018

28Arxiv: 1705.05802

29Collaborations: LHCb

30Kinematics: q2>1.1 && q2<6\.

31HLAuthor: Gal Anonim

32HLEmail<EMAIL_ADDRESS>
33HLType: HL_ProfLikelihood

The above informations are used to store the information relevant for the
bookkeeping. For instance the entries BibCite and BibEntry correspond to the
information that are used to generate a BiBtex citation file with the
measurements that have been used in the studies. The DOI corresponds to the
digital object identifier of the publication. The Decay defines the process
that has been studied. It can also be replaced by the Process entry. The Name
is a unique name of this measurement type. If the measurement gets updated
with more data or by other collaboration the Name entry in the new YAML file
should be the same as in the old one. Source entry corresponds to the source
of the measurement. This can be either a HEPData or the collaboration itself.
The SubmissionYear (PublicationYear) refers to the year of appearance
(publication) of the result. The Arxiv codes the Arxiv number, while the
Collaborations stores the information which experimental collaboration has
performed the measurement. Finally, the Kinematics stores additional
information about the kinematic region that has been measured. The HLAuthor
and HLEmail encode the information about the YAML file author and his email in
case user needs further information about the encoded measurement. Last but
not least the entry HLType contains the information about which HEPLike object
should be used to read the file.

Reading of this content in the YAML is implemented in the HL_Data class. All
other classes that construct the likelihood functions inherit from this class
its capabilities. Please note that if the information is missing in the YAML
file the program will omit reading this entry. The only exception is the
FileName, which is mandatory. If a user wants to be notified by the program
that some informations are missing the HL_debug_yaml variable has to be set to
true (default value is false).

###  2.1 Upper limits

In case where the measurement did not observe a significant access of signal
candidates the collaborations usually report an upper limit on the measured
quantity. Commonly $90\%$ or $95\%$ upper limits are quoted. Experiments use
various statistical approaches to compute this limits. It can be the $\rm
CL_{s}$ method [8], Feldman–Cousins [9] or some variation of Bayesian methods
[10]. Publication of only an upper limits does not provide enough information
to use the result in global fits. However, nowadays experiments besides the
aforementioned upper limits publish a full p-value scans. Examples of such
scans are shown in Fig. 1. The plots are usually available in digital format,
which allows the information to be extracted and used in computer program.

Figure 1: Example of p-value scans for the
$\mathup{{{B}}{}_{\scriptstyle{\mathup{{{s}}}}}^{\scriptstyle{0}}}\to\mathup{{{\tau}}^{\scriptstyle{-}}}\mathup{{{\tau}}^{\scriptstyle{+}}}$
[11] (left) and $\mathup{{{D}}}\to\mathup{{{e}}}\mathup{{{\mu}}}$ [8] (right).
Please note that the $\rm CL_{s}$ value can be interpreted as p-value as
explained in [12]. The black line corresponds to the observed $\rm
CL_{s}$/p-value.

In HEPLike a class HL_Limit is responsible for handling this type of
measurements. It reads the YAML file that contains the standard information
about the measurement (see Sec. 2 for details). The additional information of
the observed $\rm CL_{s}$/p-value is stored in the YAML file in the following
way111Please note that the besides this information the previous information
from Sec. 2 should be included.:

⬇

1Cls:

2- [0.0, 1.0]

3- [1.0e-10, 0.977091694706]

4- [2.0e-10, 0.954375824297]

5- [3.0e-10, 0.93200355343]

6- [4.0e-10, 0.910630700546]

7- [5.0e-10, 0.889382721809]

The Cls can be replaced in the YAML file by p-value as they correspond to the
same information. The first number in each array is the value of tested
hypothesis (for example branching fraction), while the second is the
corresponding $\rm CL_{s}$/p-value. These values are then interpreted using a
$\chi^{2}$ distribution with one degree of freedom:

$pdf(x)=\frac{1}{2^{1/2}\Gamma(1/2)}x^{1/2-1}e^{-x/2},$ (1)

which had the cumulative distribution function defined as:

$cdf(x)=\frac{1}{\Gamma(1/2)}\gamma(1/2,x/2).$ (2)

In the above equations the $\Gamma(x)$ and $\gamma(k,x)$ correspond to Gamma
and incomplete gamma functions. By revering the $cdf(x)$ one can obtain the
$\chi^{2}$ value:

$\chi^{2}=cdf^{-1}(1-p),$ (3)

where p corresponds to the p-value of a given x hypothesis. This $\chi^{2}$
can be then translated to the log-likelihood via Wilks theorem [13]:

$-\log(\mathcal{L})=\frac{1}{2}\chi^{2},$ (4)

where the $\mathcal{L}$ is the likelihood. The user can choose if he wants to
obtain the $\chi^{2}$, likelihood or a log-likelihood value of a given
hypothesis.

###  2.2 Single measurement with symmetric uncertainties

The simplest case of a published experimental result is a single value with a
symmetric uncertainty. This is for example a typical result of an PDG of HFLAV
average [14, 15]. The measurement is coded in the YAML file as:

⬇

1Observables:

2- [ ”Br_A2BCZ”, 0.1, 0.05, 0.01 ]

The first argument in the array ”Br_A2BCZ” corresponds to the observable name.
Then the first number corresponds to the measured central value. The 2nd and
the 3rd number are the statistical and systematic uncertainties. In case where
only one uncertainty is available the 3rd number should be omitted and it will
be automatically set to 0 in the software. We have decided to keep the plural
Observables to be more uniform in case where more observables are measured.

The module responsible for reading this YAML file is called HL_Gaussian, it
calculates the $\chi^{2}$ for an $x$ hypothesis:

$\chi^{2}=\frac{(x_{obs}-x)^{2}}{\sigma_{stat}^{2}+\sigma_{syst}^{2}},$ (5)

where the $x_{obs}$ correspond to the measured central value in the YAML file
and the $\sigma_{stat}$ and $\sigma_{syst}$ are the statistical and systematic
uncertainties respectively. This can be the again translated to the likelihood
and log-likelihood value using Eq. 4.

###  2.3 Single measurement with asymmetric uncertainties

A simple extension of the Gaussian uncertainty is when an asymmetric
uncertainty is reported. This type of measurements although less frequent
appear in the literature. The publication in this case reports the central
value and two uncertainties: $\sigma_{+}$ and $\sigma_{-}$, which correspond
to the right (for values larger than the measured central value) and left (for
values smaller than the measured central value) uncertainty. In HEPLike we
have created a HL_BifurGaussian class, which reads the following entry in the
YAML file:

⬇

1Observables:

2- [ ”Br_A2BCZ”, 0.1, 0.05, -0.06, 0.01, -0.02 ]

The first argument is again the name of the observable and the second one is
its central value. The third and fourth arguments correspond to the
statistical $\sigma_{+}$ and $\sigma_{-}$ uncertainties, while the fifth and
sixth to the systematical $\sigma_{+}$ and $\sigma_{-}$ uncertainties. It is
important to keep the minus sign before the left side uncertainties. The code
will indicate the error in case of missing sign. In some cases the
systematical uncertainty is reported to be symmetric. In such case the last
number can be omitted in the YAML entry.

In the literature there exist number of ways to interpret asymmetric
uncertainties [16]. We have chosen the most commonly used one which is the so-
called bifurcated Gaussian:

$\displaystyle\chi^{2}=\begin{cases}\frac{(x_{obs}-x)^{2}}{\sigma_{+}^{2}},&\text{if
}x\geq x_{obs}\\\ \frac{(x_{obs}-x)^{2}}{\sigma_{-}^{2}},&\text{if
}x<x_{obs},\\\ \end{cases}$ (6)

where the $\sigma_{\pm}^{2}$ is the sum of squared statistical and systematic
uncertainties for right/left case. Once $\chi^{2}$ is calculated it can be
translated to the log-likelihood using Eq. 4.

###  2.4 Multiple measurements with symmetric uncertainties

Nowadays the most common are simultaneous measurements of various quantities,
which are correlated between each other. For instance cross section
measurements in different kinematic bins, or measurements of angular
coefficients in heavy mesons decays. In HEPLike the class responsible for
handling these cases is called HL_nDimGaussian. It reads the following
information from the YAML file:

⬇

1Observables:

2- [ ”BR1”, 0.1, 0.02]

3- [ ”BR2”, 0.2, 0.01, 0.01]

4- [ ”BR3”, 0.4, 0.04]

5Correlation:

6- [ ”BR1”, ”BR2”, ”BR3”]

7- [ 1. , 0.2 , 0 ]

8- [ 0.2, 1., 0. ]

9- [ 0 , 0., 1. ]

The information in the “Observables” entry is exactly the same as in the
HL_Gaussian class. Please note that similarly to the previous class the
systematic uncertainty is not mandatory and in case if it is not provided the
code will treat it as 0. The next entry in the YAML file is the “Correlation”,
which encodes the correlation matrix. The first ”row” is the names of the
variables it is important to keep the same order of variables as in the
“Observables” entry. The HL_nDimGaussian evaluates the $\chi^{2}$ in the
following way:

$\displaystyle\chi^{2}=V^{T}{\rm Cov}^{-1}V,$ (7)

where V is a column matrix, which is the difference between the measured and
the tested i-th observable value. The $\rm Cov$ is a square matrix,
constructed from the correlation matrix (${\rm Corr}$): ${\rm Cov}_{i,j}={\rm
Corr}_{i,j}\sigma_{i}\sigma_{j}$.

Often a user does not want to use the full set of measured quantities but just
their subset. In this case a function Restrict(vector<string>) can be used. By
passing in a form of vector the list of observables to be used, the program
will create new smaller covariance matrix, which will be used to evaluate the
$\chi^{2}$. In a similar manner the $\chi^{2}$ can be translated to the
likelihood and log-likelihood value by Eq. 4.

###  2.5 Multiple measurements with asymmetric uncertainties

More complicated case is when multiple correlated measurements are reported
with asymmetric uncertainties. The case is similar to the one discussed in
Sec. 2.3 and same statistic comments apply in this case. The YAML file
encoding such a measurement will contain the following entries:

⬇

1Observables:

2- [ ”BR1”, 0.1, +0.02, -0.01, 0.02]

3- [ ”BR2”, 0.2, +0.01, -0.05, +0.03, -0.02]

4- [ ”BR3”, 0.3, +0.04, -0.03, 0.05]

5Correlation:

6- [ ”BR1”, ”BR2”, ”BR3”]

7- [ 1. , 0.1 , 0.2 ]

8- [ 0.1, 1., 0.1 ]

9- [ 0.2 , 0.1, 1. ]

The meaning of the “Observables” entry is the same as in the previous class
(cf. Sec. 2.3) and the “Correlation” encodes the same information as in the
HL_nDimGaussian class (cf.. Sec. 2.4). The rules about the minus sign and
symmetric systematic uncertainty are the same as in case of the
HL_BifurGaussian (cf. Sec. 2.3). The difference arises when one evaluates the
$\chi^{2}$, namely the $\rm cov$ matrix is constructed depending if
$\sigma_{+}$ and $\sigma_{-}$ uncertainty is relevant:

$\displaystyle{\rm Cov}_{i,j}=\begin{cases}{\rm
Corr}_{i,j}~{}\sigma^{i}_{+}\sigma^{j}_{+},&\text{if }x^{i}\geq
x^{i}_{obs}\text{ and }x^{j}\geq x^{j}_{obs}\\\ {\rm
Corr}_{i,j}~{}\sigma^{i}_{+}\sigma^{j}_{-},&\text{if }x^{i}\geq
x^{i}_{obs}\text{ and }x^{j}<x^{j}_{obs}\\\ {\rm
Corr}_{i,j}~{}\sigma^{i}_{-}\sigma^{j}_{+},&\text{if }x^{i}<x^{i}_{obs}\text{
and }x^{j}\geq x^{j}_{obs}\\\ {\rm
Corr}_{i,j}~{}\sigma^{i}_{-}\sigma^{j}_{-},&\text{if }x^{i}<x^{i}_{obs}\text{
and }x^{j}<x^{j}_{obs}\\\ \end{cases}$ (8)

The obtained $\rm Cov$ matrix is then used to calculate the $\chi^{2}$ using
Eq. 7. The rest follows the same procedure as described in Sec. 2.4.

###  2.6 One dimensional likelihood function

The best way a result can be published is by providing the (log-)likelihood
function. This type of results are more and more common in the literature. The
most easy is the one-dimensional likelihood scans as can be presented in form
of a figure, which examples are shown in Fig. 2.

Figure 2: Examples of published one-dimensional likelihoods in the Lepton
Universality Violation of the
$\mathup{{{B}}}\to\mathup{{{K}}^{\scriptstyle{\ast}}}\ell\ell$ [7] (left) and
$\mathup{{{B}}}\to\mathup{{{K}}}\ell\ell$ [17] (right).

The biggest advantage of publishing the results in this form is its
completeness. The (log-)likelihood curve contains all the information about
all the non-Gaussian effects and incorporates the systematic uncertainties.
The technical problem is how to publish such information. Usually plots are
published in the pdf or png formats which makes them hard to be used. Since
experiments are mostly using ROOT [18] framework the plots are saved also in
the C format, which contains the points in the form of arrays. This of course
makes the points accessible however it is not easy to automate retrieving this
data from the C file. The best solution is provided by the HEPData portal
[19]. It allows to download the data in a user preferred format. In HEPLike we
have chosen to use the ROOT format by default, in which the data points are
saved in the form of a TGraph object, which is also the way experimentalists
like to store this information. In the YAML file we specify the path of the
ROOT in the following way:

⬇

1ROOTData: data/HEPData-ins1599846-v1-Table_1.root

2TGraphPath: ”Table 1/Graph1D_y1”

The ROOTData encodes the location of the ROOT file, while the TGraphPath
encodes the location of the TGraph object in that ROOT file. In HEPLike the
class HL_ProfLikelihood is responsible for reading and encoding this
likelihood. The value of the log-likelihood can be ten translated again into
the $\chi^{2}$ with Eq. 4.

###  2.7 n-dimensional likelihood function

The natural extension of one dimensional likelihood is an n-dim likelihood,
where $n\geq 2$. Currently experimental collaborations publish only 2-dim
likelihood functions (cf. Fig. 3).

Figure 3: Examples of published two-dimensional likelihoods. The
$\mathcal{B}(\mathup{{{B}}{}_{\scriptstyle{\mathup{{{s}}}}}^{\scriptstyle{0}}}\to\mu\mu)$
vs
$\mathcal{B}(\mathup{{{B}}{}_{\scriptstyle{\mathup{{{d}}}}}^{\scriptstyle{0}}}\to\mu\mu)$
likelihood [20] (left) and
$\sigma(\mathup{{{t}}}\mathup{{\overline{{t}}}}\mathup{{{Z}}})$ vs
$\sigma(\mathup{{{t}}}\mathup{{\overline{{t}}}}\mathup{{{W}}})$ likelihood
[21] (right).

The natural way of encoding such information is a histogram: TH2D or TH3D and
we have chosen this way to store this information. The corresponding entry in
the YAML file looks as following:

⬇

1ROOTData: data/LHCb/RD/Bs2mumu_5fb/histB2mumu.root

2TH2Path: ”h_2DScan”

Similar to the one dimensional likelihood (Sec. 2.6) the ROOTData encodes the
location of the ROOT file, while the TH2Path(TH3Path) encodes the location of
the TH2D(TH3D) object. In the long run the community will have to address the
question how to publish higher dimensional likelihoods and this module
(HL_nDimLikelihood) will have to be extended for such use cases.

###  2.8 Fits to experimental data

It is possible that in the future experimental collaborations besides the
results will made the datasets public. The procedure and the form in which the
data should be published is not decided and there is an ongoing debate if the
published data should correspond to the raw detector data, the final selected
points used in the analysis or something between? Clearly publishing a raw
data is problematic, as people outside the collaboration do not have the
necessary knowledge about the calibration and efficiency correction procedures
or data taking conditions. The most useful way to publish the dataset is to
allow the experimentalists to perform all the selection, all the necessary
efficiency corrections and publish the final dataset that has been used for
analysis. This would allow the theory community to use the dataset directly in
their fits without knowing the technicalities about the experimental data
analysis. For this case in HEPLike we have implemented such a class
HL_ExpPoints.

The data are stored in the TTree structure located in the ROOT file. The YAML
file encodes this information in form:

⬇

1ROOTData: data/toy/data.root

2TTreePath: t

3Observables:

4- [ x ]

5- [ y ]

6- [ z ]

7Weight: w

where the ROOTData points to the ROOT file and the TTreePath stores the
information of the TTree location inside the ROOT file. It is assumed that the
experiments will provide all the corrections in form of event-by-event
weights. The name of the weight inside the TTree is encoded in the Weight
entry. In general the data points are elements of $\mathcal{R}^{n}$ vector
space, which coordinates are stored in the Observables entry.

The only thing that user needs to provide to the HL_ExpPoints object is a
pointer to the function to be fitted. The function should have a form: double
(*fun)(vector<double> par , vector<double> point), where the par vector
encodes the parameters that want to be fitted and the point corresponds to a
data point. The HL_ExpPoints will then evaluate the likelihood:

$\mathcal{L}(\omega)=f(\textbf{x}|\omega)^{w(\textbf{x})}$ (9)

for the whole dataset. In the above the x correspondents to the n-dimensional
point, $\omega$ denotes the parameters that want to be fitted par, and $f$
denotes the fitting function (fun). The HEPLike does not provide a minimalizer
or a scanner tool as it is not purpose of this type of software. It has to be
interfaced with proper scanner tool for example [1]. Again the user can decide
if he/she prefers to perform a $\chi^{2}$ or log-likelihood fit.

The biggest advantage of such format is the compatibility with the
experimental analysis. Experimentalist can in principle publish as well the
function that they have used to fit this data and therefore a theorists
reproduce the experimental result and start where the experimentalists
finished.

##  3 Code implementation

In this section we will discuss the implementation of the code used to create
likelihoods discussed in Sec. 2. The code is build in several classes:

  * 1.

HL_Data: base class from which other classes inherit their base functionality.

  * 2.

HL_Limit: class that handles the upper limit measurements.

  * 3.

HL_Gaussian: class that handles measurements with Gaussian uncertainty.

  * 4.

HL_BifurGaussian: class that handles measurements with asymmetric uncertainty.

  * 5.

HL_nDimGaussian: class that handles measurements with n-dimensional Gaussian
uncertainties.

  * 6.

HL_nDimBifurGaussian: class that handles measurements with n-dimensional
asymmetric uncertainties.

  * 7.

HL_ProfLikelihood: class that handles measurements with one-dimensional
likelihood function.

  * 8.

HL_nDimLikelihood: class that handles measurements with 2(3)-dimensional
likelihood function.

  * 9.

HL_ExpPoints: class that allows to perform the fits to experimental datasets.

In Tab. LABEL:tab:functions we present the functionality of these classes. In
addition we present the hierarchy of the structure of the class inheritance in
Fig. 4.

Table 1: Functions available in the HEPLike software.  Function |  Description
---|---
HL_Data() |  Constructor of the HL_Data class.
HL_Data(string) |  Constructor of the HL_Data class. The argument that is taken by constructor is the path for the YAML file encoding the measurement.
HL_Limit() |  Constructor of the HL_Limit class.
HL_Limit(string) |  Constructor of the HL_Limit class. The argument that is taken by constructor is the path for the YAML file encoding the measurement.
HL_Gaussian() |  Constructor of the HL_Gaussian class.
HL_Gaussian(string) |  Constructor of the HL_Gaussian class. The argument that is taken by constructor is the path for the YAML file encoding the measurement.
HL_BifurGaussian() |  Constructor of the HL_BifurGaussian class.
HL_BifurGaussian(string) |  Constructor of the HL_Gaussian class. The argument that is taken by constructor is the path for the YAML file encoding the measurement.
HL_nDimGaussian() |  Constructor of the HL_nDimGaussian class.
HL_nDimGaussian(string ) |  Constructor of the HL_nDimGaussian class. The argument that is taken by constructor is the path for the YAML file encoding the measurement.
HL_nDimBifurGaussian() |  Constructor of the HL_nDimBifurGaussian class.
HL_nDimBifurGaussian(string) |  Constructor of the HL_nDimBifurGaussian class. The argument that is taken by constructor is the path for the YAML file encoding the measurement.
HL_ProfLikelihood() |  Constructor of the HL_ProfLikelihood class.
HL_ProfLikelihood(string) |  Constructor of the HL_ProfLikelihood class. The argument that is taken by constructor is the path for the YAML file encoding the measurement.
HL_nDimLikelihood() |  Constructor of the HL_nDimLikelihood class.
HL_ProfLikelihood(string) |  Constructor of the HL_nDimLikelihood class. The argument that is taken by constructor is the path for the YAML file encoding the measurement.
HL_ExpPoints() |  Constructor of the HL_ExpPoints class.
HL_ExpPoints(string) |  Constructor of the HL_ExpPoints class. The argument that is taken by constructor is the path for the YAML file encoding the measurement.
read_standard() |  Function that reads the general information about the measurement from the YAML file.
set_debug_yaml(bool) |  Function that enables debugging the YAML file. By default the debugging is switched off and can be switched on by passing a true bool argument to this function. Debugging will print a message that for a given information in the YAML file is missing.
Read() |  Function reading the YAML file. The function
GetChi2(double) |  Function that returns the $\chi^{2}$ value for a given point (passed to the function as double). Function is available for all classes besides HL_Data.
GetLogLikelihood(double) |  Function that returns the log-likelihood value for a given point (passed to the function as double). Function is available for all classes besides HL_Data.
GetLikelihood(double) |  Function that returns the likelihood value for a given point (passed to the function as double). Function is available for all classes besides HL_Data.
GetCLs(double) |  Function that returns $\rm CL_{s}$ or p-value for a given point (passed to the function as double). The function is a member of the HL_Limit class.
Restrict(vector<string>) |  Function that restricts number of observables from the YAML file. Function is a member of the HL_nDimGaussian, HL_nDimBifurGaussian and HL_nDimLikelihood classes.
InitData() |  Function of HL_ExpPoints class that reads to the memory the data from the TTree object.
Profile() |  Function of HL_nDimLikelihood class that creates the profile log-likelihood projections.
SetFun() |  Function of HL_ExpPoints class, that sets the pointer to the function to be fitted.
Figure 4: Diagram of class inheritance of the HEPLike package.

##  4 Installation and usage

In this chapter we will present the requirements and installation for the
HEPLike package. The software is distributed via the github site:
https://github.com/mchrzasz/HEPLike.

In order to compile HEPLike the following packages (and the minimal version)
needs to be installed:

  * 1.

git

  * 2.

cmake, 2.8

  * 3.

yaml-cpp, 1.58.0

  * 4.

gsl, 2.1

  * 5.

Boost, 1.58.0

  * 6.

ROOT, 6.08

The compilation is done in the following way:

⬇

1cd <instalation dir>

2git clone https://github.com/mchrzasz/HEPLike.git

3cd HEPLike

4mkdir build

5cd build

6cmake ..

7make

In the above the make can be replaced with make -jN, where N is the number of
threads that user wants to be used for compilation. Please note that in case
of non standard installation of some packages one might have to provide cmake
with a proper path to the library. After successful compilation a libHEPLike.a
and libHEPLike.solibraries will be created in the build directory.

The HEPLike is provided with seven examples:

  * 1.

Br_example.cc: example program showing the usage of the HL_Gaussian class.

  * 2.

BrBifurGaussian_example.cc: example program showing the usage of the
HL_BifurGaussian class.

  * 3.

Data_Fit_example.cc: example program showing the usage of the HL_ExpPoints
class.

  * 4.

Limit_example.cc: example program showing the usage of the HL_Limit class.

  * 5.

Ndim_BifurGaussian_example.cc: example program showing the usage of the
HL_nDimBifurGaussian class.

  * 6.

Ndim_Gaussian.cc: example program showing the usage of the HL_nDimGaussian
class.

  * 7.

Ndim_Likelihood_example.cc: example program showing the usage of the
HL_nDimLikelihood class.

  * 8.

ProfLikelihood_example.cc: example program showing the usage of the
HL_ProfLikelihood class.

To compile them a proper variable has to be set during the cmake stage:

⬇

1 cd build

2 cmake -DEXECUTABLE=TRUE ..

3 make

After the compilation in the build directory will contain executables from the
following examples. The HEPLike package comes also with test procedures for
each of the classes. To perform the tests user has to perform the command:

⬇

ctest

or an equivalent:

⬇

make test

If the HEPLike was successfully installed the output will look as following:

⬇

Test project /storage/github/HEPLike/build

Start 1: HL_Test_YAML

1/7 Test #1: HL_Test_YAML ………………… Passed 0.01 sec

Start 2: HL_Limit

2/7 Test #2: HL_Limit ……………………. Passed 0.27 sec

Start 3: HL_Br_example

3/7 Test #3: HL_Br_example ……………….. Passed 0.02 sec

Start 4: HL_BrBifurGaussian_example

4/7 Test #4: HL_BrBifurGaussian_example ……. Passed 0.01 sec

Start 5: HL_Ndim_Gaussian

5/7 Test #5: HL_Ndim_Gaussian …………….. Passed 0.01 sec

Start 6: HL_ProfLikelihood_example

6/7 Test #6: HL_ProfLikelihood_example …….. Passed 0.25 sec

Start 7: HL_Ndim_BifurGaussian_example

7/7 Test #7: HL_Ndim_BifurGaussian_example …. Passed 0.01 sec

100% tests passed, 0 tests failed out of 7

Total Test time (real) = 0.57 sec

###  4.1 Available measurement

The YAML files that contain the stored measurements are located in a second
independent repository. The reason for this separation is that the YAML files
are expected to be updated more frequently then the code itself. It is
expected that users and experiments will contribute to this repository. By
implementing such model it is ensured that the repository will contain the
most up to date measurements.

The repository can be found at: https://github.com/mchrzasz/HEPLikeData. The
repository should be downloaded or cloned:

⬇

1cd <some new dir>

2git clone https://github.com/mchrzasz/HEPLikeData.git

Since the repository contains only YAML files there is no need for any
compilation. The repository contains a directory data, where all the YAML
files are kept. It should be linked by a symbolic link to the HEPLike package.
Inside the data the measurements are grouped by experiments (ex. LHCb, ATLAS,
CMS, etc.). Inside the experiment directory the measurements are grouped
according to type of measurement in the collaborations, for example: RD,
Semileptonic, Charmless, Exotica, etc. The names of the YAML files should be
named accordingly to publication report number. For example: CERN-
EP-2018-331.yaml. If a single publication produced more independent
measurements, user might code them in the independent files and give further
information at the end of the file, for example:CERN-PH-
EP-2015-314_q2_01_0.98.yaml.

Currently we are publishing the measurements that have been used by us in
other projects [22, 23, 24]. The list of YAML files with the context is
presented in Tab. LABEL:tab:yaml.

Table 2: Functions available in the HEPLike software.  |
---|---
File |  Description
CERN-EP-2017-100.yaml |  YAML file encoding the measurement of branching fraction of the $\mathup{{{B}}{}_{\scriptstyle{\mathup{{{d}}}}}^{\scriptstyle{0}}}\to\mu\mu$ and $\mathup{{{B}}{}_{\scriptstyle{\mathup{{{s}}}}}^{\scriptstyle{0}}}\to\mu\mu$ decays [20].
PH-EP-2015-314_q2_0.1_0.98.yaml PH-EP-2015-314_q2_11.0_12.5.yaml PH-EP-2015-314_q2_1.1_2.5.yaml PH-EP-2015-314_q2_15.0_19.yaml PH-EP-2015-314_q2_2.5_4.0.yaml PH-EP-2015-314_q2_4.0_6.0.yaml PH-EP-2015-314_q2_6.0_8.0.yaml |  YAML files encoding the measurements of the angular coefficients of $\mathup{{{B}}{}_{\scriptstyle{\mathup{{{d}}}}}^{\scriptstyle{0}}}\to\mathup{{{K}}^{\scriptstyle{\ast}}}\mu\mu$ decay in different $q^{2}$ regions [25].
CERN-EP-2016-141_q2_0.1_0.98.yaml CERN-EP-2016-141_q2_11.0_12.5.yaml CERN-EP-2016-141_q2_1.1_2.5.yaml CERN-EP-2016-141_q2_15.0_19.yaml CERN-EP-2016-141_q2_2.5_4.0.yaml CERN-EP-2016-141_q2_4.0_6.0.yaml CERN-EP-2016-141_q2_6.0_8.0.yaml |  YAML files encoding the measurements of the branching fraction of the $\mathup{{{B}}{}_{\scriptstyle{\mathup{{{d}}}}}^{\scriptstyle{0}}}\to\mathup{{{K}}^{\scriptstyle{\ast}}}\mu\mu$ decay in different $q^{2}$ regions [26].
CERN-EP-2016-215_q2_0.1_0.98.yaml CERN-EP-2016-215_q2_1.1_2.5.yaml CERN-EP-2016-215_q2_2.5_4.yaml CERN-EP-2016-215_q2_4_6.yaml CERN-EP-2016-215_q2_6_8.yaml |  YAML files encoding the measurements of the branching fraction of the $\mathup{{{B}}{}_{\scriptstyle{\mathup{{{d}}}}}^{\scriptstyle{0}}}\to\mathup{{{K}}}\mathup{{{\pi}}}\mu\mu$ decay in different $q^{2}$ regions [27].
CERN-PH-EP-2015-145_0.1_2.yaml CERN-PH-EP-2015-145_11_12.5.yaml CERN-PH-EP-2015-145_15_19.yaml CERN-PH-EP-2015-145_1_6.yaml CERN-PH-EP-2015-145_2_5.yaml CERN-PH-EP-2015-145_5_8.yaml |  YAML files encoding the measurements of the branching fraction of the $\mathup{{{B}}{}_{\scriptstyle{\mathup{{{s}}}}}^{\scriptstyle{0}}}\to\mathup{{{\phi}}}\mu\mu$ decay in different $q^{2}$ regions [27].
CERN-EP-2019-043.yaml |  YAML file encoding the measurement of the $R_{K}$ [28].
CERN-EP-2017-100_q2_0.045_1.1.yaml CERN-EP-2017-100_q2_1.1_6.yaml |  YAML file encoding the measurement of the $R_{\mathup{{{K}}^{\scriptstyle{\ast}}}}$ [7].
b2sgamma.yaml |  YAML file encoding the HFLAV average of the $\mathup{{{b}}}\to\mathup{{{s}}}\mathup{{{\gamma}}}$ [15].
RD_RDstar.yaml |  YAML file encoding the HFLAV average of the $R(\mathup{{{D}}})$ and $R(\mathup{{{D}}^{\scriptstyle{\ast}}})$ [15].
HFLAV_2016_157.yaml HFLAV_2016_160.yaml HFLAV_2016_161.yaml HFLAV_2016_162.yaml HFLAV_2016_164.yaml HFLAV_2016_165.yaml HFLAV_2016_166.yaml HFLAV_2016_167.yaml HFLAV_2016_168.yaml HFLAV_2016_169.yaml HFLAV_2016_170.yaml HFLAV_2016_171.yaml HFLAV_2016_176.yaml HFLAV_2016_177.yaml HFLAV_2016_178.yaml HFLAV_2016_179.yaml HFLAV_2016_180.yaml HFLAV_2016_181.yaml HFLAV_2016_182.yaml HFLAV_2016_183.yaml HFLAV_2016_211.yaml HFLAV_2016_212.yaml |  YAML files encoding the upper limits of $\tau$ Lepton Flavour Violation decays [27].

As already mentioned the measurements are constantly growing and there is
expected that the community will contribute to develop this repository. When a
new YAML file is wrote before merging it with the repository it should be
checked if it contains all the necessary information. It can be checked with
the Test_YAML.cc program. It can be used in the following way:

⬇

1cd HEPLike

2./build/Test_YAML <PATH_TO_YAML>

If an entry is missing the user will be notified by a printout. The
HEPLikeData repository contains also a template YAML (data/template.yaml) file
which can be used to create new measurements YAML files.

As already mentioned we provide useful utilities for the encoded measurements.
The first is the ability to create BiBtex file for the measurements that have
been used. The user should store the BiBtex items or YAML file names:

⬇

1Aaij:2017vbb

2b2mumu.yaml

To prepare the BiBtex file user should run the make_citations.py script
located in the utils directory:

⬇

1cd utils

2python make_citations.py list.txt

after this command a new file references.bib, will be created, which will
contain the full BiBtex entries. This can be directly used in preparing the
publication.

Another useful feature of HEPLike is the ability to search the measurement
database for relevant measurements. The script allowing for that utility is
also located in the utils. Currently the database can be searched for using
the year of publication, Arxiv number, author of the YAML file or the unique
name of the measurements. The syntax for running a search is the following:

⬇

1python lookup.py –Arxiv 1705.05802

2Found files:

3../data/examples/RKstar_lowq2.yaml

To see all available search options in the following script user can run it
with help option: python lookup.py -h.

##  5 Summary

We have presented a computer program HEPLike that enables to construct and
evaluate experimental likelihoods. The software is designed to handle the
interpretation of wide range of published results. It also allows to perform
direct fits to data once it is provided by the experimental collaborations.

The program can be easily interfaced with other computer programs and is aimed
to help users, who perform fits to experimental results in their scientific
work. It is especially useful for large fitting collaborations, which till now
had to implement the experimental measurements on their own. The measurement
themselves are stored in YAML files in separate repository. This allows for
easy extensions of the database without the need of compilation. Furthermore,
users and experimental collaborations can share their encoded measurements
with the community.

## Acknowledgments

This work is partly supported by the CERN FCC Design Study Program. The
research of M. Chrzaszcz is funded by the Polish National Agency for Academic
Exchange under the Bekker program. M. Chrzaszcz is also grateful to Foundation
for Polish Science (FNP) for its support.

We would like thank Mike Williams, Patrick Koppenburg, Pat Scott, Danny van
Dyk and Maria Moreno Llacer for invaluable comments about our manuscript.

## References

  * [1] P. Athron, et al., GAMBIT: The Global and Modular Beyond-the-Standard-Model Inference Tool, Eur. Phys. J. C77 (11) (2017) 784, [Addendum: Eur. Phys. J.C78,no.2,98(2018)].  arXiv:1705.07908, doi:10.1140/epjc/s10052-017-5513-2,10.1140/epjc/s10052-017-5321-8.
  * [2] J. C. Costa, et al., Likelihood Analysis of the Sub-GUT MSSM in Light of LHC 13-TeV Data, Eur. Phys. J. C78 (2) (2018) 158.  arXiv:1711.00458, doi:10.1140/epjc/s10052-018-5633-3.
  * [3] P. Bechtle, K. Desch, P. Wienemann, Fittino, a program for determining MSSM parameters from collider observables using an iterative method, Comput. Phys. Commun. 174 (2006) 47–70.  arXiv:hep-ph/0412012, doi:10.1016/j.cpc.2005.09.002.
  * [4] F. Mahmoudi, New constraints on supersymmetric models from b —¿ s gamma, JHEP 12 (2007) 026.  arXiv:0710.3791, doi:10.1088/1126-6708/2007/12/026.
  * [5] T. Feldmann, D. Van Dyk, K. K. Vos, Revisiting $B\to\pi\pi\ell\nu$ at large dipion masses, JHEP 10 (2018) 030.  arXiv:1807.01924, doi:10.1007/JHEP10(2018)030.
  * [6] J. Kumar, D. London, R. Watanabe, Combined Explanations of the $b\to s\mu^{+}\mu^{-}$ and $b\to c\tau^{-}{\bar{\nu}}$ Anomalies: a General Model Analysis, Phys. Rev. D99 (1) (2019) 015007.  arXiv:1806.07403, doi:10.1103/PhysRevD.99.015007.
  * [7] R. Aaij, et al., Test of lepton universality with $B^{0}\rightarrow K^{*0}\ell^{+}\ell^{-}$ decays, JHEP 08 (2017) 055.  arXiv:1705.05802, doi:10.1007/JHEP08(2017)055.
  * [8] R. Aaij, et al., Search for the lepton-flavour violating decay $D^{0}\to e^{\pm}\mu^{\mp}$, Phys. Lett. B754 (2016) 167–175.  arXiv:1512.00322, doi:10.1016/j.physletb.2016.01.029.
  * [9] G. J. Feldman, R. D. Cousins, A Unified approach to the classical statistical analysis of small signals, Phys. Rev. D57 (1998) 3873–3889.  arXiv:physics/9711021, doi:10.1103/PhysRevD.57.3873.
  * [10] C. Rover, C. Messenger, R. Prix, Bayesian versus frequentist upper limits, in: Proceedings, PHYSTAT 2011 Workshop on Statistical Issues Related to Discovery Claims in Search Experiments and Unfolding, CERN,Geneva, Switzerland 17-20 January 2011, CERN, CERN, Geneva, 2011, pp. 158–163.  arXiv:1103.2987, doi:10.5170/CERN-2011-006.158.
  * [11] R. Aaij, et al., Search for the decays $B_{s}^{0}\to\tau^{+}\tau^{-}$ and $B^{0}\to\tau^{+}\tau^{-}$, Phys. Rev. Lett. 118 (25) (2017) 251802.  arXiv:1703.02508, doi:10.1103/PhysRevLett.118.251802.
  * [12] A. L. Read, Modified frequentist analysis of search results (the $CL_{s}$ method) (CERN-OPEN-2000-205).
URL https://cds.cern.ch/record/451614

  * [13] S. S. Wilks, The large-sample distribution of the likelihood ratio for testing composite hypotheses, Ann. Math. Statist. 9 (1) (1938) 60–62.  doi:10.1214/aoms/1177732360.
URL https://doi.org/10.1214/aoms/1177732360

  * [14] M. Tanabashi, et al., Review of particle physics, Phys. Rev. D 98 (2018) 030001.  doi:10.1103/PhysRevD.98.030001.
URL https://link.aps.org/doi/10.1103/PhysRevD.98.030001

  * [15] Y. Amhis, et al., Averages of $b$-hadron, $c$-hadron, and $\tau$-lepton properties as of summer 2016, Eur. Phys. J. C77 (12) (2017) 895.  arXiv:1612.07233, doi:10.1140/epjc/s10052-017-5058-4.
  * [16] R. Barlow, Asymmetric systematic errorsarXiv:physics/0306138.
  * [17] R. Aaij, et al., Test of lepton universality using $B^{+}\rightarrow K^{+}\ell^{+}\ell^{-}$ decays, Phys. Rev. Lett. 113 (2014) 151601.  arXiv:1406.6482, doi:10.1103/PhysRevLett.113.151601.
  * [18] I. Antcheva, et al., ROOT: A C++ framework for petabyte data storage, statistical analysis and visualization, Comput. Phys. Commun. 180 (2009) 2499–2512.  arXiv:1508.07749, doi:10.1016/j.cpc.2009.08.005.
  * [19] E. Maguire, L. Heinrich, G. Watt, HEPData: a repository for high energy physics data, J. Phys. Conf. Ser. 898 (10) (2017) 102006.  arXiv:1704.05473, doi:10.1088/1742-6596/898/10/102006.
  * [20] R. Aaij, et al., Measurement of the $B^{0}_{s}\to\mu^{+}\mu^{-}$ branching fraction and effective lifetime and search for $B^{0}\to\mu^{+}\mu^{-}$ decays, Phys. Rev. Lett. 118 (19) (2017) 191801.  arXiv:1703.05747, doi:10.1103/PhysRevLett.118.191801.
  * [21] M. Aaboud, et al., Measurement of the $t\bar{t}Z$ and $t\bar{t}W$ cross sections in proton-proton collisions at $\sqrt{s}=13$ TeV with the ATLAS detectorarXiv:1901.03584.
  * [22] F. U. Bernlochner, et al., FlavBit: A GAMBIT module for computing flavour observables and likelihoods, Eur. Phys. J. C77 (11) (2017) 786.  arXiv:1705.07933, doi:10.1140/epjc/s10052-017-5157-2.
  * [23] P. Athron, et al., Global fits of GUT-scale SUSY models with GAMBIT, Eur. Phys. J. C77 (12) (2017) 824.  arXiv:1705.07935, doi:10.1140/epjc/s10052-017-5167-0.
  * [24] P. Athron, et al., A global fit of the MSSM with GAMBIT, Eur. Phys. J. C77 (12) (2017) 879.  arXiv:1705.07917, doi:10.1140/epjc/s10052-017-5196-8.
  * [25] R. Aaij, et al., Angular analysis of the $B^{0}\to K^{*0}\mu^{+}\mu^{-}$ decay using 3 fb-1 of integrated luminosity, JHEP 02 (2016) 104.  arXiv:1512.04442, doi:10.1007/JHEP02(2016)104.
  * [26] R. Aaij, et al., Measurements of the S-wave fraction in $B^{0}\rightarrow K^{+}\pi^{-}\mu^{+}\mu^{-}$ decays and the $B^{0}\rightarrow K^{\ast}(892)^{0}\mu^{+}\mu^{-}$ differential branching fraction, JHEP 11 (2016) 047, [Erratum: JHEP04,142(2017)].  arXiv:1606.04731, doi:10.1007/JHEP11(2016)047,10.1007/JHEP04(2017)142.
  * [27] R. Aaij, et al., Differential branching fraction and angular moments analysis of the decay $B^{0}\to K^{+}\pi^{-}\mu^{+}\mu^{-}$ in the $K^{*}_{0,2}(1430)^{0}$ region, JHEP 12 (2016) 065.  arXiv:1609.04736, doi:10.1007/JHEP12(2016)065.
  * [28] R. Aaij, et al., Search for lepton-universality violation in $B^{+}\to K^{+}\ell^{+}\ell^{-}$ decaysarXiv:1903.09252.