Unbiased determination of DVCS Compton Form Factors

TL;DR

This paper presents a comprehensive, model-independent extraction of DVCS Compton Form Factors using neural networks and advanced statistical methods, providing new insights into proton structure from experimental data.

Contribution

It introduces a neural network-based approach for unbiased CFF extraction, reducing model dependence and incorporating uncertainty propagation in DVCS analysis.

Findings

Neural network parameterizations of CFFs with reduced model bias

Successful propagation of experimental uncertainties to CFFs

Determination of the subtraction constant via dispersion relations

Abstract

The extraction of Compton Form Factors (CFFs) in a global analysis of almost all Deeply Virtual Compton Scattering (DVCS) proton data is presented. The extracted quantities are DVCS sub-amplitudes and the most basic observables which are unambiguously accessible from this process. The parameterizations of CFFs are constructed utilizing the artificial neural network technique allowing for an important reduction of model dependency. The analysis consists of such elements as feasibility studies, training of neural networks with the genetic algorithm and a careful regularization to avoid over-fitting. The propagation of experimental uncertainties to extracted quantities is done with the replica method. The resulting parameterizations of CFFs are used to determine the subtraction constant through dispersion relations. The analysis is done within the PARTONS framework.

Tables

Table 1. Table 1: Values used in the scaling of input and output variables.

$v^{\prime }$	$v_{\min }^{\prime }$	$v_{\max }^{\prime }$		$v^{\prime }$	$v_{\min }^{\prime }$	$v_{\max }^{\prime }$
${\log }_{10}\xi$	$-6$	$1$		$t$	$-1$	$0.5$
${\log }_{10}{Q}^2$	$-1$	$2$
$\xi \mathrm{Re}\mathscr{H}$	$-1$	$1.5$		$\xi \mathrm{Im}\mathscr{H}$	$-4$	$6$
$\xi \mathrm{Re}\tilde{\mathscr{H}}$	$-0.6$	$0.9$		$\xi \mathrm{Im}\tilde{\mathscr{H}}$	$-1$	$1.5$
$\xi \mathrm{Re}ℰ$	$-1.5$	$1$		$\xi \mathrm{Im}ℰ$	$-6$	$4$
$\xi \mathrm{Re}\widetilde{ℰ}$	$-120$	$180$		$\xi \mathrm{Im}\widetilde{ℰ}$	$-8$	$12$

Table 2. Table 2: DVCS data used in this analysis.

No.	Collab.	Year	Ref.	Observable		Kinematic dependence	No. of points used / all
1	HERMES	2001	Airapetian:2001yk	$A_{LU}^{+}$		$\varphi$	10 / 10
2		2006	Airapetian:2006zr	$A_C^{\cos i\varphi }$	$i=1$	$t$	4 / 4
3		2008	Airapetian:2008aa	$A_C^{\cos i\varphi }$	$i=0,1$	$x_{\mathrm{Bj}}$	18 / 24
				$A_{UT,\mathrm{DVCS}}^{\sin \left(\varphi -{\varphi }_S\right)\cos i\varphi }$	$i=0$
				$A_{UT,I}^{\sin \left(\varphi -{\varphi }_S\right)\cos i\varphi }$	$i=0,1$
				$A_{UT,I}^{\cos \left(\varphi -{\varphi }_S\right)\sin i\varphi }$	$i=1$
4		2009	Airapetian:2009aa	$A_{LU,I}^{\sin i\varphi }$	$i=1,2$	$x_{\mathrm{Bj}}$	35 / 42
				$A_{LU,\mathrm{DVCS}}^{\sin i\varphi }$	$i=1$
				$A_C^{\cos i\varphi }$	$i=0,1,2,3$
5		2010	Airapetian:2010ab	$A_{UL}^{+,\sin i\varphi }$	$i=1,2,3$	$x_{\mathrm{Bj}}$	18 / 24
				$A_{LL}^{+,\cos i\varphi }$	$i=0,1,2$
6		2011	Airapetian:2011uq	$A_{LT,\mathrm{DVCS}}^{\cos \left(\varphi -{\varphi }_S\right)\cos i\varphi }$	$i=0,1$	$x_{\mathrm{Bj}}$	24 / 32
				$A_{LT,\mathrm{DVCS}}^{\sin \left(\varphi -{\varphi }_S\right)\sin i\varphi }$	$i=1$
				$A_{LT,I}^{\cos \left(\varphi -{\varphi }_S\right)\cos i\varphi }$	$i=0,1,2$
				$A_{LT,I}^{\sin \left(\varphi -{\varphi }_S\right)\sin i\varphi }$	$i=1,2$
7		2012	Airapetian:2012mq	$A_{LU,I}^{\sin i\varphi }$	$i=1,2$	$x_{\mathrm{Bj}}$	35 / 42
				$A_{LU,\mathrm{DVCS}}^{\sin i\varphi }$	$i=1$
				$A_C^{\cos i\varphi }$	$i=0,1,2,3$
8	CLAS	2001	Stepanyan:2001sm	$A_{LU}^{-,\sin i\varphi }$	$i=1,2$	—	0 / 2
9		2006	Chen:2006na	$A_{UL}^{-,\sin i\varphi }$	$i=1,2$	—	2 / 2
10		2008	Girod:2007aa	$A_{LU}^{-}$		$\varphi$	283 / 737
11		2009	Gavalian:2008aa	$A_{LU}^{-}$		$\varphi$	22 / 33
12		2015	Pisano:2015iqa	$A_{LU}^{-}$, $A_{UL}^{-}$, $A_{LL}^{-}$		$\varphi$	311 / 497
13		2015	Jo:2015ema	$d^4{\sigma }_{UU}^{-}$		$\varphi$	1333 / 1933
14	Hall A	2015	Defurne:2015kxq	$\Delta d^4{\sigma }_{LU}^{-}$		$\varphi$	228 / 228
15		2017	Defurne:2017paw	$\Delta d^4{\sigma }_{LU}^{-}$		$\varphi$	276 / 358
16	COMPASS	2018	Akhunzyanov:2018nut	$d^3{\sigma }_{UU}^{\pm }$		t	2 / 4
17	ZEUS	2009	Chekanov:2008vy	$d^3{\sigma }_{UU}^{+}$		t	4 / 4
18	H1	2005	Aktas:2005ty	$d^3{\sigma }_{UU}^{+}$		t	7 / 8
19		2009	Aaron:2009ac	$d^3{\sigma }_{UU}^{\pm }$		t	12 / 12
						SUM:	2624 / 3996

Table 3. Table 3: Values of the χ 2 superscript 𝜒 2 \chi^{2} function per data set. For a given data set, cf . Table 2 , given are: χ 2 superscript 𝜒 2 \chi^{2} value, the number of experimental points n 𝑛 n , and the ratio between those two numbers.

