Nuclear Effects in the Deuteron and Constraints on the d/u Ratio

TL;DR

This paper investigates nuclear effects in the deuteron through various scattering experiments, focusing on how these effects influence the extraction of the d quark distribution and confirming the universality of off-shell modifications.

Contribution

It provides a comprehensive analysis of nuclear corrections in the deuteron and introduces a consistent off-shell function applicable across different nuclear environments.

Findings

Nuclear corrections significantly affect the d quark distribution extraction.

The off-shell function is universal across different nuclei.

Deuteron wave function modeling impacts nuclear correction estimates.

Abstract

We present a detailed study of nuclear corrections in the deuteron (D) by performing an analysis of data from deep inelastic scattering off proton and D, dilepton pair production in $pp$ and pD interactions, and W^+- and Z boson production in pp and p pbar collisions. In particular, we discuss the determination of the off-shell function describing the modification of the parton distribution functions in bound nucleons in the context of global QCD fits. Our results are consistent with the ones obtained independently from the study of deep inelastic scattering data off heavy nuclei with $A\geq4$, further confirming the universality of the off-shell function. We also study the sensitivity to the modeling of the deuteron wave function. As an important application we discuss the impact of nuclear corrections to the deuteron on the determination of the d quark distribution.

Tables

Table 1. Table 1: Values of the average nucleon virtuality v = ( p 2 − M 2 ) / M 2 𝑣 superscript 𝑝 2 superscript 𝑀 2 superscript 𝑀 2 v=(p^{2}-M^{2})/M^{2} , bound nucleon energy ε 𝜀 \varepsilon , and kinetic energy 𝒑 2 / 2 ​ M superscript 𝒑 2 2 𝑀 \bm{p}^{2}/2M for each deuteron wave function shown in Fig.  1 .

Wave function	Reference	$⟨v⟩$	$⟨\epsilon ⟩$ [MeV]	$⟨{𝒑}^2⟩/2M$ [MeV]
WJC1	Gross and Stadler (2008, 2010)	-0.062	-16.16	13.94
WJC2	Gross and Stadler (2008, 2010)	-0.049	-12.89	10.67
AV18	Veerasamy and Polyzou (2011)	-0.045	-11.93	9.71
Paris	Lacombe et al. (1980)	-0.043	-11.55	9.33
CD-Bonn	Machleidt (2001)	-0.037	-9.96	7.74

Table 2. Table 2: List of the various data sets used in the present analysis.

	Experiment	Reference	Beam	Target(s)	Final states	Data points
DIS collider	HERA I+II	Abramowicz et al. (2015)	$e$	$p$	$eX$	1168
					$\nu X$
	HERA I+II	Abramowicz et al. (2013)	$e$	$p$	$ecX$	52
	H1	Aaron et al. (2010)	$e$	$p$	$ebX$	12
	ZEUS	Abramowicz et al. (2014)	$e$	$p$	$ebX$	17
DIS fixed target	BCDMS	Benvenuti et al. (1989, 1990)	$\mu$	$p,D$	$\mu X$	605
	NMC	Arneodo et al. (1997)	$\mu$	$p,D$	$\mu X$	490
	SLAC E49a	Bodek et al. (1979)	$e$	$p,D$	$eX$	118
	SLAC E49b	Bodek et al. (1979)	$e$	$p,D$	$eX$	299
	SLAC E87	Bodek et al. (1979)	$e$	$p,D$	$eX$	218
	SLAC E89b	Mestayer et al. (1983)	$e$	$p,D$	$eX$	162
	SLAC E139	Gomez et al. (1994)	$e$	$D$	$eX$	17
	SLAC E140	Dasu et al. (1994)	$e$	$D$	$eX$	26
	JLab BONuS	Griffioen et al. (2015)	$e$	$D$	$eX$	5
	NOMAD	Samoylov et al. (2013)	$\nu$	Fe	${\mu }^{+}{\mu }^{-}X$	48
	CHORUS	Kayis-Topaksu et al. (2011)	$\nu$	Emul.	$\mu cX$	6
	CCFR	Goncharov et al. (2001)	$\nu$	Fe	${\mu }^{+}{\mu }^{-}X$	89
	NuTeV	Goncharov et al. (2001)	$\nu$	Fe	${\mu }^{+}{\mu }^{-}X$	89
Drell-Yan fixed target	FNAL E866	Towell et al. (2001)	$p$	$p,D$	${\mu }^{+}{\mu }^{-}$	39
	FNAL E605	Moreno et al. (1991)	$p$	Cu	${\mu }^{+}{\mu }^{-}$	119
$W,Z$ collider	D0	Abazov et al. (2013)	$p$	$\overline{p}$	$W^{+}\to {\mu }^{+}\nu$	10
					$W^{-}\to {\mu }^{-}\nu$
	D0	Abazov et al. (2015)	$p$	$\overline{p}$	$W^{+}\to e^{+}\nu$	13
					$W^{-}\to e^{-}\nu$
	ATLAS	Aad et al. (2012, 2016)	$p$	$p$	$W^{+}\to l^{+}\nu$	36
					$W^{-}\to l^{-}\nu$
					$Z\to l^{+}{l}^{-}$
	CMS	Chatrchyan et al. (2014); Khachatryan et al. (2016a)	$p$	$p$	$W^{+}\to {\mu }^{+}\nu$	33
					$W^{-}\to {\mu }^{-}\nu$
	LHCb	Aaij et al. (2015a, 2016)	$p$	$p$	$W^{+}\to {\mu }^{+}\nu$	63
					$W^{-}\to {\mu }^{-}\nu$
					$Z\to {\mu }^{+}{\mu }^{-}$
	LHCb	Aaij et al. (2015b)	$p$	$p$	$Z\to e^{+}{e}^{-}$	17

Equations

$$F_2^D=F_2^{N/D}+{\delta }_{\mathrm{MEC}}{F}_2^D+{\delta }_{\mathrm{coh}}{F}_2^D,$$

$$F_2^N(x,Q^2,p^2)=F_2^{\text{TMC}}(x,Q^2,p^2)+H_2^{(4)}(x)/Q^2+\mathcal{O}(Q^{-4})$$

$$F_2^{\text{TMC}}(x,Q^2,p^2)=\frac{x^2}{{\xi }^2{\gamma }^3}{F}_2^{\mathrm{LT}}(\xi ,Q^2,p^2)+\frac{6x^3{p}^2}{Q^2{\gamma }^4}{\int }_{\xi }^1\frac{dz}{z^2}{F}_2^{\mathrm{LT}}(z,Q^2,p^2)+\mathcal{O}(Q^{-4})$$

Full text

I Introduction

Parton distribution functions (PDFs) are universal process-independent characteristics of the target, which determine leading contributions to the cross sections of various hard processes involving leptons and hadrons Collins et al. (1989). The PDF content of both the proton and the neutron is extracted from global fits Alekhin et al. (2017); Harland-Lang et al. (2015); Ball et al. (2015); Dulat et al. (2016) to experimental data at large momentum transfer, including lepton deep inelastic scattering (DIS), lepton-pair production (Drell-Yan process), jet production, and $W$ and $Z$ boson production in hadron collisions. In order to disentangle the content of different parton flavors, global fits must include complementary data which are flavor dependent. Traditionally, the most efficient separation between $d$ and $u$ quark distributions is obtained by comparing charged-lepton DIS data for proton and deuterium targets, the latter being considered as an “effective” neutron target. Since the deuteron is a weakly bound nucleus with a binding energy of about $2.2\>\mathrm{MeV}$ – accounting only for about 0.1% of its mass – it is often assumed to be well approximated by the sum of a quasi-free proton and a quasi-free neutron in the PDF analyses.

