Progress in computing parton distribution functions from the quasi-PDF approach

We discuss the current developments by the European Twisted Mass Collaboration in extracting parton distribution functions from the quasi-PDF approach. We concentrate on the non-perturbative renormalization prescription recently developed by us, using the RI$'$ scheme. We show results for the renormalization functions of matrix elements needed for the computation of quasi-PDFs, including the conversion to the $\overline{\rm MS}$ scheme, and for renormalized matrix elements. We discuss the systematic effects present in the $Z$-factors and the possible ways of addressing them in the future.


Introduction
Parton distribution functions (PDFs) are one of the most important tools describing the rich internal structure of hadrons, resulting from the strong interactions of valence quarks with each other and with dynamically created sea quarks and antiquarks, via exchange of gluons. Due to the large value of the strong coupling constant at the relevant energy scales, PDFs are non-perturbative and thus, in principle, can be studied quantitatively only in the framework of Lattice QCD. Such an evaluation from first principles would constitute an important test of non-perturbative aspects of QCD. Moreover, it could provide important data for experiments, in particular for the transversity PDFs (and to some extent also for polarized PDFs), which are rather poorly constrained by phenomenological fits using perturbation theory. Even for the well-constrained experimentally unpolarized case, there are still some regions that are not kinematically accessible, such as the large Bjorken x region [1,2]. However, since PDFs are defined on the light cone, the standard lattice Euclidean formulation allows only to compute their low Mellin moments (cf. [3][4][5][6][7] for recent reviews). Attempts to reconstruct full PDFs from these moments are bound to fail, due to both practical problems (unfavourable signal-tonoise ratio for higher moments) and theoretical ones (inevitable power divergent mixings with lower dimensional operators).
A thoroughly investigated alternative to lattice computations of the moments is the method suggested by Ji [8]. This approach does not rely on the light cone definition, but instead one computes spatial matrix elements using fermion operators with a finite-length Wilson line that makes the object gauge invariant. The infinite nucleon momentum on the light cone is, in turn, replaced by a nucleon boost in a spatial direction, conventionally the third direction (z). The Fourier transform of such computed matrix elements defines quasi-PDFs (qPDFs), which can be matched to standard light cone PDFs using perturbation theory [9,10], presently available at one-loop. If sufficiently large nucleon momentum can be achieved, which is itself a non-trivial issue due to the signal-to-noise ratio worsening for large boosts, the matching is relatively small and the higher order corrections are small in the matching, making the extraction of PDFs from this approach robust.
There has been considerable attention devoted to Ji's approach, see e.g. Refs. [11][12][13][14][15]. A plethora of aspects has already been addressed, including effects of smearing the links in the Wilson line [11,[13][14][15], using momentum smearing [16] to improve the signal at large momenta [15], nucleon boost dependence in lattice data [11,[13][14][15] and in models [17,18] (using the relation of PDFs and transverse momentum dependent distribution functions (TMDs)) and corrections that the finite nucleon mass induces [13,14], as well as excited states effects (dependence on the source-sink separation in the computation of the three-point functions) [13,14]. Theoretical aspects have also been discussed. In particular, the role of the Euclidean signature in the computation of quasi-PDFs has been considered [19,20] and led to the conclusion that matrix elements obtained from the lattice are the same as ones from the LSZ formula in Minkowski space. Question has also been raised [21] about the presence of power divergent mixings in moments of qPDFs. However, as argued in Ref. [22], the moments of qPDFs do not exist and hence the problem of such mixings does not exist either. An alternative generalization of light cone PDFs to finite momentum has also been proposed, leading to the so-called pseudo-PDFs [23,24]. Finally, for completeness we mention that Ji's approach can also be used to compute distribution amplitudes (DAs), as shown in Refs. [25,26].
An essential aspect of PDFs computation from qPDFs has until recently been missing, which is the issue of their renormalization. Obviously, this is crucial from the point of view of applicability of the whole programme. The renormalizability of qPDFs to all orders in perturbation theory was recently addressed [27,28]. These papers used two different approaches -the former is based on the auxiliary field formalism of Dorn [29] (which can also be used to extract renormalization functions [30]), while the latter considered qPDFs defined in coordinate space. The matrix elements of qPDFs contain standard logarithmic divergences, but in addition also a power divergence is present, related to the inserted Wilson line [31]. At one-loop in perturbation theory, the divergence appears as a linear divergence, and was calculated in Ref. [32] for different fermionic and gluonic actions. A way to eliminate the linear divergence was also developed in Ref. [33], where it is made finite by using the gradient flow. An additional aspect in qPDFs renormalization is also that qPDFs operators exhibit mixing for certain Dirac structures [32]. The aim of our present work is to perform the renormalization procedure in a fully non-perturbative manner. Such a procedure was introduced, for the first time, in our paper [34] and the proposed renormalization programme was also applied in Ref. [35]. In these proceedings, we recapitulate the results of our extensive analysis that can be found in Ref. [34].