No.	Collab.	Year	Ref.	${\chi }^2$	$n$	${\chi }^2/n$
1	HERMES	2001	Airapetian:2001yk	$10.7$	$10$	$1.07$
2		2006	Airapetian:2006zr	$5.5$	$4$	$1.38$
3		2008	Airapetian:2008aa	$18.5$	$18$	$1.03$
4		2009	Airapetian:2009aa	$34.7$	$35$	$0.99$
5		2010	Airapetian:2010ab	$40.7$	$18$	$2.26$
6		2011	Airapetian:2011uq	$16.7$	$24$	$0.70$
7		2012	Airapetian:2012mq	$22.4$	$35$	$0.64$
8	CLAS	2001	Stepanyan:2001sm	—	$0$	—
9		2006	Chen:2006na	$1.0$	$2$	$0.52$
10		2008	Girod:2007aa	$376.4$	$283$	$1.33$
11		2009	Gavalian:2008aa	$28.3$	$22$	$1.29$
12		2015	Pisano:2015iqa	$306.6$	$311$	$0.99$
13		2015	Jo:2015ema	$884.7$	$1333$	$0.66$
14	Hall A	2015	Defurne:2015kxq	$231.8$	$228$	$1.02$
15		2017	Defurne:2017paw	$211.4$	$276$	$0.77$
16	COMPASS	2018	Akhunzyanov:2018nut	$3.0$	$2$	$1.50$
17	ZEUS	2009	Chekanov:2008vy	$5.49$	$4$	$1.38$
18	H1	2005	Aktas:2005ty	$22.2$	$7$	$3.17$
19		2009	Aaron:2009ac	$23.4$	$12$	$1.95$

Equations

$$l(k)+N(p)\to l(k^{\prime })+N(p^{\prime })+\gamma (q^{\prime }).$$

$$\mathcal{T}={\mathcal{T}}^{\mathrm{BH}}+{\mathcal{T}}^{\mathrm{DVCS}}.$$

$$\xi =-\frac{(q+q^{\prime }{)}^2}{2(p+p^{\prime })(q+q^{\prime })},$$

$$\eta =\frac{(q^{\prime }-q)(q+q^{\prime })}{(p+p^{\prime })(q+q^{\prime })}.$$

$$\mathrm{Re}\mathcal{G}(\xi )=C_G+$$

$$\frac{1}{\pi }{\int }_0^1\mathrm{d}{\xi }^{\prime }\mathrm{Im}\mathcal{G}({\xi }^{\prime })\left[\frac{1}{\xi -{\xi }^{\prime }}\mp \frac{1}{\xi +{\xi }^{\prime }}\right],$$

$$C_H=-C_E,$$

$$C_{\tilde{H}}=C_{\tilde{E}}=0.$$

$$v_j^i=f_{\text{includegraphics},j}^i\left(f_{\Sigma ,j}^i(\text{mathaccent}382{\beta }_j^i,\text{mathaccent}382v^{i-1})\right).$$

$$f_{\Sigma ,j}^i(\text{mathaccent}382{\beta }_j^i,\text{mathaccent}382v^{i-1})={\beta }_{j0}^i+{\beta }_{j1}^i{v}_1^{i-1}+{\beta }_{j2}^i{v}_2^{i-1}+\dots$$

$$f_{\text{includegraphics},j}^i(x)=\tanh (x),$$

$$f_{\text{includegraphics},j}^i(x)=x.$$

$$f_{\mathrm{stop}}(t_i)=\frac{f_{\mathrm{cost}}^{\mathrm{av}}(t_i)}{\underset{j<i}{\min }{f}_{\mathrm{cost}}^{\mathrm{av}}(t_j)}-1.$$

$$f_{\mathrm{cost}}^{\mathrm{av}}(t_i)=\frac{1}{2n+1}\sum _{j=-n}^n{f}_{\mathrm{cost}}(t_{i-j}),$$

$$v_{ij}\pm {\Delta }_{ij}^{\mathrm{tot}}\stackrel{\mathrm{replica}k}{\to }$$

$$\left({\mathrm{rnd}}_k(v_{ij},{\Delta }_{ij}^{\mathrm{tot}})\pm {\Delta }_{ij}^{\mathrm{tot}}\right)\times {\mathrm{rnd}}_k(1,{\Delta }_i^{\mathrm{norm}}).$$

$${\Delta }_{ij}^{\mathrm{tot}}=\sqrt{{\left({\Delta }_{ij}^{\mathrm{stat}}\right)}^2+{\left({\Delta }_{ij}^{\mathrm{sys}}\right)}^2}.$$

$${\xi }^{\prime }={\log }_{10}\xi ,$$

$$t^{\prime }=t,$$

$${Q^2}^{\prime }={\log }_{10}{Q}^2,$$

$$\mathrm{Re}{\mathcal{G}}^{\prime }=\xi \mathrm{Re}\mathcal{G},$$

$$\mathrm{Im}{\mathcal{G}}^{\prime }=\xi \mathrm{Im}\mathcal{G},$$

$$v^{\prime \prime }=-1+2\frac{v^{\prime }-v_{\min }^{\prime }}{v_{\max }^{\prime }-v_{\min }^{\prime }},$$

$$Q^2>1.5{\mathrm{GeV}}^2,$$

$$-t/Q^2<0.2.$$

$$\mathrm{Im}\mathcal{H}(\xi ,t)\stackrel{\mathrm{LO}}{=}\pi \sum _q{e}_q^2{H}^{q(+)}(x=\xi ,\xi ,t),$$

Full text

††journal: Eur. Phys. J. C

1 Introduction

An intense experimental activity has been dedicated over the last twenty years to the measurements of observables towards a 3D description of the nucleon. Most of the knowledge about the 3D quark and gluon structure is embodied in Generalized Parton Distributions (GPDs) and Transverse Momentum Dependent parton distribution functions (TMDs). The continuing efforts in understanding and determining GPDs and TMDs are partly driven by the quest for nucleon tomography in mixed position-momentum space or pure momentum space (see e.g.\xspaceRef. Diehl:2015uka ; Bacchetta:2016ccz and refs. therein), the access to the nucleon energy-momentum tensor Ji:1996ek ; Ji:1996nm and the description of the mechanical properties (radial and tangential pressures, energy density, surface tension, etc.\xspace) of the nucleon Polyakov:2002wz ; Polyakov:2002yz ; Polyakov:2018zvc ; Lorce:2018egm .

Relations between measurements and GPDs or TMDs are derived through factorization theorems, established at all order in QCD perturbation theory (pQCD). In particular, this connection can be brought under good theoretical control with a careful and systematic check of various pQCD assumptions. This task requires a sophisticated and modular computing machinery with fitting features. In the context of GPDs, the open-source PARTONS framework Berthou:2015oaw has been designed and publicly released to fulfill the needs of the experimental and theoretical GPD communities.

From a theoretical point of view, Deeply Virtual Compton Scattering (DVCS) is the cleanest channel to access GPDs. Recent status of the phenomenology of this channel can be found in Ref. Kumericki:2016ehc and of the related collected data sets in Ref. dHose:2016mda . The DVCS cross section is described in terms of Compton Form Factors (CFFs), which are convolutions of GPDs with coefficient functions computed order by order in pQCD. While GPDs are in principle renormalization and factorization scale-dependent, CFFs do not depend on any arbitrary scales, and are the off-forward analogs of the classical structure functions of Deeply Inelastic Scattering (DIS) – with GPDs being the off-forward extensions of the Parton Distribution Functions (PDFs) extracted from DIS.