However, charged-lepton DIS data from various nuclear targets demonstrate significant nuclear effects with a rate that is more than 1 order of magnitude larger than the ratio of the nuclear binding energy to the nucleon mass (for a review see Refs. Arneodo (1994); Norton (2003)). These observations indicate that the nuclear environment plays an important role even at energies and momenta much higher than those involved in typical nuclear ground state processes Arneodo (1994); Norton (2003); Geesaman et al. (1995). Similar considerations can be drawn for the DIS off the deuteron Atwood and West (1973); Landshoff and Polkinghorne (1978); Kusno and Moravcsik (1981); Bodek and Ritchie (1981); Dunne and Thomas (1986); Frankfurt and Strikman (1988); Kaptar and Umnikov (1991); Nakano and Wong (1991); Epele et al. (1992); Braun and Tokarev (1994); Gomez et al. (1994); Melnitchouk et al. (1994); Abe et al. (1995); Burov and Molochkov (1998); Alekhin et al. (2004); Burov et al. (2004); Kulagin and Petti (2006); Kulagin and Melnitchouk (2008); Arrington et al. (2009); Weinstein et al. (2011); Arrington et al. (2012). In spite of the broad range of predictions, such studies indicate that nuclear effects in the deuteron are non-negligible and rise rapidly in the region of large Bjorken $x$. A recent direct measurement of nuclear effects in the deuteron Griffioen et al. (2015) indicates a few-percent negative correction at $x\sim 0.5$-$0.6$, with a steep rise at large $x$. Therefore, if neglected or treated incorrectly, these nuclear effects can potentially introduce significant uncertainties and/or biases in the extraction of the neutron structure functions and of the $d$ quark distribution from the DIS data Accardi et al. (2011).

A microscopic model for nuclear structure functions and PDFs accounting for a number of different nuclear effects was developed in Refs.Kulagin and Petti (2006, 2007, 2010, 2014). It includes the smearing with the energy-momentum distribution of bound nucleons (Fermi motion and binding, FMB), the off-shell correction (OS) to bound nucleon structure functions, the contributions from meson exchange currents and the propagation of the hadronic component of the virtual intermediate boson in the nuclear environment. This model has been successfully used to quantitatively explain the observed $x$, $Q^{2}$ and $A$ dependence of the nuclear DIS data in a wide range of targets from 3He to 207Pb Kulagin and Petti (2006, 2007, 2010), the magnitude, the $x$ and mass dependence of the nuclear Drell-Yan (DY) data Kulagin and Petti (2014), as well as the data on the differential cross sections and asymmetries for $W^{\pm},Z$ production in $p+{\rm Pb}$ collisions at the LHC Ru et al. (2016).

A consistent description of the scattering off bound nucleons not only involves the smearing due to the nuclear momentum distribution, but also requires the knowledge of the off-shell (OS) scattering amplitudes. The model of Ref. Kulagin and Petti (2006) exploits the observation that the nucleus is a weakly bound system and thus it is sufficient to evaluate the OS correction to the bound nucleon PDFs in the vicinity of the mass shell. The shape of this correction is defined by a universal function $\delta f(x)$ of the Bjorken variable $x$, while its nuclear dependence is driven by the average virtuality of the nucleon (off-shellness) inside the nucleus. The OS function $\delta f$ can be regarded as a special nucleon structure function, which does not contribute to the cross section of the physical nucleon, but is relevant only for the bound nucleon and describes its response to the interaction with the nuclear environment. The off-shell correction proved to be an important contribution to explain nuclear effects at large $x$. The function $\delta f$ was determined from the analysis of data on ratios of DIS structure functions in various nuclei Kulagin and Petti (2006). It was also shown that in a simple single-scale model, in which the quark momentum distributions in the nucleon are functions of the nucleon core radius, the observed behavior of $\delta f$ can be interpreted in terms of an increase of the confinement radius of the bound nucleon in the nuclear environment Kulagin and Petti (2006).

The deuteron is a weakly bound state of two nucleons with peculiar attributes. Its dynamics is better understood than the dynamics of many-particle nuclei, making it an ideal benchmark tool for the study of different nuclear effects. However, it is also considerably different with respect to even a three-body nucleus like 3He. For these reasons one cannot rely on simple extrapolations of nuclear effects from heavy targets based upon nuclear density or atomic weight, as it is often assumed in phenomenological analyses. In contrast, the model of Ref. Kulagin and Petti (2006) suggests a unified treatment of the deuteron and heavier nuclei on the basis of common underlying physics mechanisms.

In this paper we discuss an independent determination of the off-shell function $\delta f$, together with the proton PDFs, from a global QCD analysis of proton and deuterium data. In Sec.II we review the model of nuclear corrections in the deuteron, while in Sec. III we discuss the details of the data analysis. In Sec. IV we compare our results with the ones obtained from heavy nuclear targets and discuss the impact on the uncertainties related to the $d/u$ ratio from global QCD fits. We summarize our results in Sec. V.

II Model of Nuclear Corrections

The nuclear corrections to the inelastic structure functions involve a number of different contributions. For the deuteron we can write (for simplicity we summarize the structure function $F_{2}$ here) Kulagin and Petti (2006, 2014):

[TABLE]

where the first term in the right-hand side stands for the incoherent scattering off the bound isoscalar nucleon $N$ including the off-shell correction, and $\delta_{\mathrm{MEC}}F_{2}^{D}$ and $\delta_{\rm coh}F_{2}^{D}$ are the corrections due to nuclear meson exchange currents (MEC) and coherent interactions of the intermediate virtual boson with nuclear target, respectively.

II.1 Incoherent Scattering off Bound Nucleons

The first term in Eq.(1) dominates at $x>0.2$ and can be written as follows Kulagin and Petti (2006): {align}

γ^2 F_2^N/D(x,Q^2) &= ∫d3p(2π)3 —Ψ_D(p)—^2 (1+γpzM) (γ’^2 +6x’2p⟂2Q2 ) F_2^N(x’,Q^2,p^2),

where the integration is over the momentum of the bound nucleon $\bm{p}$, $\Psi_{D}(\bm{p})$ is the deuteron wave function, $M=\tfrac 12(M_{p}+M_{n})$ and $F_{2}^{\smash{N}}=\tfrac 12(F_{2}^{p}+F_{2}^{n})$ are respectively the mass and the structure function of the bound nucleon with four-momentum $p=(M+\varepsilon,\bm{p})$, where $\varepsilon=\varepsilon_{D}-\bm{p}^{2}/(2M)$ and $\varepsilon_{D}=M_{D}-2M$ is the deuteron binding energy. The integration in Eq.(II.1) requires the structure function of the bound nucleon in the off-shell region and $F_{2}^{\smash{N}}$ depends on the Bjorken variable $x^{\prime}=Q^{2}/(2pq)$, the momentum transfer square $Q^{2}$, and also on the nucleon invariant mass squared $p^{2}$. In Eq.(II.1) we use a coordinate system such that the momentum transfer $\bm{q}$ is antiparallel to the $z$ axis, $\bm{p}_{\perp}$ is the transverse component of the nucleon momentum, and $\gamma^{2}=1+4x^{2}M^{2}/Q^{2}$ and ${\gamma^{\prime}}^{2}=1+4{x^{\prime}}^{2}p^{2}/Q^{2}$.

