nnpdf / nnusf

Predictions for all-energy neutrino structure functions

NNPDF: An open-source machine learning framework for global analyses of parton distributions

The NNPDF collaboration determines the structure of the proton using Machine Learning methods. This is the main repository of the fitting and analysis frameworks. In particular it contains all the necessary tools to reproduce the NNPDF4.0 PDF determinations.

Documentation

The documentation is available at https://docs.nnpdf.science/

Install

See the NNPDF installation guide for the conda package, and how to build from source.

Please note that the conda based workflow described in the documentation is the only supported one. While it may be possible to set up the code in different ways, we won't be able to provide any assistance.

We follow a rolling development model where the tip of the master branch is expected to be stable, tested and correct. For more information see our releases and compatibility policy.

Cite

This code is described in the following paper:

@article{NNPDF:2021uiq,
    author = "Ball, Richard D. and others",
    collaboration = "NNPDF",
    title = "{An open-source machine learning framework for global analyses of parton distributions}",
    eprint = "2109.02671",
    archivePrefix = "arXiv",
    primaryClass = "hep-ph",
    reportNumber = "Edinburgh 2021/13, Nikhef-2021-020, TIF-UNIMI-2021-12",
    doi = "10.1140/epjc/s10052-021-09747-9",
    journal = "Eur. Phys. J. C",
    volume = "81",
    number = "10",
    pages = "958",
    year = "2021"
}

If you use the code to produce new results in a scientific publication, please follow the Citation Policy, particularly in regards to the papers relevant for QCD NNLO and EW NLO calculations incorporated in the NNPDF dataset.

Contribute

We welcome bug reports or feature requests sent to the issue tracker. You may use the issue tracker for help and questions as well.

If you would like contribute to the code, please follow the Contribution Guidelines.

nnusf's Issues

PyPI package is broken

The relative imports from pip install nnusf is broken which means that the contents of commondata and coefficients could not be accessed. The proposed solution is the following:

Tasks

Beta Give feedback

1. Remove the commondata and coefficients from the Wheel.
2. Store these data in some server where they are publicly available.
3. Add a script in nnusf/scripts to download the data and use appdirs to both (a) store them in the platform specific directories and (b) call them in the modules.

In this case, the imports will be exactly the same for the development and pip-installation mode. Optionally, one can add versioning to the data and dump the version in the output metadata to make sure that everything is consistent.

Improve predictions in the small-x extrapolation region

With 42cb7e3, the NN function (especially for $F_2$) saturates in the small $x$ extrapolation region.

There are several possibilities to address this but the easiest check would be to use a linear activation throughout the network and use $\ln (x)$ in the input instead of $x$.

Address potential overfitting wrt matching data

It seems that the models are over-fitting on the Yadism matching data, see the following report for example. We might do something about this.

Collect various analyses plots

Given that we are more or less ready to run a full fit with our machine learning framework, it would be good to collect here various plots: data vs Yadism (in order to make sure that the coefficients are correct), kinematics, covariance matrices, etc.

Getting minimization to work

I was just thinking that one other thing we could try is to take one step/epoch per experiment instead of calculating the chi2 for all experiments and then taking the step. I don't really have high hopes, but perhaps this added "randomization" can prevent getting "stuck".

Thought I'd open this issue so I don't forget :). Other ideas are of course also welcome.

Add module that corrects isoscalarity for a given PDF set

As agreed for $A=1$ we need to use a free isoscalar nucleon with nNNPDF30_nlo_as_0118_p that contains constraints from the LHCb data. To do so, we need to create a functionality that enforces isoscalarity on a given PDF set.

Such a constraint is enforced by replacing each PDF flavour with a linear combination of the proton- and neutron-bound PDFs as follows:
$$f^{\rm iso}_i (x, Q^2) = \left( \frac{ f^{p}_i (x, Q^2) + f^{n}_i (x, Q^2) }{2} \right)$$

Assuming isospin symmetry the proton-bound $f^p_i$ and neutron-bound $f^n_i$ PDFs are related:
$$u^n (x, Q^2) = d^p (x, Q^2) , d^n (x, Q^2) = u^p (x, Q^2)~\text{(same relations hold for anti-u and anti-d)}$$
and $f^n_i (x, Q^2) = f^p_i (x, Q^2)$ for other flavours.

Quote experimental chi2 of real data in report summary

Currently, only the total $\chi^2$ values are quoted in the report summary. Since we are mainly interested in the fit onto the real data, it will be useful to add $\chi_{\rm exp}^{2, \rm real}$ to the summary.

Computation of GLS sum rules broke due to changes in EKO

