Neural network predictions of inclusive electron-nucleus cross sections

O. Al Hammal ●

IPSA-DRII, 63 Boulevard de Brandebourg, 94200 Ivry-sur-Seine, France

M. Martini ●

IPSA-DRII, 63 Boulevard de Brandebourg, 94200 Ivry-sur-Seine, France
and Sorbonne Université, Université Paris Diderot, Centre National de la Recherche Scientifique IN2P3,
Laboratoire de Physique Nucléaire et de Hautes Energies, Paris, France

J. Frontera-Pons ●

IPSA-DRII, 63 Boulevard de Brandebourg, 94200 Ivry-sur-Seine, France
and Laboratoire AIM, CEA, Centre National de la Recherche Scientifique, Université Paris-Saclay,
Université de Paris, 91191 Gif-sur-Yvette, France

T. H. Nguyen

IPSA, 63 Boulevard de Brandebourg, 94200 Ivry-sur-Seine, France

R. Pérez-Ramos

IPSA-DRII, 63 Boulevard de Brandebourg, 94200 Ivry-sur-Seine, France
and Sorbonne Université, Université Paris Diderot, Laboratoire de Physique Théorique
et de Hautes Energies, Paris, France

<figure><img src="image_1.png" /></figure>

(Received 22 September 2021; revised 21 April 2023; accepted 8 May 2023; published 5 June 2023)

We investigate whether a neural network approach can reproduce and predict the electron-nucleus cross sections in the kinematical domain of present and future accelerator-based neutrino oscillation experiments. For this purpose, we consider the large amount of data available to the community via the University of Virginia's Quasielastic Electron Nucleus Scattering Archive and use a residual, fully connected feedforward neural network. We illustrate the training performances of the neural network by comparing its results with experimental data for the electron double-differential cross section on carbon. The agreement between predictions and data is remarkable from quasielastic to deep-inelastic scattering. To test the predicting power of the neural network we consider the numerous kinematical conditions for which experimental cross sections on calcium are available. Furthermore, we show the predictions of the electron-scattering cross sections on oxygen, argon, and titanium: Nuclei of particular interest in the context of present and future accelerator-based neutrino oscillation programs. The agreement between these predictions and the data is comparable to the one of other theoretical models commonly used to calculate electron and neutrino cross sections, such as SuSAv2 and GiBUU. Results obtained with GENIE, a Monte Carlo event generator, are also discussed for comparison. The good performances obtained with our neural network suggest that neural networks could be exploited for theoretical and experimental investigations of electron- and neutrino-nucleus scattering.

DOI: 10.1103/PhysRevC.107.065501

# I. INTRODUCTION

For over 100 years, electron-scattering experiments have represented one of the most powerful and fruitful approaches for investigating the structure of physical systems. The Franck-Hertz experiment [1], which showed that the energy absorption by atoms is quantized by studying the flux of electrons through a vapor of mercury atoms, can be considered as the progenitor of an important tradition of electron-scattering experiments allowing us to move the frontiers of our knowledge of the structure of matter: From atoms, to nuclei [2], to

quarks, precisely discovered thanks to experiments [3,4] of this kind.

The possibility of using other lepton beams also appeared in parallel with the development of increasingly energetic and intense electron beams. In particular, starting from Pontecorvo's [5] idea and Schwartz's [6] project, experiments with neutrino beams have been performed, the first having led to the discovery of the muon neutrino [7].

Nowadays accelerator-based neutrino experiments, such as T2K [8] and NOvA [9], are performed for a precise determination of neutrino oscillation parameters. The next generation
experiments DUNE [10] and Hyper-Kamiokande [11] will play a central role in the neutrino oscillation program.

To ensure the success of these experiments, a reduction in systematic errors to the level of a few percent is needed. Today, one of the greatest sources of systematic error is the neutrino-nucleus cross sections. In these experiments nuclear targets (such as C, O, and Ar) are involved, and in the energy region of hundreds of MeV to a few GeV these cross sections are known to a precision not exceeding 20% [12,13]. Knowledge of these cross sections is crucial to determine the neutrino energy that enters the expression of the neutrino oscillation probability. Since neutrino beams are not monochromatic, in contrast with electron beams, the initial neutrino energy is reconstructed from the final states of the neutrino-nucleus reaction.

The current program of electron-scattering experiments continues in connection with nuclear, hadronic, and particle physics and represents a powerful tool in connection with neutrino physics. Indeed, electrons and neutrinos both being leptons, their interaction with atomic nuclei is similar: It happens via vector current in the case of electrons and via vector and vector-axial currents in the case of neutrinos. Hence the vector part of the cross section and the final state interaction of the scattered particles are identical for electron and neutrino scattering. Although the axial part of the neutrino cross section complicates the comparisons, electron-scattering experiments, which have the great advantage of employing monoenergetic electron beams, are useful in the investigation of neutrino-nucleus scattering from several points of view. First of all, the large amount of cross section data already obtained by many different electron-scattering experiments can be used to validate microscopic theoretical models employed to predict neutrino-nucleus cross sections; a nonexhaustive list of such studies includes the works of Refs. [14–17]. These world data have been recently considered in Refs. [18–20] to benchmark the neutrino Monte Carlo event generators, such as GENIE [21], commonly used in the neutrino community, as well as to check the neutrino energy reconstruction via calorimetric and quasielastic kinematics-based methods [22]. The available data also include the recently published electron cross sections on argon and titanium [23–26]. These cross sections were measured in the experiment proposed in Ref. [27] to determine the spectral function of $^{40}$Ar, the nucleus employed in DUNE, and in the detectors of the Fermi National Accelerator Laboratory short-baseline neutrino program [28].