The integrand in Eq.(II.1) factorizes into two independent terms involving the contribution from two different scales: i) the wave function $\Psi_{D}(\bm{p})$ describing the deuteron properties in momentum space, and ii) the nucleon structure function $F_{2}^{N}$ describing processes at the parton level in the nucleon. In the following we will consider several deuteron wave functions $\Psi_{D}(\bm{p})$ corresponding to different models for the nucleon-nucleon potential: Paris Lacombe et al. (1980), CD-Bonn Machleidt (2001), AV18 Veerasamy and Polyzou (2011), WJC1 and WJC2 Gross and Stadler (2008, 2010). These wave functions are constrained by high-precision fits to nucleon-nucleon scattering data at low energies. However, these models for $\Psi_{D}(\bm{p})$ can differ by more than a factor of 10 in the high momentum tail, as shown in Fig. 1. Table 1 summarizes the salient kinematic parameters associated with each deuteron wave function. To be consistent with the weak binding approximation of Ref.Kulagin and Petti (2006), we perform the integration over the nucleon momentum in Eq.(II.1) up to $|\bm{p}|<1$ GeV/$c$.

The nucleon structure function in Eq.(II.1) includes the target mass and the high-twist (HT) corrections represented as follows

[TABLE]

where $F_{2}^{\text{TMC}}$ is the leading-twist (LT) structure function corrected for the target mass effects (TMC) and $H_{2}^{(4)}$ describes the twist-4 contribution (we suppress any explicit notation to higher order terms for brevity). The LT structure function is computed using the proton and the neutron PDFs extracted from the analysis of data as described in Sec.III. The target mass correction is computed following the prescription of Ref. Georgi and Politzer (1976):

[TABLE]

where $\xi=2x/(1+\gamma)$ is the Nachtmann variable, we substitute $p^{2}$ for the mass of the bound nucleon squared $M^{2}$, and $\gamma^{2}=1+4x^{2}p^{2}/Q^{2}$. Note that the second term in Eq.(3) is suppressed as $1/Q^{2}$ and therefore can be formally considered as a kinematic HT contribution. Recent phenomenology suggests that Eq.(2) with twist-4 contributions provides a good description of data down to $Q\sim 1\>\mathrm{GeV}$ Alekhin et al. (2007). It is also worth noting that this model is consistent with duality principle and on average describes the resonance data with $W<1.8\>\mathrm{GeV}$ Alekhin et al. (2007); Malace et al. (2009).

II.2 Off-shell Correction

The structure function of the bound nucleon $F_{2}^{N}(x,Q^{2},p^{2})$ appearing in the calculation of the nuclear correction in Eq.(II.1) explicitly depends on the nucleon invariant mass squared $p^{2}$. The $p^{2}$-dependence of the structure function has two different sources Kulagin et al. (1994); Kulagin and Petti (2006): (i) the dynamic off-shell dependence of the LT structure function; (ii) the kinematic target mass correction, which generates terms of the order of $p^{2}/Q^{2}$. We evaluate the off-shell dependence of the target mass correction by replacing $M^{2}\to p^{2}$ in Eq.(3). Since the characteristic momenta of a bound nucleon are small compared to its mass, the integration in Eq.(II.1) mainly covers a region in the vicinity of the mass shell. The nucleon virtuality $v=(p^{2}-M^{2})/M^{2}$ can then be treated as a small parameter, so that we can expand the structure function in series of $v$ keeping only the leading term: {align}

F_2^LT(x,Q^2,p^2) &= F_2^LT(x,Q^2)[ 1+δf(x,Q^2) v ],

δf = ∂lnF_2^LT/∂lnp^2∣_p^2=M^2,

where the first term on the right-hand side in Eq.(II.2) is the structure function of the on-mass-shell nucleon. The off-shell (OS) function $\delta f$ can be regarded as a special nucleon structure function, which describes the relative modification of the nucleon structure functions and PDFs in the vicinity of the mass shell. This function does not contribute to the cross section of the physical nucleon, but it is relevant only for the bound nucleon and describes its response to the interaction in a nucleus.

In general, the function $\delta f$ might be flavor dependent and different for protons and neutrons. However, the study of nuclear DIS and DY data Kulagin and Petti (2006, 2010, 2014) supports the hypothesis of the OS function universality, with no significant $Q^{2}$ dependence suggested by the data, i.e. $\delta f(x,Q^{2})=\delta f(x)$. Although we assume that $\delta f$ is only a function of $x$, the overall OS correction to the nuclear structure functions also depends on $Q^{2}$, as a result of the integration of Eq.(II.1). It is important to note that this $Q^{2}$ dependence is different from the ones of both the LT and HT contributions to the structure functions in Eq.(2). This difference allows a simultaneous extraction of PDFs, HTs and $\delta f$ from global QCD fits (cf. Sec. III).

II.3 Other Nuclear Corrections

The remaining corrections appearing in Eq.(1), i.e. the nuclear meson exchange current $\delta_{\mathrm{MEC}}F_{2}^{D}$ and the nuclear shadowing (NS) $\delta_{\rm coh}F_{2}^{D}$, are relevant only at small $x$. For the details of the treatment of these terms we refer to Ref.Kulagin and Petti (2006). Here we only emphasize that nuclear effects present in different kinematical regions of $x$ are related by the DIS sum rules and normalization constraints Kulagin and Petti (2014). For example, the light-cone momentum sum rule links the FMB and the MEC corrections. We use this relation to constrain the mesonic contributions to the nuclear structure functions. Similarly, the baryon number sum rule links the shadowing and the OS corrections. In our approach, the OS effect provides the mechanism to cancel a negative nuclear-shadowing contribution to the normalization of the nuclear valence quark distributions.

II.4 Model Predictions and Phenomenology

The model described above was used to perform a detailed analysis of data on the ratio $\mathcal{R}(A/B)=F_{2}^{A}/F_{2}^{B}$ between the nuclear targets $A$ and $B$ Kulagin and Petti (2006). The data sets analyzed stem from a variety of electron and muon DIS experiments (CERN NMC, EMC, and BCDMS, SLAC E139 and E140, Fermilab E665), with targets ranging from 4He to 208Pb in a wide region of $x$ and $Q^{2}$. In this way we tested the hypothesis that the OS modification of bound nucleons is responsible for the difference between the data and all known nuclear effects in Eq.(1), including the FMB Akulinichev et al. (1985); Kulagin (1989), the nuclear shadowing and the nuclear MEC. In turn, this OS correction is controlled by the universal off-shell function $\delta f$ in Eq.(II.2), which was determined from DIS data with the corresponding uncertainty Kulagin and Petti (2006).

Such an approach leads to an excellent agreement with the available DIS data on the $x$, $Q^{2}$, and $A$ dependence of $\mathcal{R}(A/B)$. The model predictions are also in a good agreement Kulagin and Petti (2010) with the recent DIS data by the HERMES experiment at HERA Ackerstaff et al. (2000) and the E03-103 experiment at JLab Seely et al. (2009) down to 3He. Furthermore, the same model allows a calculation of nuclear PDFs, which can describe well the magnitude, the $x$ and mass dependence of the available data on Drell-Yan production off various nuclear targets Kulagin and Petti (2014), as well as the differential cross sections and asymmetries for $W^{\pm},Z$ production in $p+{\rm Pb}$ collisions at the LHC Ru et al. (2016).