The computation of the GLS sum rule is now broken due to type change in the coupling class of eko.coupling. In particular, at some point the thresholds_ratios was of type MatchingScales now it is a list. Since eko is not a direct dependency its version depends on yadism (which now has to be v0.12.3).

Add A=1 boundary conditions

A constraint that we should impose at the level of the fit is the boundary condition $F_i(x, Q^2, A=1)=F_i^p(x, Q^2)$ which can be implemented as follows (adding this here in order to not forget once the matching is done):

$$ \chi^2 \longleftarrow \chi^2 + \lambda_{\mathrm{BC}} \sum_{i} \sum_{j}^{n_x} \sum_{k}^{n_Q} \left( F_i(x_j, Q_k^2, A=1) - F_i^p(x_j, Q_k^2) \right)^2 $$

Fix coefficients for some cross sections

For some cross sections (CHORUS and NUTEV), the coefficients seem to be wrong. This would/might explain why we couldn't fit these datasets. The plots comparing the datasets (for all values of $x$) for the these two experiments with the Yadism predictions are given below:

CHORUS

• predictions-CHORUS_DXDYNUB.zip

NUTEV

• predictions-NUTEV_DXDYNUU.zip
• predictions-NUTEV_DXDYNUB.zip

For the record, the coefficients for the CDHSW experiment are fine:

• predictions-CDHSW_DXDYNUB.zip
• predictions-CDHSW_DXDYNUU.zip

Tests to be carried out for the paper once final fit is ready

In addition to the tests that are already listed/produced in the paper, the followings are the tests that need to be carried out once the final baselines (w/ and w/o Yadism pseudodata) are available:

• Study the $A$ dependence by (a) removing a given $A$ value in the fit and compare the results with the baselines, and (b) generating predictions for a given $A$ value not included in the fit at high $Q^2$ and compare with the pQCD predictions
• study the dependence on $Q^2$ in the matching region
• replace the double differential cross sections by the structure functions for the datasets for which these observables are available

Potential issue with (quoted) Validation chi2 for CDHSW_F2, CHORUS_F2, NUTEV_F2

There seems to be some issues with the (quoted) validation $\chi^2$ values for the CDHSW_F2, CHORUS_F2, and NUTEV_F2 datasets. As an example, a fit including only the NUTEV_F2 (with one single replica) yield the following results:

Dataset	Epoch	REP ID	$\rm{N}_{\rm tr}$	$\chi^2_{\rm tr}$	$\rm{N}_{\rm vl}$	$\chi^2_{\rm vl}$	$\rm{N}_{\rm tot}$	$\chi^2_{\rm exp}$
NUTEV_F2	2316	35	58	1.847	20	922.098	78	8.721

The reason why there are $\rm{N}_{\rm tot} = 78$ is because no maximum $Q^{2}$ cut is imposed. The 20 data-points forming the validation set are shown in the table below (with $Q^2_M$ the mapped values). The data vs prediction comparisons are also shown below for all the slices in $x$.

$x$	$Q^2_M~([0, 1])$
0.015	0
0.015	0.051
0.015	0.206
0.045	0.051
0.125	0.323
0.125	0.466
0.175	0.115
0.175	0.323
0.175	0.466
0.175	0.609
0.275	0.206
0.275	0.856
0.35	0.466
0.45	0.609
0.45	0.739
0.45	0.856
0.55	0.609
0.55	0.739
0.65	0.466
0.65	0.739

As shown above, there shouldn't be a reason for the validation $\chi^2$ to be large this large given that no prediction is (significantly) far way from the true validation points.

Make the Code public

Followings are a few things that should be addressed in the view of making the code public:

Tasks

Beta Give feedback

1. Merge #77
2. Finally fix Workflow #46
3. Merge #76
4. Clean the repo by removing folders that do not belong to the package. These are: pdfs_lo, validation, xsec-grids, yad_bench. These folders should instead be moved to the paper directory.
5. Rename repository and package to NNSFnu in order to be consistent with the paper ⚠️
6. Update documentation. This includes: add a couple of representative plot(s), create a page that links to the LHAPDF grids and NNSFnu model.
7. Deploy package to PyPI (write the corresponding github action).
8. Merge #80
9. Merge #82

Implement pdf+th error for yadism grids

The cli option exists, but so far only calculates the pdf error

Fix github workflows

It will be nice to have the tests up and running again, dependent on NNPDF/workflows#6. The problem right now seems that the container registry cannot be pulled.

Add theory covmat to Yadism data

