QCD analysis of $W$- and $Z$-boson production at Tevatron

Recent measurements of the $W$-boson charge asymmetry and of the $Z$-boson production cross sections, performed at the Tevatron collider in Run II by the D0 and CDF collaborations, are studied using the HERAFitter framework to assess their impact on the proton parton distribution functions (PDFs). The Tevatron measurements, together with deep-inelastic scattering data from HERA, are included in a QCD analysis performed at next-to-leading order, and compared to the predictions obtained using other PDF sets from different groups. Good agreement between measurements and theoretical predictions is observed. The Tevatron data provide significant constraints on the $d$-valence quark distribution.


Introduction
Accurate knowledge of the parton distribution functions (PDFs) is essential for predictions at hadron colliders. The primary source of information on the proton PDFs comes from deep-inelastic scattering (DIS). Measurements at fixed target experiments and at the HERA e ± p collider provide constraints on the quark and gluon densities, and discrimination of the quark flavours. The DIS proton data mostly constrain the utype quark density, due to the greater couplings to the photon at low absolute four momentum transfers, Q 2 , whereas the d-type quark densities are only constrained at high Q 2 with limited precision. Even more challenging is the separation of the d-valence quark density which relies on the HERA e + charge current data, which are statistically limited in the published HERA I combined data [1]. A better flavour separation is needed to challenge the limits of precision physics at the LHC.
Drell-Yan production of W and Z bosons in protonantiproton and proton-proton collisions can provide additional information on the d-quark PDFs. At leading order (LO) in QCD, the Drell-Yan processes probe the PDFs at energy scales Q corresponding to the boson masses, m V = m W and m V = m Z , and momentum fractions carried by the interacting partons of x 1,2 = m V / √ Se ±y , where √ S is the centre-of-mass energy and y is the boson rapidity.
At the Tevatron proton-antiproton collider, the production of W and Z bosons is dominated by valencequark interactions. The Z-boson production has similar couplings for uū and dd fusion processes, whereas W bosons are produced predominantly by ud and dū fusions for W + and W − bosons, respectively. Various measurements of Z-boson inclusive production and of W -boson charge asymmetry have been reported by the D0 and CDF collaborations [2][3][4][5][6][7]. Some of these data samples were included in previous PDF studies [8][9][10][11]. The addition of the Tevatron data resulted in improved PDFs, but some tensions were observed with global PDF fits [12][13][14][15].
In this paper the data collected at the Tevatron collider in Run II are analysed to assess their impact on the PDFs. The assumptions of the correlation model of the experimental systematic uncertainties are revised with respect to the recommendation of the experiments, leading to improved agreement with the theoretical predictions. The analysis is performed using the HERA-Fitter framework [1, [16][17][18] at next-to-leading order (NLO) QCD. The Tevatron W -and Z-boson measurements are also compared to predictions evaluated with the recent PDF sets CT10nlo [8], MMHT2014 [9] and NNPDF3.0 [10]. The impact of the Tevatron data on PDFs is studied using hessian profiling [19] and Bayesian reweighting [20][21][22] techniques. The profiling of PDF uncertainties is generalised to the case of asymmetric PDF uncertainties.
This paper is organised as follows: the data samples are introduced in Sec. 2 and the theoretical predictions are discussed in Sec. 3. The QCD analysis settings and the methods for comparing data with predictions based on existing PDFs are discussed in Sec. 4. Section 5 reports the results of the PDF analysis. The results obtained in the paper are summarised in Sec. 6.

Data Sets
The most recent measurements of W -boson charge asymmetry and Z-boson inclusive production performed in Run II of the Tevatron collider are considered in this study. They include the Z-boson differential cross section as a function of rapidity, measured by the D0 collaboration with 0.4 fb −1 of integrated luminosity in the Z → ee channel [2]; the Z-boson differential cross section as a function of rapidity, measured by the CDF collaboration with 2.1 fb −1 of integrated luminosity in the Z → ee channel [3]; the charge asymmetry of muons as a function of rapidity in W → µν decays, measured by the D0 collaboration with 7.3 fb −1 of integrated luminosity [4]; the W -boson charge asymmetry as a function of rapidity in the W → eν decay channel, measured by the CDF collaboration with 1 fb −1 of integrated luminosity [5]; the W -boson charge asymmetry as a function of rapidity in the W → eν decay channel, measured by the D0 collaboration with 9.7 fb −1 of integrated luminosity [6]. These measurements supersede the previous Run II Tevatron measurements of Wboson charge asymmetry and Z-boson inclusive produc-tion. Recently, the D0 collaboration has also released a measurement of the charge asymmetry of electrons as a function of rapidity in W → eν decays [7]. However, this measurement is performed with the same data set and event selection as the measurement of Ref.
[6], and cannot be included simultaneously in a PDF fit without provision of the correlation information. The Tevatron W -and Z-boson measurements considered in this study are summarised in Tab. 1.
Besides the Tevatron W -and Z-boson measurements, the HERA I combined measurements of the inclusive DIS neutral-and charged-current cross sections measured by the H1 and ZEUS experiments [1] are used in this study. The neutral-current measurements cover a wide range in Bjorken-x and Q 2 , which is essential for the determination of PDFs, whereas the charged-current measurements provide further information to disentangle the contributions in PDFs from u-type and d-type quarks and anti-quarks at x > 0.01. The DIS data are required to be in the kinematic region Q 2 > Q 2 min = 7.5 GeV 2 , where perturbative QCD calculations are reliable.

Experimental Uncertainties
Statistical uncertainties are considered to be uncorrelated between bins, with the exception of the D0 measurement of W -boson charge asymmetry of Ref.
[6], for which bin-to-bin statistical correlations are provided.
In general, the correlation model of the experimental uncertainties recommended by the Tevatron experiments is adapted and followed in the QCD analysis, with the exception of the experimental systematic uncertainties related to trigger and lepton identification efficiencies. These uncertainties are provided by the D0 and CDF experiments in the form of total uncertainties in each bin of the measurements. However, the trigger and lepton identification corrections are estimated from data, and they are influenced, among other effects, by statistically uncorrelated bin-to-bin fluctuations. Since the exact bin-to-bin correlation pattern of these uncertainties is not provided by the experiments, a conservative approach is followed in this study, and the uncertainties related to trigger and lepton identification efficiencies are treated as uncorrelated bin-to-bin for the nominal fit. According to this prescription, the following uncertainties are treated as uncorrelated binto-bin: the central-and forward-electron identification efficiencies of Ref.
[4], the trigger and electron identification efficiencies of Ref. [5], and the electron identification, charge misidentification and positron to electron efficiency corrections of Ref. [7].
All the other experimental systematic uncertainties are considered fully correlated bin-to-bin, with the exception of the D0 measurement of W -boson charge asymmetry of Ref.
[6], where the total experimental systematic uncertainty is treated as bin-to-bin uncorrelated, as recommended by the D0 experiment, and for the electron charge asymmetry in W → eν decays of Ref. [7], where the uncertainty of the unfolding procedure due to the limited statistics of the Monte Carlo (MC) sample is treated as uncorrelated bin-to-bin. The dependence of the measured asymmetry on the PDF set used to reconstruct the W -boson rapidity was studied in the D0 measurement of W -boson charge asymmetry of Ref.
[6]. In this paper, the W -boson charge asymmetry extracted with the CTEQ6.6 PDF set is used as the central value, and the 22 CTEQ6.6 positive and negative PDF eigenvector variations are considered as bin-to-bin correlated systematic uncertainties.
For the two measurements of Z-boson differential cross section as a function of rapidity, statistical uncertainties are scaled to the expected number of events assuming they are Poisson distributed. The experimental systematic uncertainties are treated as multiplicative, and linearly scaled to the expected cross sections, except for the background uncertainties which are treated as additive, and are not scaled. For the measurements of W -boson and lepton charge asymmetry, all the uncertainties are treated as additive, and are not scaled.
The statistical uncertainties of the HERA I data are treated as uncorrelated and scaled to the expected number of events assuming Poisson distribution, whereas the experimental systematic uncertainties are fully correlated and are scaled linearly to the expected cross sections.

Theoretical Predictions
The theoretical predictions corresponding to the Tevatron measurements of W -boson charge asymmetry and Z-boson inclusive production are included in the fits using APPLGRID [23,24] files. These predictions have been evaluated with MCFM [25,26] at NLO QCD according to the phase-space definitions of each measurement, which are as follows: the D0 and CDF measurements of the Z-boson differential cross section as a function of rapidity are defined in the full kinematic range of the decay leptons, without any requirements on the rapidity and p T of the leptons. In the D0 measurement, the invariant mass of the dielectron system is defined in the range 71 < m ee < 111 GeV, whereas in the CDF measurement, it is defined in the range 66 < m ee < 116 GeV. The charge asymmetry of muons as a function of rapidity in W → µν decays, and the charge asymmetry of electrons in W → eν decays, measured by the D0 experiment, are defined with p ℓ T > 25 GeV and p ν T > 25 GeV. The W -boson charge asymmetry as a function of rapidity in the W → eν decay channel, measured by the CDF collaboration, is defined in the full kinematic range, without any requirements on the lepton rapidity and p T . The corresponding D0 measurement of the W -boson charge asymmetry in the W → eν decay channel is defined in a kinematic region where the charged lepton and the neutrino are required to have p T > 25 GeV without further requirements on the lepton rapidity. The kinematic requirements of the Tevatron W -and Z-boson measurements are summarised in Tab. 1. Notice that the CDF and D0 measurements of W -boson charge asymmetry in the W → eν decay channel of Refs.
[5] and [6] are defined in different kinematic regions and they should not be compared without extrapolating them to a common phase space. Tables of the Tevatron measurements, with updated correlation model, and corresponding APPLGRID theoretical predictions are publicly available at herafitter.org.
The QCD predictions for the DIS cross sections are evaluated by solving the DGLAP evolution equations [27][28][29][30][31][32] at NLO in the M S scheme [33] using the QCDNUM program [34] with the renormalisation and factorisation scales set to Q 2 . The light quark coefficient functions are calculated in QCDNUM. The heavy c-and b-quark distributions are dynamically generated, and the corresponding coefficient functions for the neutral-current processes with γ * exchange are calculated in the general-mass variable-flavour-number scheme [35][36][37], with up to five active quark flavours. For the charged-current processes and the neutralcurrent processes with a Z contribution, the heavy quarks are treated as massless.

PDF Fit Settings
The QCD analysis and PDF extraction is performed with the open-source framework HERAFitter. The charm mass is set to m c = 1.38 GeV, as estimated from HERA charm production cross section [38], and the bottom mass to m b = 4.75 GeV [39]. The stronginteraction coupling constant at the Z boson mass, α s (M Z ), is set to 0.118, and two-loop order is used for the running of α s .
The PDFs for the gluon, u-valence, d-valence,ū, d quark densities are parametrised at the input scale Q 2 0 = 1.7 GeV 2 as follows: The contribution of the s-quark density is taken to be proportional to thed-quark density by setting xs(x) = r s xd(x), with r s = 1.0, as suggested in Ref.
[40]. The strange and anti-strange quark densities are taken to be equal: , is determined by the quark-counting sum rule, whereas the normalisation of the gluon density, A g , is determined by the momentum sum rule. The x → 0 limit of the uto d-sea quark densities is fixed to unity by setting Bū = Bd and Aū = Ad. A χ 2 function used for the data to theory comparison is defined as in Ref.
[1], with an additional penalty term as described in Ref. [41], and minimised with MINUIT [42] to extract the PDFs from the data.

PDF Profiling
The impact of a new data set on a given PDF set can be quantitatively estimated with a profiling procedure [19]. The profiling is performed using a χ 2 function which includes both the experimental uncertainties and the theoretical uncertainties arising from PDF variations: The correlated experimental and theoretical uncertainties are included using the nuisance parameter vectors β exp and β th , respectively. Their influence on the data and theory predictions is described by the Γ exp ij and Γ th ik matrices. The index i runs over all N data data points, whereas the index j (k) corresponds to the experimental (theoretical) uncertainty nuisance parameters. The measurements and the uncorrelated experimental uncertainties are given by σ exp i and ∆ i , respectively, and the theory predictions are σ th i . The χ 2 function of Eq. 2 can be generalised to account for asymmetric PDF uncertainties: are determined from the shifts of predictions corresponding to up (Γ th+ ik ) and down (Γ th− ik ) PDF uncertainty eigenvectors.
The minimisation of Eq. 2 in its original form leads to a system of linear equations. The generalised function, with asymmetric PDF uncertainties, is minimised iteratively: the values of Γ th+ ik are updated using β k,th from the previous iteration and following the substitution of Eq. 3. Several iterations are required to converge, and the procedure is verified using the MINUIT program which yields identical results.
The value at the minimum of the χ 2 function provides a compatibility test of the data and theory. In addition, the values at the minimum of the nuisance parameters β min k,th can be interpreted as optimisation ("profiling") of PDFs to describe the data [19]. Explicitly, the profiled central PDF set f ′ 0 is given by where f 0 is the original central PDF set and f ± k represents the eigenvector sets corresponding to up and down variations.
The shifted PDFs have reduced uncertainties. In general, the shifted eigenvectors are no longer orthogonal, but can be transformed to an orthogonal representation using a standard diagonalisation procedure, as in Ref. [43]. In this method the covariance matrix C of the PDF nuisance parameters is diagonalised as where G is an orthogonal matrix, D is a positive definite diagonal matrix, and √ D is a diagonal matrix built of √ D ii . The matrices G and D can be constructed using the eigenvectors and eigenvalues of the matrix C. The transformation G ′ can be adjusted, using orthogonal transformations, to keep the new eigenvector basis aligned along the original as much as possible. As a result of this adjustment, the transformation matrix can take a triangular form with all diagonal elements greater than zero.
The method can be extended to PDF sets with asymmetric uncertainties: the transformation matrix is determined using symmetrised uncertainties as in Eq. 5, and the orthogonal up and down PDF eigenvectors f + ′ i and f − ′ i are calculated as

Bayesian Reweighting
An alternative approach to assess the impact of new data on PDFs is the Bayesian reweighting technique, first proposed in Ref. [20] and further developed by the NNPDF collaboration [21,22]. The Bayesian reweighting can be applied to PDF sets provided in the form of MC replicas, such as the NNPDF3.0 set [10]. Recently, a variant of the method which can be used with PDFs provided in the eigenvector representation has been developed [44] and is also available in HERAFitter.
The Bayesian reweighting is based on the assumption that an ensemble of MC replicas provides a representation of the probability distribution in the space of PDFs. For a given PDF set with N rep replicas {f k }, with k = 1, 2, ..., N rep , the central value for a general observable, O({f k }), is estimated as the average of the predictions obtained from the ensemble: With the inclusion of new data, the probability distribution associated with the original PDF set is modified according to Bayes Theorem. For each replica k, a weight w k is obtained from the χ 2 function according to: where N data is the number of new data points and χ 2 k is the χ 2 value between data and predictions corresponding to the k-th PDF replica.
The prediction for a given observable, after the inclusion of the new data, is evaluated as the weighted average of predictions obtained from the ensemble: The reweighting procedure is very fast and results in a new, updated, MC PDF set. Some of the replicas of the PDF set may have very small weights (typically those which do not describe the new data), and they do not contribute to the ensemble any longer. The number of effective replicas, N eff , of a reweighted set is quantified by the Shannon entropy An un-weighting procedure can be performed on the MC set such that PDFs with small weights are suppressed and a new set is produced, which has unit weight for all PDF replicas in addition to statistically reproducing the averages from Eq. 8.

Results
The QCD fit analysis described in Sec. 4 is performed on the Tevatron W -and Z-boson data, together with the HERA I data. The fit is used to study the compatibility of the data with NLO QCD predictions, and to assess the impact of the Tevatron data on PDFs.
The profiling and the reweighting techniques are used to asses the impact of the Tevatron data on various PDF sets. The optimal parametrisation for the PDF fit is found through a parametrisation scan, a procedure first introduced in Ref.
[1]. The scan is performed by starting from a parametrisation with a basic polynomial form, where D f , E f , and the exponential parameters F f of Eq. 1 are set to zero. After application of the quarkcounting and momentum sum rules, and of the x → 0 constraints onū andd, the initial PDFs parametrisation has 10 free parameters. The 15 D f , E f and F f additional parameters are allowed to vary, one parameter at a time, and the parameter which induces the largest reduction of χ 2 min is added as a free parameter for the next iteration of the scan. The PDF fits which lead to solutions with negative high-x PDFs are discarded. For each PDF, the exponential term e F f x and the polynomial term (1+D f x+E f x 2 ) are considered as mutually exclusive, that is, when the exponential term is preferred, the polynomial term is no longer considered in the scan, and vice versa. The procedure is stopped when the reduction in the χ 2 min value, ∆χ 2 min , is less than unity. Table 2 shows the results of the parametrisation study: the parameters which induce the largest ∆χ 2 min are, in order, F dv , F uv , D g , Dd, and Dū. The optimal parametrisation found with this procedure has 15 free parameters, and the PDFs are expressed as: The parametrisation of Eqs. (10-14) is used for a fit to the HERA I data, and for a combined fit to the HERA I and Tevatron W -and Z-boson data. Table 3 shows the χ 2 min per degrees of freedom (dof) of the two fits. The contribution to the total χ 2 min of each data set, henceforth referred to as partial χ 2 , is also shown. The inclusion of the Tevatron W -and Z-boson data in the fit, which corresponds to 93 additional points, results in an increase of about 110 in the overall χ 2 min of the fit, and the partial χ 2 per number of points of each of the Tevatron and HERA I data set is close to unity. Figures 1 and 2 show the Tevatron Z-and Wboson measurements, respectively, compared to the theoretical predictions evaluated with the PDFs extracted from the combined fit to the HERA I and Tevatron data.
The central value and the uncertainties of the PDFs are evaluated with MC replicas [45]: the data points are smeared using Gaussian distributions, according to their experimental uncertainties, and the PDF fit is repeated 1000 times, using different random seeds for the smearing. The central PDFs are calculated as the average of the replicas and the PDF uncertainties are calculated from their standard deviation. Figure 3 shows the comparison of the PDFs extracted with the MC-replica method by fitting the HERA I data, and by fitting the HERA I and Tevatron W -and Z-boson data. Figure 4 shows the comparison of the relative uncertainty of the two PDFs. A significant reduction of the PDF uncertainties is observed in the fit which includes the Tevatron W -and Z-boson measurements, in particular for the valence quarks andd quarks.
A fit of the HERA I and Tevatron W -and Z-boson data with the same settings, but with a correlation model in which trigger and identification uncertainties are treated as correlated bin-to-bin, yields very similar central PDFs and PDF uncertainties [not shown]. The χ 2 of all the data sets are also very similar, except for the χ 2 of the CDF W -boson asymmetry measurements, which is about twice as large.
The W -boson charge asymmetries rely on the reconstruction of the W -boson rapidity, which is measured assuming a fixed W -boson mass, and inferring the unmeasured longitudinal momentum of the neutrino on a statistical basis [5, 6]. The reconstruction of the Wboson rapidity introduces a model dependence in the measurement. To study the possible bias due to the W -boson rapidity reconstruction, an alternative fit is performed in which the W -boson charge asymmetries measured by CDF and D0 are excluded, and the latest D0 measurement of the electron asymmetry is included. The χ 2 min / dof and the partial χ 2 of the fit are shown in Tab. 5. Also for this fit the partial χ 2 of each of the Tevatron and HERA I data set is close to unity. The d v PDF determined from the fit is shown in Fig. 5, and compared to the nominal fit. The fit to the lepton asymmetry data yields very compatible results, but the uncertainties on the d v PDF are up to twice as large.
The impact of posterior inclusion of the Tevatron W -and Z-boson measurements on the PDF uncertainties as estimated by CT10nlo, MMHT2014, and NNPDF3.0 is assessed by profiling and reweighting. For consistency with the other PDFs, the uncertainties of the CT10nlo PDFs are scaled to 68% confidence interval by applying a factor of 1.645. The three PDF sets already include the CDF and D0 Z-boson differential cross sections as a function of rapidity, and the MMHT2014 fit also includes the D0 muon charge asymmetry in W → µν decays and the CDF W charge asymmetry in the W → eν decay channel. Only the measurements that are not included in each of the PDF sets are considered for the corresponding profiling or reweighting. The compatibility of the Tevatron data with the CT10nlo, MMHT2014 and NNPDF3.0 sets is tested by evaluating the χ 2 function of Eq. (2), accounting for asymmetric PDF uncertainties according to Eq. (3). To perform this calculation for the NNPDF3.0 set, the covariance matrix for the predictions is decomposed using the eigenvector representation. Table 4 shows the compatibility between the Tevatron measurements and the above PDF sets, together with the partial χ 2 of each data set. The partial χ 2 per number of points of each of the Tevatron data set, and the total χ 2 / dof, are close to unity for all the PDFs, when the χ 2 evaluation includes the PDF uncertainties. The quality of the agreement significantly deteriorates if the χ 2 evaluation neglects the PDF uncertainty. This effect is more pronounced for the CT10nlo and NNPDF3.0 sets which include fewer data from the Tevatron. This indicates the significant constraining power of the Tevatron data.
The CT10nlo and MMHT2014 PDFs are profiled according to Eq. (4). The results of the profiling on the d-valence PDFs, and on their relative uncertainties, are shown in Fig. 6. The profiling affects the shape of the distribution more for the CT10nlo when compared to MMHT2014 set. Significant reduction of the uncertainties is observed for both sets, in particular in the low-and medium-x range. The NNPDF3.0 PDFs are reweighted to the Tevatron data. The number of effect-ive replicas remaining after reweighting, N eff , is only 1 and hence the resulting PDFs are not shown.
The original and profiled d-valence PDFs, and the result of the fit to the HERA I and Tevatron W -and Zboson data, are compared in Figs. 7 and 8, respectively. The profiling using Tevatron data improves agreement of the d-valence distribution between the MMHT2014 and CT10nlo PDF sets.

Summary
The HERAFitter framework is used to perform a QCD analysis of the DIS data from HERA, together with Wand Z-boson production measurements performed at the Tevatron collider in Run II. The correlation model of the systematic uncertainties of the Tevatron data is investigated, and a modification is proposed which accounts for the statistical nature of some of the systematic uncertainties. The Tevatron and HERA data are well described by a NLO fit with 15 free parameters, with a new parametrisation of PDFs which adds to the basic form a combination of linear and exponential terms. The impact of the Tevatron W -and Z-boson measurements is assessed by comparing PDF uncertainties from a fit to the HERA data alone, and a fit to the HERA and Tevatron data. A significant reduction of the uncertainties is observed in the latter case, for the valence quarks andd quarks in particular.
The Tevatron measurements are also compared to predictions evaluated with modern PDF sets, and the impact of the data on the PDFs is assessed using profiling and reweighting techniques. The profiling techniques take into account asymmetric PDF uncertainties. A good agreement between measurements and predictions is observed, if the PDF uncertainties of the predictions are taken into account. After the inclusion of the Tevatron data, the PDF uncertainties on the dvalence quarks are significantly reduced, especially for the PDF sets which include only the Z-boson data from the Tevatron, and the agreement between the various PDF sets is improved. These findings highlight the importance of the Tevatron W -and Z-boson production data to constrain d-quark and valence PDFs, and suggest that the data should be used in the future global PDF analyses. All the supporting material to allow fits of the Tevatron data, including the updated correlation model and the grid files for fast theory calculations, are publicly available on the web page of the HERAFitter project.   Table 2 Results of the parametrisation study. For each additional free parameter D, E, and F , of the d v , u v , gluon,ū, and d PDF, the reduction of χ 2 min of a fit to the Tevatron and HERA I data, ∆χ 2 min , is shown. For each of the fit with n free parameters, with n = 10, 11, 12, 13, 14, the largest ∆χ 2 min is shown in bold, and the corresponding parameter is added as a free parameter for the n + 1-parameters fit. The fits which leads to negative high-x PDFs are shown in parenthesis, and are not considered in the parametrisation study.      Figure 1 Theoretical predictions evaluated with the PDFs extracted from a fit to the HERA I and Tevatron data are compared to (a) Z-boson differential cross section as a function of rapidity, measured by the D0 collaboration and (b) Z-boson differential cross section as a function of rapidity, measured by the CDF collaboration. The red continuous lines correspond to the theoretical predictions, the red dashed lines are the theoretical predictions shifted by the experimental shift terms j Γ exp ij β j,exp of Eq.
(2). The yellow bands show the total experimental uncertainty, the black vertical bars show the quadratic sum of statistical and systematic uncorrelated uncertainties. Note, the theoretical predictions and the theoretical predictions shifted by the experimental shift terms are nearly identical, and virtually indistinguishable in these plots.