The aim of this paper is to investigate the use of neural networks to reproduce and predict the electron-nucleus cross sections in the kinematical domain of present and future accelerator-based neutrino oscillation experiments.

In the last decade, applications using neural networks have become ubiquitous and are found in many tasks beyond fundamental research. In the high energy physics community, deep learning provides a faster alternative to both standard data analysis techniques and Monte Carlo approaches to simulate the detector outputs [29–33]. Deep learning for event reconstruction in accelerator-based neutrino experiments is used for example in Refs. [34–40]. Recent studies [41,42] applied generative adversarial networks to construct an AI-based Monte Carlo event generator for deep-inelastic electron-proton

<table><caption>TABLE I. Percent contribution of the different nuclear targets electron-scattering data [50] considered in our paper.</caption><tr><td><sup>4</sup>He</td><td><sup>6</sup>Li</td><td><sup>9</sup>Be</td><td><sup>12</sup>C</td><td><sup>16</sup>O</td><td><sup>24</sup>Mg</td><td><sup>27</sup>Al</td><td><sup>40</sup>Ar</td><td><sup>40</sup>Ca</td><td><sup>48</sup>Ca</td><td><sup>48</sup>Ti</td><td><sup>56</sup>Fe</td><td><sup>59</sup>Ni</td></tr><tr><td>21.5</td><td>1.3</td><td>3.1</td><td>24.0</td><td>1.0</td><td>0.3</td><td>6.3</td><td>1.4</td><td>10.5</td><td>9.6</td><td>1.4</td><td>19.3</td><td>0.3</td></tr></table>

scattering, free of theoretical assumptions about the underlying particle dynamics. In the context of electron and neutrino scattering, theoretical studies have used neural networks to obtain information on the two main ingredients of the lepton-nucleus cross sections: nuclear responses and nucleon form factors. More precisely, in Ref. [43] a physics-informed artificial neural network is employed to reconstruct the electromagnetic response functions. In Ref. [44] a Bayesian approach for feed-forward neural networks has been applied to extract from the neutrino-deuteron scattering data the nucleon axial form factor, a quantity which has been widely debated in the last ten years, following the MiniBoNE measurement of the quasielastic-like neutrino cross section on carbon [45]. The possibility of using machine learning algorithms in reconstructing neutrino energy has been explored in Ref. [46]. For a recent review on the current trends and perspectives of artificial intelligence in nuclear physics we refer to Ref. [47].

## II. DATA

A large amount of data is needed to predict electron-nucleus cross sections via a neural network approach. These data, accumulated by many different electron-scattering experiments which began in the mid-1970s, have been assembled by the authors of Refs. [48,49] and made available in the University of Virginia’s Quasielastic Electron Nucleus Scattering Archive [50].

This archive contains about 600 different combinations of targets, energies, and angles consisting of some 20 000 data points. These data cover the energy region of interest for present and future long-baseline neutrino oscillation experiments, from giant resonance excitations up to deep-inelastic scattering, with a predominant contribution of quasielastic and $\Delta$ resonance excitations.

The nuclei included in the archive [50] vary from hydrogen to uranium. Driven by the idea of having a relatively homogeneous dataset, containing the nuclei employed in neutrino detectors, we decided to discard the very light nuclei and the heavy nuclei characterized by a large neutron-proton asymmetry. The final subset of nuclei considered varies between $^{4}$He and $^{59}$Ni. They are specified in Table I together with their relative contribution to our dataset.

In Fig. 1 we plot some examples of inclusive double-differential cross sections $d^2\sigma/d\omega d\Omega$ for different values of incoming electron energy $E$ and electron-scattering angle $\theta^1$