Global GPD fits over the whole DVCS kinematic domain, from the glue to the valence region, have not been achieved yet. However, a great deal is already known about CFFs, in particular in the valence region. This is indeed precious knowledge: due to the scale dependence, GPDs can only be promoted to the status of quasi-observables thanks to universality (the same GPDs should allow the simultaneous analysis of several independent exclusive processes), while CFFs are genuine observables. Hence measuring CFFs is an interesting task per se. Moreover, having a unified description of all existing DVCS measurements in terms of CFFs is important in view of the definition of future experimental programmes: one has to select an observable to be measured, specify the kinematic region, and evaluate the required uncertainty in order to make substantial progress in the understanding of the DVCS process. This is particularly relevant for future experimental campaigns at Jefferson Lab or for Electron-Ion Collider (EIC) design studies.

As mentioned in Ref. Kumericki:2016ehc , the attempts to determine CFFs from experimental data mostly fall into two categories:

Local fits

CFFs are independently determined from measurements between different kinematic bins. This amounts to a sampling of CFFs over the probed kinematic domain. The model-dependence of the result is low, but most of the time, the problem is ill-posed by lack of uniqueness.

Global fits

The free coefficients of a CFF parameterization are matched to experimental data. Kinematic bins are no more treated independently. Interpolating between measurements of the same observable on neighboring kinematic bins is feasible. Extrapolating out of the domain constrained by measurements becomes possible, paving the way for impact studies. However the estimation of the systematic uncertainty associated to the choice of a parameterization is an extremely difficult task, and this is a limit to potential impact studies.

A possible solution to the tension between the advantages and drawbacks of these two approaches may be found in Artificial Neural Networks (ANNs) Hassoun:1995:FAN:526717 . While being of common use for PDFs today, their role in GPD-related physics has been assessed (to the best of our knowledge) only in the pioneering work of Kumerički, Müller and Schäfer Kumericki:2011rz . ANNs offer the joint promise of a great flexibility and of a common description of data sets in different kinematic regions. This difficulty of this last point should not be underestimated. In our previous fit of CFFs from almost all existing DVCS data Moutarde:2018kwr , our physically-motivated choice of CFF parameterizations restricted the scope of our study to the quark sector. There is a need for flexible parameterizations while escaping at the same time the curse of dimensionality, and in this respect, ANNs are quite appealing solutions.

We present here the first global fit of most existing DVCS measurements in terms of CFFs described by neural networks. The paper is organized as follows. In Sec. 2 we briefly introduce the DVCS channel, define CFFs and stress their role in dispersion relations. Sec. 3 is a reminder of ANN technique and of the genetic algorithm. Our implementation of those techniques and the main elements of our analysis are detailed in Sec. 4. We point out the experimental data being used in this fit in Sec. 5, and outline our results with a focus on the DVCS subtraction constant in Sec. 6. This discussion is of major contemporary interest since it is related to the measurability of the distribution of pressure forces in the proton Polyakov:2018zvc ; Lorce:2018egm ; Burkert:2018bqq ; Kumericki:2019ddg . At last we summarize in Sec. 7.

2 Theory framework

The golden channel in the GPD extraction programme is the leptoproduction of a single real photon off a nucleon:

[TABLE]

Here, the symbols between parentheses denote the four-momentum vectors of the lepton, $l$, the nucleon, $N$, and of the real photon, $\gamma$. The amplitude for this process, $\mathcal{T}$, is given by a sum of amplitudes for two processes having the same initial and final states: $\mathcal{T}^{\mathrm{BH}}$ for the purely electromagnetic Bethe-Heitler process and $\mathcal{T}^{\mathrm{DVCS}}$ for the hadronic DVCS process, such as:

[TABLE]

The Bethe-Heitler part can be expressed with a great precision in terms of the nucleon electromagnetic form factors. The DVCS part is generally parameterized by twelve helicity amplitudes Kroll:1995pv , or equivalently twelve CFFs. However, in this analysis we restrict ourselves only to four of them, which can be related to the leading-twist, chiral-even GPDs. Those CFFs are usually denoted by $\mathcal{H}$, $\mathcal{E}$, $\widetilde{\mathcal{H}}$ and $\widetilde{\mathcal{E}}$, and nowadays they are the most extensively studied ones in the context of GPD phenomenology. The exploration of other CFFs, which can be related to either higher-twist or chiral-odd GPDs, suffers from the sparsity of data collected in the kinematic domain where the factorization theorem applies.

The cross section for single photon production is usually expressed in terms of four variables: (i) the usual Bjorken variable, $x_{\mathrm{Bj}}=Q^{2}/\left(2p\cdot\left(k-k^{\prime}\right)\right)$, (ii) the square of the four-momentum transfer to the target nucleon, $t=(p^{\prime}-p)^{2}$, (iii) the virtuality of the exchanged photon, $Q^{2}=-q^{2}=-(k^{\prime}-k)^{2}$, and (iv) the azimuthal angle between the leptonic (spanned by incoming and outgoing lepton momenta) and production (spanned by virtual and outgoing photon momenta) planes. In the case of a transversely polarized target one also introduces the angle $\phi_{S}$ between the leptonic plane and the nucleon polarization vector. It is also convenient to exchange the usual Bjorken variable for the generalized one:

[TABLE]

which in the case of DVCS is approximately equal to the skewness variable:

[TABLE]

In this paper a set of lengthy formulae relating CFFs with cross sections is omitted for brevity, but those can be easily found e.g.\xspacein Ref. Belitsky:2012ch . The unpublished analytic expressions of BH and DVCS amplitudes by Guichon and Vanderhaeghen are used in this analysis, as implemented and publicly available in the PARTONS framework Berthou:2015oaw . Let us note at this moment that the leptoproduction of a single real photon is sensitive to the beam charge, but also to the beam and target polarization states. This gives rise to exclusive measurements performed with variety of experimental setups, allowing for a combination of cross sections to probe specific sub-processes (BH, DVCS or the interference between them) or combinations of CFFs. Experimental data used in this analysis will be introduced in Sec. 5.

Our goal here is a global analysis of DVCS data avoiding any model-dependency whatsoever. We separately extract the real and imaginary parts of the four aforementioned CFFs in the three-dimensional phase-space of ($\xi$, $t$, $Q^{2}$). During the extraction we do not exploit the fact that the real and imaginary parts of a single CFF can be related together by a fixed-$t$ dispersion relation (see Ref. Diehl:2007jb for a study at all orders in pQCD):

[TABLE]

Here, $\mathcal{G}\in\{\mathcal{H},\mathcal{E},\widetilde{\mathcal{H}},\widetilde{\mathcal{E}}\}$ denotes a single CFF with $t$ and $Q^{2}$ dependencies suppressed for brevity, one has minus sign for $\mathcal{G}\in\{\mathcal{H},\mathcal{E}\}$ and plus sign for $\mathcal{G}\in\{\widetilde{\mathcal{H}},\widetilde{\mathcal{E}}\}$ in the square brackets, and where $C_{G}$ is a subtraction constant associated to the corresponding CFF (and GPD), where:

[TABLE]

The dispersion relation is only used afterwards to determine the subtraction constant from the extracted CFF parameterizations. This quantity has on its own important interpretation in terms of strong forces acting on partons inside the nucleon as shown in the recent review Ref. Polyakov:2018zvc . We also use the requirement of having $C_{G}$ independent on $\xi$ as a consistency check of our analysis.

Although hadron tomography Burkardt:2000za ; Burkardt:2002hr ; Burkardt:2004bv is a very important motivation for the phenomenology of the DVCS process, in the present analysis we restrict ourselves to (quasi) model-independent statements. Therefore we will only sketch what we can foresee from the obtained results and leave detailed quantitative tomographic interpretation for future studies.

3 Methodology

To make this paper self-contained, we remind here some elements of the theory and practice of ANNs and genetic algorithms, and introduce the associated terminology.

3.1 Artificial neural networks

ANNs Hassoun:1995:FAN:526717 are systems processing information. What distinguishes those systems from other information techniques is a unique inspiration from Nature, namely by biological neural networks making human and animal brains. All neural networks have a similar structure made out of simple, but highly connected elements called neurons. In this work we will exclusively focus on feed-forward neural networks, but many other types of ANNs exist. Their usefulness in the study of hard exclusive reactions has not been assessed yet.

An example of a typical feed-forward neural network structure is shown in Fig. 1a. The data are processed layer-by-layer. One can distinguish three type of layers: (i) one input layer, storing and distributing input information, being the first layer of the network, (ii) a number of hidden layers processing information and (iii) one output layer, aggregating and returning output information, being the last layer of the network. The number of hidden layers is a parameter of the network architecture.

In machine learning neurons share a similar structure shown in Fig. 1b. A single neuron processes information by evaluating two functions: (i) combination function, $f_{\Sigma,j}^{i}(\mathaccent 382{\beta}_{j}^{i},\mathaccent 382{v}^{i-1})$, and (ii) activation function, $f_{\includegraphics{plots/symbol/step.pdf},j}^{i}(x)$, so that the overall neuron’s output, $v_{j}^{i}$, is given by:

[TABLE]

Here, the index $i$ identifies a specific layer of the neural network, while $j$ denotes a single neuron in that layer. The vector $\mathaccent 382{\beta}_{j}^{i}$ contains a so-called bias, $\beta_{j0}^{i}$, and all weights associated to the input nodes, $\beta_{jk}^{i}$, with $k=1,2,\dots$ The vector $\mathaccent 382{v}^{i-1}$ contains values returned by neurons in the previous layer, $v_{k}^{i-1}$, where again $k=1,2,\dots$ The role of the combination function is the aggregation of input information to a single value to be processed by the activation function. A typical combination function is based on the dot product of two vectors:

[TABLE]

The activation function gives an answer returned by the neuron. It can be a numeric value, but it can be also a binary answer (“yes” or “no”), which is typically used in classification applications. Many types of activation functions exist. In this work we use the hyperbolic tangent function for neurons in the hidden layers,

[TABLE]

and the linear function (identity) for neurons in the output layer,

[TABLE]

We illustrate both functions with Fig. 2. For the hyperbolic tangent function one can spot an active part near $-1<x<1$ where the signal is processed almost linearly. Outside this range a saturation appears, so that the input signal is blocked, i.e.\xspaceit is not processed by a neuron in an effective way. Those features put extra constraints on input data to be processed by the network. In particular, those data must be normalized to the range of $(-1,1)$, so they can be processed by the first layer of the network.

Biases and weights are free parameters of neural networks. Instead of predefining those parameters a training is done, and only during that training biases and weights are fixed to accommodate the network’s answer to the stated problem. Such training is done with a training sample, which should be as representative as possible for the problem under consideration. Both input and output information are known for the training sample.

In a wide variety of applications neural networks are trained with the back-propagation algorithm goodfellow2016deep . This algorithm minimizes the so-called cost function (e.g.\xspaceroot mean square error, RMSE) against the training sample. Such a minimization usually is straightforward, as typically the cost function is analytically differentiable with respect to each parameter of the network, so that its gradient is known in the space of those parameters. However, there are cases where the back-propagation algorithm (or similar) cannot be used. Those are in particular problems with the differentiation of the cost function being practically difficult. In those cases one can evaluate the gradient numerically. Eventually, one can use the genetic algorithm Mitchell:1998:IGA:522098 , which we will describe in the next paragraph. It is a popular replacement for the back-propagation algorithm in the context of machine learning. In particular, the genetic algorithm is used in this analysis, as a differentiation of cross sections and other observables with respect to the free parameters of CFF parameterizations is not straightforward.

3.2 Genetic algorithm

The genetic algorithm Mitchell:1998:IGA:522098 is a heuristic technique used for search and optimization purposes. It is applicable in particular whenever standard minimization techniques fail because of the complexity of the problems under consideration. Because of that, but also because of a unique ability to avoid local minima, the genetic algorithm is commonly used in physics, mathematics, engineering and manufacturing. It is frequently used for the learning of neural networks, where it supervises the search of biases and weights against the training sample.

The genetic algorithm mimics mechanisms of reproduction, natural selection and mutation, all being the cornerstones of the evolution process known from Nature. It utilizes nomenclature inspired by the natural evolution process. A population is a set of possible solutions to the stated problem. A single solution is referred to as a candidate and it is characterized by a set of genes making up a genotype. In the following, we will thus speak equivalently of candidates and genotypes. The usability of single candidates in solving problems is examined by a fitness function, whose argument is the genotype. To help the reader understand this nomenclature we explain it with an example of a multidimensional fit to experimental data. In such a fit a single fitted parameter may be identified as one gene, while a complete set of those parameters as the genotype. The ability to describe experimental data by a given set of parameters is typically examined with a $\chi^{2}$ function, which in this case we will refer to as the fitness function. A set of possible solutions makes the population. We point out that in general the implementation of introduced terms strongly depends on the problem under consideration.

The genetic algorithm is iterative, that is, it consists of several steps being consecutively executed in a loop, until one of well defined conditions, referred to as stopping criteria, is fulfilled. The algorithm is described with the help of Fig. 3. (i) In the first step an initial population is generated. A number of candidates making up this initial population is a parameter of this algorithm. The creation of a new candidate consists for instance of the random generation of its genotype with gene values laying in defined ranges. (ii) Each candidate is evaluated by acting on its genotype with the fitness function. The resulting fitness score indicates how good the candidate is accommodated to the stated problem. Then, all the candidates are sorted according to those scores. The one having the most preferable score is considered as the best one. (iii) Stopping criteria are checked and the algorithm is stopped if one of them is fulfilled. In such a case the genotype of the best candidate is considered as the solution. A possible stopping criterion is for instance a threshold on the fitness score of the best candidate or a threshold on the number of iteration epochs that has passed already. (iv) A fraction of candidates being identified by the worst fitness scores is removed from the population and forgotten. Their place in the population is taken by new candidates carrying genes of the remaining candidates. That is, two or more remaining candidates pass a fraction of their genotypes to a new candidate in a hope of creating even a better adapted individual. This step is known as the reproduction or cross-over. (v) A random fraction of genes in the population is modified. This modification may consist of the random generation of a new gene value within a specified range (as during the initialization) or it can be a small shift with respect to the original gene value. The purpose of this step, which is known as the mutation, is to keep the population diverse. The mutation probability cannot be too high, as in such a case the genetic algorithm becomes a random search. The mutation ends a single epoch, so that the algorithm proceeds with the evaluation as the next step.