In this paper we perform an independent analysis of deuterium and proton data in the context of global QCD fit. In Fig.2 we show the predictions of Ref. Kulagin and Petti (2006) for the ratio $\mathcal{R}(D/N)=F_{2}^{D}/F_{2}^{N}$ of the deuteron and the isoscalar nucleon structure functions at $Q^{2}=20\>\mathrm{GeV}^{2}$. The region of $0.35\leq x\leq 0.55$ is characterized by an almost linear dependence from $x$, with a slope $\mathrm{d}{\mathcal{R}}(D/N)/\mathrm{d}x=-0.099\pm 0.006$, including model uncertainties. This slope is often used in the analysis of experimental data Seely et al. (2009) since it is less affected by the experimental uncertainties (especially the overall normalization) than the absolute value of the nuclear correction. At large $x>0.1$ nuclear corrections are dominated by the FMB and OS effects, as shown in Fig. 2. In particular, the off-shell correction is a crucial contribution in this kinematic region, which is studied in more detail in the present analysis.

III Off-shell Correction from Global QCD Fit

In this paper we discuss the impact of nuclear effects in the deuteron in the context of global QCD fits. Our goals are twofold: (i) an independent determination of the off-shell correction preferred by the deuteron data; (ii) an estimate of the PDF uncertainties (in particular for the $d/u$ ratio at large $x$) introduced by the nuclear corrections to deuterium data. The analysis framework and the main data sets used are common to the ABMP16 fit Alekhin et al. (2017).

III.1 Analysis Framework

In our analysis we use the next-to-next-to-leading-order (NNLO) approximation in the QCD perturbation theory to calculate the partonic cross sections entering the LT terms for the hard interaction processes considered. We set the renormalization and factorization scales to $\mu_{r}=\mu_{f}=\mu$ and we identify this scale $\mu$ with the relevant kinematics of each process, e.g. $\mu=Q$ for DIS. The individual PDFs are parametrized as in Ref. Alekhin et al. (2017) at the starting scale $\mu^{2}=Q^{2}_{0}=9$ GeV2. The PDFs are subject to sum rule constraints due to the conservation of the quark number and the momentum in the nucleon.

The splitting functions controlling the scale dependence of the PDFs in the evolution equations are evaluated at NNLO in perturbation theory Moch et al. (2004); Vogt et al. (2004). The Wilson coefficients entering the massless DIS structure functions are calculated at NNLO van Neerven and Zijlstra (1991); Zijlstra and van Neerven (1991, 1992a, 1992b); Moch and Vermaseren (2000); Moch et al. (2005); Vermaseren et al. (2005). Similarly, we use NNLO calculations for the partonic cross sections of the Drell-Yan process and the hadronic $W$ and $Z$ boson production Hamberg et al. (1991); Harlander and Kilgore (2002); Anastasiou et al. (2003, 2004); Catani et al. (2009).

In our PDF analysis we use a fixed flavor number scheme with $n_{f}=3$ light flavors from which heavy quark PDFs are generated. The heavy quark masses $m_{q}$ are defined in the $\overline{\rm MS}$ renormalization scheme as running masses $m_{q}(\mu)$ depending upon the scale $\mu$ of the hard scattering in analogy to the running coupling $\alpha_{s}(\mu)$. As discussed in Ref. Alekhin and Moch (2011), the use of the $\overline{\rm MS}$-mass allows better convergence properties and greater perturbative stability at higher orders. The heavy quark Wilson coefficients entering the DIS structure functions for heavy quark production are known exactly only to the next-to-leading-order (NLO) for both the charged current (CC) Gottschalk (1981); Gluck et al. (1996) and neutral current (NC) Laenen et al. (1993) processes. For the NC case approximate NNLO coefficients are employed Alekhin et al. (2017).

In the kinematic range of our analysis the twist-6 terms give negligible contributions to the structure functions Alekhin et al. (2007). Therefore, in addition to the leading twist we only include two twist-4 contributions $H_{2}$ and $H_{T}$ – as defined in Eq.(2) – to the structure functions $F_{2}$ and $F_{T}$, respectively. We also considered the target dependence of the HT parametrization. The isospin asymmetry in $H_{T}$ is consistent with zero within uncertainties Alekhin et al. (2004) and therefore is neglected in our analysis. The isospin asymmetry in $H_{2}$ is also small Alekhin et al. (2004). Although the values of this latter have a better statistical significance, we set it to zero as well in order to avoid a potential bias in the nuclear corrections extracted from the global QCD fits. In summary, we fit two twist-4 coefficients for the isoscalar nucleon, $H_{2}^{N}$ and $H_{T}^{N}$. These power corrections are parametrized as cubic spline functions of $x$.

The nuclear corrections for the deuteron are calculated according to the model described in Sec. II. We do not include meson exchange currents and coherent nuclear effects (shadowing) for the deuteron, since their impact is negligible in the kinematic coverage of our analysis (see Fig. 2) and we are mainly focused on the study of the off-shell correction. The only free parameters entering the nuclear corrections are the ones describing the off-shell function $\delta f(x)$, which are extracted simultaneously with the PDFs and HT terms. To this end, we use a parametrization with generic second- and third-order polynomials for $\delta f(x)$. We verified that there is no statistically significant difference between these two options within the accuracy of the data samples used in our analysis.

III.2 Data Samples

In our analysis most of the information about the deuteron is provided by the inclusive DIS data off deuterium from the SLAC E49, E87, E89, E139, E140 Bodek et al. (1979); Atwood et al. (1976); Mestayer et al. (1983); Gomez et al. (1994); Dasu et al. (1994) and the CERN BCDMS Benvenuti et al. (1990), NMC Arneodo et al. (1997) experiments, as well as by the ratio of Drell-Yan production in $p{\rm D}$ and $pp$ collisions from the Fermilab E866 experiment Towell et al. (2001). In addition, the recent direct measurement of the ratio $F_{2}^{D}/F_{2}^{N}$ Griffioen et al. (2015) by the BONuS experiment Tkachenko et al. (2014) at Jefferson laboratory allows a better disentanglement of the nuclear corrections in the deuteron from possible variations of the $d/u$ ratio in the nucleon. Since most of the BONuS data either have low values of $Q^{2}$ or are in the resonance region, we only include the BONuS points with $Q^{2}>1.5$ GeV2 and $W>1.6$ GeV. Although these cuts are less stringent than the ones we apply for the other data sets, they are justified by a partial cancellation of HT effects in the ratio $F_{2}^{D}/F_{2}^{N}$ and by the relevance of the direct BONuS measurement for our study.

For consistency with the ABMP16 fit Alekhin et al. (2017) and to better constrain the sea quark distributions we include the Drell-Yan data in $p{\rm Cu}$ interactions by the E605 experiment, as well as charm production data in (anti)neutrino interactions off heavy targets by the CCFR Goncharov et al. (2001), NuTeV Goncharov et al. (2001), NOMAD Samoylov et al. (2013), and CHORUS Kayis-Topaksu et al. (2011) experiments. We verified that these data using nuclear targets do not affect our results on the deuteron by performing dedicated fits with and without such data. It is also worth noting that NOMAD and CHORUS measured the ratios of charm to inclusive charged-current cross sections and it was shown that the corresponding nuclear corrections cancel out at the subpercent level in such ratios Alekhin, S. and Blümlein, J. and Caminadac, L. and Lipka, K. and Lohwasser, K. and Moch, S. and Petti, R. and Placakyte, R. (2015).

Nuclear corrections in the deuteron can be determined by comparing the available data for deuterium targets with the ones originated from interactions on the free nucleons. However, inclusive DIS data off protons do not allow to disentangle the $d$ and $u$ quark distributions because of the lack of a corresponding free neutron target.