<footnote n="1">¹We recall that the relation between the differential solid angle $d\Omega$ in the direction specified by the scattered electron momentum and the electron-scattering angle $\theta$ is $d\Omega = 2\pi d \cos\theta$.</footnote>
<figure><img src="image_2.png" /><figcaption>FIG. 1. Inclusive (e, e') $^{12}$C double-differential cross section data as a function of the transferred energy in four different kinematical conditions written above each panel. The vertical lines correspond to the kinematics of quasielastic and Δ excitations given by Eqs. (2) and (4), respectively.</figcaption></figure>

as a function of the energy transferred to the nucleus ω (also called energy loss, ω = E − E ′ being the difference between the incoming and outgoing electron energy E ′ ). For this il- lustration we have selected $^{12}$C, a nucleus for which much data exist. Two characteristics that would complicate the neu- ral network prediction task can already be seen from these few examples which refer to different kinematical conditions. First, the cross section may span many orders of magnitude: Five orders in the four panels of Fig. 1, which refer to the same nucleus, and 12 orders of magnitude (from 4 × 10 − 6 to 10 $^{6}$) when considering the whole dataset. Secondly, the shapes of the cross section can vary greatly, reflecting the different reaction mechanisms that can be induced by different kinematics. Starting from the top left panel of Fig. 1, the first recognizable excitations are the nuclear giant resonances which correspond to the sharp peaks in the cross sections for ω < 20 MeV. Moreover, the quasielastic bump appears in the first three panels. It corresponds to one nucleon knockout and is peaked around

where M$_{N}$ is the nucleon mass and q the momentum transfer to the nucleus, given by the difference between the incoming and scattered electron momentum, q = k − k $^{′}$. In terms of the electron kinematics variable, when electron mass is neglected, Eq. (1) can be written as

$$
\omega_{QE} = \frac{E^2(1 - \cos \theta)}{M_N + E(1 - \cos \theta)}. \qquad (2)
$$

A vertical line corresponding to this value is plotted in each panel of Fig. 1. The shift of the position of the real quasielastic peak with respect to the value of Eqs. (1) and (2) is due to the absence of nucleon binding and nuclear collective effects in these formulas. The broadening of the quasielastic bump is due to nucleon Fermi motion. The second bump at a higher ω corresponds to the Δ resonance excitation. In the case of scattering with a free nucleon at rest, it would peak at

$$
\omega_{QE} = \frac{Q^2}{2M_N} = \frac{q^2 - \omega^2}{2M_N} = \sqrt{q^2 + M_N^2} - M_N, \quad (1)
$$

$$
\omega_{\Delta} = \sqrt{\mathbf{q}^2 + M_{\Delta}^2} - M_N = \frac{Q^2}{2M_N} + \Delta M, \quad (3)
$$
with $\Delta M = (M_{\Delta}^2 - M_N^2)/2M_N = 338$ MeV. Equivalently, Eq. (3) can be written as

$$ \omega_{\Delta} = \frac{M_N \Delta M + E^2(1 - \cos \theta)}{M_N + E(1 - \cos \theta)}, \quad (4) $$

also shown in Fig. 1. Beyond the $\Delta$ excitations, higher energy nucleon resonances contribute to the cross section up to the onset of the deep inelastic scattering. On the other hand, in the $\omega$ region between the quasielastic and $\Delta$ peaks, the so-called dip region, a large part of the cross section is due to multinucleon excitations, arising from nucleon-nucleon correlations and meson exchange processes. These complicated many-body mechanisms have attracted a lot of attention in the neutrino community in recent years, starting from the suggestion [51] of their inclusion as a natural explanation of the MiniBooNE cross sections [45]. These mechanisms play a crucial role in the neutrino energy reconstruction problem [52–56].

In this paper, we remove the small number of data points corresponding to the low energy giant resonances. Their presence worsens the neural network predictions for essentially two reasons: They are characterized by the largest (and spiked) cross sections, and they represent a small proportion of the whole dataset, which focuses on quasielastic excitations and beyond. We postpone the inclusion of these excitations to future work.

Several microscopic quantum-mechanical and phenomenological approaches exist to describe the different reaction mechanisms that contribute to the inclusive electron-nucleus cross section. While some approaches allow a description of all the nuclear excitations mentioned earlier, others focus on specific excitations only. For a review of these models we refer the reader to Refs. [49,57], where different approaches are discussed in connection with both electron and neutrino scattering. In the following we analyze to what extent neural networks can reproduce the electron-nucleus cross sections.

### III. NEURAL NETWORKS

In regression tasks, machine learning models are designed to learn mappings between an input domain and a predefined set of outputs for the target variables. The dataset considered in this paper is composed of different features: The atomic number Z, the nucleus mass number A, the electron beam energy E, the electron-scattering angle $\theta$, the energy loss $\omega$, and the inclusive double-differential cross section $d^2\sigma/d\omega d\Omega$. This last quantity is the value to be predicted by the supervised learning algorithm we present. Beyond the five data features (Z, A, E, $\theta$, and $\omega$), we add four additional complementary variables, obtained as combinations of the five original ones. These added variables are $\cos \theta$, $\omega_{QE}$ defined in Eq. (2), $\omega_{\Delta}$ defined in Eq. (4), and $Q^2 = q^2 - \omega^2$, here defined as

$$ Q^2 = 2E(E - \omega)(1 - \cos \theta), \quad (5) $$

a relation that is only valid when the electron mass is set to zero.

The addition of handcrafted features, such as the four-momentum transfer in the reaction and the approximate position of the quasielastic and $\Delta$ peaks [despite a possible offset between real and approximated values due to the lack of removal energy parameters in Eqs. (2) and (4)], is found to be useful to drive the network.

A larger neural network is expected to perform equally well without adding handcrafted features but our choice to include these additional variables makes the optimization process easier and allows us to use a smaller network, that is faster to train and less prone to overfitting. However, exactly quantifying how these four variables actually improve the prediction is not possible since there is no guarantee that the fine tuning is equally well performed when we add the four additional variables and when we do not add them.

The model is built upon deep neural networks (DNNs) [58]. DNNs are designed to learn hierarchical and increasingly abstract representations of the data. A DNN is a composition of L parametric functions named layers. The output of each layer, $f_l$, $l \in \{0, \dots, L-1\}$, is understood as a representation of the input samples. More specifically, the layer is composed of neurons, which are the building blocks of a layer. Hence, the layer $f_l$ takes the output of the previous layer $f_{l-1}$ and applies a nonlinear transformation to compute its output. These transformations use the model parameters, $W_l$ for each layer, commonly called weights and conveniently represented by a rectangular matrix. They relate the neurons of a layer to the previous layer and contain the information extracted by the model from the training data. Thus, given an input $x$, a neural network $f$ performs the following computation to infer its output:

$$ f(W, x) = f_{L-1}(W_{L-1}, f_{L-2}(W_{L-2}, \dots, f_0(W_0, x))) \dots. \quad (6) $$

In Eq. (6), the state of layer $l$ can be denoted $x^l$ and is called a representation. In a standard feedforward network, $x^{l-1}$ is employed to compute $x^l$ as follows:

$$ x^l = f_{l-1}(W_{l-1}, x^{l-1}) = g(W_{l-1}x^{l-1} + b_{l-1}), \quad (7) $$

where $W_{l-1}x^{l-1}$ is a matrix product, $x^{l-1}$ can be seen as a column vector of $n_{l-1}$ components, $W_{l-1}$ is then a matrix of $n_{l-1}$ columns and $n_l$ lines, and $b_{l-1}$ is a column vector of $n_l$ lines. g is an elementwise nonlinear function. We used the standard rectified linear unit (ReLU) function defined as $ReLU(z) = max\{0, z\}$. Other choices might work equally well or even better, and we make no claim of optimality in this paper. The result of Eq. (7) is the input of the next layer of the neural network, layer $l$.

During the initial exploratory phase of this paper we used fully connected neural networks as described by Eq. (7). We found that the richness of the data was such that deep neural networks with less than 35 layers were unable to fit the training set. We therefore used deeper networks that then became difficult to train for the reasons explained in Ref. [59]. The solution proposed by the authors of Ref. [59] is to use what they named a residual network. The architecture of our neural network is a fully connected residual network [59] which means that we perform an addition that shortcuts the architecture to help the optimization
<figure><img src="image_3.png" /><figcaption>FIG. 2. Illustration of the chosen neural network, composed of ten fully connected residual blocks. Details of one of these residual blocks are shown on the right of the figure. The activation function used is a rectified linear unit.</figcaption></figure>

procedure converge. The skipped connection is performed as follows:

$$x^l = f_{l-1}(W_{l-1}, x^{l-1}, x^{l-5}) = g[W_{l-1}g(x^{l-1} + x^{l-5}) + b_{l-1}]. \quad (8)$$

The skipped connections or residual blocks are used in deep learning to help training deep neural networks [58,59]. Our neural network, implemented with TENSORFLOW [60], is composed of ten residual blocks of five fully connected layers with increasing and then decreasing sizes ranging from input size (nine input data) up to size 50 and then down to size 5 before connecting to the output of size 1. Our neural network architecture is represented in Fig. 2, and the resulting model has 50 connected layers. To help information flow

in the backpropagation stage we have chosen a layer width that varies progressively. Abrupt changes in widths made the optimization process harder and local minima were difficult to escape from. The constraint we have is that we start with nine inputs and have an output of dimension 1.

An essential element of neural networks training is the choice of the cost function. For the neural network to be able to accurately predict the cross section structure over its wide range and different shapes, we chose to use the relative absolute error given by

$$C(\hat{y}) = \frac{1}{n} \sum_{i=1}^{n} \frac{|y_i - \hat{y}_i|}{y_i}. \quad (9)$$
In Eq. (9), $\hat{y}_i$ is a model prediction and y is the cross section divided by A, $y = (d^2\sigma/d\omega d\Omega)/A$ obtained from the dataset. In our case,

$$ \hat{y}_i = f(W, x_i), \qquad (10) $$

where the input data vector is the nine-dimensional vector of normalized components: $x_i = (Z, A, E, \theta, \omega, \cos\theta, \omega_{QE}, \omega_\Delta, Q^2)$ and $(x_i, y_i)$ is an input-output pair. Another choice for the loss function could be the $\chi^2$ that would allow to take into account the experimental errors. For this first exploratory study, driven by a major simplicity and rapidity of calculations and by the fact that in the present dataset the experimental errors are in general small, we prefer to choose the relative absolute error as cost function. We leave the use of $\chi^2$ to future works.

During the training phase, the model is presented with a large number of input-output pairs. The weights, W, are initialized randomly with a normal distribution [58] and an initial prediction of the output is computed using Eq. (6). During the forward pass, the models prediction error is estimated by computing the value of a cost function that quantifies the discrepancy between the current prediction $\hat{y}$ given by Eq. (9) and the true target y. This cost function is computed over a batch of data and differentiated with respect to the model weights, W, in the backward pass. Model weights are then updated with the computed differential to improve the predictions of the neural network. The general idea is to update W in the following way:

$$ W = W - \epsilon \frac{\partial C}{\partial W}. \qquad (11) $$

The variable $\epsilon$ is called the learning rate and the value we chose was standard, $\epsilon = 0.001$. In practice, we use the ADAM optimizer [61] which provides an empirical improvement over Eq. (11) with a slowly decreasing learning rate, reminiscent of annealing, to help the optimization process. The values of the neural network weights that minimize the cost function are obtained during the training phase by iteratively taking forward and backward passes on the training dataset. Once a predetermined number of epochs (a forward and a backward pass) is computed, the training phase is finished and the model is deployed on new samples, unseen during the training, to make predictions. Training takes about 30 min for 5000 epochs on a 12GB NVIDIA Tesla K80 GPU. We include standard $L_2$ regularization in order to reduce the variance of the model and hopefully improve generalization: The predictions on data not used during the training phase. Weights are stored during training and for inference we roll back to the epoch that yielded the lowest validation error during the whole training. In the inference phase, the output of the model is computed for the test dataset. One hopes that the model will generalize well on the test data and infer target values close to the actual values thanks to the patterns it has learned from its training data.

## IV. RESULTS

The training set is composed of the available samples excluding one or more chosen isotopes ($Z', A'$), or examples corresponding to some particular kinematical conditions.

<figure><img src="image_4.png" /><figcaption>FIG. 3. Evolution of the cost function during the 4000 epochs training process for a successful training run. This illustrates a slight overfitting, controlled by $L_2$ regularization.</figcaption></figure>

These excluded data will be examined in the testing stage. In other words, we used the trained neural network to make predictions for the isotopes ($Z', A'$) or for some specific kinematical conditions. This allows us to assess the ability of the trained neural network to infer the desired target for cases unseen during the training stage.

### A. Training results

To further illustrate our procedure let us consider some specific cases. The whole dataset, representing our starting point, is the subset of the electron-scattering data [50] corresponding to nuclei from $^4$He to $^{59}$Ni. Their relative contribution is illustrated in Table I. As a first experiment, we remove from the entire dataset only the isotope ($Z = 20, A = 40$), i.e., the $^{40}$Ca, our test nucleus, for which we intend to predict the cross section. The rest of the dataset is randomly split: 90% for the training set and the remaining 10% for the validation set that allows us to monitor overfitting.

Stochastic gradient descent on minibatches is used for parameters optimization [58]. The minibatches of 512 instances of the training set are constituted randomly. This promotes regularization and speeds up the training. As a standard practice, at the end of the training process over 4000 epochs, we plot the evolution of the cost function, displayed in Fig. 3, to check that the training has converged as expected and to monitor overfitting. In order to reduce the variance of the model and to control overfitting we used $L_2$ regularization on the first and largest layers. It is standard practice to regularize large layers, that are more prone to overfitting, more than small layers. To limit overfitting we use for our prediction the parameters that gave the best validation value.

In order to illustrate the overall performance of the neural network on the training set, we consider its predictions for $^{12}$C, a nucleus currently employed in the T2K and NOvA neutrino detectors. For this nucleus, many electron-scattering cross section data are available, for different kinematical conditions, inducing different nuclear excitations, from nuclear giant resonances up to deep inelastic scattering, as discussed in Sec. II.
The neural network double-differential cross section results for $^{12}$C are illustrated in Figs. 4 and 5 for a selected set of the values of incoming electron energy $E$ and electron-scattering angle $\theta$ present in the database. The overall agreement with data is remarkable for all the incoming energies and electron-scattering angles.<fnref n="2" />

These good results suggest that neural networks can be considered as a tool for representing data that have this kind of behavior (changing shape depending on the type of induced excitation). In this context, we recall that for $^{12}$C the values of $d^2\sigma/d\omega d\Omega$ vary between a minimum of $7.7 \times 10^{-5}$ nb/(GeV sr) and a maximum of $9.6 \times 10^5$ nb/(GeV sr). Despite this variation of ten orders of magnitude, the choice of the mean relative absolute error of Eq. (9) as cost function allows for efficient training of the neural network.

## B. Testing results

Following the results obtained for the training set, in this section we investigate the predictive power of the neural networks for the test set. First, we consider $^{40}$Ca as test nucleus. This choice is motivated by the existence of many available cross section data for this nucleus. The neural network predictions for this nucleus, which, we reiterate, was not included in the training set, are shown in Figs. 6 and 7. These figures show the result of an average over $N = 100$ different neural networks trained on the same dataset. After initializing the same neural network architecture with $N$ different realizations of a random variable, we train the neural networks. For each training run we keep the weights that performed best on the validation set. These different neural networks give different predictions on the test set. For each scattering angle $\theta$, each energy $E$, and each transferred energy $\omega$, we obtain $N$ independent values for $\hat{y}_i$, the prediction for $d^2\sigma/d\omega d\Omega$. We can then estimate the average based on these $N$ different predictions. Assuming a Gaussian distribution of these predictions, average and 95% confidence intervals for the true average $\mu$ are estimated as follows:

$$ \bar{y} = \frac{1}{N} \sum_{i=1}^{N} \hat{y}_i, \qquad (12) $$

$$ \sigma^2 = \frac{1}{N-1} \sum_{i=1}^{N} (\hat{y}_i - \bar{y})^2, \qquad (13) $$

$$ |\bar{y} - \mu| < 1.96 \times \frac{\sigma}{\sqrt{N}}. \qquad (14) $$

The estimated average of the predictions is plotted in Figs. 6 and 7 together with the estimated 95% confidence interval, bounded by two red lines and shaded in blue. One can observe that the overall agreement of our neural network predictions with the data for $^{40}$Ca is satisfactory from quasielastic up to the $\Delta$ resonance. The prediction results could be further

improved by adding more examples into the training set, with different choices of the neural-network hyperparameters or including more handcrafted features. One could also use alternative approaches, like Bayesian neural networks, to better take into account the uncertainty. However, for the exploratory purpose of the present paper, we do not claim optimality in any sense, and our focus is on presenting the interest of the method.

As an additional illustration of the predictive performances of the neural network approach we consider some tests on nuclei of particular interest for the neutrino-oscillation program, namely, oxygen, employed in the present and future water Cherenkov detectors Super-Kamiokande and Hyper-Kamiokande, and argon, the nucleus employed in three detectors of the present short-baseline Fermi National Accelerator Laboratory neutrino programs as well as by the future DUNE experiment.

In order to investigate the predictive power of our neural network for some specific kinematical conditions, in the case of oxygen (for which only few data are available, as shown in Table I) we remove from the training dataset all the events corresponding to electron scattering on $^{16}$O with a scattering angle of $\theta = 37.1^\circ$ and an incoming electron energy of $E = 0.737$ GeV. The predictions are shown in the left-hand panel of Fig. 8. In the same figure, in its middle and right-hand panel, we show the neural network predictions for argon and titanium. For these two nuclei, the inclusive electron-scattering data have been made available only recently [23–25]. Due to the increasing interest in argon detectors, the results of titanium ($Z = 22$) represent an important source of information for the charged current interaction of neutrinos with a neutron of $^{40}$Ar ($Z = 18, N = 22$). For $^{40}$Ar and $^{48}$Ti, experimental data only exist for one particular experimental setup (the same for both nuclei), corresponding to $E = 2.222$ GeV and $\theta = 15.541^\circ$. We removed from the usual dataset data points corresponding to these two nuclei in order to avoid the risk of overfitting due to the presence in the training dataset of other similar nuclei with the same kinematical conditions. For all three panels of Fig. 8, the represented 95% confidence interval is computed by considering $N = 40$ independent neural network predictions, following the methodology previously explained for $^{40}$Ca (where it was taken to be $N = 100$). In the case of $^{16}$O the band is very narrow because the confidence intervals are small but in the case of $^{40}$Ar and $^{48}$Ti it is more visible in the figure. We explain this because of the complete absence of $^{40}$Ar and $^{48}$Ti in the training set.

For all three cases shown in Fig. 8, the neural network predictions are excellent for the whole dataset, from the quasielastic bump to the $\Delta$ resonance, passing by the dip region. Similar agreement has also been obtained by other theoretical approaches, such as SuSAv2 [57] and Giessen Boltzmann-Uehling-Uhlenbeck (GiBUU) [64], largely employed in neutrino cross section studies. We briefly remind the reader that SuSAv2, an updated version of the super scaling analysis [65] approach, is based on the observation of the superscaling behavior (i.e., on the simultaneous $q$ and Fermi momentum independence) of the scaling function (a sort of nuclear response function in terms of a scaling

<footnote n="2">²There are some additional kinematical conditions that are not shown in Figs. 4 and 5 to avoid overwhelming figures but the agreement between the data and the neural network fit on the training set remains remarkable also in the cases not shown.</footnote>
<figure><img src="image_5.png" /><figcaption>FIG. 4. Comparisons between the experimental inclusive (e, e') double differential cross section on $^{12}$C, in black squares with experimental error bars, and the neural network's predictions, red circles, for the different values of incoming electron energy *E* and scattering angle *θ* (indicated above each panel), as a function of the transferred energy. To obtain these predictions, the $^{12}$C data occurrences were used in the training set as well as all the other nuclei of the dataset except $^{40}$Ca.</figcaption></figure>
<figure><img src="image_6.png" /><figcaption>FIG. 5. The same as Fig. 4 but for other kinematical conditions.</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>FIG. 6. Comparisons between the experimental inclusive (e, e') double-differential cross section on $^{40}$Ca, in black squares with experimental error bars, and the neural network's predictions, red circles, for the different values of incoming electron energy *E* and scattering angle *θ* written above each panel, as a function of the transferred energy. These are the predictions for the test data for which the neural network has not been trained. A 95% confidence interval is shown by the shaded blue area, bounded by red lines.</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>FIG. 7. The same as Fig. 6 but for other kinematical conditions.</figcaption></figure>
<figure><img src="image_9.png" /><figcaption>FIG. 8. Comparison between the inclusive (e, e') double differential cross section data and the neural network's predictions for ¹⁶O (left panel), ⁴⁰Ar (middle panel), and ⁴⁸Ti (right panel). The different values of incoming electron energy and scattering angle are written above each panel. For the left panel, the neural network was trained without the data represented on the figure. For the middle and right panels, the neural network was trained without both ⁴⁰Ar and ⁴⁸Ti data. Theoretical predictions of SuSAv2 and GiBUU approaches are also shown. They are taken respectively from Refs. [14,57] for ¹⁶O and from Refs. [62,63] for ⁴⁰Ar and ⁴⁸Ti. GENIE Monte Carlo event generator results, taken from Ref. [20] in the case of ⁴⁰Ar and from Ref. [18] in the case of ⁴⁸Ti, are also displayed. The red lines are the boundary of the 95% confidence interval shaded in blue.</figcaption></figure>

variable) extracted from electron-scattering data. Similarly to the neural network, it can be considered as a phenomenological data-driven approach, obviously with a much lower number of parameters and with a sound physical basis which microscopically justify it [66]. The GiBUU framework is an implementation of transport theory which allows a description of many nuclear reactions. It allows consistent treatment of the initial reaction vertex and of the final state processes. It takes into account various nuclear effects via the local density approximation for the nuclear ground state, mean-field potentials, and in-medium spectral functions. In the implementation of Ref. [14], from which we take the electron cross section results of oxygen shown in Fig. 8, multinucleon excitations were not included, which explains the GiBUU's underestimation of data in the dip region in this case. The underestimation was reduced in the GiBUU calculations for argon and titanium [63] by using an updated version of the model that takes into account the multinucleon contributions. We note that the dip region, which is particularly difficult to microscopically describe and reproduce with theoretical approaches, is well predicted by the neural network.

For completeness, in Fig. 8 we also plot the results obtained with two different versions of GENIE (one of the most commonly used Monte Carlo event generators in the neutrino community) for ⁴⁰Ar and ⁴⁸Ti, taken from Refs. [20] and [18], respectively. One can observe that starting from the dip region GENIE results disagree with the data. This disagreement is analyzed in Refs. [18,20], where several results are also shown and discussed for inclusive electron-scattering on carbon.

## V. SUMMARY AND PERSPECTIVES

We have deployed a neural network model to predict the electron-scattering inclusive double-differential cross sections on nuclei. The neural network predictions have been compared with a large amount of electron-scattering data collected in the past as well as with recent data on argon and

titanium that are of interest to the neutrino community. Our results show that neural networks can reproduce and predict the electron-scattering cross section with an accuracy comparable to that provided by microscopic approaches developed in the past and nowadays generalized to investigate the neutrino-nucleus scattering.

Neural network approaches can be impaired by a limited amount of data. Therefore, it was not a given that the data available were sufficient to train the model to make inference for unseen data. The many orders of magnitude spanned by the cross sections and their varying shape due to the different types of nuclear excitations (quasielastic, multinucleon, resonance excitations, and deep inelastic) induced in different scattering kinematics make inference challenging. The neural network performance suggests that it could be used as an additional tool in the studies of electron and neutrino scattering on nuclei. Moreover, this can be used to predict the electron cross sections for nuclei and/or for kinematical conditions where experimental data are absent and employ these predictions to validate Monte Carlo simulations. It could also be employed as a support to drive and speed up microscopic evaluation of cross sections and/or response functions, in some sense by generalizing what has been performed in Ref. [43] beyond the quasielastic excitation considered in that work.

Another perspective involves the use of neural networks to directly predict the neutrino-nucleus cross sections. This task is not trivial since the neutrino beams are not monochromatic. So in this case the measured quantity is the flux-integrated double differential cross section in terms of the final state measurable variables. For the charged-current scattering process ($v_l, l$) on nuclei, where $l$ is the charged lepton, this neutrino cross section reads

$$ \frac{d^2\sigma_v}{dE_l d\Omega} = \frac{1}{\int \Phi(E_{v_l}) dE_{v_l}} \int dE_{v_l} \left[ \frac{d^2\sigma_v}{d\omega d\Omega} \right]_{\omega=E_{v_l}-E_l} \Phi(E_{v_l}), \quad (15) $$
where $E_l$ is the charged lepton energy and $d\Omega$ is the differential solid angle in the direction specified by the charged lepton momentum. They represent the two final state measured variables. This expression formally reduces to that analyzed in this paper, when the flux of incoming particles, $\Phi$, reduces to a delta distribution. This is the case of monochromatic electron beams. For neutrino scattering, the problem is more complex since for a given set of the measured variables $E_l$ and $\Omega$ one explores the full energy spectrum of neutrinos above the charged lepton energy, being $E_{\nu_l} = E_l + \omega$. As a consequence, all the reaction channels (giant resonances, quasielastic, multinucleon excitations, pion production arising from nucleon resonances decay, and deep-inelastic scattering) are entangled and isolating a primary vertex process from the measurement of the neutrino flux-integrated differential cross section is much more difficult. When employing a deep learning approach for neutrino scattering, the technical challenge would be to design a network sufficiently rich to encode the

complexity of the cross section for different primary vertex processes over the phase space relevant to the signal process. A similar approach is currently under study in the high energy collider community [30] in connection with the so-called matrix element method. This possibility should be investigated with the perspective of the development of an AI-based Monte Carlo event generator for neutrino-nucleus scattering.

## ACKNOWLEDGMENTS

We would like to thank Fabrice Couderc for his suggestions on neural network optimization. We warmly thank Stefano Fortunati for interesting conversations on statistical topics. We also acknowledge Christophe Coudé, Raphaël Lasseri, Gabriel Perdue, Pierre-Etienne Pion, David Regnier, Guillaume Scamps, and Kazuhiro Terao for fruitful discussions and exchanges on neural networks in connection with nuclear and neutrino physics.

[1] J. Franck and G. Hertz, Verh. Dtsch. Phys. Ges. **16**, 457 (1914).

[2] R. Hofstadter, Rev. Mod. Phys. **28**, 214 (1956).

[3] E. D. Bloom, D. H. Coward, H. DeStaebler, J. Drees, G. Miller, L. W. Mo, R. E. Taylor, M. Breidenbach, J. I. Friedman, G. C. Hartmann et al., Phys. Rev. Lett. **23**, 930 (1969).

[4] M. Breidenbach, J. I. Friedman, H. W. Kendall, E. D. Bloom, D. H. Coward, H. C. DeStaebler, J. Drees, L. W. Mo, and R. E. Taylor, Phys. Rev. Lett. **23**, 935 (1969).

[5] B. Pontecorvo, Zh. Eksp. Teor. Fiz. **37**, 1751 (1959).

[6] M. Schwartz, Phys. Rev. Lett. **4**, 306 (1960).

[7] G. Danby, J. M. Gaillard, K. A. Goulianos, L. M. Lederman, N. B. Mistry, M. Schwartz, and J. Steinberger, Phys. Rev. Lett. **9**, 36 (1962).

[8] K. Abe et al. (T2K), Nucl. Instrum. Methods Phys. Res. Sect. A **659**, 106 (2011).

[9] D. S. Ayres et al. (NOvA Collaboration), The NOvA Technical Design Report, Report No. FERMILAB-DESIGN-2007-01, Fermi National Accelerator Laboratory, 2007.

[10] R. Acciarri et al. (DUNE), arXiv:1601.05471.

[11] K. Abe et al. (Hyper-Kamiokande Proto-Collaboration), Prog. Theor. Exp. Phys. **2015**, 053C02 (2015).

[12] T. Katori and M. Martini, J. Phys. G **45**, 013001 (2018).

[13] L. Alvarez-Ruso et al. (NuSTEC), Prog. Part. Nucl. Phys. **100**, 1 (2018).

[14] T. Leitner, O. Buss, L. Alvarez-Ruso, and U. Mosel, Phys. Rev. C **79**, 034601 (2009).

[15] A. M. Ankowski, O. Benhar, and M. Sakuda, Phys. Rev. D **91**, 033005 (2015).

[16] V. Pandey, N. Jachowicz, T. Van Cuyck, J. Ryckebusch, and M. Martini, Phys. Rev. C **92**, 024606 (2015).

[17] G. D. Megias, J. E. Amaro, M. B. Barbaro, J. A. Caballero, and T. W. Donnelly, Phys. Rev. D **94**, 013012 (2016).

[18] A. M. Ankowski and A. Friedland, Phys. Rev. D **102**, 053001 (2020).

[19] J. L. Barrow, S. Gardiner, S. Pastore, M. Betancourt, and J. Carlson, Phys. Rev. D **103**, 052001 (2021).

[20] A. Papadopoulou et al. (e4v), Phys. Rev. D **103**, 113003 (2021).

[21] C. Andreopoulos, A. Bell, D. Bhattacharya, F. Cavanna, J. Dobson et al., Nucl. Instrum. Methods Phys. Res. Sect. A **614**, 87 (2010).

[22] M. Khachatryan et al. (CLAS, e4v), Nature (London) **599**, 565 (2021).

[23] H. Dai et al. (Jefferson Lab Hall A), Phys. Rev. C **98**, 014617 (2018).

[24] H. Dai et al., Phys. Rev. C **99**, 054608 (2019).

[25] M. Murphy et al., Phys. Rev. C **100**, 054606 (2019).

[26] L. Gu et al. (Jefferson Lab Hall A), Phys. Rev. C **103**, 034604 (2021).

[27] O. Benhar et al., arXiv:1406.4080.

[28] P. A. Machado, O. Palamara, and D. W. Schmitz, Annu. Rev. Nucl. Part. Sci. **69**, 363 (2019).

[29] P. Baldi, P. Sadowski, and D. Whiteson, Nat. Commun. **5**, 4308 (2014).

[30] K. Albertsson et al., J. Phys.: Conf. Ser. **1085**, 022008 (2018).

[31] A. Radovic, M. Williams, D. Rousseau, M. Kagan, D. Bonacorsi, A. Himmel, A. Aurisano, K. Terao, and T. Wongjirad, Nature (London) **560**, 41 (2018).

[32] X. Ju et al., arXiv:2003.11603.

[33] P. Calafiura, D. Rousseau, and K. Terao, Artificial Intelligence for High Energy Physics (World Scientific, Singapore, 2022).

[34] R. Acciarri et al. (MicroBooNE), J. Inst. P03011 (2017).

[35] G. N. Perdue et al. (MINERvA), J. Inst. P11020 (2018).

[36] P. Baldi, J. Bian, L. Hertel, and L. Li, Phys. Rev. D **99**, 012011 (2019).

[37] B. Abi et al. (DUNE), Phys. Rev. D **102**, 092003 (2020).

[38] S. Alonso-Monsalve, D. Douqa, C. Jesús-Valls, T. Lux, S. Pina-Otey, F. Sánchez, D. Sgalaberna, and L. H. Whitehead, Phys. Rev. D **103**, 032005 (2021).

[39] P. Abratenko et al. (MicroBooNE), Phys. Rev. D **103**, 052012 (2021).

[40] A. Ghosh et al. (MINERvA), arXiv:2103.06992.

[41] Y. Alanazi, N. Sato, T. Liu, W. Melnitchouk, M. P. Kuchera, E. Pritchard, M. Robertson, R. Strauss, L. Velasco, and Y. Li, arXiv:2001.11103.

[42] Y. Alanazi et al., arXiv:2008.03151.

[43] K. Raghavan, P. Balaprakash, A. Lovato, N. Rocco, and S. M. Wild, Phys. Rev. C **103**, 035502 (2021).

[44] L. Alvarez-Ruso, K. M. Graczyk, and E. Saul-Sala, Phys. Rev. C **99**, 025204 (2019).
[45] A. Aguilar-Arevalo et al. (MiniBooNE), Phys. Rev. D **81**, 092005 (2010).

[46] S. Nagu, J. Singh, J. Singh, and R. B. Singh, arXiv:2105.12733.

[47] P. Bedaque et al., Eur. Phys. J. A **57**, 100 (2021).

[48] O. Benhar, D. Day, and I. Sick, arXiv:nucl-ex/0603032.

[49] O. Benhar, D. day, and I. Sick, Rev. Mod. Phys. **80**, 189 (2008).

[50] http://discovery.phys.virginia.edu/research/groups/qes-archive/.

[51] M. Martini, M. Ericson, G. Chanfray, and J. Marteau, Phys. Rev. C **80**, 065501 (2009).

[52] M. Martini, M. Ericson, and G. Chanfray, Phys. Rev. D **85**, 093012 (2012).

[53] M. Martini, M. Ericson, and G. Chanfray, Phys. Rev. D **87**, 013009 (2013).

[54] J. Nieves, F. Sanchez, I. R. Simo, and M. J. Vicente Vacas, Phys. Rev. D **85**, 113008 (2012).

[55] O. Lalakulich, U. Mosel, and K. Gallmeister, Phys. Rev. C **86**, 054606 (2012).

[56] A. M. Ankowski and C. Mariani, J. Phys. G **44**, 054001 (2017).

[57] J. E. Amaro, M. B. Barbaro, J. A. Caballero, R. González-Jiménez, G. D. Megias, and I. Ruiz Simo, J. Phys. G **47**, 124001 (2020).

[58] I. Goodfellow, Y. Bengio, and A. Courville, Deep Learning (MIT Press, Cambridge, MA, 2016).

[59] K. He, X. Zhang, S. Ren, and J. Sun, arXiv:1512.03385.

[60] M. Abadi, A. Agarwal, P. Barham, E. Brevdo, Z. Chen, C. Citro, G. S. Corrado, A. Davis, J. Dean, M. Devin et al., arXiv:1603.04467 (2015); software available at https://www.tensorflow.org/.

[61] D. P. Kingma and J. Ba, arXiv:1412.6980.

[62] M. B. Barbaro, J. A. Caballero, A. De Pace, T. W. Donnelly, R. González-Jiménez, and G. D. Megias, Phys. Rev. C **99**, 042501(R) (2019).

[63] U. Mosel and K. Gallmeister, Phys. Rev. C **99**, 064605 (2019).

[64] O. Buss, T. Gaitanos, K. Gallmeister, H. van Hees, M. Kaskulov et al., Phys. Rep. **512**, 1 (2012).

[65] J. E. Amaro, M. B. Barbaro, J. A. Caballero, T. W. Donnelly, A. Molinari, and I. Sick, Phys. Rev. C **71**, 015501 (2005).

[66] R. González-Jiménez, G. D. Megias, M. B. Barbaro, J. A. Caballero, and T. W. Donnelly, Phys. Rev. C **90**, 035501 (2014).