3.3 Regularization

Two extreme cases of failure in the learning procedure are known as under-fitting and over-fitting. A compact explanation of those effects goes as follows. Under-fitting results in a lack of success in describing the training sample, as well as any other sample, while over-fitting is having a model being too flexible, which is not able to describe data other than the training sample. Over-fitting usually happens either because the model memorizes a noise present in the training sample, or because it becomes too flexible in a range between training points. The extrapolation does not count here, as neural networks typically become pliant outside ranges covered by the training sample. Both effects, together with the case of proper training, are graphically illustrated in Fig. 4. At this point we should stress that both under-fitting and over-fitting are not specific features of neural network training, but are known and common problems met in data fitting whenever the number of free parameters becomes large.

Regularization is the introduction of additional information to improve the resolution of badly behaved problems. In the case of neural networks, the regularization helps to avoid over-fitting. In the following we describe the early stopping regularization that is used in this analysis, however many other types of regularization exist. We note that in this analysis we have found regularization techniques using a penalty term in the cost function impractical. Such techniques, like Lasso or Ridge regressions goodfellow2016deep , make the values of neural network weights small and as a consequence replicas similar to each other (the replication method will be introduced in Sec. 3.4). While this may be seen as a feature in some analyses, in this one it does not lead to a reliable estimation of uncertainties, as replicas are not spread significantly in a domain that is not sufficiently constrained by data.

Our procedure of dealing with over-fitting makes the training supervised. It means that we randomly divide our initial sample into two independent sub-samples, which in the following we will refer to as training sample and validation sample. While the purpose of the training sample remains the same, i.e.\xspaceit is used for the training, the validation sample is exclusively used for the detection of over-fitting.

As a principle, too short a learning leads to under-fitting, while a too long one is responsible for over-fitting, unless the network’s architecture prevents it. Therefore, between those two extremes there must be an optimal epoch in the training corresponding to a most favorable solution. The search of such an optimal epoch is the goal of the early stopping technique. The principle of this method is to monitor the validation sample with the cost function, $f_{\mathrm{cost}}$, and to stop the training whenever the value of this function starts to grow. This idea is illustrated in Fig. 5. In practice, such a growth is difficult to be detected due to fluctuations of the cost function values. Therefore, various quantities are used for its detection. In this analysis we monitor the following quantity evaluated for the test sample, which is inspired by Ref. Prechelt2012 :

[TABLE]

Here, $t_{i}$ and $t_{j}$ denote different epochs in the training. In the numerator of Eq. (12) one has the value of the cost function for the epoch under consideration, and in the denominator the minimal value of this function obtained in all previous epochs. To minimize the impact of random fluctuations on $f_{\mathrm{stop}}(t_{i})$ the cost function is averaged using the running average technique:

[TABLE]

where $n$ is the length of a strip of numbers that is used for the averaging. In this analysis $n=25$. With those definitions $f_{\mathrm{stop}}(t_{i})>0$ indicates a need of stopping the training. However, as fluctuations of the cost function cannot be completely avoided we stop the training only when $f_{\mathrm{stop}}(t_{i})>0$ for more than one hundred consecutive iterations. We take the set of parameters obtained in the first iteration of such a sequence as the result of our training procedure. In rare cases the stopping criterion is not fulfilled during the whole training procedure. In those cases the result is obtained from the last allowed epoch of this procedure, $t_{\max}=1000$.

3.4 Replica method

As in our last analysis Moutarde:2018kwr , the replica method is used to propagate uncertainties coming from experimental data to CFF parameterizations. In addition to our nominal extraction of CFFs, the fit is repeated one hundred times, each time independently modifying the central values of experimental points using the following prescription:

[TABLE]

Here, $v_{ij}$ is the measured value associated to the experimental point $j$ coming from the data set $i$. It is linked to statistical, ${\Delta}_{ij}^{\mathrm{stat}}$, systematic, ${\Delta}_{ij}^{\mathrm{sys}}$, and normalization, ${\Delta}_{i}^{\mathrm{norm}}$, uncertainties. The latter one is correlated across bins of data, like for instance uncertainties related to beam and target polarization measurements. The total uncertainty, which is used to evaluate $\chi^{2}$ function utilized in this analysis to perform the fit, is given by:

[TABLE]

The generator of random numbers following a specified normal distribution, $f(x|\mu,\sigma)$, is denoted by $\mathrm{rnd}_{k}(\mu,\sigma)$, where $k$ both identifies the replica and is a unique random seed.

4 Implementation

4.1 Neural network architecture

In this subsection we describe the architecture of the used neural networks. In particular a determination of the number of hidden neurons in those networks is described. We start with a general statement: any neural network is a container to store information and as such it cannot be too small to serve that purpose. It can be larger than needed though, however such neural networks require a careful regularization to avoid over-fitting and their training is more time-consuming. We note that usually the optimal network’s architecture is selected prior the training, however there exist regularization techniques like dropout Srivastava:2014:DSW:2627435.2670313 , which alter an initial architecture during the training to obtain best results. In those techniques a number of neurons is either added or dropped as the training proceeds.

We extract four CFFs, each one being a complex quantity having both real and imaginary parts. Although the real and imaginary parts of the a given CFF can be connected together via a dispersion relation (with the subtraction constant being involved), in the extraction presented in this paper we keep them fully independent quantities. The dispersion relation is only used for an interpretation of the obtained results, which is presented in Sec. 6. Therefore eight independent neural networks are used.

The architecture of a single network used in this analysis is shown in Fig. 6. It consists of three input neurons, one hidden layer containing six neurons, and one output neuron. Both input and output variables of the network are linearized and normalized (as much as possible). The linearization is achieved with a logarithmic projection and $\xi$ pre-factors:

[TABLE]

where the prime symbol is used to distinguish between original and linearized variables. For the normalization we utilize the min-max method:

[TABLE]