Renormalization prescription
We consider the following matrix elements with Dirac structure Γ (γ µ for unpolarized, γ µ γ 5 for helicity and σ µν for transversity PDFs): taken between boosted nucleon states |P with P = (P 0 , 0, 0, P 3 ), i.e. with the spatial momentum along the third direction. W 3 (z) represents the Wilson line of length z in the direction of the boost. The qPDFs can be obtained by Fourier transforming these matrix elements. The adopted renormalization prescription relies on the non-perturbative RI ′ scheme [36] and proceeds along the lines of the renormalization programme developed for local operators [37]. In absence of mixing, it consists in imposing the following momentum space renormalization condition, for each Wilson line length z: where Z q is the quark field renormalization function that satisfies: In both of the above equations, the momentum p is set to the RI ′ renormalization scaleμ 0 , which we choose in such a way that its third spatial component is P 3 . We test two choices for the other two spatial components -we leave them to be 0 ("parallel" choice) or equal to P 3 ("diagonal"). Previous studies of local matrix elements renormalization [38] show that the latter choice is expected to have much smaller discretization effects, as it has a much smaller value of the ratioP≡ ρμ Further, V(p, z) denotes the amputated vertex function of the operator, with V Born its tree-level value, e.g. iγ 3 γ 5 e ipz for the helicity operator. S (p) is the quark propagator, which equals S Born (p) at treelevel.
The above RI ′ renormalization conditions remove all the divergences present in the matrix elements. The linear divergence resums into a multiplicative exponential factor, e −δm|z|/a+c|z| , where the operator-independent δm is the strength of the divergence, while c is an arbitrary scale [39] that is fixed by the renormalization condition. The presence of this divergence is the major difference with respect to the case of local operators. As in the case of the latter, the standard logarithmic divergence is also removed and finite renormalization is fixed with RI ′ equations.
Analogous conditions can be written for the case when mixing is present, e.g. for the vector operator with Dirac structure parallel to the direction of the Wilson line that mixes with the scalar operator. In such case, one constructs a 2×2 mixing matrix [34].
Having obtained the RI ′ scheme renormalization functions, it is desirable to convert them to the conventional renormalization scheme for phenomenology, i.e. the MS scheme. Moreover, change of the scale fromμ 0 to an MS scale ofμ, typically 2 GeV, can also be peformed. We use the oneloop conversion factor from Ref. [32]. The evolution to the MS renormalization scaleμ=2 GeV is performed using the intermediate Renormalization Group Invariant scheme (RGI), see Refs. [38,40] for more details.
We illustrate the conversion factor for the unpolarized and helicity operators from two RI ′ scales, aμ 0 = 2π 32 ( 7 2 + 1 4 , 3, 3, 3) ("diagonal" one) and aμ 0 = 2π 32 ( 4 2 + 1 4 , 0, 0, 3) ("parallel"), to the MS scheme at the same scale in Fig. 1. Each value of the conversion factor is calculated numerically (which includes the evaluation of integrals of modified Bessel functions). The real part is an order of magnitude larger than the imaginary part and increases with increasing z, whereas the imaginary part stabilizes for |z|/a 2. At z=0, the conversion equals 1, since the corresponding Z-factors are just ones of the local vector or axial vector currents. For other values ofμ 0 , the behaviour of the conversion factor is qualitatively similar.

Lattice setup
We use the setup of Ref. [15] for calculating qPDFs matrix elements, i.e. we take one ensemble of N f =2+1+1 maximally twisted mass fermions [41][42][43][44] produced by the ETM Collaboration [45,46]. The parameters of this 32 3 ×64 ensemble correspond to a lattice spacing of a≈0.082 fm [47] and a pion mass of around 375 MeV. The nucleon boost that we use, P 3 = 6π L , is around 1.4 GeV in physical units. We performed 30000 measurements on 1000 configurations, using the standard Gaussian smearing of quark fields.
For the computation of the renormalization functions in the RI ′ scheme, we employ the momentum smearing technique [38,48]. This technique allows to obtain good statistical precision already with 10-20 gauge field configurations and this is the number of configurations that we have used in practice.