The limitations of DIS data can be partially overcome with the addition of Drell-Yan, $W^{\pm}$, and $Z$ production at Tevatron and the LHC Alekhin et al. (2017); Abazov et al. (2013, 2015); Aad et al. (2012, 2016); Chatrchyan et al. (2014); Khachatryan et al. (2016a); Aaij et al. (2015a, b, 2016). In particular, the data on $W^{+}$ and $W^{-}$ production allow a $d/u$ separation independent from the deuterium data. We note that the same $W^{\pm}$ production data sets collected at Tevatron and the LHC may result in two distinct (but correlated) measurements: (i) the $l^{\pm}$ lepton asymmetry from the $W^{\pm}$ decays; (ii) the actual $W^{\pm}$ asymmetry. The former is more closely related to the experimental observables, while the latter requires model-dependent acceptance corrections to account for the kinematics of the $W^{\pm}$ decay. As discussed in Ref. Alekhin et al. (2016), some inconsistencies between the $l^{\pm}$ lepton asymmetries and the $W^{\pm}$ asymmetries obtained from the same experimental data sets are observed. For this reason whenever both measurements are available, we only consider the $l^{\pm}$ lepton asymmetry data in our analysis.

Table 2 summarizes all the data sets used in the present analysis. In order to exclude the region of resonance production and to reduce the impact of HT corrections we require $Q^{2}>2.5$ GeV2 and $W>1.8$ GeV for DIS data.

IV Results and Discussion

The general features of the PDFs extracted from the global fits described above as well as a detailed discussion of the individual data samples considered were presented elsewhere Alekhin et al. (2017). In this paper we focus on the nuclear correction extracted from deuterium data and on the corresponding impact for the $d/u$ ratio.

IV.1 Off-shell Function $\delta f$

The results of our determination of the off-shell function $\delta f(x)$ from the global QCD fits described in Sec. III are shown in Fig. 3. A simultaneous extraction of the off-shell function with both PDFs and HTs is possible because of the different $Q^{2}$ dependence of these three contributions and the wide $Q^{2}$ coverage of the data sets listed in Table 2. In general, nuclear corrections to the deuteron data are partially correlated to the $d$-quark distribution. In order to reduce this correlation, the role of Drell-Yan and $W^{\pm}$ production at $pp$ and $p\bar{p}$ colliders is crucial. In particular, the recent combined D0 data and the LHC data from LHCb, reaching values of $x\sim 0.8$ due to the wide rapidity coverage, offer precisions comparable to the ones of older fixed-target DIS experiments.

Figure 3 illustrates the dependence of the fitted $\delta f(x)$ function upon the choice of the deuteron wave function $\Psi_{D}(\bm{p})$ among the models listed in Sec. II.1. The main differences are related to the high momentum component of the wave function, as shown in Fig. 1. Since this high momentum tail controls the region of large nucleon virtuality $v$, the off-shell correction in the large $x$ region is in principle sensitive to the corresponding nuclear smearing in Eq.(II.1), which modifies the $x$ and $Q^{2}$ dependence of the structure functions. A general trend can be observed from Fig. 3, with the harder wave function resulting in a slightly higher off-shell function at large $x$. Since the overall off-shell correction has opposite sign with respect to $\delta f(x)$ in Eq.(II.2), this trend implies an anticorrelation between $\Psi_{D}(\bm{p})$ and $v\delta f(x)$ in global QCD fits. From Fig. 3 we note that our results obtained with the Paris, CD-Bonn, AV18, WJC1, and WJC2 wave functions are all consistent within the corresponding uncertainties and indicate a relatively limited spread. This robustness against the modeling of the deuteron wave function can be explained by the use of data samples which can reduce the correlation between the nuclear correction and the $d$-quark distribution. In this context the recent direct measurement of the ratio $F_{2}^{D}/F_{2}^{N}$ from the BONuS experiment contributes to constrain the overall normalization of the nuclear corrections in our fits.

A more precise determination of the off-shell function $\delta f(x)$ was obtained in Ref. Kulagin and Petti (2006) from heavy nuclei with $A\geq 4$, as described in Sec. II.4. In order to further study the sensitivity to the nuclear smearing in Eq.(II.1), we repeat the standalone extraction of $\delta f(x)$ in Ref Kulagin and Petti (2006) after rescaling the nuclear spectral function describing the properties of heavy nuclei Kulagin and Petti (2006) by the ratios of the various deuteron wave functions shown in Fig. 1. The results summarized in Fig. 4 demonstrate a small sensitivity to the choice of the nuclear spectral function and/or of the deuteron wave function, as well as a dramatic reduction of the uncertainties with respect to Fig. 3. This reduction can be explained by the different observables considered in the two independent extractions. In the global QCD fits we use the absolute DIS cross sections off the deuteron, while in the standalone determination of Ref Kulagin and Petti (2006) we consider only ratios $\mathcal{R}(A/B)=F_{2}^{A}/F_{2}^{B}$ between two nuclear targets $A$ and $B$. Many model uncertainties largely cancel out in such ratios. For the same reason the data sets used in Ref Kulagin and Petti (2006) are more accurate than the deuteron ones, making them an excellent tool to study the off-shell function $\delta f(x)$. The $\pm 1\sigma$ uncertainty band shown in Fig. 4 includes model systematics due to the spectral and wave functions, the functional form, the PDFs, as well as corrections due to meson exchange currents and nuclear shadowing.

A comparison between the two independent determinations of the off-shell function $\delta f$ is given in Fig. 3. Since the five individual determinations of the off-shell function $\delta f$ from deuterium data using different wave functions are characterized by a comparable fit quality, we combine them by taking an average of both the central values and the corresponding uncertainties. The resulting $\pm 1\sigma$ uncertainty band is shown in the right panel of Fig. 3. This band summarizes our determination of $\delta f$ from deuteron data and is consistent with the more precise determination of $\delta f$ from the analysis of the ratios of nuclear structure functions for $A\geq 4$ Kulagin and Petti (2006). Since we are using a generic polynomial to parametrize $\delta f$ (Sec. III.1), no functional form bias is present in this comparison. The agreement between the two independent determinations supports the interpretation of the off-shell function $\delta f$ as a universal structure function of the nucleon, validating the unified treatment of the deuteron and heavier nuclei developed in Ref Kulagin and Petti (2006).

IV.2 Nuclear Corrections to $F_{2}^{D}/F_{2}^{N}$

The nuclear correction stemming from the FMB and OS effects on the ratio $F_{2}^{D}/F_{2}^{N}$ is given in Fig. 5. This ratio is particularly interesting because it represents the overall nuclear corrections in the deuteron. The variation due to the choice of the deuteron wave function in the global QCD fit appears to be even smaller in the ratio $F_{2}^{D}/F_{2}^{N}$ than in the off-shell functions $\delta f$ given in Fig. 3. This behavior can be explained by the anticorrelation between $\Psi_{D}(\bm{p})$ and $v\delta f(x)$ discussed in Sec. IV.1: a larger off-shell function partially compensates a reduced strength of the high momentum component of the wave function so that the observable structure function remains consistent with the fitted data. In addition, the recent BONuS measurement significantly constrains the ratio $F_{2}^{D}/F_{2}^{N}$, as mentioned in Sec. IV.1. The results obtained from deuteron data with such a constraint agree with the predictions from Ref Kulagin and Petti (2006) based upon a standalone analysis of heavy nuclei with $A\geq 4$, as shown in Fig. 5. 111 The result with the off-shell function $\delta f$ from Ref. Kulagin and Petti (2006) shown in Fig. 5 is slightly different with respect to the calculation in Fig. 2. The differences are mainly at large $x$ values and appear due to the fact that the results shown in Fig. 5 are obtained with the PDFs and HT terms extracted from our global QCD fit.