where $v^{\prime}$ and $v^{\prime\prime}$ is a given linearized variable before and after the normalization, respectively, and where $v_{\mathrm{min}}^{\prime}$ and $v_{\mathrm{max}}^{\prime}$ are normalization parameters specified in Table 1. The min-max values used for the normalization of $\xi^{\prime}$, $t^{\prime}$ and ${Q^{2}}^{\prime}$ have been selected to cover the target phase-space, while those for $\mathrm{Re}\mathcal{G}^{\prime}$ and $\mathrm{Im}\mathcal{G}^{\prime}$ roughly correspond to $v_{\mathrm{min,GK}}^{\prime}-\Delta v_{\mathrm{GK}}^{\prime}$ and $v_{\mathrm{max,GK}}^{\prime}+\Delta v_{\mathrm{GK}}^{\prime}$ values, where $\Delta v_{\mathrm{GK}}^{\prime}=v_{\mathrm{max,GK}}^{\prime}-v_{\mathrm{min,GK}}^{\prime}$. Here, $v_{\mathrm{min,GK}}^{\prime}$ and $v_{\mathrm{max,GK}}^{\prime}$ are min-max values found in a poll of CFFs evaluated from the Goloskokov-Kroll (GK) GPD model Goloskokov:2005sd ; Goloskokov:2007nt ; Goloskokov:2009ia for the experimental data kinematics that are used in this analysis, see Sec. 5. Both linearization and normalization significantly improve the performance of ANNs. We point out that the ranges specified in Table 1 are not absolute, i.e.\xspacethe networks may still reasonably describe data covering the exterior of domains defined by $v_{\mathrm{min}}^{\prime}$ and $v_{\mathrm{max}}^{\prime}$ values. In particular in this analysis the activation function of output neurons is set to the identity and thus does not show any saturation effects.

The number of neurons in the hidden layer is determined with a benchmark sample made out of one thousand CFF points evaluated with the GK GPD model. Those CFF points are randomly generated in a broad range of $10^{-6}<\xi<1$, $0~{}\mathrm{GeV}^{2}<|t|<1~{}\mathrm{GeV}^{2}$ and $1~{}\mathrm{GeV}^{2}<Q^{2}<100~{}\mathrm{GeV}^{2}$. We have checked how many neurons in the hidden layer are needed to provide a good description of the benchmark sample. This test is performed with the FANN library nissen03 for the neural network implementation and with the quickprop training algorithm fahlman:faster available in that library. The post-training RMSE for CFF $\mathcal{H}$ is shown in Fig. 7 as a function of the number of hidden neurons. From this figure one may conclude that the addition of a new neuron to the network made already out of six neurons does not significantly improve the performance. We note that a high number of neurons slows down the training and makes the regularization more difficult. That is why in this analysis the number of hidden neurons in each network is set to six, which we consider to be a sufficient number.

4.2 Initial fit

Without any improvement the fit either does not converge or converges very slowly. This is due to a large number of parameters to be constrained in the training and a widely opened phase-space for CFF values, cf.\xspaceTable 1. To overcome this difficulty the network is initially trained on values obtained in a local extraction of CFFs. That is, CFF values multiplied by $\xi$ are independently extracted in each kinematic bin of data, analogously as in Ref. Dupre:2016mai ; Dupre:2017hfs ; Burkert:2018bqq , and then used in a training based on the back-propagation algorithm.

The local extraction procedure includes all four CFFs, so effectively eight values per each kinematic bin of $(x_{\mathrm{Bj}},t,Q^{2})$ are fitted. This may make a single fit underconstrained, in particular if a given kinematic bin is populated only by a single experimental point. To overcome this difficulty it is imperative to set sensible ranges in which the fitted values are allowed to vary. In this analysis those ranges are set to ($5\times v_{\mathrm{min}}^{\prime},5\times v_{\mathrm{max}}^{\prime}$), where the values of $v_{\mathrm{min}}^{\prime}$ and $v_{\mathrm{max}}^{\prime}$ are specified in Table 1. The ranges are wider than those used in the construction of our neural networks, which prevents from introducing a bias on the final extraction of CFF parameterizations.

The initial fit reduces the $\chi^{2}$ value per single data point from $\mathcal{O}(10^{6})$ to $\mathcal{O}(10)$. The further minimization, which reduces the same quantity to $\mathcal{O}(1)$, is done by our genetic algorithm minimizer, where values of weights and biases obtained in the initial fit are used as starting parameters.

4.3 Genetic algorithm

The parameters of our genetic algorithm minimizer are as follows. The gene values (weights and biases) may vary between $-10$ and $10$. The population consists of $1000$ candidates and the fraction of candidates left after the selection process (survivors) is $30\%$. We consider two types of mutation and each of them occurs with the same probability of $0.1\%$ per gene. The mutation type-$A$ is a random generation of gene values from the range of $-10$ and $10$. This type of mutation keeps the population diverse. The mutation type-$B$ is a small shift of gene values. This shift is randomly generated from a distribution of gene values in the population. This distribution is made for the modified gene. The mutation type-$B$ provides a fine-tuning of fitted parameters. The fraction of $20\%$ and $10\%$ best candidates in the population is resistant to the mutation type-$A$ and type-$B$, respectively, so the best results are not destroyed by the mutation process.

We demonstrate the performance of our genetic algorithm minimizer with Fig. 8, where values of a single gene are plotted against the training epoch. Intensity in population equal $1$ means that all candidates in the population share the same gene value. Low values of this quantity indicate that only few or none candidates share the same gene value. One can see few features in this plot being typical for a performance of genetic algorithms: (i) the diversity of values is large at the beginning of minimization, which means that the whole phase-space available is equally scanned to find a cost function minimum, (ii) usually, few such local minima are found at the beginning of minimization and simultaneously explored by the algorithm for a number of iterations, (iii) close to the actual minimum the algorithm concentrates on its neighborhood, trying to perform a fine-tuning of fitted parameters, (iv) even at this stage of minimization, the whole phase-space is homogeneously scanned due to the mutation process to eventually find a new minimum and to avoid the convergence to a local one.

4.4 Regularization

A typical distribution of the cost function as a function of the training epoch is shown in Fig. 9. In addition, we show the distribution of the $f_{\mathrm{stop}}$ quantity that we use to detect over-fitting, see Eq. (12). With our stopping criterion: $f_{\mathrm{stop}}>0$ for more than one hundred iterations, for this example we consider the solution obtained in the epoch $382$ as the valid one and we use this solution in the further analysis.

4.5 Feasibility test

In order to check our procedure a feasibility test was performed. In this test the extraction of CFF parameterizations is done on a sample of pseudo-data generated in the leading-order formalism with the GK GPD model Goloskokov:2005sd ; Goloskokov:2007nt ; Goloskokov:2009ia . The generation of pseudo-data is straightforward: for each experimental data point, a single point of pseudo-data is generated. This single point of pseudo-data is generated with the same type (differential cross section, or beam spin asymmetry, or beam charge asymmetry, etc.\xspace), and at the same kinematics as the original point. The obtained observable value is then smeared according to the uncertainties of the original point. As a result of this procedure, pseudo-data are faithfully generated, corresponding to the experimental data; however “input” CFFs are known for this sample. These pseudo-data are used in the same procedure of CFF extraction as the genuine experimental data.

The outcome of this feasibility test is in general positive. The value of the $\chi^{2}$ function per single data point is $1.13$, or $1.06$ for data other than HERA and COMPASS points (see Sec. 6 for an explanation of this behavior). The agreement between the CFF parameterizations obtained from the pseudo-data and the GK GPD model used to generate those data is demonstrated with Fig. 10, where the imaginary part of CFF $\mathcal{H}$ and the model are compared as a function of $\xi$ for example kinematics of $t=-0.3~{}\mathrm{GeV}^{2}$ and $Q^{2}=2~{}\mathrm{GeV}^{2}$. The curve representing the model stays within the uncertainty band of the extracted parameterization. The only exception is the region of $\xi\approx 1$, which is not covered by pseudo-data and where the model rapidly goes to zero. However, in this region the difference between the model and the central value of the extracted parameterization is still reasonable (e.g.\xspace$1.23\sigma$ for $\xi=0.95$).

