N3LO Anomalous Dimensions

The N3LO QCD anomalous dimensions \(\gamma^{(3)}\) are not yet fully known, since they rely on the calculation of 4-loop DIS integrals. Moreover, the analytical structure of these function is already known to be complicated since in Mellin space they include harmonics sum up to weight 7, for which an analytical expression is not available.

Here, we describe the various assumptions and limits used in order to reconstruct a parametrization that can approximate their contribution. In particular we take advantage of some known physical constrains, such as the large-x limit, the small-x limit, and sum rules, in order to make our reconstruction reasonable.

Generally, we remark that the large-x limit correspond to large-N in Mellin space where the leading contribution comes from the harmonics \(S_1(N)\), while the small-x region corresponds to poles at \(N=0,1\) depending on the type of divergence.

In any case N3LO DGLAP evolution at small-x, especially for singlet-like PDFs, is not reliable until the splitting function resummation is available up to NNLL.

Non-singlet sector

In the non-singlet sector we construct a parametrization for \(\gamma_{ns,-}^{(3)},\gamma_{ns,-}^{(3)},\gamma_{ns,s}^{(3)}\) where:

\[\gamma_{ns,s}^{(3)} = \gamma_{ns,v}^{(3)} - \gamma_{ns,-}^{(3)}\]

In particular, making explicitly the dependence on \(n_f\), the non-singlet anomalous dimensions include the following terms:

non-singlet 4-loop Anomalous Dimensions








Some of these parts are known analytically exactly (\(\propto n_f^2,n_f^3\)), while others are available only in the large \(N_c\) limit (\(\propto n_f^0,n_f^1\)). In EKO they are implemented as follows:

  • the part proportional to \(n_f^3\) is common for \(\gamma_{ns,+}^{(3)},\gamma_{ns,-}^{(3)},\gamma_{ns,v}^{(3)}\) and is exact [DVR+17] (Eq. 3.6).

  • In \(\gamma_{ns,s}^{(3)}\) the part proportional to \(n_f^2\) is exact [DVR+17] (Eq. 3.5).

  • In \(\gamma_{ns,s}^{(3)}\) the part proportional to \(n_f^1\) is parametrized in x-space and copied from [MRU+17] (Eq. 4.19, 4.20).

  • The remaining contributions include the following constrains:

    • The small-x limit, given in the large \(N_c\) approximation by [DKMV22] (see Eq. 3.3, 3.8, 3.9, 3.10) and coming from small-x resummation. This part contains the so called double logarithms:

      \[\gamma_{ns} \approx \sum_{k=1}^{6} c^{(k)} \ln^k(x) \quad \text{with:} \quad \mathcal{M}[\ln^k(x)] = \frac{1}{N^{k+1}}\]

      Note the expressions are evaluated with the exact values of the QCD Casimir invariants, to better agree with the [MRU+17] parametrization.

    • The large-N limit [MRU+17], which reads (Eq. 2.17):

      \[\gamma_{ns} \approx A^{(f)}_4 S_1(N) - B^{(f)}_4 + C^{(f)}_4 \frac{S_1(N)}{N} - D^{(f)}_4 \frac{1}{N}\]

      This limit is common for all \(\gamma_{ns,+}^{(3)},\gamma_{ns,-}^{(3)},\gamma_{ns,v}^{(3)}\). The coefficient \(A^{(f)}_4\), being related to the twist-2 spin-N operators, can be obtained from the QCD cusp calculation [HKM20], while the \(B^{(f)}_4\) is fixed by the integral of the 4-loop splitting function and has been firstly computed in [MRU+17] in the large \(n_c\) limit. More recently [DMV22], it has been determined in the full color expansion by computing various N3LO cross sections in the soft limit. \(C^{(f)}_4,D^{(f)}_4\) instead can be computed directly from lower order splitting functions. From large-x resummation [DVR+17], it is possible to infer further constrains on sub-leading terms \(\frac{\ln^k(N)}{N^2}\), since the non-singlet splitting functions contain only terms \((1-x)^a\ln^k(1-x)\) with \(a \ge 1\).

    • The 8 lowest odd or even N moments provided in [MRU+17], where from quark number conservation we can trivially obtain: \(\gamma_{ns,s}(1)=\gamma_{ns,-}(1)=0\).

    • The difference between the known moments and the known limits is parametrized in Mellin space. The basis includes:

      \(\gamma_{ns,\pm}^{(3)}\) parametrization basis











      \(- \rm{Li_2}(x) + \zeta_2\)






      \(x^{2}, x^{3}\)


      The first five functions model the sub-leading differences in the \(N\to \infty\) limit, while the last three help the convergence in the small-N region. Finally, we add a polynomial part \(x^{2}\) or \(x^{3}\) respectively for \(\gamma_{ns,+},\gamma_{ns,-}\). For large-N we have the limit:

      \[\mathcal{M}[(1-x)\ln^k(1-x)] \approx \frac{S_1^k(N)}{N^2}\]

      Note that the constant coefficient is included in the fit, following the procedure done in [MRU+17] (section 4), to achieve a better accuracy. It is checked that this contribution is much more smaller than the values of \(B_4\).