Although the off-shell function $\delta f$ is extracted in our analysis as a generic polynomial, we are still calculating the nuclear correction to the structure functions using the nuclear convolution in Eq.(II.1), following the prescriptions of the model of Ref. Kulagin and Petti (2006). In order to verify whether this procedure introduces any indirect model dependence in our results, we perform a separate fit in which we parametrize the overall off-shell correction to the structure function $F_{2}^{D}$ as a generic polynomial added to the standard FMB correction. In this approach the fitted off-shell correction is model independent as it is not part of the nuclear convolution in Eq.(II.1). The results obtained with such a parametrization shown in Fig. 6 are in good agreement with the corresponding fits based upon the nuclear convolution with the off-shell function $\delta f$. We can thus conclude that the functional form we are using in our fits for $\delta f$ is flexible enough to reproduce the data and that our modeling of nuclear effects does not introduce any significant bias.

IV.3 Systematic Studies

As we discuss in Sec. IV.1, the uncertainty on the off-shell function $\delta f(x)$ related to the modeling of the deuteron wave function turns out to be negligible compared to the statistical uncertainty from our global QCD fits. Our final $\pm 1\sigma$ band on the ratio $F_{2}^{D}/F_{2}^{N}$, computed by averaging the results from the five individual fits with different wave functions in Fig. 5, is given in Fig. 7, including the statistical uncertainty from the fit and the one related to the choice of the deuteron wave function added in quadrature.

In Ref. Alekhin et al. (2007) it was shown that there is some tension between the DIS data sets from the BCDMS and SLAC experiments, resulting in significant modifications of the extracted HT terms and PDFs. In order to mitigate the impact of this tension on our studies we allow the overall normalization of both the BCDMS proton and deuteron data sets to vary freely in our fits. This approach is justified by the use of separate normalizations for the deuterium and proton data sets in the BCDMS measurements Benvenuti et al. (1989, 1990). The normalization of the BCDMS proton data is essentially defined by the precise HERA data in the overlap region, resulting in an overall factor consistent with the corresponding normalization uncertainties quoted by the experiment, up to 3%. The only constraint on the normalization of the BCDMS deuteron data comes from the SLAC experiments. However, the partial correlation between this normalization and the determination of the deuteron nuclear correction can potentially introduce an additional uncertainty in the global fits. The recent direct measurement of the deuteron nuclear correction by the BONuS experiment substantially reduces this uncertainty by constraining the normalization of the overall nuclear corrections. As a result, the normalization factor for the BCDMS deuterium data obtained from our fits is stable against variations of the deuteron wave function and is very close to the one obtained for the BCDMS proton data.

As an additional test of the robustness of our analysis, we perform separate fits to different subsets of the data listed in Table 2. These variants of our analysis do not indicate any anomalous tension related to individual data samples, but rather suggest that our results on the off-shell correction are originated from the combined fit of all deuteron data sets.

IV.4 Discussion

The results of our analysis discussed in Sec. IV.1 support the predictions of Ref. Kulagin and Petti (2006) for the nuclear effects in the deuteron and the unified treatment of all nuclei. We can then exploit the higher precision offered by the existing DIS data off heavier nuclear targets ($A\geq 4$) to fix the value of the off-shell function $\delta f$ used in global QCD fits following Ref. Kulagin and Petti (2006). The corresponding reduction of the overall uncertainties on the deuteron nuclear corrections is illustrated in Fig. 7. Within a simple single-scale model, relating the quark momentum distributions in the nucleon to the nucleon radius Kulagin and Petti (2006), the off-shell function $\delta f$ obtained from nuclei with $A\geq 4$ suggests an increase of the nucleon core radius by about 2% in the deuteron, taking an average virtuality of $-0.045$ from Table 1. This value is comparable to estimates obtained with a different model Close et al. (1985) in relation to an increase of the overlap of nucleons in nuclei with the nuclear density.

The predictions from both the present analysis and the ones from Ref. Kulagin and Petti (2006) are compared with the recent direct measurement of the ratio $F_{2}^{D}/F_{2}^{N}$ by the BONuS experiment Griffioen et al. (2015) in Fig. 8. A good agreement with BONuS data is observed, although the current experimental accuracy is somewhat limited and many data points fall in the resonance region with $W<1.6$ GeV or correspond to low $Q^{2}<1.5$ GeV2. The model-dependent extrapolation of the $F_{2}^{A}$ measurements with $A\geq 4$ performed by the SLAC E139 experiment Gomez et al. (1994) is also compared with our results and with the ones from Ref. Kulagin and Petti (2006) in Fig. 8. Although the basic assumption of Ref. Gomez et al. (1994) about the scaling of the magnitude of nuclear effects with the nuclear density Frankfurt and Strikman (1988) was excluded by the recent measurements on 9Be target by the E03-103 experiment Seely et al. (2009), the results of Ref. Gomez et al. (1994) with $x<0.7$ are consistent with our predictions. It is worth noting that the typical $Q^{2}$ values of the E139 data are substantially larger than in the BONuS sample (2-15 GeV2 vs. 1-4 GeV2). This difference allows to demonstrate the $Q^{2}$ dependence in the nuclear corrections, which appears mainly at large $x$ due to the combined effect of the TMC and the off-shell corrections, cf. Fig. 8.

Since nuclear corrections are almost linearly dependent from $x$ in the region $0.35<x<0.55$, they are often quantified by the corresponding linear slope. The main advantage of this slope is that it can be measured more accurately since it is not affected by the normalization uncertainties. The empirical model-independent determination of the slope $\mathrm{d}{\mathcal{R}}(D/N)/\mathrm{d}x=-0.100\pm 0.050$ Griffioen et al. (2015) of the BONuS data agrees well with the value $-0.099\pm 0.006$ predicted by Ref. Kulagin and Petti (2006). To this end, the model-dependent extrapolation of the SLAC E139 data Gomez et al. (1994) gives a consistent value of $-0.098\pm 0.005$, while the empirical extrapolation using the short range correlation scale factors from Ref. Weinstein et al. (2011) results in a somewhat smaller slope $-0.079\pm 0.006$. However, while being useful for the analysis of experimental data, the slope $\mathrm{d}{\mathcal{R}}(D/N)/\mathrm{d}x$ describes the behavior of the nuclear corrections in a limited region only. Meanwhile, the microscopic model of Ref. Kulagin and Petti (2006) reproduces not only the measured slopes, but also the shape and magnitude of the nuclear corrections in the entire kinematic range covered by existing data.

Our results on the off-shell correction differ from the ones obtained using a similar formalism in Refs. Accardi et al. (2011, 2016a). The analysis of Ref. Accardi et al. (2011) is based on a modified model of Ref. Kulagin and Petti (2006) to relate the off-shell function $\delta f$ to an increase in the nucleon confinement radius in the nuclear medium. The analysis of Ref. Accardi et al. (2016a) follows more closely the model of Ref. Kulagin and Petti (2006), with the off-shell function $\delta f$ being detrmined from a global QCD fit to the deuteron and proton data. The differences in the results of those analyses with respect to our study can be attributed to the implementation of the deuteron model and the details of the calculations Accardi et al. , as discussed in the Appendix.

As discussed in Sec. IV.1, the ratio of DIS structure functions for two different nuclear targets, $\mathcal{R}(A/B)=F_{2}^{A}/F_{2}^{B}$, offers a good tool to study the off-shell function $\delta f(x)$, due to a large cancellation of both experimental and model uncertainties. In this respect we note that the addition of data on nuclear ratios $\mathcal{R}(A/B)$ to global QCD fits should improve the determination of the nucleon off-shell function $\delta f$ and of the $d$ quark distribution in the proton.