5 Experimental data

Table 2 summarizes DVCS data used in this analysis. For the explanation of symbols used to distinguish between observable types see Ref. Moutarde:2018kwr . Only proton data are used. Available neutron data are sparse ones, but in principle they may be used to attempt a flavor separation. This is however beyond the scope of the presented analysis.

As in our previous analysis Moutarde:2018kwr , recent Hall A data Defurne:2015kxq ; Defurne:2017paw for cross sections measured with unpolarized beam and target, $d^{4}\sigma_{UU}^{-}$, are not used in the final extraction of CFF information. Again, only the corresponding differences of cross sections probing longitudinally polarized beam, $\Delta d^{4}\sigma_{LU}^{-}$, are used. We will elaborate on the inclusion of Hall A data in Sec. 6.

The difference with respect to our previous analysis Moutarde:2018kwr in terms of used data comes from the inclusion of low-$x_{\mathrm{Bj}}$ sets. For COMPASS, instead of a single point for the measurement of the $t$-slope $b$, four points for the cross section measurement, $d^{3}\sigma_{UU}^{\pm}$, are used. The measurement of both slope and cross sections is reported in Ref. Akhunzyanov:2018nut . The cross sections provide richer information than the slope itself. In addition, the slope measurement relies on the assumption of the $t$-dependence being exponential, which can be avoided analyzing the cross section points. However, the COMPASS cross sections are measured in broad kinematic bins and because of that they cannot be compared to models evaluated at the average kinematics of events corresponding to those bins. That is why in this analysis, the COMPASS points are compared to respective integrals of cross sections. Those multidimensional integrals are time consuming and significantly extend the computing time of our fits.

HERA data for DVCS Chekanov:2008vy ; Aktas:2005ty ; Aaron:2009ac are included in this analysis. As we will demonstrate in Sec. 6, those data provide important constraints on CFF parameterizations in the low-$x_{\mathrm{Bj}}$ region. However, the sparsity of data covering this region makes our extraction of CFFs difficult, as will be explained in Sec. 6.

We apply two kinematic cuts on experimental data111There is a typo in Eq. (86) of our last publication Moutarde:2018kwr . Instead of $-t/Q^{2}<0.25$ there should be $-t/Q^{2}<0.2$ , as in Eq. (23) here.:

[TABLE]

The purpose of those cuts is the restriction of the phase-space covered by experimental data to the deeply virtual region, where one can rely on the factorization between GPDs and the hard scattering kernel. In principle, those cuts can be avoided in the extractions of amplitudes. However, we keep them here to allow an interpretation of the extracted CFFs in terms of GPDs. In addition, the cuts are applied to keep a correspondence with our previous analysis Moutarde:2018kwr , allowing a straightforward comparison. Let us note at this point, that the neural network approach developed in this work provides an easy way for the addition of CFFs identified with higher-twist contributions. With those contributions included one could relax the cuts (22) and (23), allowing more data to be included.

In total, in this analysis $2624$ points are used out of $3996$ available in all considered data sets. The coverage of phase-space by those data in the sections of $(x_{\mathrm{Bj}},Q^{2})$ and $(x_{\mathrm{Bj}},-t/Q^{2})$ is shown in Fig. 11. With those plots one can easily conclude about the actual coverage by various experiments, but also identify parts of phase-space not probed at all by existing data.

6 Results

6.1 Performance

The value of the $\chi^{2}$ function obtained after the training of our neural network system with $2624$ experimental points is $2243.5$ for the central replica. It gives the average deviation of that system’s answer from experimental data $2243.5/2624\approx 0.85$ per single data point. This value per experimental data set is given in Table 3. The number of free parameters, that is the number of weights and biases in all networks, is $248$. There is not much sense of studying the goodness of fit with this information taken into account, as we know a priori that the size of a single network may be too large, and that is why we are using a regularization. In other words, we are able to increase the number of neurons in our network, and therefore the number of weight and biases, keeping at the same time the same precision of data description.

Figures 12 - 15 provide a straightforward comparison between our fit and the selected data sets. As indicated by the $\chi^{2}$ values summarized in Table 3, we are able to describe the data well within a single phenomenological framework based on the neural network approach. This includes data ranging from HERA to JLab kinematics. Beyond the results of this analysis, predictions coming from the GK Goloskokov:2005sd ; Goloskokov:2007nt ; Goloskokov:2009ia and VGG Vanderhaeghen:1998uc ; Vanderhaeghen:1999xj ; Goeke:2001tz ; Guidal:2004nd GPD models are also shown. Those two models originate from the exploratory phase of GPD studies and are able to describe only a general behavior followed by experimental data. This confirms the need for new GPD models constrained in global analyses, preferably multi-channel ones.

As has been already mentioned in Sec. 5, Hall A cross sections are not included in the nominal data set that is used in this analysis. We point out, that in our last analysis Moutarde:2018kwr we were not able to describe those cross sections with the proposed Ansatz and those data were excluded there as well. The question is: can those data be described in this analysis, where a flexible CFF parameterization based on the neural network approach is used? The answer consists of two parts: (i) the parameterizations of CFFs obtained in this analysis from the nominal data set do not have the predictive power to describe Hall A cross sections. The value of the $\chi^{2}$ function for those cross sections is $5916.6$, which for $594$ points gives the reduced value of $9.96$. The bad description is also seen by eye in Fig. 13. (ii) the inclusion of Hall A cross sections in the extraction of CFFs significantly improves the $\chi^{2}$ evaluated for this subset, from aforementioned $9.96$ to $1.12$ per single data point. However, such an inclusion makes the $\chi^{2}$ for the CLAS cross sections worse, from $0.66$ reported in Table 3 to $0.83$ per single data point. Taking both observations into account we found the situation unclear at this moment. Hall A cross sections will be investigated further in a future analysis with higher-twist contributions taken into account.

The inclusion of low-$x_{\mathrm{Bj}}$ experimental data has significantly extended the coverage of phase-space. However, the sparsity of those data creates problems in the supervised training that we utilize in this analysis, see Sec. 3.3. Namely, the random division of the available experimental data into training and test subsets causes an insufficient coverage of some parts of the phase-space by either the training or the test points. This may lead to over-fitting, if a given part of the phase-space is only covered by training points, or under-fitting, if kinematics is only covered by test points. The effect increases the spread of replicas and because of that gives an additional contribution to the estimated uncertainties. The goodness of fit for those data is also worse as shown in Table 3. Foreseen data coming from electron-ion facilities will improve this situation.

6.2 Compton form factors

The parameterizations of CFFs are shown in Figs. 16 - 19 for example kinematics of $t=-0.3~{}\mathrm{GeV}^{2}$ and $Q^{2}=2~{}\mathrm{GeV}^{2}$ as a function of $\xi$. As expected the data provide the best constraints on $\mathrm{Im}\mathcal{H}$, and some on $\mathrm{Re}\mathcal{H}$, $\mathrm{Im}\widetilde{\mathcal{H}}$ and $\mathrm{Re}\widetilde{\mathcal{E}}$. Other CFFs are poorly constrained by the available data, in particular $\mathcal{E}$ related to GPD $E$, being of a great importance for the study of parton densities in a transversely polarized proton and the determination of the orbital angular momentum through Ji’s sum rule.

