{'fraction_non_alphanumeric': 0.03494529795158287, 'fraction_numerical': 0.011027700186219739, 'mean_word_length': 4.607603197911568, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Next-to-leading order evolution of twist-two conformal operators: The abelian case

We present the method, based on the use of the broken conformal Ward identities, for the calculation of the anomalous dimensions of conformal operators beyond the leading order of perturbation theory. By means of this technique we find the complete set of two-loop corrections to the off-diagonal elements of the anomalous dimension matrix (the diagonal ones are known from deep inelastic scattering calculations) of the parity odd and even singlet local operators in the Abelian gauge theory. We reconstruct corresponding exclusive two-loop evolution kernels and give the reduction formulae for the restoration of the diagonal part of the ER-BL kernels from the DGLAP splitting functions.


Introduction.
Conformal invariance [1,2,3] has been known for a long time as a powerful tool for the derivation of the constraints on the form of the Green functions in massless field theory. However, the real world could not be conformal invariant even approximately, since otherwise the mass spectrum of the elementary particles might not be discrete. When the interaction is switched on the necessity to construct well-defined theory forces us to introduce a UV cut-off, which inevitably breaks the symmetry. It is reflected in the fact that the field operators acquire the anomalous dimensions and the coupling constant becomes running. In the description of the physical quantities new objectslocal composite operators -appear which define the dynamical content of the former. The computation of the spectrum of the anomalous dimensions of these composite operators is one of the most important problems of the theory, since these quantities determine the scaling behaviour of experimental observables, i.e. the moments of the parton distributions or wave functions depending on the particular processes in question. These non-perturbative entries satisfy corresponding evolution equations, known as Dokshitzer-Gribov-Lipatov-Altarelli-Parisi (DGLAP) [4,5,6] one for inclusive reactions and Efremov-Radyushkin-Brodsky-Lepage (ER-BL) equation [7,8] for the exclusive ones [9]. In order to solve them one has to determine the eigenfunctions and eigenvalues of renormalization group kernels. The main problem arises in the study of exclusive processes [7] where contrary to the forward scattering one needs to handle the whole tower of local operators with total derivatives which mix with each other under renormalization. In the leading order (LO) this problem was analyzed in a number of papers [7,8,10,11,12,13,14,15]. It was recognized that the operators which diagonalize the anomalous dimension matrix to the lowest order have a conformally covariant transformation behaviour at tree level.
Nowadays it becomes a common wisdom to perform the analyses in the next-to-leading order (NLO) of perturbation theory. The self-consistent approach to this problem requires the knowledge of the two-loop exclusive evolution kernels. Because of the complications spelled out above, the direct calculation of the anomalous dimension matrices or the evolution kernels beyond the leading order turns out to be tremendous task and was resolved up to now in the non-singlet sector [16,17,18,19,20,21,22]. Since the conformal operators diagonalize corresponding ER-BL evolution kernels in the leading order of perturbation theory and the off-diagonal part appears only in the NLO approximation, the symmetry breaking starts at the order O(α) and might be expected to be evaluated from the Green functions of some quasi-one-loop objects. Therefore, the first step toward the accomplishment of the program of the evaluation of the O(α 2 )-corrections to the evolution kernels, is to specify the particular pattern of the conformal symmetry breaking in gauge theories to leading order. To this end true conformal Ward identities are required.
Due to the UV divergences of the theory there exist anomalies, but only part of them could be related to those which appear in the violation of scale invariance [23]. Note, that this is already true in scalar theories, where there is one-to-one correspondence between dilatation and special conformal anomalies induced by the renormalization of the fields [24] only. In gauge theories, a priori, one expects in addition a breaking due to the gauge fixing and ghost terms in the Lagrangian, but it should not affect the final result provided we are interested in some physically relevant quantities. Thus, the limit of non-running coupling constant corresponds to the exact scale invariance; however, it is by no means true that the special conformal symmetry also respects this limit because there exists a new source of symmetry breaking which can be attributed to the renormalization of the product of composite operators and trace anomaly of the improved energy-momentum tensor. However, in the study of gauge invariant quantities these effects could be always embedded into the redefined conformal representation of the algebra, so that the symmetry will be broken explicitly only by the terms proportional to the β-function of the theory. Thus, conformal covariance is restored in the non-trivial fixed point of the renormalization group function. This fact gives us the opportunity to apply the powerful machinery of conformal symmetry. It has been exploited by us quite recently to give several QCD predictions coming from the use of the conformally covariant operator product expansion [25,26,27].
In the present analysis we will construct the formalism and calculate the two-loop results for the singlet ER-BL kernels, restricting ourselves to Abelian gauge theory. The main tools for these purposes will be the renormalized Ward identities for the dilatation and special conformal transformations [28,29,30,31,32] and the commutator constraint coming from the algebra of the conformal group. The most important feature of the former is that they are derived without any assumptions about the possible form of the symmetry breaking and, thus, turn out to be true anomalous Ward identities in gauge theory.
Consequent presentation will be organized as follows. Next section is mainly introductory. We summarize there some basic definitions of the main objects of our analysis: singlet conformal operators, renormalized Lagrangian of the Abelian theory, as well as other generalities related to the transformation of the latter under the action of generators of the conformal algebra. In section 3 we construct renormalized conformal Ward identities for the Green functions with conformal operator insertions. In section 4 they are used to derive the commutator constraint on the off-diagonal elements of the anomalous dimension matrix. The most important feature of the conformal consistency relation is that it allows to find corresponding entries in the n-th order from the calculation of the, so-called, special conformal anomaly matrix in the (n − 1)-loop approximation. This is used by us to derive the two-loop results for the anomalous dimensions in terms of the one-loop quantities which are evaluated in section 5. These results are also transformed in the language of ER-BL evolution kernels and checked in section 6 for the particular set of two-loop diagrams with fermion bubble insertions in the GQ-kernels. Next, we present the reduction formulae which allows to reconstruct the diagonal part of the ER-BL kernels starting from DGLAP analogues.
The final section is devoted to the summary. Trying to make the paper self-consistent as much as possible we add two appendices. In the first one we give the Feynman rules used in the main body of the paper. In the second one we list the extended ER-BL evolution kernels responsible for the evolution of the generalized non-forward distribution functions which attract lots of attention recently.

Preliminaries.
To start with, we review few basic properties of the theory invariant under conformal algebra of space-time. The conformal algebra [1,2,3] is a 15-dimensional Lie algebra, isomorphic to O(4, 2). Its generators are the Poincaré generators P µ , M µν together with dilatation D and special conformal transformation K µ . The latter could be defined as product K µ = RP µ R −1 , where R produces inversion of space-time points x µ → x µ /x 2 . Apart from commutators which manifest the scalar and vector nature of D and K µ , respectively, there are extra non-trivial commutation The transformations of the field operators under dilatation and special conformal transformations 2 are expressed by the formulae where the parameters of the infinitesimal transformations have been set equal to unity. Using the method of induced representation [2] one can construct from the, so-called, stability subgroup the differential representationĜ for the generators D and K µ of the conformal algebra. Here we deal with the irreducible standard representation, namely: where d φ is a scale dimension of the field φ. The particular choice for them will be discussed below. The action of the spin operator Σ µν on the fermion and vector fields is defined as usual: In this paper we intend to construct and study the Ward identities which originates from the dilatation and special conformal transformations. The simplest way to derive them is to use the path integral formalism [29,30,31,32] and dimensional regularization in d = 4 − 2ǫ dimensions.
Their generic form looks like where δ G -symbol means some infinitesimal transformation from the algebra and . . . corresponds to the path integral averaging with the weight function exp iS. S = d d xL stands for the action defined in terms of the renormalized Lagrangian of Abelian gauge fields coupled to fermions where G µν = ∂ µ A ν − ∂ ν A µ is the Abelian field strength tensor. Here we have also added a free ghost Lagrangian to make L invariant under BRST transformation. In Eq. (5) X = k φ(x k ) is a product of elementary fields (entering in the Lagrangian); its variation is expressed in terms of the generators G of the conformal group acting on the field operators: In Eq. (5) we have considered Green functions with a conformal operator insertion which is defined as two-dimensional vector An explicit form of the conformal covariant tree operators is 3 We have used the following relation ǫ µν−+ G µ+ . The renormalized operators are expressed through the non-renormalized ones in the following way where the matrix Z jk (ǫ, g) is triangular and independent on l due to Lorentz invariance of the theory. It is given in the minimal subtraction scheme bŷ Since we are restricted to the leading twist operators (9) defined on the light-cone, it is sufficient to deal consequently with the collinear conformal subgroup SU(1, 1) of the whole group which simplify considerably the analysis. For that we project the Lorentz indices with two light-like vectors n and n * , which are defined along opposite tangents to the light cone, with n 2 = n * 2 = 0 and nn * = 1.
Beyond the tree level we are not limited to any unique choice of the scale dimensions d φ , rather we can define them as we wish. It turns out, that the most optimal and convenient set is where a jk (l) = δ jk a(j, l) = δ jk 2{(j + 1)(j + 2) − (l + 1)(l + 2)}.
In the consequent discussion and all practical calculations we will always employ the lightcone position formalism for the evaluation of the renormalization matrices of different composite operators. Therefore, they should be replaced for these purposes by the non-local string operators.
Thus, we introduce the quark and gluon gauge invariant operators where we have omitted for brevity the path ordered exponential in the first entry.
To simplify the derivation of the Ward identities it is instructive to use the following variations of the Lagrangian They lead to the non-renormalized scale and special conformal variations of the action expressed in terms of the following equation-of-motion (EOM) operators (Ωψ ψ = Ω ψ + Ωψ), and the operators of class A and B [35] O Because of the factor ǫ appearing in Eqs. (17,18), the conformal variation of the action can be different from zero only due to the UV divergencies. The appearing operators can be renormalized along the line of Refs. [35]. However, since these are exactly the operators which enter in the Lagrangian density, their renormalization could be gained from the study of the differential vertex operator insertion g ∂ ∂g S and ξ ∂ ∂ξ S. Namely, a little algebra leads to Here 5 From the finiteness of the differential vertex operator insertions the counterterms of can be fixed: in complete agreement with general renormalization theorems [35]. These equalities hold up to counterterms which are total derivatives of some other operators R µ i . By considering the Green functions of the differential operator insertion ξ ∂ ∂ξ S and from the finiteness of the special conformal Ward identity, we can conclude that they are proportional to the gauge fixing parameter ξ and is finite by itself, since there is no other dimension three operator with the same quantum numbers which it could mix with due to the renormalization. Now the final step can be easily done with the result where the coefficients in front of the operator insertions are the well-known renormalization group (σ = −2γ G ) and the anomalous dimensions of the elementary fields The first expansion coefficients of Aiming consequent generalization of the present considerations to the non-Abelian theory we keep dependence on the colour factors explicit. In QED we have to set C F = T F = 1.

Conformal Ward identities.
Combining the results, we have obtained so far in Eq. (5), we will come to conformal Ward identities, which, however, contain unrenormalized products of operators. The aim of this section is to renormalize these products [36]. Since the coefficients appearing in the Ward identities depend on the regularization parameter ǫ, this renormalization procedure leads to the conformal anomalies.

Scale Ward identity.
The renormalization of the dilatation Ward identity is quite straightforward and does not require any additional information as one which comes from the consideration of the differential vertex insertions. By using the action principle [37] When substituted into the Ward identity, the terms proportional to the functional derivative acting on the conformal operator will cancel upon similar terms coming from the renormalization For quark equation of motion we have where the projectors on the quark and gluon sectors are defined bŷ It should be noted that the definition of the renormalization constant we have used above does not correspond to the minimal subtraction of UV divergencies, rather it comes from the mere integration by parts in the functional integral for the Green functions.
Putting all these equations together, we obtain the renormalized dilatation Ward identity Note, that only the whole anomalous dimension matrix of the composite operators,γ ≡γ Z + 2γ GPG + 2γ ψPQ , is gauge independent. Consequently, the off-diagonal elements ofγ Z are gauge independent too, while its diagonal entries have a very simple dependence on the gauge parameter.
Taking into account that the improved energy-momentum tensor [38,39,41] in the massless gauge theory admits the anomalous trace [40,41,42] (35) and inserting it in the above equality (33), making use of the relation where N φ is a number of φ-particles operator, we come to the familiar Callan-Symanzik equation [43] for the Green function, albeit with a conformal operator insertion, where the canonical scale dimensions of the field operators are shifted by the anomalous ones [24].

Special conformal Ward identity.
The derivation of the renormalized special conformal Ward identity requires new inputs. We define the renormalization constants for the product of the renormalized conformal operator and finite For the product with EOM operator we have Simple calculations lead to the resultŝ where we have introduced ab-matrix with elements An additional renormalization constant in Eq. (38) is extracted from the 1 ǫ -pole part of the operator insertion In the non-local representation, to be used below, the left hand side of this equality is given by in terms of the non-local string operators defined above.
Finally, the renormalized special conformal Ward identity takes the form Since we have expressed everything in terms of the renormalized operator products, we can safely put the regularization parameter ǫ = 0. Moreover, because the LHS of this equality is finite the same property is fulfilled by the RHS. We have defined above the special conformal anomaly where the second equality holds becauseγ c is finite.  (33) and (42).
For this purpose we need the following commutators Since σ = −2γ G is gauge independent in QED, the RHS of Eq. (45) vanishes identically. One should be careful in treatment of these equalities acting on the differential vertex operator insertions, since they are only valid when both differentials act on the same objects. Naive application could lead to an error. An explicit derivation, which makes use of these commutators, gives the following renormalization group equations for the operator insertions [∆ i − ] in QED: The following steps are quite straightforward. Now we will act repeatedly with the operator D on every term in the special conformal Ward identity without taking into account the action of the latter on the field monomial X , since at the end these terms would anyway die out by the use of the conformal Ward identity. First where we have taken into account that (physical part of) the special conformal anomaly matrix is gauge independent quantity. Next The result of Eq. In the NLO approximation the constraint (53) fixes unambiguously the off-diagonal elements of the anomalous dimension matrices: Adding the two-loop results for the anomalous dimensions γ total derivatives 6 , which are available for about two decades from the calculations of unpolarized [44] and, recently, from polarized [45] deep inelastic scattering splitting functions, we obtain the complete entryγ which governs the evolution of composite operators in Abelian gauge theory to NLO.
5 Evaluation of the conformal anomaly matrix.

Local anomalous dimensions.
In this section we evaluate all one-loop renormalization constants which enter in the special conformal anomaly matrix and determine the anomalous dimensions AB γ First, in the leading order 2Ẑ [1] = α 2πγ (0) Z . Second, the calculation of the diagrams depicted in Fig. 1

and the extraction of the pole part from the renormalization constantẐ
gives: Here we took into account that according to the Feynman rules (A.4) these diagrams are given by the derivatives with respect to the external momenta of the familiar one-loop graphs; and that 6 Due to the particular normalization of the Gegenbauer polynomials the QG and GQ results of Refs. [44,45] should be multiplied by appropriate factors, i.e. 6/j and j/6, respectively. the differentiation of the Gegenbauer polynomials leads to theb-matrix: Defining also the one-loop approximation of the constantẐ ⋆ eff = 1 ǫẐ ⋆ [1] eff + · · · in terms of thê w-matrix, we finally obtain the result for the special conformal anomaly to leading order: The first entries, which are the usual anomalous dimensions of conformal operators, are well known and read × j 2 + 3j + 4, for even parity j(j + 3), for odd parity (64) GQ γ (0) j = −C F 3(j + 1)(j + 2) × j 2 + 3j + 4, for even parity j(j + 3), for odd parity (65) The calculation of theŵ-matrix for the QQ-sector can be found in Ref. [32]. For the GQ-channel we employ the non-local representation for the evaluation of the corresponding kernel. Simple calculation provides the result with the same kernel in the parity even and odd channels Note that this kernel is gauge independent as it should be (the ξ-dependent piece we have found in the calculation is symmetrical with respect to the permutation y ↔ z, and vanishes because of the symmetry properties of the operators). The evaluation of the conformal moments of this simple kernel is quite straightforward and we get QQ w jk = −C F θ j−2,k 1 + (−1) j−k (2k + 3)a(j, k) (k + 1)(k + 2) × A jk − ψ(j + 2) + ψ(1) + 4(k + 1)(k + 2) A jk a(j, k) , (69) In addition QG w = GG w = 0. Defining new matrices d jk = b jk /a(j, k) and g jk = w jk /a(j, k), we can rewrite the off-diagonal elements:

Evolution kernels.
The representations (71-74) allow us to write the off-diagonal part of the evolution kernels as convolution of known kernels. The LO kernels AB V (0) in the expansion V (x, y|α) = α 2π V (0) (x, y) + These results are already known or can be derived from the DGLAP kernels by the procedure described below. In the restricted region |x, y| ≤ 1 these kernels have the general structure GQ F (x, y) = C F x y 2y − x, for even parity x, for odd parity , Here one comment is in order to clarify the differences of the even parity QG and GQ kernels available in the literature. Indeed these kernels are not uniquely defined and different expressions could be obtained from various calculation techniques [11,12,46,47,48]. However, apparently different functional dependence on the momentum fraction x, y leads to the same eigenvalues in the physical sectors of the corresponding channels. The only difference arises in the space of the non-physical moments. Namely, the kernels derived by Chase [11] (and listed above) possess the off-diagonal matrix elements in the basis of the Gegenbauer polynomials in the non-physical sector, while the kernels derived below are exactly diagonal. The difference can be attributed to the following additional pieces ∆ QG F (x, y) which should be added to the results of Chase to get the equations given behind. Note also that the kernels we have listed (taken without colour factors) respect the SUSY relation [49]: Now we are able to reconstruct the off-diagonal parts of the singlet two-loop ER-BL kernels from the local anomalous dimensions we have derived above, i.e. g jk and ( AB γ In the QQ-channel the g-kernel is given in Ref. [32], while for the GQ-channel we have been able to perform the summation of the Gegenbauer polynomials by the procedure discussed below in section 7: Intuitively, one could anticipate this last result. For that one should note that in Eq. (70) the factor 1/(k + 1)(k + 2) gives just the eigenvalues of the ER-BL kernel in scalar φ 3 (6) -theory. Next, the off-diagonal matrix elements of theb-matrix are (2k +3). According to Eq. (60), the latter can be induced by derivative of the Gegenbauer polynomials C ν k (2x − 1) with respect to the argument x times the factor (2x − 1). Since the derivative increase the index and reduce the order by one unity, we come to the equation we have just presented.
Note that the kernels (83,84) possess also diagonal conformal moments. They can be subtracted by applying the projection operator (I−D), where D extracts the diagonal part of any test function τ (x, y): The kernel AB g(x, y) that possesses the AB g ij -moments is defined as follows The factor ( AB γ k )d jk corresponds to the expected breaking of the diagonality of the NLO exclusive kernels [15,20] and it is related to the one-loop renormalization of the composite operators. These contributions can be embedded into a redefined conformal representation, which differs from the original one by the shift of the scale dimensions of local operators by the anomalous ones γ j . This means that new eigenfunctions should be orthogonal to each other with the weight factor (yȳ) n+γ j , with n = 1, 2 corresponding to quarks and gluons, respectively. This gives us a guide how to construct new kernels V δ (x, y) from the leading order ones. For that we take into account that V δ (x, y) should be symmetrical after multiplication by the factor (yȳ) n+δ : with n = 1, 2 for A = Q, G-channels, respectively, and for mixed ones. These equations tell us that the two-loop off-diagonal termsV (x, y) should be proportional to the first term in the Taylor expansion of V δ (y, x) with respect to δ. This recipe leads to the following dotted kernels Note, however, that simple logarithmic modification of the leading order result does not work for the parity even mixing channels. The main difference from the other cases comes from the fact that Chase's kernels possess off-diagonal elements in the non-physical space of conformal moments, which immediately manifest themselves when these kernels are multiplied by factor ln(x/y). However, the above symmetry consideration could not fix unambiguously the modification of extra pieces ∆F and, therefore, it requires an explicit analysis. We find the following results for ∆ ABḞ : Combining the equations (83,84,88,89), the off-diagonal parts of the NLO evolution kernels read: Here we use a shorthand ⊗ = 1 0 dz.

Fermion bubble insertions in the singlet kernels.
In this section we present the check of the self-consistency of our approach by explicit evaluation of the fermion vacuum insertions in the GQ-evolution kernels. From the conformal constraints we know that corresponding contribution to the off-diagonal part of the evolution kernel is of the Therefore, evaluating the diagrams in Fig. 2 we should immediately obtain the desired form of the dotted and g-kernels. To this end we have exploited the formalism developed in Refs. [18,20], which is the most appropriate for the multi-loop calculations of the evolution kernels and makes use of the, so-called, α-representation for the Feynman propagators. Corresponding calculations run along the same line as for the one-loop diagrams but now we insert the fermion vacuum polarization iΠ µν (k) into the gluon lines (see Fig. 2) (with µ * 2 = 4πµ 2 ), which results to the mere shift of the power of the denominator of the gluon propagator by the magnitude ǫ. Since we are interested only in the off-diagonal part of the corresponding diagrams we can safely omit, after subtraction of sub-divergences, the contribution coming from the surviving overall constants, i.e. Π(ǫ) etc., since it would be proportional to the leading order evolution kernel. The remainder of the graphs leads to the integrals where the functions v i read for the parity even and odd channels To perform the final integrations we need the following formulae [50]: for the first one, and for the second. Putting ǫ = 0 in the square brackets in Eq. (96), we immediately obtain twice the leading order kernels 7 . Evaluating the integral for arbitrary ǫ and differentiating it with respect to ǫ, we finally obtain the contribution to the two-loop kernel with F (x, y) = 2C F x y 1 − 2y, for even parity 1 − 2x, for odd parity , 7 To be precise, in the parity even channel, we get Chase's kernel.
while the other functions are given by Eqs. (78,79). Obviously, these expressions consist of both diagonal and off-diagonal parts of the evolution kernel. Subtracting the sum of the dotted and the g-kernels derived in the preceding section it is a trivial task to convince oneself that the remainder is diagonal in the basis of the Gegenbauer polynomials. Namely, this provides a new diagonal and an equivalent representation for the g-kernel in the GQ-channel (84) which holds true for the off-diagonal matrix elements. Note that it is most close to the functional form of QQ g.
Finally, we can conclude that an explicit calculation supports our results derived from the use of the conformal consistency relation (53).

DGLAP → ER-BL reduction.
Aiming an explicit numerical integration of the two-loop evolution equations one is interested in the momentum fraction splitting kernels. While the off-diagonal parts of these kernels are known explicitly, it is not the case for the diagonal entries. In this section we outline the formalism which fills this gap. We have been able to restore the one-loop ER-BL evolution kernels analytically. The reduction formulae we have derived here are sufficient to find numerically the diagonal part (in the basis of the Gegenbauer polynomials) of the exclusive singlet kernels from the corresponding DGLAP analogues.

QQ channel.
The formal expression for the diagonal part of the ER-BL evolution kernel can be given by the where N j (ν) = 2 −4ν+1 Γ 2 ( 1 b 1 1 u Figure 3: The integration contour for the evaluation of the integral (109). There is a n-th order pole in the point u = b and the cut on the real axis from −1 to 1.
Therefore, the main problem is to perform the infinite summation and construct the generating function for the exclusive evolution kernels. Using the Gegenbauer's summation theorem for C ν jpolynomials [51] this goal can be easily achieved. We have evaluated the following conversion function from the forward to non-forward limits: In order to find the ER-BL evolution kernel we have to convolute the above equality with the DGLAP splitting function and evaluate the integral on the RHS of this equation. For this purpose we use the Euler's integral representation of the hypergeometric function, so that the above formula Performing the external integration first leads to the simple elementary function, so that the final integral looks like: It can be easily evaluated by considering its analytical continuation to the whole complex plane and choosing the integration contour as displayed in Fig. 3. It reduces to the evaluation of the residues in the n-th order pole in the point u = b ≡ (1 − (2x − 1)(2y − 1))/4 √ xxyȳ and the point u = ∞: The first residue provides the contribution proportional to difference of θ-functions, while the second one -to the sum. Namely, for the n = 1, 2 which appear in the transformation of the one-loop splitting functions we get for n = 1 and Thus, we finally obtain For instance, for the QQ- and after little algebra we come to the known non-singlet kernel.

QG and GQ channels.
We can make use of the same summation formula Eq. (107) for the mixed channels too. For that we have to differentiate both sides of this identity with respect to the first x or the second y arguments. Since the anomalous dimensions of the conformal operators differ from the one of the usual local operators by the factor 6/j and j/6 for QG and GQ channels, respectively, we will introduce these factors by representing the forward evolution kernels as a convolution with appropriate functions. The main reason for that is to use the same Eq. (107), since otherwise the evaluation of the infinite series may turn out to be quite complicated. For the QG and GQ sectors we define the kernels and corresponding summation formulae Since the parity odd channel does not introduce any new features as compared with results available in the literature we turn to the vector sector. The integrals appeared in the calculation are just the same as before, namely, we get for the QG P V (z) = (z 2 + (1 − z) 2 ) splitting function and for the GQ P V (z) = (1 + (1 − z) 2 )/z kernel Discarding the constant term J 0 = π 2 we reproduce the Chase's kernels, however, if the latter is taken into account we come to the expressions given in Eqs. (81).

Summary.
In the present paper we have developed a formalism, which allows to evaluate the n-th order nondiagonal part of the exclusive twist-2 evolution kernels starting form the (n−1)-loop approximation for the special conformal anomaly matrix, for the flavour singlet channel in Abelian gauge theory.
The main tools of our analysis are the renormalized anomalous conformal Ward identities and commutator constraint obtained from the algebra of the conformal group. Using this approach we have derived O(α 2 ) corrections to singlet anomalous dimensions of the conformal operators. This is sufficient to reconstruct the complete ER-BL splitting functions AB V (1) (x, y) in the NLO approximation. We have partially checked our calculations by explicit evaluation of the simple two-loop diagrams with fermion vacuum polarization insertions. The formalism we have described allows to write in a simple way the solution of the corresponding evolution equations, since the corrections to the conformal partial wave expansion are defined by known one-loop evolution kernels and g-functions [32,26]. We have also derived, for the first time, the reduction formulae which allow to find the diagonal part of the ER-BL kernels from the corresponding DGLAP analogues. This procedure fix unambiguously the form of the diagonal part of the exclusive kernels which are exactly diagonal in the physical as well as in non-physical spaces of moments.
We wish to thank B. Geyer, A.V. Radyushkin

A Feynman rules.
In this appendix we summarize Feynman rules which have been used in the calculation of the renormalization constants in the main text. Everywhere it is implied that all momenta on the lines are incoming.
For quark and gluon non-local string operators we have In the quark sector we have omitted the contributions coming form the phase factor, since they are not of relevance for our present consideration.
The only modified Feynman rule coming from the differential operator insertion [∆ g − ], required to the order we are interested in, is It should be compared to the usual one for the quark-gluon vertex V µ = igγ µ (2π) 4 δ (4) (k 1 + k 2 + k 3 ). (A.5) The resulting expression, composed from propagators and vertices we have just derived, should be integrated with respect to the momentum of every internal line, i.e. multiplied by the factor ℓ d d k ℓ (2π) ℓ·d .
Recently, there is renewed interest in the deeply virtual Compton scattering [22,52,46,48]. The evolution of the corresponding non-forward parton distribution functions [46] h where the non-local string operators are defined by Eqs. (13) and the upper (lower) signs corresponds to the parity-even (odd) channels, is governed by generalized ER-BL evolution equations We perform the calculations and end up with ER-BL type kernels which look like for the parityeven sector and for parity odd one We have used here the generalized Θ-functions defined by the integral To cure this problem in appropriate way, so that the forward limit would not be affected, we should add the terms (B.11) In the forward limit Θ 0 2 (x) = δ(x), so that these extra terms die out. As usual these kernels possess the following DGLAP and ER-BL limits: P (z) ≡ −K(z, 1, 0), V (x, y) ≡ −K(x, y, 1). (B.12)