IV.5 Constraints on $d/u$ and $F_{2}^{n}/F_{2}^{p}$

Correlations between the deuteron nuclear corrections and the $d$-quark distribution can substantially limit the PDF accuracy achievable in the PDF fits based on the proton and deuterium DIS data. In this context the data from flavor sensitive processes like $W^{\pm}$ production in $pp(\bar{p})$ collisions play a major role in reducing such correlations. A possible approach to avoid the effects of the deuteron nuclear corrections is to avoid any DIS data off the deuteron in global QCD fits, as in the recent ABMP16 analysis Alekhin et al. (2017). The corresponding results for the $d/u$ ratio shown in Fig. 9 indicate that the recent precision data on $W^{\pm}$ boson production from D0 and the LHC experiments (cf. Table 2) provide a good sensitivity to the $d$-quark distribution. In particular, the $d/u$ ratio at large $x>0.7$ is well constrained, mainly due to the large rapidity data from the recent LHCb measurement of $W^{\pm}$ boson production Aaij et al. (2015a, 2016). This sample indeed probes values of $x$ up to 0.8 and its accuracy is comparable to the one of DIS experiments.

The universality of $\delta f$ allows a further improvement of the accuracy in the determination of the $d/u$ ratio, by using the deuteron DIS data in combination with the more precise off-shell function obtained from the analysis of the nuclear targets with $A\geq 4$ (see Sec. IV.1). In Fig. 9 we show the $d/u$ ratio obtained in such a way in the ABMP16 fit. The impact of the DIS deuteron data on the $d/u$ ratio is more evident in the region of $x>0.4$, where the uncertainties are substantially reduced, as compared to the ABMP16 results obtained without the DIS deuteron data.

An interesting observation is that the $d/u$ ratio tends to vanish as $x\to 1$ (see Fig. 9). In order to verify that this behavior is not biased by the functional form of the PDF parametrization, we multiply the $d$- and $u$-quark parametrizations by an additional free polynomial. We do not find any significant impact on the corresponding large-$x$ behavior of the $d/u$ ratio, thus confirming its stability. Furthermore, the value of $d/u$ obtained in the present analysis is consistent with the one of the ABMP16 fit performed without using the deuteron data, cf. Fig. 9. These results indicate that our $d$- and $u$-quark parametrizations are flexible enough to be driven by the data sets, rather than by the functional form used (cf. Appendix).

The $d/u$ ratio is related to the neutron to proton structure function ratio, $F_{2}^{n}/F_{2}^{p}$, displayed in Fig. 10. The impact of the deuteron nuclear correction on $F_{2}^{n}/F_{2}^{p}$ is somewhat larger than on the $d/u$ ratio. Note, however, that the behavior of $F_{2}^{n}/F_{2}^{p}$ at $x\to 1$ is dominated by the HT contributions, which introduce a significant uncertainty on this ratio, as demonstrated in Fig. 10.

It is instructive to compare the nuclear corrections applied to the deuteron data in various PDF analyses. The CJ15 analysis Accardi et al. (2016a) is based upon a formalism similar to the model of Ref. Kulagin and Petti (2006), with the corresponding off-shell correction in the deuteron determined from the global QCD fit (see discussion in Sec. IV.4). The MMHT14 analysis Harland-Lang et al. (2015) is based on an empirical parametrization of the nuclear correction to the ratio $F_{2}^{D}/F_{2}^{N}$, which is extracted from data. The CT14 Dulat et al. (2016) and NNPDF3.0 Ball et al. (2015) do not use any nuclear correction arguing that nuclear corrections would introduce additional uncertainties in the analysis Ball et al. (2015). Furthermore, it was claimed that nuclear corrections can be neglected when using more stringent cuts in $Q^{2}$ and $W^{2}$ Dulat et al. (2016). The nuclear corrections used by CJ15 and MMHT14, while mutually consistent, are characterized by relatively large uncertainties (cf. Fig. 11). Since these two determinations are driven by the available deuterium data and are largely correlated with the $d/u$ ratio, the possibility to improve such uncertainties in the context of global QCD fits appears to be limited. As discussed above, the use of the microscopic model of Ref. Kulagin and Petti (2006) and of the corresponding results on the universal off-shell function $\delta f$, allows a substantial reduction of the uncertainties related to the nuclear correction in the deuteron. Such a reduction is illustrated by the comparison with the model of Ref. Kulagin and Petti (2006) given in Fig. 11. In turn, the uncertainty on the deuteron nuclear correction is a bottleneck in the overall accuracy of the $d/u$ ratio obtained from global QCD fits. The differences in the treatment of the nuclear corrections in the deuteron illustrated in Fig. 11 would translate into a corresponding spread in the $d$-quark distribution. It is worth noting that these systematic effects cannot be mitigated by more stringent $Q^{2}$ and $W^{2}$ cuts, since nuclear effects survive even at very large energy and momentum, as demonstrated in DIS experiments (for a review see Refs. Arneodo (1994); Norton (2003)) and by recent observations of nuclear modifications in $p+{\rm Pb}$ and ${\rm Pb}+{\rm Pb}$ collisions at the LHC Khachatryan et al. (2015, 2016b); Aad et al. (2015); Aaij et al. (2014); Senosi (2015); Aad et al. (2013); Chatrchyan et al. (2015).

The largest deviations between the CJ15 and MMHT14 analyses and the model of Ref. Kulagin and Petti (2006) are observed in the intermediate region $x\sim 0.15$ and at large $x>0.6$ (cf. Fig. 11), although their significance is limited by the current uncertainties. Our analysis of deuteron data is consistent with Ref. Kulagin and Petti (2006), with an uncertainty band extending close to the CJ15 and MMHT14 ones (cf. Fig. 7). The small enhancement present in the model of Ref. Kulagin and Petti (2006) for $0.05<x<0.3$ is the result of an interplay of the off-shell correction and the meson exchange currents (cf. Fig. 2). Nuclear corrections at $x>0.6$ are instead dominated by the FMB and OS corrections. This kinematic region is very sensitive to the treatment of both the bound nucleon momentum distribution and the target mass corrections to the nucleon structure functions (Sec. II). It is worth noting that the prescription of Ref. Georgi and Politzer (1976) for the TMC is known to have an incorrect behavior for $x\to 1$, which can affect the calculations at very large $x$ values. More detailed comparisons with the CJ15 and MMHT14 results can be found in the Appendix. Future DIS measurements from the BONuS experiment with the 12 GeV JLab upgrade Bueltmann et al. (2012), from neutrino and antineutrino scattering off free proton in the DUNE experiment Acciarri et al. (2015), and from the electron-ion collider Accardi et al. (2016b) can further improve our understanding of the ratio $F_{2}^{D}/F_{2}^{N}$.

V Summary

We performed a study of nuclear effects in the deuteron using the data from DIS off proton and deuterium, Drell-Yan production in $pp$ and $p{\rm D}$ interactions, and $W^{\pm}$- and $Z$-boson production in $pp$ and $p\bar{p}$ collisions in the context of global QCD fits. We found that it is possible to determine simultaneously PDFs, high twist terms, and the off-shell function describing the modification of PDFs in bound nucleons due to their different $Q^{2}$ dependence and the wide kinematic coverage of existing data. Flavor sensitive processes like $W^{\pm}$ production in $pp(\bar{p})$ collisions play an important role in disentangling the impact of the nuclear corrections in the deuteron from the $d$ quark distribution function, allowing a more accurate determination of both contributions. We also evaluated the sensitivity of our results to various models of the deuteron wave function and found that the corresponding model dependence is reduced by the recent BONuS measurement of the ratio $F_{2}^{D}/F_{2}^{N}$.