Singlet sector

In the singlet sector we construct a parametrization for \(\gamma_{gg}^{(3)},\gamma_{gq}^{(3)},\gamma_{qg}^{(3)},\gamma_{qq}^{(3)}\) where:

\[\gamma_{qq}^{(3)} = \gamma_{ns,+}^{(3)} + \gamma_{qq,ps}^{(3)}\]

In particular, making explicitly the dependence on \(n_f\), the singlet anomalous dimensions include the following terms:

singlet 4-loop Anomalous Dimensions









The parts proportional to \(n_f^3\) are known analytically [DVR+17] and have been included so far. For \(\gamma_{qq,ps}\) and \(\gamma_{gq}\) also the component proportional to \(n_f^2\) has been computed in [GvMSY24] and [FHM+24] respectively and it’s used in our code through an approximations obtained with 30 moments.

The other parts are approximated using some known limits:

  • The small-x limit, given in the large \(N_c\) approximation by [DKMV22] (see Eq. 5.9, 5.10, 5.11, 5.12) and coming from small-x resummation of double-logarithms which fix the leading terms for the pole at \(N=0\):

    \[\gamma_{ij} \approx c^{(6)}_{ij} \ln^6(x) + c^{(5)}_{ij} \ln^5(x) + c^{(4)}_{ij} \ln^5(x) + \dots \quad \text{with:} \quad \mathcal{M}[\ln^k(x)] = \frac{1}{N^{k+1}}\]
  • The small-x limit, coming from BFKL resummation [BM18] (see Eq. 2.32, 2.20b, 2.21a, 2.21b) which fix the leading terms (LL, NLL) for the pole at \(N=1\):

    \[\gamma_{ij} \approx d^{(3)}_{ij} \frac{\ln^3(x)}{x} + d^{(2)}_{ij} \frac{\ln^2(x)}{x} + \dots \quad \text{with:} \quad \mathcal{M}[\frac{\ln^k(x)}{x}] = \frac{1}{(N-1)^{k+1}}\]

    Note that in principle also the term \(\frac{\ln^6(x)}{x}\) could be present at N3LO, but they are vanishing. These terms are way larger than the previous ones in the small-x limit and are effectively determining the raise of the splitting functions at small-x. In particular only the expansion for \(\gamma_{gg}^{(3)}\) is known at NLL. LL terms respect the representation symmetry :

    \[\begin{split}\gamma_{gq} & \approx \frac{C_F}{C_A} \gamma_{gg} \\ \gamma_{qq,ps} & \approx \frac{C_F}{C_A} \gamma_{qg} \\\end{split}\]
  • The large-x limit of the singlet splitting function is different for the diagonal part and the off-diagonal. It is known that [AB01, MRU+22] the diagonal terms diverge in N-space as:

    \[\gamma_{kk} \approx A^{(r)}_4 S_1(N) + B^{(r)}_4 + C^{(r)}_4 \frac{S_1(N)}{N} - D^{(r)}_4 \frac{1}{N}\]

    Where again the coefficient \(A^{(r)}_4\) is the QCD cusp anomalous dimension for the adjoint or fundamental representation, the coefficient \(B^{(r)}_4\) has been extracted from soft anomalous dimensions [DMV22]. and \(C^{(r)}_4,D^{(r)}_4`can be estimate from lower orders :cite:`Dokshitzer:2005bf\). However, \(\gamma_{qq,ps}^{(3)}\) do not constrain any divergence at large-x or constant term so its expansion starts as \(\mathcal{O}(\frac{1}{N^2})\). The off-diagonal do not contain any +-distributions or delta distributions but can include divergent logarithms of the type [SMVV10]:

    \[\ln^k(1-x) \quad k=1,..,6\]

    where also in this case the term \(k=6\) vanish. The values of the coefficient for \(k=4,5\) can be guessed from the lower order splitting functions. These logarithms are not present in the diagonal splitting function, which can include at most terms \((1-x)\ln^4(1-x)\). While for \(\gamma_{gg}\) these contributions are beyond the accuracy of our implementation, they are relevant for \(\gamma_{qq,ps}\). At large-x we have [SMVV10]:

    \[\gamma_{qq,ps} \approx (1-x)[c_{4} \ln^4(1-x) + c_{3} \ln^3(1-x)] + \mathcal{O}((1-x)\ln^2(1-x))\]
  • The 5 lowest even N moments provided in [MRU+24, MRU+22], where momentum conservation fixes:

    \[\begin{split}& \gamma_{qg}(2) + \gamma_{gg}(2) = 0 \\ & \gamma_{qq}(2) + \gamma_{gq}(2) = 0 \\\end{split}\]

    For \(\gamma_{qq,ps}, \gamma_{qg}\) other 5 additional moments are available [FHMV23a, FHMV23b]. making the parametrization of this splitting function much more accurate.