Results
We now show (Fig. 2) the renormalization function Z ∆h for the helicity case, i.e. the Zfactors that renormalize the matrix element ∆h(P 3 , z). The left panel shows the "parallel" case, aμ 0 = 2π 32 ( 4 2 + 1 4 , 0, 0, 3), while the right panel corresponds to the "diagonal" choice of the renormalization scale, aμ 0 = 2π 32 ( 7 2 + 1 4 , 3, 3, 3). In both cases, the real part of the Z-factors increases significantly for increasing Wilson length lines, in accordance with the presence of the power divergence. The conversion to the MS scheme further enhances the values. The "parallel" case leads to slightly larger values of the real parts. We can expect that these values are more contaminated by lattice artifacts in the "parallel" case. This can be argued on general grounds from the higher value ofP, which is around 0.54, as compared to approx. 0.26 for "diagonal". Moreover, we can compare the local value, for z = 0, when the Z-factor is real and scheme-and scale-independent. For the "parallel" case, we find Z ∆h (0) = 0.8620 (15), as compared to Z ∆h (0) = 0.7727(2) for the "diagonal" choice. Both values should be compared to the value of Z A =0.7556 (5) found in Ref. [38], where a robust calculation has been performed, including usage of several scales, subtraction of lattice artifacts computed in lattice perturbation theory and extrapolation to the limit (ap) 2 → 0. It is clear that the value that we find for the "diagonal" case is much closer to this reference value and hence, this choice of the renormalization scale aμ 0 is subject to smaller discretization effects.
For the imaginary parts of Z ∆h in the RI ′ scheme, we find smaller values in the "diagonal" case as compared to the "parallel" one and they are further reduced after conversion to the MS scheme. Again, the difference between the two choices suggests that the "diagonal" choice is subject to smaller lattice artifacts, since the perturbative Z-factor in dimensional regularization is extracted only from the poles and hence it is real to all orders in perturbation theory. As such, the non-zero value of the imaginary part of Z ∆h is an indication of both lattice artifacts and truncation effects resulting from the usage of only one-loop conversion between the schemes.
The above discussion makes it clear what needs to be done to make the extraction of Z-factors more robust for our case of interest. The subtraction of lattice artifacts in lattice perturbation theory has been shown to be very efficient for the local Z-factors case [38] and is expected to work in a similar fashion here. In this way, we expect that after performing this subtraction, different scales will lead to similar values of the Z-factors, i.e. the slope of the (ap) 2 -dependence of Z-factors will be significantly reduced. The other unwanted effect in the MS renormalization functions comes from truncating the conversion factor to one-loop in perturbation theory. We plan a two-loop computation to investigate this aspect.
Finally, we employ the above discussed Z-factors to renormalize the helicity matrix elements. In Fig. 3 left/right, we show the real/imaginary part of bare matrix elements (blue circles), together with MS (μ = 2 GeV) renormalized matrix elements for the two choices of the renormalization scale, the "parallel" one (green triangles) and the "diagonal" one (red squares). The differences between the latter are a measure of lattice artifacts present in the Z-factors. Note that these differences can be treated as an upper bound of these artifacts, since they involve two very different choices of renormalization scales. In typical extractions of the renormalization functions, one typically restricts oneself to scales with corresponding values ofP below 0.35-0.40 [38]. In such case, lattice artifacts are much smaller and moreover, they can be reliably subtracted with a one-loop lattice perturbative calculation. As we argued above, another source of uncertainty in the Z-factors and hence also in the renormalized matrix elements is the truncation of the conversion factor between the RI ′ and MS schemes to one-loop. This manifests itself, for example, in the fact that the real part of renormalized matrix elements be- comes negative for large Wilson line lengths. This behaviour can occur even when the bare real part has already decayed to zero (as happens around |z|/a = 9), but when the imaginary part is still nonzero and it is multiplied by the non-zero imaginary part of the MS Z-factor. Obviously, both of the unwanted effects would propagate to the extracted quasi-PDFs and would contaminate the extracted light-front PDFs. Hence, addressing them is of utmost importance for the robustness of the whole research programme. After the LATTICE 2017 conference, we have indeed addressed the systematic effects shortly described here and we refer to our paper [34] for more details. In short, we have investigated there several "diagonal" scales and estimated the magnitude of lattice artifacts as well as of truncation effects. The explicit computation of these effects, i.e. the one-loop evaluation of lattice artifacts in lattice perturbation theory and calculation of the two-loop RI ′ -MS conversion formulae is the necessary next step.

Conclusions and future directions
In this proceeding, we have summarized our prescription for the non-perturbative renormalization of matrix elements involved in the computation of quasi-PDFs. The employed RI ′ scheme correctly handles both kinds of divergences that are present, i.e. the standard logarithmic divergence as well as the Wilson line related power divergence. We have shown our results for the RI ′ and MS renormalization functions for the relevant Wilson line lengths and we have discussed the systematic effects present in these Z-factors. Addressing them is a necessary step to obtain robust results for the appropriate matrix elements and, thus, also for quasi-PDFs and the light-front PDFs as the final step. In Ref. [34], we have made first quantitative assessment of such effects, i.e. lattice artifacts resulting from the breaking of the rotational symmetry and truncation effects in the perturbative conversion to the phenomenologically widely used MS scheme. However, the explicit computation of these effects remains essential.
Obviously, the systematic effects that need to be addressed for a robust extraction of light-front PDFs are not restricted to the ones present in the renormalization functions. Apart from this, several other effects have to be taken into account and this direction is also currently pursued by us. Arguably the most important of these effects is the influence of the pion mass. So far, our computations [13,15,34] have focused on an ensemble of N f = 2 + 1 + 1 twisted mass fermions at a non-physical pion mass of around 375 MeV. However, as shown e.g. in Ref. [49], the moments of PDFs with such a pion mass are considerably different than the experimentally measured ones. Moreover, the moments extracted from the 375 MeV ensemble using the quasi-PDF approach [15] (including matching, but with renormalization proxied only by HYP smearing and the usage of the local currents Z-factor) agree with the extraction of Ref. [49] within uncertainties. This indicates that the pion mass can indeed considerably influence the shape of the quasi-distributions. To investigate these effects, we are currently performing computations for an ensemble of N f = 2 twisted mass fermions with a clover term, at the physical pion mass, see Ref. [50] in these proceedings.