The results on the off-shell function $\delta f$ reported in this paper are in good agreement with the earlier determination obtained in the analysis of the ratio of nuclear structure functions with mass number $A\geq 4$ Kulagin and Petti (2006). This result confirms the universality of $\delta f$, which can be regarded as a special structure function of the nucleon describing the modification of the bound nucleons in the nuclear medium. This study supports the unified treatment of the deuteron and heavier nuclei developed in Ref. Kulagin and Petti (2006).

We also studied the impact of deuteron nuclear corrections on the $d/u$ ratio within global QCD fits. We found that the recent precision data on $W^{\pm}$ boson production from D0 and the LHC experiments allow a reduction of the uncertainties in the $d/u$ ratio at large $x$. Our results indicate that the accuracy in the determination of the $d/u$ ratio can be further substantially improved by including the DIS data off a deuterium target corrected for nuclear effects using the model of Ref. Kulagin and Petti (2006) with the universal off-shell function $\delta f$.

Acknowledgements.

We thank F. Gross and W. Polyzou for providing the parametrizations of the WJC and AV18 deuteron wave functions, respectively. We thank R. Machleidt for the tables of the CD-Bonn deuteron wave function. We thank A. Accardi, W. Melnitchouk, S. Moch, and J. Owens for fruitful discussions. R.P. thanks the II Institute for Theoretical Physics at the University of Hamburg for hospitality during the manuscript preparation. S.A. was supported by Bundesministerium für Bildung und Forschung (Contract No. 05H15GUCC1). R.P. was supported by Grant No. DE-SC0010073 from the Department of Energy, USA. S.K. was supported by the Russian Science Foundation Grant No. 14-22-00161.

Appendix A Results from different analyses

A few phenomenological studies aimed to extract the nuclear corrections in the deuteron from global QCD fits are available in literature Accardi et al. (2011); Harland-Lang et al. (2015); Accardi et al. (2016a), including the analyses of Refs. Accardi et al. (2011, 2016a) based on a deuteron model similar to that used in our studies. In this Appendix we investigate the differences with respect to those studies Accardi et al. and benchmark the recent CJ15 analysis Accardi et al. (2016a).

Clear differences with the CJ15 fit appear in the data sets used and in the proton PDFs correspondingly obtained. From Table I of Ref. Accardi et al. (2016a) we conclude that most of the sensitivity to the nuclear corrections in the CJ15 fit manifests for the D0 data on the $W^{\pm}$-boson production asymmetry and the DIS deuteron data by the SLAC experiments (cf. reduction of $\chi^{2}$ obtained by adding the nuclear corrections). While we include the SLAC DIS data (Table 2), we use significantly different $W,Z$-boson collider data in our analysis (see Sec. III.2). More specifically, we include the D0 data on the lepton asymmetry from the $W$-boson decays rather than the actual $W$-boson asymmetry data, as well as the recent LHC DY data Aad et al. (2012, 2016); Chatrchyan et al. (2014); Khachatryan et al. (2016a); Aaij et al. (2015a, b, 2016). In order to test the impact of such differences, we perform a variant of our fit in which we drop all the $W,Z$-boson collider data, replacing them with the D0 $W$-asymmetry data Abazov et al. (2014). The corresponding results for the function $\delta f$ are consistent with the ones presented in Sec. IV.1 (cf. Fig. 12). We also perform a separate NLO fit to all deuteron data sets using fixed CJ15 proton PDFs from the LHAPDF library Buckley, Andy and Ferrando, James and Lloyd, Stephen and Nordström, Karl and Page, Ben and Rüfenacht, Martin and Schönherr, Marek and Watt, Graeme (2015). The off-shell function $\delta f$ obtained in this way is consistent with our results in Sec. IV.1 and Ref. Kulagin and Petti (2006) (cf. Fig. 12). Therefore, we can conclude that the differences in the data samples and proton PDFs cannot explain the different off-shell correction obtained in Ref. Accardi et al. (2016a).

Meanwhile, it is instructive to compare the values of the $d/u$ ratio obtained in the various fits. The $d/u$ ratio obtained from the variant of our fit with the D0 $W$-asymmetry data is consistent within uncertainties with the results presented in Sec. IV.5 (cf. Fig. 13). The differences between the central values suggest that our PDF parametrizations are flexible enough to describe different data sets without limitations from the functional form used. In particular, we do not explicitly constrain the $d/u$ ratio to vanish for $x\to 1$. The flexibility of our PDF parametrization is confirmed by the fact that we obtain similar results by multiplying the $d$- and $u$-quark parametrizations by an additional free polynomial (cf. Sec. IV.5). The uncertainty of the $d/u$ ratio determined in the MMHT14 analysis Harland-Lang et al. (2015) at large $x$ does not allow quantitative comparisons (cf. Fig. 13), due to the lack of relevant experimental constraints. It is worth noting that the $d/u$ ratio of the CJ15 analysis Accardi et al. (2016a) displays a substantially different behavior at large $x$. Since this kinematics is largely controlled by the data on $W$-boson asymmetry at large rapidity, we compute this quantity at the NLO approximation using the FEWZ package Li and Petriello (2012); Gavin et al. (2013) with CJ15 proton PDFs from the LHAPDF library Buckley, Andy and Ferrando, James and Lloyd, Stephen and Nordström, Karl and Page, Ben and Rüfenacht, Martin and Schönherr, Marek and Watt, Graeme (2015). A comparison with the D0 data Abazov et al. (2014) (c.f. Fig. 13) shows deviations from these predictions at large values of the $W$-boson rapidity, in contrast with the corresponding results of the CJ15 analysis (cf. Fig. 13 of Ref. Accardi et al. (2016a)).

A possible source of differences between this analysis and Ref.Accardi et al. (2016a) may stem from the implementation of the convolution model, as well as from the treatment of the TMC and/or the HT contributions. In order to check the sensitivity to these effects we compute the ratio $F_{2}^{D}/F_{2}^{N}$ using the $\delta f(x)$ function from Ref.Accardi et al. (2016a) and the AV18 wave function Veerasamy and Polyzou (2011). The result shown in Fig. 14 (left panel) deviates significantly from both the one presented in Sec.IV.1 and the one of Ref. Accardi et al. (2016a) using the same deuteron wave function, suggesting a different implementation of the convolution model in Ref.Accardi et al. (2016a). Figure 14 (right panel) also illustrates the sensitivity of the ratio $F_{2}^{D}/F_{2}^{N}$ to various implementations of the TMC and HT corrections. 222Both TMC and HT corrections are characterized by a strong $Q^{2}$ dependence, increasing their impact at lower $Q^{2}$ values. In particular, we compare the standard TMC scheme of Ref. Georgi and Politzer (1976) (on-shell TMC) with its off-shell continuation by Eq.(3). The TMC and HT corrections to $F_{2}^{D}/F_{2}^{N}$ are treated differently in the CJ15 and MMHT14 analyses, and could contribute to the disagreement in the region of large $x$ (cf. left panel of Figure 14).

In summary, the present studies indicate that we cannot reproduce the CJ15 results of Ref. Accardi et al. (2016a) on the function $\delta f$. All our systematic checks are consistent with the determination presented in Sec. IV.1. The differences with the results of Ref.Accardi et al. (2016a) cannot be explained by the different data samples nor by the PDFs used. Instead, we find indications pointing towards the implementation of the deuteron model, TMC and HT corrections in the CJ15 fit.