The inclusion of HERA and COMPASS data in a global extraction of CFFs is not trivial, but those data provide important constraints in the low and intermediate range of $x_{\mathrm{Bj}}$. We demonstrate it with the example of Fig. 20, where the imaginary part of the CFF $\mathcal{H}$ is shown without the inclusion of low-$x_{\mathrm{Bj}}$ data. We point out that the gap between the collider and fixed target experiments seen in the coverage of the $(x_{\mathrm{Bj}},Q^{2})$ phase-space in Fig. 11 is expected to be filled by future experiments at expected electron-ion collider facilities. The precision of foreseen data should allow for a precise phenomenology in that domain.

The $Q^{2}$ evolution of the CFF $\mathcal{H}$ for $\xi=0.002$ and $\xi=0.2$ is shown in Fig. 21. One can note a rather mild $Q^{2}$-dependence followed by the extracted CFF parameterizations, consistent with the expected logarithmic behavior.

6.3 Nucleon tomography

In Fig. 22 we show our results for the slope $b$ of the DVCS cross section described by a single exponential function, $d^{3}\sigma_{UU}\propto\exp(bt)$. This slope is evaluated as indicated in our previous analysis Moutarde:2018kwr . It can be converted into a transverse extension of partons under assumptions that are also specified in Ref. Moutarde:2018kwr . The uncertainties on the extracted values of $b$ are larger than expected, mainly because of the aforementioned problems with the supervised training caused by the sparsity of the low-$x_{\mathrm{Bj}}$ data. In addition, without an explicit assumption about the exponential $t$-behavior of $\mathrm{Im}\mathcal{H}$, the estimation of $b$ from pliant replicas gets additional uncertainties.

The extraction of tomography information from CFFs is also possible in the high-$x_{\mathrm{Bj}}$ range, see for instance Refs. Dupre:2016mai ; Dupre:2017hfs . However, this requires a “de-skewing” of $\mathrm{Im}\mathcal{H}$, i.e.\xspaceone needs to evaluate the GPD $H^{q}(x,0,t)$ from:

[TABLE]

where the sum runs over quark flavors $q$, $e_{q}$ is the electric charge of a specific quark flavor in units of the positron charge $e$ and where $H^{q(+)}(x,\xi,t)=H^{q}(x,\xi,t)-H^{q}(-x,\xi,t)$. Because of the model-dependency of this procedure we refrain from doing it here. We point out that a straightforward, however still model-dependent, access to the nucleon tomography is a feature of the CFF Ansätze proposed in our previous analysis Moutarde:2018kwr .

6.4 Subtraction constant

The subtraction constant $C_{H}$ is evaluated with our results on the CFF $\mathcal{H}$ and the dispersion relation introduced in Sec. 2. More precisely, the imaginary part of the CFF $\mathcal{H}$ is integrated according to Eq. (5) in the range of $\epsilon<\xi^{\prime}<1$ and then subtracted from the corresponding real part. By studying $C_{H}$ as a function of $\epsilon$ the latter has been chosen to be $10^{-6}$, which introduces a little bias on $C_{H}$ comparing to $\epsilon=0$. This bias is estimated to be smaller than $1\%$.

The subtraction constant for a given $Q^{2}$ and $t$ should be the same independently on the value of $\xi$ used in its extraction. This is demonstrated in Fig. 23 (top), where $C_{H}$ is shown as a function of $\xi$ for example kinematics of $t=-0.3~{}\mathrm{GeV}^{2}$ and $Q^{2}=2~{}\mathrm{GeV}^{2}$. We consider this test as a proof of consistency between the parameterizations of the real and imaginary parts of the CFF $\mathcal{H}$.

The subtraction constant as a function of $t$ and $Q^{2}$ is shown in Fig. 23 (bottom) for $\xi=0.2$. As expected, the uncertainties are large in domains sparsely covered by data, in particular for $t\rightarrow 0$ and $Q^{2}\rightarrow\infty$, which are important for a direct interpretation of $C_{H}$ in terms of the energy-momentum tensor. However, the uncertainties obtained in domains covered by data are encouraging, so our result coming from a model-independent extraction can provide important constraints on models of $C_{H}$. Our findings are consistent with the recent observation of Ref. Kumericki:2019ddg .

7 Summary

In this paper we report the extraction of CFF parameterizations from proton DVCS data. The extracted quantities are the most basic observables as one can unambiguously access by exploring the DVCS process. We analyze the data in the region in which the interpretation in the language of GPDs is applicable. In this analysis, a dispersion relation is used to access the DVCS subtraction constant, which is related to the mechanical forces acting on partons in the nucleon. Also, the nucleon tomography that describes a spatial distribution of partons inside the nucleon is discussed.

The extraction of CFFs is done with the help of the artificial neural networks technique, allowing for an essential reduction of model dependency. The presented analysis includes such elements as the training of neural networks with the genetic algorithm, the careful regularization to avoid over-fitting and the propagation of experimental uncertainties with the replica method. The work is done within the PARTONS framework Berthou:2015oaw .

The results of this analysis include in particular unbiased CFF parameterizations extracted in the three-dimensional phase-space of $(x_{\mathrm{Bj}},t,Q^{2})$, with a reliable estimation of uncertainties. In addition, a direct extraction of the subtraction constant from the experimental data is presented.

The analysis is complementary to our previous one Moutarde:2018kwr , where CFF parameterizations were constructed with the basic GPD properties acting as Ansatz building blocks. Although a physically motivated Ansatz gives more insight into GPD physics, its inherent model-dependency introduces an extra theoretical uncertainty, which usually can be only roughly estimated. The situation is opposite for the analysis presented in this paper – the extraction of CFFs is unbiased, but the link with GPD physics is not as straightforward as in the case of previous analysis. Therefore, both analyses, which are performed on a similar data set, provide a complete picture of the DVCS process.

This work provides a benchmark for a powerful tool to be used in future GPD analyses. This tool allows in particular for a nearly model-independent estimation of the impact of future experiments, in particular those to be performed in the foreseen electron-ion facilities. CFF parameterizations that are presented in this paper can be used for a fast generation of DVCS cross sections, which may be useful e.g.\xspacefor Monte Carlo generation. Studies of Timelike Compton Scattering and higher twist contributions are also possible within the proposed framework.

Acknowledgements.

The authors would like to thank C. Weiss and K. Kumerički for fruitful discussions. This work was supported by the Grant No. 2017/26/M/ST2/01074 of the National Science Centre, Poland and was done in part under CRADA agreement No. 2017S002 between Jefferson Lab, USA and National Centre for Nuclear Research, Poland. The project is co-financed by the Polish National Agency for Academic Exchange and by the COPIN-IN2P3 Agreement. The computing resources of Świerk Computing Centre, Poland are greatly acknowledged.