The difference between the known moments and the known limits is parametrized in Mellin space using different basis, in order to estimate the uncertainties of our determination.

Uncertainties estimation

Since the available constrains on the singlet anomalous dimension are not sufficient to determine their behavior exactly, for instance the poles at \(N=1\) and \(N=0\) are not fully known, we need to account for a possible source of uncertainties arising during the approximation. This uncertainty is neglected in the non-singlet case.

The procedure is performed in two steps for each different anomalous dimension separately. First, we solve the system associated to the 5 (10) known moments, minus the known limits, using different functional bases. Any possible candidate contains 5 elements and is obtained with the following prescription:

  1. one function is leading small-N unknown contribution, which correspond to the highest power unknown for the pole at \(N=1\),

  2. one function is the leading large-N unknown contribution,

  3. the remaining functions are chosen from of a batch of functions describing sub-leading unknown terms both for the small-N and large-N limit.

This way we generate a large set of independent candidates, roughly 70 for each anomalous dimension, and by taking the standard deviation of the solutions we get as an estimate of the parametrization uncertainties. When looking at the x-space results we must invert/perform the evolution with each solution and then compute the statical estimators on the final ensemble. The “best” result is always taken as the average on all the possible variations.

In the second stage we apply some “post fit” selection criteria to reduce the number of candidates (to \(\approx 20\)) selecting the most representative elements and discarding clearly unwanted solutions. This way we can achieve a smoother result and improve the speed of the calculation.

  • Among the functions selected at point 3 we cherry pick candidates containing at least one of the leading sub-leading small-N (poles N=0,1) or large-N unknown contributions, such that the spread of the reduced ensemble is not smaller than the full one.

  • By looking at the x-space line integral, we discard any possible outlier that can be generated by numerical cancellations.

The following tables summarize all the considered input functions in the final reduced sets of candidates.

\(\gamma_{gg}^{(3)}\) parametrization basis








\(\frac{1}{N^4},\ \frac{1}{N^3},\ \frac{1}{N^2},\ \frac{1}{(N+1)},\ \frac{1}{(N+2)},\ \mathcal{M}[(1-x)\ln^2(1-x)],\ \mathcal{M}[(1-x)\ln(1-x)]\)

\(\gamma_{gq}^{(3)}\) parametrization basis








\(\frac{1}{N^4},\ \frac{1}{N^3},\ \frac{1}{N^2},\ \frac{1}{(N+1)},\ \frac{1}{(N+2)},\ \mathcal{M}[\ln^2(1-x)],\ \mathcal{M}[\ln(1-x)]\)

Following [MRU+24] we have assumed no violation of the scaling with \(\gamma_{gg}\) also for the NLL small-x term, to help the convergence. We expect that any possible deviation can be parametrized as a shift in the NNLL terms which are free to vary independently.

Slightly different choices are performed for \(\gamma_{gq}^{(3)}\) and \(\gamma_{qq,ps}^{(3)}\) where 10 moments are known. In this case we can select a larger number of functions in group 3 and following [FHMV23a, FHMV23b] we use:

\(\gamma_{qg}^{(3)}\) parametrization basis






\(\frac{1}{N^4},\ \frac{1}{N^3},\ \frac{1}{N^2},\ \frac{1}{N},\frac{1}{N-1}-\frac{1}{N},\ \mathcal{M}[\ln^2(1-x)]\)


\(\mathcal{M}[\ln(x) \ln(1-x)],\ \mathcal{M}[\ln(1-x)],\ \mathcal{M}[(1-x)\ln^3(1-x)],\ \mathcal{M}[(1-x)\ln^2(1-x)],\ \mathcal{M}[(1-x)\ln(1-x)],\ \frac{1}{1+N}\)

\(\gamma_{qq,ps}^{(3)}\) parametrization basis


\(-\frac{1}{(N-1)^2} + \frac{1}{N^2}\)




\(\frac{1}{N^4},\ \frac{1}{N^3},\ \mathcal{M}[(1-x)\ln(1-x)],\ \mathcal{M}[(1-x)^2\ln^2(1-x)],\ \frac{1}{N-1}-\frac{1}{N},\ \mathcal{M}[(1-x)\ln(x)]\)


\(\mathcal{M}[(1-x)(1+2x)],\ \mathcal{M}[(1-x)x^2],\ \mathcal{M}[(1-x) x (1+x)],\ \mathcal{M}[(1-x)]\)

Note that for \(\gamma_{qq,ps},\gamma_{qg}\) the parts proportional to \(n_f^0\) are not present.