Since the grids with theory uncertainties have now been generated and parset (#47) we should now add a function that construct them and add them to the already present matching covmats.

Add the possibility to include a Docker container

It will be a nice feature to have the possibility to include a container in the on.workflow_call.inputs in case some actions require specific programs. A status check (as illustrated here) will be required in case no container is passed.

Missing datasets

Dataset that needs to be still implemented:

• CDHSW tables 30-40 F_W, A
• CHARM table 2, Qbar
• CDHS all tables

Re-compute Yadism theory with FFNS5 @ NLO

In order to account for top-quark mass effects, we need to re-compute the Yadism pseudodata using the FFNS5. For this we need two classes of grids:

• Grids used to compute the matching pseudodata at high $Q^2$ which typically probe the $Q^2$ range between $30~\mathrm{GeV}^2$ and $500~\mathrm{GeV}^2$.
• Dense grids used for switching the Yadism predictions when computing the total inclusive cross sections.

Matching for high-Q2 and A=1 does not properly work

So far, fits with #29 and #35 (which attempt to perform the matching) lead to extremely bad $\chi^2$.

Some improvements we may want to make to the code

• Make the name of the fit output folder the same as the name of the runcard instead of "fits"
• Make code more user-friendly by storing info in dataclass objects instead of dicts
• Make a default folder to store the fits, so the path doesn't need to be passed each time plots are made
• Make a single interface for all plots
• Add the option to either use patience_epochs or patience_fraction to the runcard. Current stopping_patience is in epochs, which is different from n3fit where it's a fraction, this might lead to confusion
• add a generate_pseudodata boolean flag option to the fit runcard
• Make names in callbacks.py more sensible
• do feature scaling outside model, this will also affect how models are saves/loaded
• Stop rich from messing with ipdb
• related to above, just print chi2 info line-by-line so the logs contain a history
• make gradient descent stochastic to reduce computational costs

General plan & roadmap

Below are a few things that we should give some thoughts beforehand. These concern some of the technical difficulties I foresee mainly in terms of implementations. There exist various solutions to approach all of them, we just need to find the most optimal ones. Hence, this thread will not only serve as a place to collect ideas, but most importantly to converge on final decisions concerning the various aspects of the implementation (design choices, etc).

Fitting parametrization, masks, observables

Ideally, we would like to parametrize the three structure functions ($F_2$, $F_3$, and $F_L$) as the outputs of the NN. In this way, the object that we could compute at the level of the $\chi^2$ is the double-differential cross section $d^2 \sigma / (dxdQ^2)$ (or equivalently $d^2 \sigma / (dxdy)$). However, for a given $(x,Q^2)$-pair, most of the datasets provide measurements for either $d^2 \sigma / (dxdy)$ or both $F_2$ and $F_3$ at the same time. This means that for some of the measurements, one needs to leave out $F_L$ in the computations. Here, we could do the same as in n3fit, meaning that for each measurement one specifies which structure functions are active. Better suggestion on what should be compared at the $\chi^2$-level and how to implement it is more than welcome.

Covariance matrix, correlations

Given that we will implement our own module(s) for the reading & parsing of the experimental datasets, we will also have to compute the covariance matrix by ourselves as we will no longer be able to rely on validphys.

Custom/Early stopping

One of the technical parts we need to think about in terms of the fitting code is an easy implementation of the stopping. n3fit has a module that could perform a custom callbacks on a PDF-like fitting model. However, this module is very complicated and contains a huge amount of things we do not want. Thus, we need something simple enough that is able to check all the various features we would like on top of tracking the history of the $\chi^2$.

Theoretical constraints: sum rules, etc.

An important part of the fitting is the constraints that one can impose on the structure functions both along the momentum fraction $x$ and the energy scale $Q^2$. While no literature (as far as I have looked) provides constraints on the $Q^2$ (apart from the obvious fact that $F_i (1,Q^2)=0 \text{ } \forall \text{ } Q^2$), below are a couple of sum rules that we could potentially use:

Gross Llewellyn-Smith:

$$ \int_{0}^{1} dx F_3^{\nu N} \sim \int_{0}^{1} \left( u_v(x) + d_v(x) \right) = 3 $$

Adler:

$$ \int_{0}^{1} \frac{dx}{x} \left( F_{2}^{\bar{\nu} P} - F_{2}^{\nu P} \right) \sim 2 \int_{0}^{1} dx \left( u_v(x) - d_v(x) \right) =2$$

From momentum sum rules:

$$ x \Sigma(x)= F_{2}^{\nu N} \Longrightarrow \int_{0}^{1} dx F_{2}^{\nu N} \sim \frac{1}{2} $$

Generate missing tables

The tables below are missing:

• CHARM: table 2 (Obs: Qbar)
• CDHS: all tables
• CDHSW: tables 30-40 (Obs: A, F_W)

Add option for log frequency of chi2 history to runcard