#  A description of charmonia decays $J/\psi\rightarrow B\bar{B}$ within the
QCD factorisation framework

Nikolay Kivel
Physik-Department, Technische Universität München,
James-Franck-Str. 1, 85748 Garching, Germany

###### Abstract

Decays $J/\psi\to B\bar{B}$ into octet baryon-antibaryons pairs are studied
within the QCD factorisation framework. The next-to-leading power correction
to the leading-order decay amplitude is calculated for the first time. The
phenomenological analysis also includes the interference with the
electromagnetic decay $J/\psi\to\gamma^{*}\to B\bar{B}$. The branching
fractions and the angular coefficients $\alpha_{B}$ are obtained with an
accuracy of about $20\%$. The obtained results allow us to conclude that the
calculated amplitudes provide a dominant contribution and therefore QCD
factorisation provides a suitable basis for understanding of these decays.

##  1 Introduction

A systematic description of exclusive charmonia decays can be performed using
QCD factorisation framework, see e.g. Refs.[1, 2] and references therein.
Within this approach one constructs the expansion of an amplitude in inverse
powers of $1/m_{c}$ combining non-relativistic QCD and hard-collinear
effective field theory. As a rule, most of these calculations are performed
with leading power accuracy. In many cases this gives a reliable description,
but in many cases this approximation is insufficient even for a good
qualitative description of the data. Potentially such discrepancies can be
associated with various large corrections and with a complicated underlying
QCD dynamics. Therefore, a description of exclusive charmonia decays still
includes many open questions, see, for example, reviews in Refs.[3, 4].

Decays of $S$-wave charmonia into proton-antiproton pairs have been discussed
already long time ago within the QCD factorisation framework in Refs[1, 5]. It
was found that with the appropriate models of baryon matrix elements one can
obtain a reasonable description of the branching ratio for $J/\psi\to
p\bar{p}$. This observation was further used in order to constrain the nucleon
DAs in Refs. [6, 7]

In recent years, the BESII and BESIII collaborations have provided a lot of
new and accurate data in Refs. [8, 9, 10, 11, 12, 13]. The description of
these data is a non-trivial task, since requires calculations beyond the
leading-power approximation, which are often are very problematic. In
addition, the phenomenological consideration suggested in Refs.[14, 15]
indicates potential problems for the description based on QCD factorisation.

The new data provide not only the values of the cross sections, but also
interesting information about the angular behaviour of the $e^{+}e^{-}\to
J/\psi\to B\bar{B}$ cross section. This behaviour is sensitive to the
amplitude, which is power suppressed within the systematic QCD description.
The contribution of this amplitude to the decay width is suppressed by
additional factor $1/m_{c}^{2}$ and therefore is expected to be a small
correction. On the other hand, new data show that the effect of this amplitude
on the description of the angular behaviour is very significant.

The various attempts to describe the angular coefficient $\alpha_{B}$ have
been considered in Refs. [16, 17, 18] using different phenomenological models
and crude approximations. For the first time the power suppressed amplitude
was calculated in the QCD factorisation framework in Refs.[19, 20]. The
qualitative phenomenological analysis carried out in these works indicates
that the corresponding description works sufficiently well. However, in order
to better understand the role of the power corrections some theoretical
uncertainties need to be clarified.

In this paper we continue the investigation of the decays into baryon-
antibaryon pairs within the QCD factorisation framework. Two additional
effects are included into phenomenological consideration, which was ignored in
Ref.[20]. First, we calculate the next-to-leading power correction
$\sim\Lambda^{2}/m^{2}_{c}$ to the dominant leading-power amplitude, where
$\Lambda$ is a typical hadronic scale. This correction is interesting because
the mass of charm quark is not very large, and such contribution could provide
rather large numerical effect. In addition, this correction is described by
the same non-perturbative parameters as the subleading amplitude, which gives
them better phenomenological constraints.

This power correction is closely associated with the higher Fock states of the
baryon wave function and do not involve relativistic corrections with respect
to the small heavy quark velocity in NRQCD. The decay mechanism associated
with annihilation into the three hard gluons suggests that such contributions
can also be described within the collinear factorisation framework that gives
us an opportunity of a more detailed understanding of the systematic expansion
in the inverse mass $1/m_{c}$. The higher twist baryon light-cone distribution
amplitudes (DAs), associated with the non-perturbative collinear matrix
elements, have already been investigated in Refs. [21, 22, 23, 24, 25, 26, 27,
28]. These results allows one to reduce theoretical uncertainties associated
with the non-perturbative QCD input.

In this work we also discuss the effect that occurs due to the contribution of
the electromagnetic amplitude $J/\psi\to\gamma^{*}\to B\bar{B}$. The given
amplitude is sensitive to the baryon time-like electromagnetic form factors.
For a given energy scale these form factors can not be computed within the QCD
framework because of relatively large effect from the soft-overlap
rescattering (or Feynman) mechanism. Time-like baryonic FFs are also poorly
studied experimentally. Different cross sections for the $e^{+}e^{-}\to
B\bar{B}$ channels have been measured by BESIII recent years, and the so-
called effective FFs are obtained from data fitting, see Refs. [29, 30, 31,
32, 33, 34, 35, 36]. We use these results in order to estimate the numerical
influence of the electromagnetic admixture in the angular distribution
$\alpha_{B}$.

Our paper is organised as follows. In Sec. 2 we discuss the higher twist
baryon DAs and describe the models, which are used in our calculations. The
calculation and the analytical expressions for the hard kernels and for the
total amplitudes are discussed in the Sec. 3. The phenomenological analysis is
presented in Sec. 4. In Sec. 5 we summarise our findings. In Appendix we
provide additional useful details and a more detailed discussion of some
technical issues.

##  2  Baryon light-cone distribution amplitudes

###  2.1 The basic set of the twist-5 DAs and their properties

In this section we describe the properties of twist-5 DAs and define the
twist-5 auxiliary matrix elements, which are used in our calculation of the
power corrections. We consider contributions with only three quark operators,
which greatly simplifies the calculations. The total power correction is given
by the sum of various contributions with the different combinations of DAs.
The leading-power (lp) contribution to amplitude $A_{1}$ only includes the
integral with the twist-3 DAs $\varphi_{3}$, schematically:
$A_{1}^{(\text{lp})}\sim\varphi_{3}\ast H\ast\varphi_{3}$, where the asterisk
is used as a short notation for the convolution integral with respect to quark
momentum fractions. The next-to-leading power (nlp) correction is more
complicated, it depends on twist-4 $\varphi_{4}$ and twist-5 $\varphi_{5}$ DAs
and can be represented as

$A_{1}^{(\text{nlp})}\sim\varphi_{3}\ast H_{35}\ast\varphi_{5}+\varphi_{4}\ast
H_{44}\ast\varphi_{4}+\frac{m_{B}^{2}}{m_{Q}^{2}}\varphi_{3}\ast
H_{33}\ast\varphi_{3},$ (1)

where the last term can be associated with the kinematical contribution. The
various higher twist DAs and their properties have been discussed in Refs.[23,
25, 26, 27] and we use these results in the present work.

Let us define the light-cone matrix elements, which define the light-cone DAs.
The analysis of the twist-4 auxiliary DAs are given in Refs.[19, 20, 25].
Therefore, here we will focus only on the discussion of the twist-5 matrix
elements.

We use the same notations as in Refs.[19, 20]. The light-like vectors
$n,\bar{n}$ ( $n^{2}=\bar{n}^{2}=0$, $(n\bar{n})=2$ ) define the longitudinal
and transverse components of a Lorentz vector $V^{\mu}$

$V_{+}=(kn),\quad V_{-}=(V\bar{n}),\ \ \
V^{\mu}=V_{+}\frac{\bar{n}^{\mu}}{2}+V_{-}\frac{n^{\mu}}{2}+V_{\bot}^{\mu}.$
(2)

In the following we also assume that the baryon momentum $k$ is directed along
$z$-axis and can be expanded as

$k=k_{-}\frac{{n}}{2}+k_{+}\frac{\bar{n}}{2},\ \ k_{+}\gg
k_{-}=\frac{m_{B}^{2}}{k_{+}},~{}~{}k^{2}=m_{B}^{2}.$ (3)

For the baryon spinor $N(k,s)$ we define the large and small components
$N_{\bar{n}}$ and $N_{{n}}$, respectively:

$N_{\bar{n}}=\frac{\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\hbox to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n}{4}N(k,s),~{}\ ~{}N_{{n}}=\frac{\hbox
to0.0pt{\hbox to6.00235pt{\hfil/\hfil}\hss}n\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/}{4}N(k,s)=\frac{m_{B}}{k_{+}}\frac{\hbox
to0.0pt{\hbox to6.00235pt{\hfil/\hfil}\hss}n}{2}N_{\bar{n}}.~{}$ (4)

The similar relations are also used for the antibaryon spinors.

The relevant three-quark light-cone operators are constructed from the QCD
quark fields $q=\\{{u,d,s\\}}$ and from the light-cone the Wilson lines

$W_{\bar{n}}[x_{-},z_{-}]=\text{P}\exp\left\\{ig\int_{(z_{-}-x_{-})/2}^{0}ds\
A_{+}(x_{-}n/2+sn)\right\\}.$ (5)

The basic operator can be written as

$\mathcal{O}_{\alpha_{1}\alpha_{2}\alpha_{3}}(x_{-},y_{-},z_{-})=\varepsilon^{ijk}\
q_{\alpha_{1}}^{i^{\prime}}(x_{-})W_{\bar{n}}[x_{-},v_{-}]_{i^{\prime}i}\
q_{\alpha_{2}}^{j^{\prime}}(y_{-})W_{\bar{n}}[y_{-},v_{-}]_{j^{\prime}j}\
q_{\alpha_{3}}^{k^{\prime}}(z_{-})W_{\bar{n}}[z_{-},v_{-}]_{k^{\prime}k},$ (6)

where we use the short notation $q(x_{-})\equiv q(x_{-}n/2)$. Below, for
brevity, we also simplify the notation of the Dirac indices setting
$\\{\alpha_{1},\alpha_{2},\alpha_{3}\\}\equiv\\{1,2,3\\}$ and we also will not
show the colour structure. Following to Ref. [26] we assume the following
flavour content of the operators $\langle
0|q_{\alpha_{1}}q_{\alpha_{2}}q_{\alpha_{3}}|B\rangle\equiv\langle
0|q_{1}q_{2}q_{3}|B\rangle$:

$\displaystyle\langle 0|u_{1}u_{2}d_{3}|p\rangle,\,\langle
0|u_{1}d_{2}s_{3}|\Lambda\rangle,$ (7) $\displaystyle\langle
0|u_{1}d_{2}s_{3}|\Sigma^{0}\rangle,\,\langle
0|u_{1}u_{2}s_{3}|\Sigma^{+}\rangle,$ (8) $\displaystyle\langle
0|s_{1}s_{2}u_{3}|\Xi^{0}\rangle,\,\langle 0|s_{1}s_{2}d_{3}|\Xi^{-}\rangle.$
(9)

For simplicity, this will not be shown explicitly.

Let us write the light-cone matrix element for the baryon state $B$ as

$\displaystyle\left\langle
0\right|\mathcal{O}_{123}(x_{-},y_{-},z_{-})\left|B(k)\right\rangle$
$\displaystyle=\left\langle
0\right|\mathcal{O}_{123}\left|B(k)\right\rangle_{\text{tw3}}+\left\langle
0\right|\mathcal{O}_{123}\left|B(k)\right\rangle_{\text{tw4}}$
$\displaystyle+\left\langle
0\right|\mathcal{O}_{123}\left|B(k)\right\rangle_{\text{tw5}}+\ldots,$ (10)

where the subscript “tw3” indicate that we only count in this term the
contributions of twist-3 and so on. The leading twist-3 contribution read

$\displaystyle\left\langle
0\right|\mathcal{O}_{123}\left|B(k)\right\rangle_{\text{tw3}}$
$\displaystyle=~{}\frac{1}{8}k_{+}\left[\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/C\right]_{12}~{}\left[\gamma_{5}N_{\bar{n}}\right]_{3}\text{FT}\left[V_{1}^{B}\right]+~{}\frac{1}{8}k_{+}\left[\hbox
to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\gamma_{5}C\right]_{12}~{}\left[N_{\bar{n}}\right]_{3}\text{FT}\left[A_{1}^{B}\right]$
$\displaystyle+\frac{1}{8}k_{+}\left[\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\gamma_{\bot}^{\alpha}C\right]_{12}~{}\left[\gamma_{\bot}^{\alpha}\gamma_{5}N_{\bar{n}}\right]_{3}\text{FT}\left[T_{1}^{B}\right].$
(11)

The twist-4 contribution is defined as

$\displaystyle\left\langle
0\right|\mathcal{O}_{123}\left|B(k)\right\rangle_{\text{tw4}}=~{}\frac{m_{B}}{16}\left[\hbox
to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/C\right]_{12}~{}\left[\gamma_{5}\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}nN_{\bar{n}}\right]_{3}\text{FT}\left[V_{2}^{B}\right]+\frac{m_{B}}{8}\left[\gamma_{\bot}C\right]_{12}~{}\left[\gamma_{\bot}\gamma_{5}N_{\bar{n}}\right]_{3}~{}\text{FT}\left[V_{3}^{B}\right]$
$\displaystyle+~{}\frac{m_{B}}{16}\left[\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\gamma_{5}C\right]_{12}~{}\left[\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}nN_{\bar{n}}\right]_{3}\text{FT}\left[A_{2}^{B}\right]+\frac{m_{B}}{8}\left[\gamma_{\bot}\gamma_{5}C\right]_{12}~{}\left[\gamma_{\bot}N_{\bar{n}}\right]_{3}~{}\text{FT}\left[A_{3}^{B}\right]$
$\displaystyle+\frac{m_{B}}{4}\left[C\right]_{12}~{}\left[\gamma_{5}N_{\bar{n}}\right]_{3}\text{FT}\left[S_{1}^{B}\right]+\frac{m_{B}}{4}\left[\gamma_{5}C\right]_{12}~{}\left[N_{\bar{n}}\right]_{3}\text{FT}\left[P_{1}^{B}\right]$
$\displaystyle+\frac{m_{B}}{16}\left[\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\gamma_{\bot}^{\alpha}C\right]_{12}~{}\left[\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\gamma_{\bot\alpha}\gamma_{5}N_{\bar{n}}\right]_{3}\text{FT}\left[T_{2}^{B}\right]+\frac{m_{B}}{8}k_{+}\left[i\sigma^{\bar{n}n}C\right]_{12}~{}\left[\gamma_{5}N_{\bar{n}}\right]_{3}\text{FT}\left[T_{3}^{B}\right]$
$\displaystyle+\frac{m_{B}}{8}\left[\sigma_{\bot\bot}^{\alpha\beta}C\right]_{12}~{}\left[\sigma_{\alpha\beta}\gamma_{5}N_{\bar{n}}\right]_{3}\text{FT}\left[T_{7}^{B}\right].$
(12)

The twist-5 contribution read

$\displaystyle\left\langle
0\right|\mathcal{O}_{123}\left|B(k)\right\rangle_{\text{tw5}}=~{}\frac{1}{4}\frac{m_{B}^{2}}{4k_{+}}\left[\gamma_{\bot}^{\alpha}C\right]_{12}~{}\left[\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\gamma_{\bot\alpha}\gamma_{5}N_{\bar{n}}\right]_{3}~{}\text{FT}\left[V_{4}^{B}\right]+\frac{1}{4}~{}\frac{m_{B}^{2}}{4k_{+}}\left[\gamma_{\bot}^{\alpha}\gamma_{5}C\right]_{12}~{}\left[\gamma_{\bot\alpha}\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}nN_{\bar{n}}\right]_{3}\text{FT}\left[A_{4}^{B}\right]$
$\displaystyle+~{}~{}\frac{1}{4}\frac{m_{B}^{2}}{2k_{+}}\left[\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}nC\right]_{12}~{}\left[\gamma_{5}N_{\bar{n}}\right]_{3}\text{FT}\left[V_{5}^{B}\right]+~{}\frac{1}{4}\frac{m_{B}^{2}}{2k_{+}}\left[\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\gamma_{5}C\right]_{12}~{}\left[N_{\bar{n}}\right]_{3}\text{FT}\left[A_{5}^{B}\right]$
$\displaystyle-\frac{1}{4}\frac{m_{B}^{2}}{2k_{+}}~{}\left[C\right]_{12}\left[\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\gamma_{5}N_{\bar{n}}\right]_{3}~{}\text{FT}\left[S_{2}^{B}\right]+\frac{1}{4}\frac{m_{B}^{2}}{2k_{+}}~{}\left[\gamma_{5}C\right]_{12}\left[\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}nN_{\bar{n}}\right]_{3}~{}\text{FT}\left[P_{2}^{B}\right]$
$\displaystyle+\frac{1}{4}\frac{m_{B}^{2}}{4k_{+}}\left[i\sigma_{\bar{n}n}C\right]_{12}\left[\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\gamma_{5}N_{\bar{n}}\right]_{3}~{}\text{FT}\left[T_{4}^{B}\right]+\frac{1}{4}\frac{m_{B}^{2}}{2k_{+}}\left[i\sigma_{\alpha
n}C\right]_{12}\left[\gamma_{\bot}^{\alpha}\gamma_{5}N_{\bar{n}}\right]_{3}~{}\text{FT}\left[T_{5}^{B}\right]$
$\displaystyle-\frac{1}{4}\frac{m_{B}^{2}}{4k_{+}}\left[\sigma_{\bot\bot}^{\alpha\beta}C\right]_{12}\left[\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\sigma_{\alpha\beta}\gamma_{5}N_{\bar{n}}\right]_{3}~{}\text{FT}\left[T_{8}^{B}\right],$
(13)

where we always assume that
$\sigma_{\bar{n}n}=\sigma_{\alpha\beta}\bar{n}^{\alpha}n^{\beta}$ . The symbol
“FT” denotes the Fourier transformation with respect to the quark momentum
fractions

$\text{FT}[f]=\int
Du_{i}~{}e^{-i(k_{1}x)-i(k_{2}y)-i(k_{3}z)}f(u_{1},u_{2},u_{3}),~{}\ $ (14)

with the measure $Du_{i}=du_{1}du_{2}du_{3}\delta(1-u_{1}-u_{2}-u_{3})$, the
quark momenta in (14) are defined as

$k_{i}=u_{i}k_{+}\frac{\bar{n}}{2}.$ (15)

The DAs defined in Eqs.(11)-(13) are symmetric or antisymmetric with respect
to interchange $u_{1}\leftrightarrow u_{2}$

$F_{i}^{B}(x_{2},x_{1},x_{3})=-(-1)_{B}F_{i}^{B}(x_{1},x_{2},x_{3}),\
\text{for }\ \ F=\\{S,P,A\\},$ (16)
$F_{i}^{B}(x_{2},x_{1},x_{3})=+(-1)_{B}F_{i}^{B}(x_{1},x_{2},x_{3}),\
\text{for }\ \ F=\\{V,T\\},$ (17)

with

$(-1)_{B}=\left\\{\begin{array}[c]{c}+1\ \ \ B\neq\Lambda\\\ -1\ \ \
B=\Lambda\end{array}\right..\text{ }$ (18)

The DAs in Eqs(11)-(13) can be rewritten in terms of the basis light-cone DAs
as [26, 27] ( below we use short notation
$(x_{i})\equiv(x_{123})\equiv(x_{1},x_{2},x_{3})$)

$f_{B}\varphi_{3\pm}^{B}(x_{123})=\frac{c_{B}^{\pm}}{2}\left\\{\left(V_{1}-A_{1}\right)^{B}(x_{123})\pm\left(V_{1}-A_{1}\right)^{B}(x_{321})\right\\},$
(19) $f_{B}\Pi_{3}^{B}(x_{123})=c_{B}^{-}(-1)_{B}T_{1}^{B}(x_{132}),$ (20)
$\Phi_{4\pm}^{B}(x_{123})=c_{B}^{\pm}\left\\{\left(V_{2}-A_{2}\right)^{B}(x_{123})\pm(-1)_{B}\left(V_{3}-A_{3}\right)^{B}(x_{231})\right\\},$
(21)
$\Xi_{4\pm}^{B}(x_{123})=(-1)_{B}3c_{B}^{\pm}\left[~{}\left(T_{3}-T_{7}+S_{1}+P_{1}\right)^{B}(x_{123})\pm\left(T_{3}-T_{7}+S_{1}+P_{1}\right)^{B}(x_{132})\right],$
(22)
$\Pi_{4}^{B}(x_{123})=c_{B}^{-}\left(T_{3}+T_{7}+S_{1}-P_{1}\right)^{B}(x_{312}),$
(23) $\Upsilon_{4}^{B}(x_{123})=6c_{B}^{-}T_{2}^{B}(x_{321}),$ (24)
$\Phi_{5\pm}^{B}(x_{123})=c_{B}^{\pm}\left\\{\left(V_{2}-A_{2}\right)^{B}(x_{123})\pm(-1)_{B}\left(V_{3}-A_{3}\right)^{B}(x_{231})\right\\},$
(25)
$\Xi_{5\pm}^{B}(x_{123})=(-1)_{B}3c_{B}^{\pm}\left[~{}\left(T_{4}-T_{8}+S_{2}+P_{2}\right)^{B}(x_{123})\pm\left(T_{4}-T_{8}+S_{2}+P_{2}\right)^{B}(x_{132})\right],$
(26)
$\Pi_{5}^{B}(x_{123})=c_{B}^{-}\left(T_{4}+T_{8}+S_{2}-P_{2}\right)^{B}(x_{312}),$
(27) $\Upsilon_{5}^{B}(x_{123})=6c_{B}^{-}T_{5}^{B}(x_{321}),$ (28)

where

$c_{B}^{+}=\left\\{\begin{array}[c]{c}1,~{}B\neq\Lambda\\\
\sqrt{\frac{2}{3}},B=\Lambda\end{array}\right.,~{}~{}\
c_{B}^{-}=\left\\{\begin{array}[c]{c}1,~{}B\neq\Lambda\\\
-\sqrt{6},B=\Lambda\end{array}\right.,$ (29)

The coefficients $c_{B}^{\pm}$ are chosen in such a way that these DAs satisfy
the simple relations in the $SU(3)$-flavour symmetry limit, see e.g.
Ref.[26].

For the nucleon state, the basic set of DAs can be further simplified by the
isospin relationships, namely, it is sufficient to determine

$f_{N~{}}\varphi_{3}(x_{123})=\phi_{3+}^{N}(x_{123})+\phi_{3-}^{N}(x_{123})=\left(V_{1}-A_{1}\right)^{N}(x_{123}),$
(30)
$\frac{1}{2}\Phi_{4}(x_{123})=\frac{1}{2}\left(\Phi_{4+}^{N}+\Phi_{4-}^{N}\right)(x_{123})=\left(V_{2}-A_{2}\right)^{N}(x_{123}),$
(31)
$\frac{1}{2}\Psi_{4}(x_{123})=\frac{1}{2}\left(\Phi_{4+}^{N}-\Phi_{4-}^{N}\right)(x_{312})=\left(V_{3}-A_{3}\right)^{N}(x_{231}),$
(32)
$\frac{\lambda_{2}}{6}\Xi_{4}(x_{123})=\left(\Xi_{4+}^{N}+\Xi_{4-}^{N}\right)(x_{123})=\left(T_{3}-T_{7}+S_{1}+P_{1}\right)^{N}(x_{123}).$
(33)
$\frac{1}{2}\Phi_{5}(x_{123})=\frac{1}{2}\left(\Phi_{5+}^{N}+\Phi_{5-}^{N}\right)(x_{123})=\left(V_{5}-A_{5}\right)^{N}(x_{123}),$
(34)
$\frac{1}{2}\Psi_{5}(x_{123})=\frac{1}{2}\left(\Phi_{5+}^{N}-\Phi_{5-}^{N}\right)(x_{312})=\left(V_{4}-A_{4}\right)^{N}(x_{231}),$
(35)
$\frac{1}{6}\Xi_{5}(x_{123})=\frac{1}{6}\left(\Xi_{5+}^{N}+\Xi_{5-}^{N}\right)(x_{123})=\left(T_{4}-T_{8}+S_{2}+P_{2}\right)^{N}(x_{123}).$
(36)

The twist-4 and twist-5 DAs can be decomposed into components corresponding to
the contributions of the light-cone operators of geometrical twist-3, twist-4
and twist-5. Let us write such decompositions in the following form

$\Phi_{4\pm}^{B}(x_{123})=f_{B}\Phi_{4\pm}^{B(3)}(x_{123})+\lambda_{1}^{B}\bar{\Phi}_{4\pm}^{B}(x_{123}),$
(37)
$\Pi_{4}^{B}(x_{i})=f_{\bot}^{B}\Pi_{4}^{B(3)}(x_{123})+\lambda^{B}_{1}\bar{\Pi}_{4}^{B}(x_{i}),\
B\neq\Lambda,$ (38)
$\Pi_{4}^{\Lambda}(x_{i})=f_{\Lambda}\Pi_{4}^{\Lambda(3)}(x_{123})+\lambda^{\Lambda}_{\bot}\bar{\Pi}_{4}^{\Lambda}(x_{i}),\
$ (39) $\Xi_{4\pm}^{B}(x_{i})=\lambda_{2}^{B}\bar{\Xi}_{4}^{B}(x_{i}),$ (40)
$\Upsilon_{4}^{B}(x_{i})=\lambda_{2}^{B}\bar{\Upsilon}_{4}^{B}(x_{i}).$ (41)
$\Phi_{5\pm}^{B}(x_{i})=f_{B}~{}\Phi_{5\pm}^{B(3)}(x_{i})+\lambda_{1}^{B}\Phi_{5\pm}^{B(4)}(x_{i})+\bar{\Phi}_{5\pm}^{B}(x_{i}),$
(42)
$\Xi_{5\pm}^{B}(x_{i})=\lambda_{2}^{B}\Xi_{5}^{B(4)}(x_{i})+\bar{\Xi}_{5\pm}^{B}(x_{i}),$
(43)
$\Pi_{5}^{B}(x_{i})=f_{\bot}^{B}\Pi_{5}^{B(3)}(x_{i})+f_{\bot}^{B}\Pi_{5}^{B(4)}(x_{i})+\bar{\Pi}_{5}^{B}(x_{i}),$
(44)
$\Pi_{5}^{\Lambda}(x_{i})=\lambda_{\bot}^{\Lambda}\Pi_{5}^{\Lambda(3)}(x_{i})+\lambda_{\bot}^{\Lambda}\Pi_{5}^{\Lambda(4)}(x_{i})+\bar{\Pi}_{5}^{\Lambda}(x_{i}),$
(45)
$\Upsilon_{5}^{B}(x_{i})=\lambda_{2}^{B}\Upsilon_{5}^{B(4)}(x_{i})+\bar{\Upsilon}_{5}^{B}(x_{i}).$
(46)

Here the superscript $(i)$ indicates that the corresponding function is
completely defined in terms of the moments of the light-cone operators of
geometrical twist-$i$. The functions with the bar denote the genuine twist-4/5
contributions. The moments of genuine twist-5 DAs are not known. In the
present work we assume that the corresponding matrix elements are sufficiently
small and therefore be neglected that yields

$\bar{\Phi}_{5\pm}^{B}(x_{i})=\bar{\Xi}_{5\pm}^{B}(x_{i})=\bar{\Pi}_{5}^{B}(x_{i})=\bar{\Pi}_{5}^{\Lambda}(x_{i})=\bar{\Upsilon}_{5}^{B}(x_{i})=0.$
(47)

The normalisation couplings $f_{B},\ f_{\bot}^{B},\ \lambda_{1,2}^{B}$ and
$\lambda_{\bot}^{\Lambda}$ have dimension of GeV2 and their values will be
discussed later.

The conformal expansion and evolution of the moments for different three-quark
DAs was discussed in Refs. [22, 23, 24]. The explicit analytical expressions
for the different terms in rhs in Eqs.(42)-(46) can be found in Refs.[23, 26,
27]. For convenience, in Appendix A we provide the explicit expressions for
the models of basic DAs, which are used in present work.

###  2.2 Auxiliary DAs and their properties

Let us also consider the auxiliary  set of the matrix elements, which is more
convenient for the actual calculation of the hard kernels. More details about
this calculation will be given in the next section. Here we only discuss the
formal definition of the auxiliary matrix elements and their relations with
the basic set of DAs.

The calculation of the decay amplitudes implies that the QCD fields can be
splitted into sum of hard and collinear fields. The letter describe the non-
perturbative overlap with the baryon-antibaryon state. Therefore performing a
matching onto collinear operators one implies the effective field theory
constructed from the hard ($p_{h}\sim m_{Q}$, $x\sim 1/m_{Q}$) and collinear
fields. Such theory provides a definite power counting, which is used for the
ranging of the operators. The power counting of the collinear fields is
associated with their momentum, if we assume the collinear sector associated
with the momentum $k$ in Eq.(3) then ($\lambda^{2}\sim\Lambda/m_{Q}$)

$\displaystyle\frac{\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\hbox to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n}{4}q(x_{-})\sim\mathcal{O}(\lambda^{2}),\ \
\frac{\hbox to0.0pt{\hbox to6.00235pt{\hfil/\hfil}\hss}n\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/}{4}q(x_{-})\sim\mathcal{O}(\lambda^{4}),\
\ $ (48) $\displaystyle\ A_{+}(x_{-})\sim\mathcal{O}(1),\
A_{\bot}(x_{-})\sim\mathcal{O}(\lambda^{2}),\
A_{-}(x_{-})\sim\mathcal{O}(\lambda^{4}).$ (49)

Performing the multipole expansions of collinear fields one also obtains
derivatives of the collinear fields, which counts as

$\displaystyle(n\partial)q(x_{-})$ $\displaystyle\sim\mathcal{O}(1)q(x_{-}),\
\ \partial_{\bot}q(x_{-})\sim\mathcal{O}(\lambda)q(x_{-}),\ \ \
(\bar{n}\partial)q(x_{-})\sim\mathcal{O}(\lambda^{4})q(x_{-}),\ $ (50)
$\displaystyle(n\partial)A_{+}(x_{-})$
$\displaystyle\sim\mathcal{O}(1)A_{+}(x_{-}),\partial_{\bot}A_{+}(x_{-})\sim\mathcal{O}(\lambda^{2})A_{+}(x_{-}),(\bar{n}\partial)A_{+}(x_{-})\sim\mathcal{O}(\lambda^{4})A_{+}(x_{-}).$
(51)

Remind, that we do not consider the contributions with the quark-gluon
operators, therefore terms like $A_{\bot},$ $A_{-}$ , $\partial_{\bot}A_{+}$
and $(\bar{n}\partial)A_{+}$ will be excluded from our consideration.

The auxiliary DAs can be defined in terms of the matrix elements of 3-quark
twist-5 operators with the transverse derivatives. These operators are defined
using the gauge invariant collinear fields, which have a definite scaling
behaviour within the effective field theory framework

$\chi_{1}(x_{-})=W^{{\dagger}}(x_{-})\frac{\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\hbox to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n}{4}q_{1}(x_{-}),\ \ \ \ \
W(x_{-})=\text{P}\exp\left\\{ig\int_{-\infty}^{0}ds\
A_{+}(x_{-}n/2+sn)\right\\}.$ (52)

Obviously

$\chi(x_{-})\sim\mathcal{O}(\lambda^{2}),\
\partial_{\bot}\chi(x_{-})\sim\mathcal{O}(\lambda^{4}),\ \
(\bar{n}\partial)\chi(x_{-})\sim\mathcal{O}(\lambda^{6}).$ (53)

The most general 3-quark twist-3 operator is

$\chi_{1}(x_{-})\chi_{2}(y_{-})\chi_{3}(z_{-})\sim\mathcal{O}(\lambda^{6}),$
(54)

where for simplicity we do not show the colour and flavour structure.

The 3-quark twist-4 operators include one transverse derivative

$\left[\partial_{\bot}\chi_{1}\right]\chi_{2}\chi_{3}\sim\chi_{1}\left[\partial_{\bot}\chi_{2}\right]\chi_{3}\sim\chi_{1}\chi_{2}\left[\partial_{\bot}\chi_{3}\right]\sim\mathcal{O}(\lambda^{8}),$
(55)

where the derivative is applied inside the brackets only. Notice that we do
not need to write the covariant derivative because the field $\chi(x_{-})$ is
gauge invariant. The matrix elements of the twist-4 operators in (55) are
considered in Ref.[19]. The results read

$~{}\left\langle
0\right|\left[i\partial_{\bot}^{\alpha}\chi(x_{-})\right](x_{+})C\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\chi(y_{-})\chi_{3}(z_{-})\left|B(k)\right\rangle=k_{+}m_{B}\left[\gamma_{\bot}^{\alpha}\gamma_{5}N_{\bar{n}}\right]_{3}\text{FT}\left[\mathcal{V}_{1}^{B}\right],$
(56) $\left\langle 0\right|~{}\chi(x_{-})C\hbox to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\left[i\partial_{\bot}^{\alpha}\chi(y_{-})\right]\chi_{3}(z_{-})\left|B(k)\right\rangle=k_{+}m_{B}~{}\left[\gamma_{\bot}^{\alpha}\gamma_{5}N_{\bar{n}}\right]_{3}\text{FT}\left[\mathcal{V}_{2}^{B}\right].$
(57) $\left\langle
0\right|\left[i\partial_{\bot}^{\alpha}\chi(x_{-})\right]C\hbox to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\gamma_{5}\chi(y_{-})\chi_{3}(z_{-})\left|B(k)\right\rangle=k_{+}m_{B}\left[~{}\gamma_{\bot}^{\alpha}N_{\bar{n}}\right]_{3}\text{FT}\left[\mathcal{A}_{1}^{B}\right],$
(58) $\left\langle 0\right|\chi(x_{-})C\hbox to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\gamma_{5}\left[i\partial_{\bot}^{\alpha}\chi(y_{-})\right]\chi_{3}(z_{-})\left|B(k)\right\rangle=k_{+}m_{B}\left[\gamma_{\bot}^{\alpha}N_{\bar{n}}\right]_{3}\text{FT}\left[\mathcal{A}_{2}^{B}\right].$
(59) $\displaystyle\left\langle
0\right|\left[i\partial_{\bot}^{\alpha}\chi(x_{-})\right]C\hbox to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\gamma_{\bot}^{\nu}\chi(y_{-})\chi_{3}(z_{-})\left|B(k)\right\rangle$
$\displaystyle=~{}k_{+}m_{B}~{}g_{\bot}^{\alpha\nu}\left[\gamma_{5}N_{\bar{n}}\right]_{3}\text{FT}\left[\mathcal{T}_{21}^{B}\right]$
(60)
$\displaystyle+k_{+}m_{B}\left[i\sigma_{\bot\bot}^{\nu\alpha}\gamma_{5}N_{\bar{n}}\right]_{3}\text{FT}\left[\mathcal{T}_{41}^{B}\right],$
(61) $\displaystyle\left\langle 0\right|\chi(x_{-})C\hbox to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\gamma_{\bot}^{\nu}\left[i\partial_{\bot}^{\alpha}\chi(y_{-})\right]\chi_{3}(z_{-})\left|B(k)\right\rangle$
$\displaystyle=k_{+}m_{B}~{}g_{\bot}^{\alpha\nu}\left[\gamma_{5}N_{\bar{n}}\right]_{3}\text{FT}\left[\mathcal{T}_{22}^{B}\right]$
$\displaystyle+k_{+}m_{B}\left[i\sigma_{\bot\bot}^{\nu\alpha}\gamma_{5}N_{\bar{n}}\right]_{3}\text{FT}\left[\mathcal{T}_{42}^{B}\right],$
(62)

where the color and the flavour structure of the operators are the same as
described before. The DAs in the rhs of these equations can be expressed in
terms of DAs $\\{V_{i},A_{i},T_{i},S_{i},P_{i}\\}^{B}$ using the Lorentz
symmetry and QCD equations of motion. This gives111 If the all function in an
equation have the same standard arguments $(x_{i})\equiv(x_{1},x_{2},x_{3})$
we use the simplification $F(x_{i})\equiv F$ to simplify the equations.  [25,
19]

$\displaystyle 4\mathcal{V}_{i}^{B}$
$\displaystyle=x_{3}\left(V_{2}+A_{2}\right)^{B}+(-1)^{i}\left[(x_{1}-x_{2})V_{3}^{B}+\bar{x}_{3}A_{3}^{B}\right]+\frac{m_{3}}{m_{B}}\left(V_{1}^{B}+(-1)^{i}A_{1}^{B}\right)$
$\displaystyle-2(-1)^{i}\frac{m_{1}-m_{2}}{m_{B}}T_{1}^{B},$ (63)
$\displaystyle 4\mathcal{A}_{i}^{B}$
$\displaystyle=-x_{3}\left(V_{2}+A_{2}\right)^{B}+(-1)^{i}\left[(x_{1}-x_{2})A_{3}^{B}+\bar{x}_{3}V_{3}^{B}\right]+\frac{m_{3}}{m_{B}}\left(A_{1}^{B}+(-1)^{i}V_{1}^{B}\right)$
$\displaystyle-2(-1)^{i}\frac{m_{1}+m_{2}}{m_{B}}T_{1}^{B},$ (64)
$\displaystyle\left(\mathcal{T}_{21}+\mathcal{T}_{41}\right)^{B}$
$\displaystyle=\frac{x_{1}}{2}\left(T_{3}-T_{7}+S_{1}+P_{1}\right)^{B}-\frac{1}{2}\frac{m_{1}}{m_{B}}\left(V_{1}^{B}+A_{1}^{B}\right),\
$ (65) $\displaystyle\ \ \ \left(\mathcal{T}_{22}+\mathcal{T}_{42}\right)^{B}$
$\displaystyle=\frac{x_{2}}{2}\left(T_{3}-T_{7}-S_{1}-P_{1}\right)^{B}-\frac{1}{2}\frac{m_{2}}{m_{B}}\left(V_{1}^{B}-A_{1}^{B}\right),$
(66) $\displaystyle\left(\mathcal{T}_{21}-\mathcal{T}_{41}\right)^{B}$
$\displaystyle=\frac{x_{1}}{2}\left(T_{3}+T_{7}+S_{1}-P_{1}\right)^{B}+\frac{1}{2}\frac{m_{1}}{m_{B}}\left(A_{1}^{B}-V_{1}^{B}\right),\
\ \ $ (67) $\displaystyle\ \left(\mathcal{T}_{22}-\mathcal{T}_{42}\right)^{B}$
$\displaystyle=\frac{x_{2}}{2}\left(T_{3}+T_{7}-S_{1}+P_{1}\right)^{B}-\frac{1}{2}\frac{m_{2}}{m_{B}}\left(V_{1}^{B}+A_{1}^{B}\right),$
(68)

where the quark masses $m_{i}$ correspond to quarks $q_{1}q_{2}q_{3}$ in the
matrix elements in Eqs.(7)-(9), respectively.

The additional relation is derived in Appendix B

$T_{2}^{B}=\frac{1}{x_{3}}(\mathcal{T}_{21}+\mathcal{T}_{22}+\mathcal{T}_{41}+\mathcal{T}_{42})^{B}+\frac{~{}m_{3}}{m_{B}}T_{1}^{B}.$
(69)

Combaining Eq.(69) with Eqs.(65), (66) and definitions in Eqs. (22) and (24)
one finds for $\Upsilon_{4}^{B}$ the following equation

$\displaystyle\Upsilon_{4}(x_{123})$
$\displaystyle=\frac{c_{B}^{-}}{c_{B}^{+}}\frac{1}{x_{1}}\left\\{(-1)_{B}\
x_{3}\Xi_{4+}^{B}(x_{321})+x_{2}\Xi_{4+}^{B}(x_{231})\right\\}$
$\displaystyle+\frac{1}{x_{1}}\left\\{(-1)_{B}\
x_{3}\Xi_{4-}^{B}(x_{321})+x_{2}\Xi_{4-}^{B}(x_{231})\right\\}$
$\displaystyle+\frac{3c_{B}^{-}}{x_{1}}\left\\{\frac{m_{2}-m_{1}}{m_{B}}A_{1}^{B}(x_{321})-\frac{m_{1}+m_{2}}{m_{B}}V_{1}^{B}(x_{321})+\frac{~{}2m_{3}}{m_{B}}T_{1}^{B}(x_{321})\right\\}.$
(70)

Therefore the DA $\Upsilon_{4}^{B}$ is not an independent function in the
3-quark approximation.

The 3-quark twist-5 operators include two transverse derivatives

$\partial_{\bot}^{2}\chi_{1}\chi_{2}\chi_{3}\sim\mathcal{O}(\lambda^{10}).$
(71)

The longitudinal derivatives $(\bar{n}\partial)\chi$ can be excluded with the
help of QCD equation of motion

$(i\bar{n}\partial)\chi=(in\partial)^{-1}\partial_{\bot}^{2}\chi.$ (72)

Therefore the basis set of the 3-quarks operators can be defined as

$\displaystyle\left[\partial_{\bot}^{\alpha}\partial_{\bot}^{\beta}\chi_{1}\right]\chi_{2}\chi_{3},\quad\chi_{1}\left[\partial_{\bot}^{\alpha}\partial_{\bot}^{\beta}\chi_{2}\right]\chi_{3},\quad\left[\partial_{\bot}^{\alpha}\chi_{1}\right]\left[\partial_{\bot}^{\beta}\chi_{2}\right]\chi_{3}.$
(73)

The operators with $\left[\partial_{\bot}\chi_{3}\right]$ can be excluded
because

$\chi_{1}\chi_{2}\left[\partial_{\bot}\chi_{3}\right]=\partial_{\bot}\left[\chi_{1}\chi_{2}\chi_{3}\right]-\left[\partial_{\bot}\chi_{1}\right]\chi_{2}\chi_{3}-\chi_{1}\left[\partial_{\bot}\chi_{2}\right]\chi_{3}\simeq-\left[\partial_{\bot}\chi_{1}\right]\chi_{2}\chi_{3}-\chi_{1}\left[\partial_{\bot}\chi_{2}\right]\chi_{3},$
(74)

where it is used that the matrix element from the operator with the total
transverse derivative vanishes because the baryon momenta have longitudinal
components only. The matrix elements of the three operators in Eq.(73) define
the auxiliary set of the twist-5 DAs. It is convenient to write their
definitions for the operators with projected Dirac indices $\\{12\\}$. We
introduce the following set of the matrix elements

$\left\langle 0\right|\left[\partial_{\bot}^{2}\chi(x_{-})\right]C\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\chi(y_{-})\chi_{3}(z_{-})\left|B(k)\right\rangle=-m_{B}^{2}k_{+}~{}\left[\gamma_{5}N_{\bar{n}}\right]_{3}~{}\text{FT}[V_{1xx}^{B}],$
(75) $\left\langle 0\right|\chi(x_{-})C\hbox to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\left[\partial_{\bot}^{2}\chi(y_{-})\right]\chi_{3}(z_{-})\left|B(k)\right\rangle=-m_{B}^{2}k_{+}~{}\left[\gamma_{5}N_{\bar{n}}\right]_{3}~{}\text{FT}\left[V_{1yy}^{B}\right],$
(76) $\left\langle 0\right|\left[\partial_{\bot}^{2}\chi(x_{-})\right]C\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\gamma_{5}\chi(y_{-})\chi_{3}(z_{-})\left|B(k)\right\rangle=-m_{B}^{2}k_{+}~{}\left[N_{\bar{n}}\right]_{3}~{}\text{FT}\left[A_{1xx}^{B}\right],$
(77) $\left\langle 0\right|\chi(x_{-})C\hbox to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\gamma_{5}\left[\partial_{\bot}^{2}\chi(y_{-})\right]\chi_{3}(z_{-})\left|B(k)\right\rangle=-m_{B}^{2}k_{+}\left[N_{\bar{n}}\right]_{3}~{}\text{FT}\left[A_{1yy}^{B}\right],$
(78) $\displaystyle\left\langle
0\right|\left[\partial_{\bot}^{\alpha}\chi(x_{-})\right]C\hbox to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\left[\partial_{\bot}^{\beta}\chi(y_{-})\right]\chi_{3}(z_{-})\left|B(k)\right\rangle$
$\displaystyle=-\frac{m_{B}^{2}}{2}k_{+}~{}g_{\bot}^{\alpha\beta}\left[\gamma_{5}N_{\bar{n}}\right]_{3}~{}\text{FT}\left[V_{1xy}^{B}\right]$
$\displaystyle-\frac{m_{B}^{2}}{2}k_{+}\left[i\sigma_{\bot\bot}^{\alpha\beta}\gamma_{5}N_{\bar{n}}\right]_{3}\text{FT}\left[T_{xy}^{B}\right],$
(79) $\displaystyle\left\langle
0\right|\left[\partial_{\bot}^{\alpha}\chi(x_{-})\right]C\hbox to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\gamma_{5}\left[\partial_{\bot\beta}\chi(y_{-})\right]\chi_{3}(z_{-})\left|B(k)\right\rangle$
$\displaystyle=-\frac{m_{B}^{2}}{2}k_{+}g_{\bot}^{\alpha\beta}\left[N_{\bar{n}}\right]_{3}~{}\text{FT}\left[A_{1xy}^{B}\right]~{}$
$\displaystyle-\frac{m_{B}^{2}}{2}k_{+}\left[i\sigma_{\bot\bot}^{\alpha\beta}N_{\bar{n}}\right]_{3}\text{FT}\left[G_{xy}^{B}\right],$
(80) $\displaystyle\left\langle
0\right|\left[\partial_{\bot}^{\alpha}\partial_{\bot}^{\beta}\chi(x_{-})\right]C\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\gamma_{\bot}^{\sigma}\chi(y_{-})\chi_{3}(z_{-})\left|B(k)\right\rangle=\frac{m_{B}^{2}}{4}~{}k_{+}\left[\gamma_{\bot}^{\rho}\gamma_{5}N_{\bar{n}}\right]_{3}$
$\displaystyle\times\left(g_{\bot}^{\alpha\beta}g_{\bot}^{\sigma\rho}\text{FT}\left[T_{0xx}^{B}\right]+~{}G_{\bot}^{\alpha\beta\sigma\rho}~{}\text{FT}\left[T_{2xx}^{B}\right]\right),$
(81) $\displaystyle\left\langle
0\right|\chi(x_{-})C\bar{n}\gamma_{\bot}^{\sigma}\left[\partial_{\bot}^{\alpha}\partial_{\bot}^{\beta}\chi(y_{-})\right]\chi_{3}(z_{-})\left|B(k)\right\rangle=\frac{m_{B}^{2}}{4}~{}k_{+}\left[\gamma_{\bot}^{\rho}\gamma_{5}N_{\bar{n}}\right]_{3}~{}$
$\displaystyle\left(g_{\bot}^{\alpha\beta}g_{\bot}^{\sigma\rho}\text{FT}\left[T_{0yy}^{B}\right]+G_{\bot}^{\alpha\beta\sigma\rho}~{}\text{FT}\left[T_{2yy}^{B}\right]\right),$
(82) $\left\langle
0\right|\left[\partial_{\bot}^{\alpha}\chi(x_{-})\right]C\hbox to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\gamma_{\bot}^{\sigma}\left[\partial_{\bot}^{\beta}\chi(y_{-})\right]\chi_{3}(z_{-})\left|B(k)\right\rangle=\frac{m_{B}^{2}}{4}k_{+}\left[\gamma_{\bot}^{\rho}\gamma_{5}N_{\bar{n}}\right]_{3}~{}T_{xy}^{\alpha\beta\sigma\rho},$
(83)

with

$T_{xy}^{\alpha\beta\sigma\rho}=~{}g_{\bot}^{\sigma\rho}g_{\bot}^{\alpha\beta}~{}\text{FT}\left[T_{0xy}^{B}\right]+G_{\bot}^{\alpha\beta\sigma\rho}~{}\text{FT}\left[T_{2xy}^{B}\right]+\left(g_{\bot}^{\alpha\rho}g_{\bot}^{\beta\sigma}-g_{\bot}^{\alpha\sigma}g_{\bot}^{\beta\rho}\right)\text{FT}\left[\tilde{T}_{2xy}^{B}\right],$
(84)

where we use the short notaion for the following combination

$G_{\bot}^{\alpha\beta\sigma\rho}=g_{\bot}^{\alpha\sigma}g_{\bot}^{\beta\rho}+g_{\bot}^{\alpha\rho}g_{\bot}^{\beta\sigma}~{}-g_{\bot}^{\alpha\beta}g_{\bot}^{\sigma\rho}.$
(85)

This set of DAs can be splitted into two groups: chiral even and chiral-odd
DAs defined in Eqs.(75)-(80) and Eqs.(81)-(83), respectively. These groups can
be considered independently.

In order to establish relations of the auxiliary DAs with the basic ones, one
needs to relate the operator $\mathcal{O}_{123}$ defined in Eq.(6) with the
operators in Eq.(73). This can be done with the help of the QCD equations of
motion. Let us write the quark field as a sum

$q(x_{-})=\xi(x_{-})+\eta(x_{-}),$ (86)

where

$\xi(x_{-})=\frac{\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\hbox to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n}{4}q(x_{-}),\ \ \eta(x_{-})=\frac{\hbox
to0.0pt{\hbox to6.00235pt{\hfil/\hfil}\hss}n\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/}{4}q(x_{-}).$ (87)

The small field $\eta(x_{-})$ can be written as, see e.g. Ref.[43]

$\eta(x_{-})=-\frac{\hbox to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n}{2}\left(in\cdot D\right)^{-1}i{\hbox
to0.0pt{\hbox to8.55695pt{\hfil/\hfil}\hss}D}_{\bot}\xi(x_{-})=-\frac{\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n}{2}W(x_{-})\left(in\partial\right)^{-1}i\hbox
to0.0pt{\hbox to5.30904pt{\hfil/\hfil}\hss}\partial_{\bot}\chi(x_{-})+\ldots,$
(88)

where dots denote the terms with insertions of the transverse gluon components
$\left[W^{{\dagger}}(x_{-})i{\hbox to0.0pt{\hbox
to8.55695pt{\hfil/\hfil}\hss}D}_{\bot}W(x_{-})\right]$, which give the quark-
gluon operators and therefore can be omitted. Hence

$\ W^{{\dagger}}(x_{-})\eta(x_{-})\simeq-\frac{\hbox to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n}{2}\left(in\partial\right)^{-1}i\hbox
to0.0pt{\hbox to5.30904pt{\hfil/\hfil}\hss}\partial_{\bot}\chi(x_{-}).$ (89)

It is clear that each field $\eta$ gives the one transverse drivative,
therefore the twist-5 operators in the expansion of $\mathcal{O}_{123}$ in
Eq.(6) have two fields $\eta$

$\displaystyle\left[\mathcal{O}_{123}(x_{-},y_{-},z_{-})\right]_{\text{tw5}}$
$\displaystyle=W^{{\dagger}}(x_{-})\eta_{1}(x_{-})W^{{\dagger}}(y_{-})\eta_{2}(y_{-})W^{{\dagger}}(z_{-})\xi_{3}(z_{-})$
$\displaystyle+W^{{\dagger}}(x_{-})\eta_{1}(x_{-})W^{{\dagger}}(y_{-})\xi_{2}(y_{-})W^{{\dagger}}(z_{-})\eta_{3}(z_{-})$
$\displaystyle+W^{{\dagger}}(x_{-})\xi_{1}(x_{-})W^{{\dagger}}(y_{-})\eta_{2}(y_{-})W^{{\dagger}}(z_{-})\eta_{3}(z_{-}).$
(90)

The main idea is the following: we rewrite the terms in Eq.(90) using Eq.(89)
in terms of the basis operators (73) and compute the matrix element using
Eqs.(75)-(28). On the other hand, using the definitions (87) and one can
express the $\left[\mathcal{O}_{123}\right]_{\text{tw5}}$ in terms of
different projections of the three quark operator defined in Eq.(6) and use
Eq.(10) in order to get the matrix element. Comparing the two results one
finds the relations, which allows one to establish a connection of the
auxiliary DAs with the basis ones. It is convenient to perform such
calculations for the chiral-even and chiral-odd operators separately. Let us
first consider the chiral-even sector, which includes 8 DAs:

$V_{1xx}^{B},V_{1yy}^{B},\ V_{1xy}^{B},\ A_{1xx}^{B},\ A_{1yy}^{B},\
A_{1xy}^{B},T_{xy}^{B},G_{xy}^{B}.$ (91)

Let us consider as example the following operator

$\mathcal{O}_{3}^{V}(x_{-},y_{-},z_{-})=W^{{\dagger}}(x_{-})\eta(x_{-})C\hbox
to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/W^{{\dagger}}(y_{-})\eta(y_{-})W^{{\dagger}}(z_{-})\xi_{3}(z_{-}),$
(92)

where the Dirac indices $\\{12\\}$ are contracted with the matrix $C\hbox
to0.0pt{\hbox to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/$. Using Eq. (89) one
finds

$\displaystyle\mathcal{O}_{3}^{V}(x_{-},y_{-},z_{-})=\left[-\left(in\partial\right)^{-1}\frac{\hbox
to0.0pt{\hbox to6.00235pt{\hfil/\hfil}\hss}n}{2}i\hbox to0.0pt{\hbox
to5.30904pt{\hfil/\hfil}\hss}\partial_{\bot}\chi(x_{-})\right]C\hbox
to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\left[-\left(in\partial\right)^{-1}\frac{\hbox
to0.0pt{\hbox to6.00235pt{\hfil/\hfil}\hss}n}{2}i\hbox to0.0pt{\hbox
to5.30904pt{\hfil/\hfil}\hss}\partial_{\bot}\chi(y_{-})\right]\chi_{3}(z_{-})$
(93)
$\displaystyle=\left[\left(in\partial\right)^{-1}\partial_{\bot}^{\alpha}\chi(x_{-})\right]C\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\left[\left(in\partial\right)^{-1}\partial_{\bot}^{\alpha}\chi(y_{-})\right]\chi_{3}(z_{-})$
(94)
$\displaystyle+\frac{1}{2}i\varepsilon[\alpha\beta]\left[\left(in\partial\right)^{-1}\partial_{\bot}^{\alpha}\chi(x_{-})\right]C\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\gamma_{5}\left[\left(in\partial\right)^{-1}\partial_{\bot}^{\beta}\chi(y_{-})\right]\chi_{3}(z_{-}),$
(95)

where we use the short notation

$i\varepsilon[\alpha\beta]=\frac{1}{2}i\varepsilon_{\alpha\beta\sigma\rho}n^{\sigma}\bar{n}^{\rho}.$
(96)

The matrix element of this operator can be computed with the help of Eqs.(2.2)
and (80). One finds

$\left\langle
0\right|\mathcal{O}_{3}^{V}(x_{-},y_{-},z_{-})\left|B(k)\right\rangle=-\frac{m_{B}^{2}}{k_{+}}\left[\gamma_{5}N_{\bar{n}}\right]_{3}~{}\text{FT}\left[\frac{1}{x_{1}x_{2}}\left(V_{1xy}^{B}+G_{xy}^{B}\right)\right],$
(97)

where the factor $1/(x_{1}x_{2})$ appears due to the inverse derivatives
$\left(in\partial\right)^{-1}$ acting on the first and second fields. Remind,
that the inverse derivatives $\left(in\partial\right)^{-1}$ must be undestood
as nonlocal operator

$\left(in\partial\right)^{-1}\phi(x)\equiv\left(in\partial+i\varepsilon\right)^{-1}\phi(x)=-i\int_{-\infty}^{0}ds\
\phi(x+sn).$ (98)

On the other hand, using Eqs.(87) one finds

$\displaystyle\mathcal{O}_{3}^{V}(x_{-},y_{-},z_{-})$
$\displaystyle=W^{{\dagger}}(x_{-})q(x_{-})\left(\frac{n\bar{n}}{4}\right)^{\top}C\hbox
to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\frac{n\bar{n}}{4}W^{{\dagger}}(y_{-})q(y_{-})\
W^{{\dagger}}(z_{-})\left[\frac{\bar{n}n}{4}q(z_{-})\right]_{3}$ (99)
$\displaystyle=\left[\frac{\bar{n}n}{4}\right]_{33^{\prime}}\
W^{{\dagger}}(x_{-})q(x_{-})C\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/W^{{\dagger}}(y_{-})q(y_{-})W^{{\dagger}}(z_{-})q_{3^{\prime}}(z_{-}).$
(100)

Computing the matrix element of this operator with the help of Eq.(13) one
gets

$\left\langle
0\right|\mathcal{O}_{3}^{V}(x_{-},y_{-},z_{-})\left|B(k)\right\rangle=-\frac{m_{B}^{2}}{k_{+}}\left[\gamma_{5}N_{\bar{n}}\right]_{3}\text{FT}\left[V^{B}_{5}\right].$
(101)

Comparing Eqs.(97) and (101) one obtains

$\frac{1}{x_{1}x_{2}}\left(V_{1xy}^{B}+G_{xy}^{B}\right)=V_{5}^{B},$ (102)

where, for simplicity, the arguments are not shown: $F(x_{123})\equiv F$.

The similar consideration for the operators $\eta C\gamma_{\bot}\xi\eta_{3}$
and $\xi C\gamma_{\bot}\eta\eta_{3}$222The Wilson lines are not shown for
simplicity. yields

$V^{B}_{1xx}+V^{B}_{1xy}+A^{B}_{1xx}+A^{B}_{1xy}+T^{B}_{xy}+G^{B}_{xy}=-x_{1}x_{3}\left(V^{B}_{4}-A^{B}_{4}\right),$
(103)
$V^{B}_{1yy}+V^{B}_{1xy}-A^{B}_{1yy}-A^{B}_{1xy}-T^{B}_{xy}+G^{B}_{xy}=-x_{2}x_{3}\left(V^{B}_{4}+A^{B}_{4}\right),$
(104)

respectively. The analysis for the axial operators $\eta
C\bar{n}\gamma_{5}\eta\xi_{3}$ gives

$A^{B}_{1xy}+T^{B}_{xy}=x_{1}x_{2}A^{B}_{5}.$ (105)

Two more equations can be obtained using the QCD EOM, which can be rewritten
as

$i(\bar{n}\partial)\xi(x_{-})\simeq(in\partial)^{-1}\partial_{\bot}^{2}\xi(x_{-}),$
(106)

and describes the free collinear quark. Consider the following operator with
the total derivative

$\displaystyle i(\bar{n}\partial)\left[\chi(x_{-})C\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\chi(y_{-})\chi_{3}(z_{-})\right]\simeq\left[(in\partial)^{-1}\partial_{\bot}^{2}\chi(x_{-})\right]C\hbox
to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\chi(y_{-})\chi_{3}(z_{-})\phantom{i(\bar{n}\partial)\chi(x_{-})C\hbox
to0.0pt{\hbox to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\chi(y_{-})}$
$\displaystyle+\chi(x_{-})C\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\left[(in\partial)^{-1}\partial_{\bot}^{2}\chi(y_{-})\right]\chi_{3}(z_{-})+\chi(x_{-})C\hbox
to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\chi(y_{-})\left[(in\partial)^{-1}\partial_{\bot}^{2}\chi_{3}(z_{-})\right].$
(107)

Taking the matrix elements from the both sides one finds

$\frac{x_{1}+x_{3}}{x_{1}x_{3}}V^{B}_{1xx}+\frac{x_{2}+x_{3}}{x_{2}x_{3}}V^{B}_{1yy}+\frac{2}{x_{3}}V^{B}_{1xy}=V^{B}_{1}.$
(108)

The same consideration for the axial operator $\chi C\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\gamma_{5}\chi\chi_{3}$ gives

$\frac{x_{1}+x_{3}}{x_{1}x_{3}}A^{B}_{1xx}+\frac{x_{2}+x_{3}}{x_{2}x_{3}}A^{B}_{1yy}+\frac{2}{x_{3}}A^{B}_{1xy}=A^{B}_{1}.$
(109)

Therefore we obtain six linear algebraic relations. The investigation of other
possible 3-quark light-cone operators using the same technique does not
provide any new information. We did not find any method in order to get two
more pure algebraic equations. However one needs at least two more equations
in order to express the eight auxiliary DAs (91) in terms of the basis ones.
In order to derive these equations we introduce off light-cone correlation
function (CF). Performing the expansion of this CF around the light-cone we
derive two more equations, which allow us to solve the problem. However these
relations are already not algebraic. The derivation of the corresponding
relations is more involved and the details are given in Appendix B. The
obtained equations read

$2i(kx)\text{FT}\left[T^{B}_{xy}-G^{B}_{xy}\right]=\text{FT}\left[\left(x_{1}+x_{2}\right)\left(A^{B}_{5}-V^{B}_{5}\right)-\left(x_{2}+x_{3}\right)\left(V^{B}_{4}+A^{B}_{4}\right)-2\left(\mathcal{V}^{B}_{2}+\mathcal{A}^{B}_{2}\right)\right],$
(110)
$2i(ky)\text{FT}\left[T^{B}_{xy}+G^{B}_{xy}\right]=\text{FT}\left[(x_{1}+x_{2})\left(V^{B}_{5}+A^{B}_{5}\right)+\left(x_{1}+x_{3}\right)\left(V^{B}_{4}-A^{B}_{4}\right)+2\left(\mathcal{V}^{B}_{1}-\mathcal{A}^{B}_{1}\right)\right],$
(111)

where the factors $2(kx)=k_{+}x_{-}$, $2(ky)=k_{+}y_{-}$ and the coordinates
$x$ and $y$ enter in the definition of the FT in Eq.(14)333 We also assume
here that the argument of the third quark field is trivial: $z=0$ . Remind,
that the DAs $\mathcal{V}_{i}$ and $\mathcal{A}_{i}$ on the rhs denote the
auxiliary twist-4 DAs from Eqs.(63)-(64). For instance, for the nucleon,
taking the asymptotic contributions of twist-3 and twist-4 DAs one finds

$\left.\left(\mathcal{A}_{1}-\mathcal{V}_{1}\right)^{N}(x_{i})\right|_{as}=-20\
f_{N}~{}x_{1}^{2}x_{2}x_{3},\ \ \
\left.\left(\mathcal{V}_{2}+\mathcal{A}_{2}\right)^{N}(x_{i})\right|_{as}=-20\
f_{N}\ x_{1}x_{2}^{2}x_{3}.$ (112)

Solving two equations (110) and (111) one finds expressions for $T_{xy}$ and
$G_{xy}$. Then one can easily solve the system of the linear equations.

The properties of Eqs. (110) and (111) imply that in case of $x=0$ or $y=0$
the one lhs vanishes and therefore must vanish the corresponding rhs. This
gives two non-trivial conditions ($x_{3}=1-x_{1}-x_{2}$)

$\int_{0}^{\bar{x}_{2}}dx_{1}\left[\left(x_{1}+x_{2}\right)\left(A^{B}_{5}-V^{B}_{5}\right)-\left(x_{2}+x_{3}\right)\left(V^{B}_{4}+A^{B}_{4}\right)-2\left(\mathcal{V}^{B}_{2}+\mathcal{A}^{B}_{2}\right)\right](x_{1},x_{2},x_{3})=0.$
(113)
$\int_{0}^{\bar{x}_{1}}dx_{2}\left[(x_{1}+x_{2})\left(V^{B}_{5}+A^{B}_{5}\right)+\left(x_{1}+x_{3}\right)\left(V^{B}_{4}-A^{B}_{4}\right)+2\left(\mathcal{V}^{B}_{1}-\mathcal{A}^{B}_{1}\right)\right](x_{1},x_{2},x_{3})=0.$
(114)

These conditions allows one to perform the integrations by parts on rhs of
Eqs.(110) and (111), that reduce these equations to a simple algebraic form.
Notice, that the DAs on the  rhs of Eqs.(113) and (114) depends on the twist-3
and twist-4 moments only. We have found that Eqs.(113) and (114) are satisfied
if one uses for description of the DAs in Eqs.(113) and (114) the first
moments (or equivalently the asymptotic contributions) for twist-3 and twist-4
DAs only. We believe that this observation is closely related to the disregard
of the quark-gluon contributions and the genuine twist-5 moments Eq.(47). Such
simple recipe does not work for the description of the twist-5 contributions
associated with higher moments of the DAs and this makes such analysis to be
more complicated in comparison with the twist-4 DAs. In order to illustrate
this point we discuss a simple example, which is considered in the Appendix C.

Therefore in this work we only use the asymptotic expressions for the twist-3
and twist-4 DAs for the description of auxiliary twist-5 DAs. Let us also
mention that obtained description of the twist-5 DAs ensures their proper
normalisations, which are consistent with the Lorentz structure for the local
matrix elements. This is automatically provided by the algebraic equations
derived above.

It is convenient to solve the equations for twist-5 auxiliary DAs using the
explicit expressions for the asymptotic DAs. Then the solutions are given by
the simple analytical expressions. Assuming below $B\neq\Lambda$ we write the
final result as

$\displaystyle
T_{xy}^{B}=x_{1}x_{2}x_{3}(x_{1}-x_{2})(5f_{B}-\lambda_{1}^{B}),\ \
T_{xy}^{\Lambda}=-x_{1}x_{2}x_{3}\left\\{\frac{5f_{\Lambda}}{c_{\Lambda}^{+}}(3x_{3}-1)+\frac{\lambda_{1}^{\Lambda}}{c_{\Lambda}^{-}}\left(1+x_{3}\right)\right\\},$
(115) $\displaystyle
G_{xy}^{B}=x_{1}x_{2}x_{3}\left\\{5f_{B}(3x_{3}-1)+(1+x_{3})\lambda_{1}^{B}\right\\},\
G_{xy}^{\Lambda}=-~{}x_{1}x_{2}x_{3}(x_{1}-x_{2})\left\\{\frac{5f_{\Lambda}}{c_{\Lambda}^{+}}-\frac{\lambda_{1}^{\Lambda}}{c_{\Lambda}^{-}}\right\\},$
(116) $\displaystyle A_{1xy}^{B}$
$\displaystyle=x_{1}x_{2}x_{3}(x_{1}-x_{2})(5f_{B}-\lambda_{1}^{B}),$ (117)
$\displaystyle A_{1xy}^{\Lambda}$
$\displaystyle=-~{}x_{1}x_{2}x_{3}\left\\{~{}\frac{5f_{\Lambda}}{c_{\Lambda}^{+}}(3x_{3}-1-8x_{1}x_{2})+\frac{\lambda_{1}^{\Lambda}}{c_{-}}\left(1+x_{3}\right)\right\\},$
(118) $\displaystyle V_{1xy}^{B}$
$\displaystyle=x_{1}x_{2}x_{3}(~{}5f_{B}(3x_{3}-1-8x_{1}x_{2})+\lambda_{1}^{B}(1+x_{3})),$
(119) $\displaystyle V_{1xy}^{\Lambda}$
$\displaystyle=-~{}x_{1}x_{2}x_{3}(x_{1}-x_{2})\left\\{\frac{5f_{\Lambda}}{c_{\Lambda}^{+}}-\frac{\lambda_{1}^{\Lambda}}{c_{\Lambda}^{-}}\right\\},$
(120) $\displaystyle
V_{1xx}^{B}=2x_{1}^{2}x_{2}x_{3}\left(5f_{B}(5-4x_{1})-\lambda_{1}^{B}\right),\
\
V_{1xx}^{\Lambda}=2x_{1}^{2}x_{2}x_{3}\left(\frac{5f_{\Lambda}}{c_{\Lambda}^{+}}-\frac{\lambda_{1}^{\Lambda}}{c_{\Lambda}^{-}}\right),$
(121) $\displaystyle
V_{1yy}^{B}=2x_{1}x_{2}^{2}x_{3}\left(5f_{B}(5-4x_{2})-\lambda_{1}^{B}\right),\
\
V_{1yy}^{\Lambda}=-2x_{1}x_{2}^{2}x_{3}\left(\frac{5f_{\Lambda}}{c_{\Lambda}^{+}}-\frac{\lambda_{1}^{\Lambda}}{c_{\Lambda}^{-}}\right),$
(122) $\displaystyle
A_{1xx}^{B}=2x_{1}^{2}x_{2}x_{3}\left(-5f_{B}+\lambda_{1}^{B}\right),\ \
A_{1xx}^{\Lambda}=2x_{1}^{2}x_{2}x_{3}\left(\frac{5f_{\Lambda}}{c_{\Lambda}^{+}}(4x_{1}-5)+\frac{\lambda_{1}^{\Lambda}}{c_{\Lambda}^{-}}\right),$
(123) $\displaystyle
A_{1yy}^{B}=2x_{1}x_{2}^{2}x_{3}\left(5f_{B}-\lambda_{1}^{B}\right),\ \
A_{1yy}^{\Lambda}=2x_{1}x_{2}^{2}x_{3}\left(\frac{5f_{\Lambda}}{c_{\Lambda}^{+}}(4x_{2}-5)+\frac{\lambda_{1}^{\Lambda}}{c_{\Lambda}^{-}}\right).$
(124)

The similar discussion also holds for the chiral-odd set of the auxiliary DAs.
The details are given in Appendix B. The results read

$T_{0xx}^{B}=40f_{B}^{~{}\bot}~{}x_{1}^{2}x_{2}x_{3}\left(3-2x_{1}\right),\
T_{0xx}^{\Lambda}=8\frac{\lambda_{\bot}^{\Lambda}}{c_{\Lambda}^{-}}\vskip
6.0pt plus 2.0pt minus 2.0ptx_{1}^{2}x_{2}x_{3},$ (125)
$T_{0yy}^{B}=40f_{B}^{~{}\bot}~{}x_{1}x_{2}^{2}x_{3}\left(3-2x_{2}\right),\
T_{0yy}^{\Lambda}=-8\frac{\lambda_{\bot}^{\Lambda}}{c_{\Lambda}^{-}}x_{1}x_{2}^{2}x_{3},$
(126)
$T_{0xy}^{B}=-20f_{B}^{~{}\bot}~{}x_{1}x_{2}x_{3}\left(1-x_{3}+4x_{1}x_{2}\right),\
~{}T_{0xy}^{\Lambda}=-4\frac{\lambda_{\bot}^{\Lambda}}{c_{\Lambda}^{-}}x_{1}x_{2}x_{3}(x_{1}-x_{2}),$
(127) $\displaystyle~{}T_{2xy}^{B}$
$\displaystyle=-\lambda_{2}^{B}\frac{2}{3}x_{1}x_{2}x_{3}(1+x_{3}),~{}~{}T_{2xy}^{\Lambda}=\frac{\lambda_{2}^{\Lambda}}{c_{\Lambda}^{+}}\frac{2}{3}x_{1}x_{2}x_{3}(x_{1}-x_{2}),$
(128) $\displaystyle\tilde{T}_{2xy}^{B}(x_{i})$
$\displaystyle=-20f_{B}^{~{}\bot}~{}x_{1}x_{2}x_{3}(x_{1}-x_{2}),\ \
\tilde{T}_{2xy}^{\Lambda}=-4\frac{\lambda_{\bot}^{\Lambda}}{c_{\Lambda}^{-}}x_{1}x_{2}x_{3}(1-x_{3}),$
(129) $\displaystyle T_{2xx}^{B}$
$\displaystyle=\lambda_{2}^{B}\frac{4}{3}x_{1}^{2}x_{2}x_{3},~{}~{}T_{2xx}^{\Lambda}=-\frac{\lambda_{2}^{\Lambda}}{c_{\Lambda}^{+}}\frac{4}{3}x_{1}^{2}x_{2}x_{3},$
(130) $\displaystyle T_{2yy}^{B}$
$\displaystyle=\lambda_{2}^{B}\frac{4}{3}x_{1}x_{2}^{2}x_{3},\ \ \
T_{2yy}^{\Lambda}=\frac{\lambda_{2}^{\Lambda}}{c_{\Lambda}^{+}}\frac{4}{3}x_{1}x_{2}^{2}x_{3}.$
(131)

For the nucleon set one has to substitute in the formulas
$f_{B}^{\bot}\rightarrow f_{N}$. The coefficients $c_{\Lambda}^{\pm}$ are
defined in Eq.(29). It is useful to note that all auxiliary DAs are
proportional to the factor $x_{1}x_{2}x_{3}$ that simplifies an investigation
of analytic properties of the collinear integrals.

In addition to the identities, which follow from the Lorentz symmetry, one can
also derive symmetry relations that follow from the condition that a baryon
state has a certain isospin. For the auxiliary matrix elements this can be
done in the same way as for the basis one, see e.g. Ref.[21]. For instance,
for the nucleon ($I_{N}=1/2$) this gives 15 different relations for the
auxiliary DAs. We provide these relations in Appendix D. These relations can
also be used for the finding of the relation between the auxiliary and basic
DAs but different baryons has different isospin, which complicates the
analysis. At the same time, the relations that follow from the Lorentz
symmetry are universal and valid for all baryons. Nevertheless, the isospin
relations provide a powerful check of the obtained result. We have verified
that the obtained expressions for the auxiliary DAs in Eqs.(115)-(124) and
(125)-(131) satisfy all 15 isotopic relations given in Appendix D.

##  3 Calculation of hard kernels

The decay amplitude can be written as [19]

$M_{B}=\bar{N}(k)\left\\{A_{1}^{B}~{}\hbox to0.0pt{\hbox
to4.05904pt{\hfil/\hfil}\hss}\epsilon_{\psi}+A_{2}^{B}\left(\epsilon_{\psi}\right)_{\mu}(k^{\prime}+k)_{\nu}\frac{i\sigma^{\mu\nu}}{2m_{N}}\right\\}V(k^{\prime})~{},$
(132)

where $\epsilon_{\psi}$ is the charmonium polarisation vector, $\bar{N}(k)$
and $V(k^{\prime})$ denote baryon and antibaryon spinors, respectively. The
scalar amplitudes $A_{i}^{B}$ can be expanded within the effective field
theory framework in small heavy quark velocity $v^{2}$ (NRQCD) and in powers
$\lambda^{2}\sim\Lambda/m_{c}$ (collinear factorisation), where $\Lambda$ is
the typical hadronic scale. In this work we consider the leading-oder
contribution with respect to small velocity $v$ and compute the next-to-
leading power (NLP) correction to the amplitude $A_{1}^{B}$, which is provided
by the diagrams in Fig. 1$a)$.

Figure 1:  Hadronic $a)$ and electromagnetic $b)$ contributions to the decay
amplitudes.

The expansion of the amplitudes $A_{i}^{B}$ in powers of $\lambda$ is
associated with the power counting of the collinear matrix elements, which
describe the overlap with outgoing baryons. The structure of the collinear
expansion allows one to write the following expansions

$A_{1}^{B}=A_{1}^{\text{lp}}+A_{1}^{\text{nlp}}+\mathcal{O}(\lambda^{14}),\ \
A_{2}^{B}=A_{2}^{\text{lp}}+\mathcal{O}(\lambda^{14}).$ (133)

The leading-power (LP) contribution scales as
$A_{1}^{\text{lp}}\sim\mathcal{O}(\lambda^{8})$ because it is described by the
two light-cone matrix elements, each of which is of order
$\mathcal{O}(\lambda^{4})$. The NLP correction $A_{1}^{\text{nlp}}$ includes
collinear operators, which are suppressed by the overall factor $\lambda^{4}$
and therefore $A_{1}^{\text{nlo}}$ $\sim\mathcal{O}(\lambda^{12}).$ For the
second amplitude one finds
$A_{2}^{\text{lp}}\sim\mathcal{O}(\lambda^{12})$444The additional power
$\lambda^{2}$ appears due to the factor $1/m_{B}$ in the definition in
Eq.(132)..

The expansion in powers of $1/m_{c}$ mixes up power corrections from the
expansions of the scalar amplitudes $A_{i}^{B}$ and from the baryon spinors.
The calculation of the amplitude in the factorisation framework yields the
result in the following form

$M_{B}=\left\\{\delta A_{1}^{\text{lp}}+\delta
A_{1}^{\text{nlp}}\right\\}\bar{N}_{\bar{n}}\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\epsilon$\hfil}\hss}/_{\psi}^{\bot}V_{n}+\bar{N}_{\bar{n}}\boldsymbol{1}V_{n}~{}\delta
A_{2}^{\text{lp}},$ (134)

where the coefficients $\delta A_{i}$ scale as

$\delta A_{1}^{\text{lp}}\sim\mathcal{O}(\lambda^{8}),\ \delta
A_{1}^{\text{nlp}}\ \sim\delta
A_{2}^{\text{lp}}\sim\mathcal{O}(\lambda^{12}).$ (135)

The expressions for $\delta A_{i}$ are given by the convolution integrals of
the hard kernels with the different baryon DAs.

In order to perform the interpretation of $\delta A_{i}$ in terms of the
scalar amplitudes $A_{1,2}^{B}$ one needs to perform the expansion of the rhs
of Eq.(132). This gives

$\displaystyle M_{B}$ $\displaystyle\simeq\bar{N}_{\bar{n}}\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\epsilon$\hfil}\hss}/_{\psi}^{\bot}V_{n}\left(A_{1}^{\text{lp}}+A_{1}^{\text{nlp}}+A_{2}^{\text{lp}}-A_{1}^{\text{lp}}\frac{m_{B}^{2}}{k_{+}k_{-}^{\prime}}\right)$
$\displaystyle+\frac{1}{4m_{B}}\left\\{\left(\epsilon_{\psi}\cdot
n\right)k_{-}^{\prime}~{}-k_{+}\left(\epsilon_{\psi}\cdot\bar{n}\right)\right\\}~{}\bar{N}_{\bar{n}}\boldsymbol{1}V_{n}~{}\left(A_{2}^{\text{lp}}+\frac{4m_{B}^{2}}{k_{+}k_{-}^{\prime}}A_{1}^{\text{lp}}\right),$
(136)

where we use that $k_{+}\sim k_{-}^{\prime}\sim m_{c}$. Comparing Eqs.(134)
and (136) one finds

$A_{1}^{\text{lp}}=\delta A_{1}^{\text{lp}},\ \ A_{2}^{\text{lp}}=\ \delta
A_{2}^{\text{lp}}-\frac{4m_{B}^{2}}{k_{+}k_{-}^{\prime}}A_{1}^{\text{lp}}$
(137) $\ A_{1}^{\text{nlp}}=\delta A_{1}^{\text{nlp}}-\ \delta
A_{2}^{\text{lp}}+\frac{5m_{B}^{2}}{k_{+}k_{-}^{\prime}}\delta
A_{1}^{\text{lp}}$ (138)

The expressions for observables have a compact analytical form if we rewrite
them in terms of the linear combinations $\mathcal{G}_{M}$ and
$\mathcal{G}_{E}$, which read

$\mathcal{G}_{M}^{B}=A_{1}^{B}+A_{2}^{B}\simeq\delta
A_{1}^{\text{lp}}\left(1+\frac{m_{B}^{2}}{M_{\psi}^{2}}\right)+\delta
A_{1}^{\text{nlp}},\ \ \
\mathcal{G}_{E}^{B}=A_{1}^{B}+\frac{M_{\psi}^{2}}{4m_{B}^{2}}A_{2}^{B}\simeq\frac{M_{\psi}^{2}}{4m_{B}^{2}}\delta
A_{2}^{\text{lp}}.$ (139)

Notice that the total power correction to $A_{1}^{B}$ also includes the simple
kinematical term $\delta A_{1}^{\text{lo}}\times m_{B}^{2}/M_{\psi}^{2}$.

The LP term $\delta A_{1}^{\text{lo}}$ is calculated long time ago, see e.g.
Refs.[1, 5]. The contribution $\delta A_{2}^{\text{lo}}$ is obtained in
Ref.[19]. The correction $\delta A_{1}^{\text{nlo}}$ will be considered for
the first time. In this calculation we only consider the contributions of
light-cone 3-quark operators and do not include the various quark-gluon matrix
elements.

We use the same technique as in Ref.[19]. Here we just provide a schematic
description and discuss some important points. The different contributions in
$\delta A_{1}^{\text{nlo}}$ can be divided into two groups, which are given by
the product of the two twist-4 operators
$\left[\chi_{1}\chi_{2}\chi_{3}\right]_{\text{tw4}}\otimes\left[\bar{\chi}_{1^{\prime}}\bar{\chi}_{2^{\prime}}\bar{\chi}_{3^{\prime}}\right]_{\text{tw4}}$
or by the product twist-3 and twist-5 operators
$\left[\chi_{1}\chi_{2}\chi_{3}\right]_{\text{tw3}}\otimes\left[\bar{\chi}_{1^{\prime}}\bar{\chi}_{2^{\prime}}\bar{\chi}_{3^{\prime}}\right]_{\text{tw5}}$.
Remind, $\chi_{1}$ and $\bar{\chi}_{1^{\prime}}$ denote the collinear quark
fields, see Eq.(52). Let us also remind that $\chi_{1}\equiv\chi_{\alpha_{1}}$
and $\bar{\chi}_{1^{\prime}}\equiv\bar{\chi}_{\alpha_{1}^{\prime}}$. Each
group can be further divided into two subgroups with respect to Dirac
projections: $[\chi\Gamma_{even}\chi]\chi_{3}$ and
$[\chi\Gamma_{odd}\chi]\chi_{3}$, where $\Gamma_{even/odd}$ denote the chiral-
even or chiral-odd Dirac matrices. The explicit definitions of the various
operators and their matrix elements are discussed in Sec.2.2.

The contribution of the sum of diagrams in Fig. 1$\,a)$ to $\delta
A_{1}^{\text{nlo}}$ can be written as

$\delta A_{1}^{\text{nlo}}=f_{\psi}\int Dx_{i}\int Dy_{i}\
\left[P(x_{i})\right]_{123}\left[\bar{P}(y_{i})\right]_{1^{\prime}2^{\prime}3^{\prime}}\frac{1}{4}\text{Tr}\left[(1-\hbox
to0.0pt{\hbox to6.22453pt{\hfil/\hfil}\hss}\omega)\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\epsilon$\hfil}\hss}/_{\psi}\
\sum_{i}D^{i}_{\left(123\right)\left(1^{\prime}2^{\prime}3^{\prime}\right)}\right],$
(140)

where $D^{i}_{\left(123\right)\left(1^{\prime}2^{\prime}3^{\prime}\right)}$
denotes the analytical expression of the diagram $i$ in momentum space. It can
be written as

$\ D_{\left(123\right)\left(1^{\prime}2^{\prime}3^{\prime}\right)}^{i}=\
D_{Q}^{(i)\alpha\beta\gamma}(k_{i},k_{i}^{\prime},k,k^{\prime})\
(ig_{s})^{3}\left[\gamma^{\alpha}\right]_{11^{\prime}}[\gamma^{\beta}]_{22^{\prime}}\left[\gamma^{\gamma}\right]_{33^{\prime}}.$
(141)

Here the matrices $\gamma^{\alpha}\otimes\gamma^{\beta}\otimes\gamma^{\gamma}$
corresponds to the light quark-antiquark vertices,
$D_{Q}^{(i)\alpha\beta\gamma}$ denotes the remnant part of the diagram, the
trace in Eq.(140) is taken with respect of Dirac indices of the heavy quark
fermion line. For simplicity, the colour structure of diagrams is not shown.
The coupling $f_{\psi}$ describes the matrix element in NRQCD and can be
related with the radial wave function at the origin

$f_{\psi}=\sqrt{2M_{\psi}}\sqrt{\frac{3}{2\pi}}~{}R_{10}(0).$ (142)

The definition of the NRQCD matrix element is the same as in Ref.[19]. The
expression $(1-\hbox to0.0pt{\hbox to6.22453pt{\hfil/\hfil}\hss}\omega)\hbox
to0.0pt{\hbox to5.00002pt{\hfil$\varepsilon$\hfil}\hss}/_{\psi}$ in the trace
in Eq.(140) represent the Dirac part of the projector on the NRQCD matrix
element, $\omega=(1,0,0,0)$ and $\varepsilon_{\psi}^{\mu}$ denote the velocity
and the polarisation vector of charmonium, respectively.

The operators $P(x_{i})$ and $\bar{P}(y_{i})$ describe the projections on the
matrix elements of the light-cone 3-quark operators for baryon and antibaryon,
respectively. For the twist-3 case these operators depend on the nucleon
spinors, Dirac matrices and twist-3 DAs and can be obtained just from the
definition in Eq.(11). For instance, for the matrix element $\left\langle
0\right|\mathcal{O}_{123}\left|B(k)\right\rangle_{\text{tw3}}$ one easily
obtains

$\displaystyle P_{\text{tw3}}(x_{i})$
$\displaystyle=~{}\frac{1}{8}k_{+}\left[\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/C\right]_{12}~{}\left[\gamma_{5}N_{\bar{n}}\right]_{3}V_{1}^{B}(x_{i})+\frac{1}{8}k_{+}\left[\hbox
to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\gamma_{5}C\right]_{12}~{}\left[N_{\bar{n}}\right]_{3}\
A_{1}^{B}(x_{i})$ (143) $\displaystyle+\frac{1}{8}k_{+}\left[\hbox
to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\gamma_{\bot}^{\alpha}C\right]_{12}~{}\left[\gamma_{\bot}^{\alpha}\gamma_{5}N_{\bar{n}}\right]_{3}T_{1}^{B}(x_{i}).$
(144)

For the higher twist matrix elements these projections also include the
derivatives with respect to quark and antiquark momenta $k_{i}$ and
$k_{i}^{\prime}$, which are assumed to be off-shell $k_{i}^{2}\neq 0$,
$k_{i}^{\prime}\neq 0$. Only after the differentiations one can neglect the
small momentum components and set $k_{i}\simeq x_{i}k_{+}\bar{n}/2$ and
$k_{i}^{\prime}\simeq y_{i}k_{-}^{\prime}n/2$. The derivation of the twist-4
projections have been considered in Ref.[19] and here we follow the same
technique. The analytical results for the different projectors are given in
Appendix E. Substituting the various projectors in Eq.(140) one obtains the
expression, which allows one to calculate the hard kernels. The analytical
calculations were performed using the package FeynCalc [37, 38].

The obtained result can be written as

$\delta
A_{1}^{\text{nlo}}=A_{0}^{B}\frac{m_{B}^{2}}{m_{c}^{2}}\left\\{J_{44}^{B}+J_{35}^{B}\right\\},$
(145)

where the normalisation coefficient $A_{0}$ read

$A_{0}^{B}=\frac{f_{\psi}}{m_{c}^{2}}\frac{f_{B}^{2}}{m_{c}^{4}}(\pi\alpha_{s})^{3}\frac{10}{81}.$
(146)

The dimensionless integrals $J_{ij}^{B}$ have subscripts “$44$” and “$35$” ,
which are related with the structure of the collinear matrix elements:
twist-4$\times$twist-4 and twist-3$\times$twist-5. Each integral in Eq.(145)
is given by the sum of simpler integrals $J[X,Y]$ where $X$ and $Y$ denote the
DAs, which are included in the integrand

$J[X,Y]=\frac{1}{f_{B}^{2}}\int Dx_{i}\ X(x_{i})\int Dy_{i}Y(y_{i})\
K_{XY}(x_{i},y_{i}),$ (147)

where $K_{XY}(x_{i},y_{i})$ is the universal hard kernel. The calculation of
the traces in Eq.(140) gives the sum of the collinear integrals, which can be
simplified using the symmetry properties of the DAs and the hard kernels with
respect to interchange $x_{1}\leftrightarrow x_{2}$ and $x_{i}\leftrightarrow
y_{i}$. Below we provide the analytical expressions for the integrals in
Eq.(145), where such simplifications are taken into account. Let us write each
term in Eq.(145) as

$J_{ij}^{B}=\left.J_{ij}^{B}\right|_{\text{even}}+\left.J_{ij}^{B}\right|_{\text{odd}}.$
(148)

The subscript “even/odd” corresponds to the chiral-even and chiral-odd subsets
of the DAs. Then

$\left.J_{44}^{B}\right|_{\text{even}}=J[\mathcal{V}^{B}_{1},\mathcal{V}^{B}_{1}]+J[\mathcal{V}^{B}_{2},\mathcal{V}^{B}_{1}]+J[\mathcal{A}^{B}_{1},\mathcal{A}^{B}_{1}]+J[\mathcal{A}^{B}_{2},\mathcal{A}^{B}_{1}]+J[\mathcal{A}^{B}_{1},\mathcal{V}^{B}_{1}]+J[\mathcal{A}^{B}_{2},\mathcal{V}^{B}_{1}].$
(149)
$\left.J_{44}^{B}\right|_{\text{odd}}=J[\mathcal{T}^{B}_{21}+\mathcal{T}^{B}_{41},\mathcal{T}^{B}_{21}+\mathcal{T}^{B}_{41}]+J[\mathcal{T}^{B}_{21}-\mathcal{T}^{B}_{41},\mathcal{T}^{B}_{21}-\mathcal{T}^{B}_{41}]+J[\mathcal{T}^{B}_{22}+\mathcal{T}^{B}_{42},\mathcal{T}^{B}_{21}+\mathcal{T}^{B}_{41}].$
(150)
$\displaystyle\left.J_{35}^{B}\right|_{\text{even}}=J[V^{B}_{5},V^{B}_{1}]+J[V^{B}_{5},A^{B}_{1}]+J[A^{B}_{5},A^{B}_{1}]+J[A^{B}_{5},V^{B}_{1}]$
(151)
$\displaystyle\phantom{J[V^{B}_{5},V^{B}_{1}]}+J_{Z}[V^{B}_{Z},V^{B}_{1}]+J_{Z}[A^{B}_{Z},V^{B}_{1}]+J_{Z}[A^{B}_{Z},A^{B}_{1}]+J_{Z}[V^{B}_{Z},A^{B}_{1}].$
(152) $\left.J_{35}\right|_{\text{odd}}=J_{Z}[T^{B}_{Z},T^{B}_{1}].$ (153)

Each integral on the rhs of Eqs.(149)-(153) is IR-finite, i.e. does not have
endpoint singularities. The integrals $J_{Z}[X,Y]$ are given by the sum of
singular integrals and read

$\displaystyle J_{Z}[V^{B}_{Z},V^{B}_{1}]$ $\displaystyle=$ $\displaystyle
J[V^{B}_{1xx},V^{B}_{1}]+J[V^{B}_{1yy},V^{B}_{1}]+J[V^{B}_{1xy},V^{B}_{1}]+(1-\delta_{B\Lambda})J[V^{B}_{1},V^{B}_{1}],$
(154) $\displaystyle J_{Z}[V^{B}_{Z},A^{B}_{1}]$ $\displaystyle=$
$\displaystyle
J[V^{B}_{1xx},A^{B}_{1}]+J[V^{B}_{1yy},A^{B}_{1}]+J[V^{B}_{1xy},A^{B}_{1}]+(1-\delta_{B\Lambda})J[V^{B}_{1},A^{B}_{1}],$
(155) $\displaystyle J_{Z}[A^{B}_{Z},V^{B}_{1}]$ $\displaystyle=$
$\displaystyle
J[A^{B}_{1xx},V^{B}_{1}]+J[A^{B}_{1yy},V^{B}_{1}]+J[A^{B}_{1xy},V^{B}_{1}]+\delta_{B\Lambda}J[A^{B}_{1},V^{B}_{1}],$
(156) $\displaystyle J_{Z}[A_{Z},A_{1}]$ $\displaystyle=$ $\displaystyle
J[A^{B}_{1xx},A^{B}_{1}]+J[A^{B}_{1yy},V^{B}_{1}]+J[A^{B}_{1xy},V^{B}_{1}]+\delta_{B\Lambda}J[A^{B}_{1},A^{B}_{1}],$
(157) $\displaystyle J_{Z}[T^{B}_{Z},T^{B}_{1}]$ $\displaystyle=$
$\displaystyle
J[T^{B}_{0xx},T^{B}_{1}]+J[T^{B}_{0yy},T^{B}_{1}]+J[T^{B}_{0xy},T^{B}_{1}]+J[\tilde{T}^{B}_{2xy},T^{B}_{1}]$
(158)
$\displaystyle\phantom{J[T^{B}_{0xx},T^{B}_{1}]+J[T^{B}_{0yy},T^{B}_{1}]+J[T^{B}_{0xy},T^{B}_{1}}+(1-\delta_{B\Lambda})J[T^{B}_{1},T^{B}_{1}].$

where $\delta_{B\Lambda}=1$ if $B=\Lambda$ and $0$ otherwise. Notice that
expressions in Eqs.(154) - (158) include the integrals, which have twist-3 DAs
only. These are contributions, which appears from twist-5 matrix elements with
the total derivatives, see discussion in Appendix E. The Kronecker symbol
$\delta_{B\Lambda}$ appears in Eqs.(154) - (157) because we only use
asymptotic twist-3 and twist-4 DAs for the evaluation of the twist-5 DAs. The
integrals with the twist-3 DAs in Eqs.(154) - (158) are singular and they are
included in the given expressions because they are required for the
cancellation of the endpoint singularities in the integrals $J_{Z}[X,Y]$.

The analytical expressions for the hard kernels $K_{XY}(x_{i},y_{i})$, as
defined in Eq.(147), are given in the Appendix F. It has been analytically
verified that these kernels correctly describe the relations between the
amplitudes in the $SU(3)$ limit. The corresponding convolution integrals can
be relatively easily calculated. Substituting the explicit expressions for DAs
with arbitrary moments and summing the singular integrands gives analytical
expressions with the well defined four-dimensional integrals. These integrals
can be easily computed numerically using the standard mathematica packages for
numerical integrations.

##  4 Phenomenology

A systematic numerical analysis of existing data is still very challenging due
to various theoretical and experimental uncertainties. Therefore, the main
goal of the following discussion is to understand how well the existing data
can be described using plausible assumptions about various nonperturbative
parameters.

In this analysis we take into account the electromagnetic amplitudes, which
are shown by the diagram in Fig.1$\,(b)$. With these contributions the decay
amplitudes read

$M[J/\psi\rightarrow B\bar{B}]=\bar{N}\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\epsilon$\hfil}\hss}/_{\psi}V~{}\left(A_{1}^{B}+F_{1}^{B}\right)+\frac{1}{2m_{N}}\epsilon_{\psi}^{\mu}(k+k^{\prime})^{\nu}~{}\bar{N}i\sigma^{\mu\nu}V~{}\left(A_{2}^{B}+F_{2}^{B}\right),$
(159)

where $F_{1,2}^{B}$ are the time-like e.m. form factors. Expression for the
width reads

$\Gamma[J/\psi\rightarrow
B\bar{B}]=\frac{M_{\psi}\beta}{12\pi}|\mathcal{\tilde{G}}_{M}^{B}|^{2}\left(1+\frac{2m_{B}^{2}}{M_{\psi}^{2}}\gamma_{B}^{2}\right),\
\ $ (160)

where $\beta=\sqrt{1-4m_{B}^{2}/M_{\psi}^{2}}$,

$~{}\mathcal{\tilde{G}}_{M}^{B}=\mathcal{G}_{M}^{B}+\zeta G_{M}^{B},~{}~{}\
\mathcal{\tilde{G}}_{E}^{B}=\mathcal{G}_{E}^{B}+\zeta G_{E}^{B},\ \
\gamma_{B}\equiv\left|\mathcal{\tilde{G}}_{E}^{B}\right|/|\mathcal{\tilde{G}}_{M}^{B}|.$
(161)

The amplitudes $\mathcal{G}_{M}$ and $\mathcal{G}_{E}$ are defined in
Eq.(139), the time-like Sachs form factors $G^{B}_{E,M}$ read

$G_{M}^{B}=F_{1}^{B}+F_{2}^{B},\ \
G_{E}^{B}=F_{1}^{B}+\frac{M_{\psi}^{2}}{4m_{B}^{2}}F_{2}^{B}.$ (162)

The coupling $\zeta$ in Eq.(161) is defined as

$\zeta=4\pi\alpha_{em}\ e_{c}\frac{f_{\psi}}{M_{\psi}^{2}}.$ (163)

Therefore, the observables depend on the time-like form factors of baryons,
which are still not well known. For the given value of the
$q^{2}=M_{\psi}^{2}$ their values are still depend on the large non-
perturbative corrections and can not be sufficiently well estimated within the
QCD factorisation framework. This introduces additional uncertainty into the
phenomenological analysis.

The hadronic amplitudes calculated to this accuracy within the factorisation
framework are real. The e.m. FFs have real and imaginary parts, therefore to
this accuracy one finds

$\left|~{}\mathcal{\tilde{G}}_{M}^{B}\right|^{2}/\left|\mathcal{G}_{M}^{B}\right|^{2}\simeq
1+2\zeta\left|G_{M}^{B}\right|\cos\varphi_{M}/\left|\mathcal{G}_{M}^{B}\right|+\zeta^{2}\left|G_{M}^{B}\right|^{2}/\left|\mathcal{G}_{M}^{B}\right|^{2},$
(164)
$\left|~{}\mathcal{\tilde{G}}_{E}^{B}\right|^{2}/\left|\mathcal{G}_{E}^{B}\right|^{2}\simeq
1+2\zeta\left|G_{E}^{B}\right|\cos\varphi_{E}/\left|\mathcal{G}_{E}^{B}\right|+\zeta^{2}\left|G_{E}^{B}\right|^{2}/\left|\mathcal{G}_{E}^{B}\right|^{2},$
(165)

where $\varphi_{M,E}$ are phases of the e.m. form factors. Let us rewrite the
ratio as

$\gamma_{B}=\gamma_{B}^{g}\mathcal{R}_{em}^{B},$ (166)

with

$\gamma_{B}^{g}=\frac{\left|\mathcal{G}_{E}^{B}\right|}{\left|\mathcal{G}_{M}^{B}\right|},\
\ \
\mathcal{R}_{em}^{B}=\sqrt{\frac{1+2\zeta\left|G_{E}^{B}\right|/\left|\mathcal{G}_{E}^{B}\right|\cos\varphi_{E}+\zeta^{2}\left|G_{E}^{B}\right|^{2}/\left|\mathcal{G}_{E}^{B}\right|^{2}}{1+2\zeta\left|G_{M}^{B}\right|/\left|\mathcal{G}_{M}^{B}\right|\cos\varphi_{M}+\zeta^{2}\left|G_{M}^{B}\right|^{2}/\left|\mathcal{G}_{M}^{B}\right|^{2}}},$
(167)

where the factor $\mathcal{R}_{em}^{B}$ defines how the measured value
$\gamma_{B}$ depends on the interference with e.m. amplitude.

At present, the time-like baryon FFs are little known. The best experimental
data for the proton can be found in Ref.[31]. For other baryons one can only
find the data for the effective FF obtained from the analysis of the e.m.
cross section under assumption $G_{\text{eff}}^{B}=|G_{M}^{B}|=|G_{E}^{B}|$,
see Refs. [29, 30, 32, 33, 34, 35, 36]. Therefore an estimation of
$\mathcal{R}_{em}^{B}$ can not be done without assumptions. In the following
we assume that $0\ll\cos\varphi_{M,E}<1$555The small value of
$\cos\varphi_{M,E}\simeq 0$ implies that $G_{M,E}$ dominates by imaginary
part. We assume that such a scenario is unlikely.. The value of $|G_{M}^{p}|$
we take from data in Ref.[31], the maximal value $|G_{E}^{p}|$ is taken to be
$|G_{E}^{p}|\leq 3.29\times 10^{-2}$ as it follows from analysis of data in
Ref.[31]. In other cases we assume that $|G_{M}^{B}|\approx$
$G_{\text{eff}}^{B}$. The value of $\cos\varphi_{M}$ can be estimated from the
data for the width (160), where $\gamma_{B}$ is fixed by the experimental
value. An estimate of $G_{E}^{B}$ is the most challenging point. We assume
that $\left|G_{E}^{B}\right|\leq$1.5 $G_{\text{eff}}^{B}$, where the numerical
coefficient $1.5$ is our conservative estimate. Numerical value of
$\zeta^{2}\left|G_{E}\right|^{2}/\left|\mathcal{G}_{E}\right|^{2}$ is quite
small and can be neglected, which yields

$\mathcal{R}_{em}^{B}\approx\sqrt{\frac{1+2\zeta\left|G_{E}^{B}\right|/\left|\mathcal{G}_{E}^{B}\right|\cos\varphi_{E}}{1+2\zeta\left|G_{M}^{B}\right|/\left|\mathcal{G}_{M}^{B}\right|\cos\varphi_{M}+\zeta^{2}\left|G_{M}^{B}\right|^{2}/\left|\mathcal{G}_{M}^{B}\right|^{2}}}.$
(168)

Therefore, assuming the estimate

$0<\left|G_{E}^{B}\right|\cos\varphi_{E}\leq\ 1.5G_{\text{eff}}^{B},$ (169)

we get an interval, which is considered as an estimate of the possible
numerical effect provided by $\mathcal{R}_{em}^{B}$.

Let us now discuss the non-perturbative input, which is used for calculation
of the hadronic amplitudes. Recall that we consider only three-quark DAs. It
is convenient to describe the non-perturbative input in terms of basic DAs as
discussed in the Sec.2.1. Analytical expressions for the corresponding models
are described in Appendix A. The values of DAs parameters are summarised in
Table 1.

Table 1: The parameters, which define the twist-3 and twist-4 models of the
baryon DAs (upper and bottom tables, respectively). All values are given in
units of $10^{-2}$ at the scale $\mu^{2}=1.5~{}$GeV2.

$B$ |  $f_{B},$ GeV2 | $\phi_{10}$ | $\phi_{11}$ | $\phi_{20}$ | $\phi_{21}$ | $\phi_{22}$ |  $f_{\bot}^{B},$ GeV2 | $\pi_{10}^{B}$ | $\pi_{11}^{B}$
---|---|---|---|---|---|---|---|---|---
$N$ | $0.54$ | $5.1$ | $5.15$ | $7.82$ | $-2.84$ | $17.95$ | $-$ | $-$ | $-$
$\Lambda$ | $0.60$ | $10.9$ | $4.42$ | $0$ | $0$ | $0$ | $-$ | $3.83$ | $-$
$\Sigma$ | $0.51$ | $1.90$ | $4.18$ | $0$ | $0$ | $0$ | $0.53$ | $-$ | $-1.90$
$\Xi$ | $0.62$ | $6.10$ | $-0.25$ | $0$ | $0$ | $0$ | $0.61$ | $-$ | $7.10$

$B$ |  $\lambda_{1}^{B}$, GeV2 | $\eta_{10}^{B}$ | $\eta_{11}^{B}$ |  $\lambda_{\bot}^{B}$, GeV2 | $\zeta_{10}^{B}$ | $\zeta_{11}^{B}$ |  $\lambda_{2}$, GeV2 | $\xi_{10}$
---|---|---|---|---|---|---|---|---
$N$ | $-2.81$ | $-4.00$ | $16.34$ | $-$ | $-$ | $-$ | $5.29$ | -27.24
$\Lambda$ | $-3.93$ | $-4.00$ | $18.86$ | $-4.87$ | $-4.00$ | $-$ | $9.21$ | $-$
$\Sigma$ | $-4.29$ | $-4.62$ | $21.99$ | $-$ | $-$ | $20.95$ | $7.93$ | $-37.56$
$\Xi$ | $-4.64$ | $-2.05$ | $13.62$ | $-$ | $-$ | $15.71$ | $9.27$ | $-27.24$

The parameters of nucleon chiral-even DAs correspond to the model ABO-I from
Ref.[25]. This set of parameters have beed obtained from the fit of the space-
like electromagnetic FF data using the light-cone sum rules. The moments
$\lambda_{2}$ and $\xi_{10}$ are evaluated in Ref.[21] using the QCD sum rule
approach. Some twist-3 moments for different baryons are also obtained using
QCD sum rules in Ref.[5]. Recently, the different twist-3 and twist-4 moments
have also been computed on the lattice, see Ref.[28]. We use these results in
order to fix parameter values in DA models. In particular, the values for the
moments $\phi_{10}^{B},\phi_{11}^{B}$, $\lambda_{1}^{B}$,$\lambda_{\bot}^{B}$,
$\lambda_{2}^{B},\pi_{10}^{B}$ and $\pi_{11}^{B}$ ($B\neq N$) are fixed
according to data in Ref.[28].

For the nucleon twist-3 normalisation coupling we take the value $f_{N}(1$
GeV$)=5.5\times 10^{-3}\ $GeV2, which is in agreement with the sum rule
estimate $f_{N}(1$ GeV$)=(5.3\pm 0.5)\times 10^{-3}$GeV2 obtained in Ref.[5].
The other twist-3 normalisation couplings $f_{B}$ are fixed from the fit of
data. The values obtained in Ref.[28] were used as the first estimate and then
changed if this improves the description of the data. Therefore the resulting
value $f_{\Lambda}$ in Tab.1 is larger than the estimates in Refs.[28, 5],
which read $0.50$ and $0.48$, respectively.666For a convenience, these values
are given in the same notation as in the Tab.1 .  We also find a relatively
small difference for values of $f_{\Sigma}$ ($0.54$ [28] and $0.48$ [5]) and
$f_{\bot}^{\Xi}$ ($0.64$ [28] and $0.50$ [5]). The values for
$f_{\bot}^{\Sigma}$ and $f_{\Xi}$ are the same as in Ref.[28].

The values of the twist-4 moments $\eta_{10}^{B},\eta_{11}^{B}$,
$\zeta_{10}^{B},\zeta_{11}^{B}$ and $\xi_{10}^{B}$ are fixed from the data,
assuming that their values are not significantly different from the nucleon
ones because of approximate $SU(3)$ symmetry.

In our calculation we fix relatively low renormalisation scale
$\mu^{2}=1.5\,$GeV2. This gives for the strong coupling
$\alpha_{s}(\mu^{2})=0.35$. The twist-4 DAs include terms with the quark
masses, see Eqs.(63)-(68). For these terms we set $m_{u}\approx m_{d}\approx
0$, $m_{s}(1.5$GeV${}^{2})=105\,$MeV. The value of $J/\psi$ total width
$\Gamma_{tot}=93\,$keV is taken from Ref.[39] .

For the value $R_{10}(0)$ in Eq.(142) we use the estimate obtained for the
Buchmuller-Tye potential [40]

$\left|R_{10}(0)\right|^{2}\simeq 0.81\text{GeV}^{3},$ (170)

which implies for the value of charm quark mass $m_{c}=1.48\,$GeV.

One of the elements of the current analysis is the estimation of the power-law
correction to $\mathcal{G}_{M}^{B}$, see Eq.(139). In Table 2 we show the
numerical results obtained for these corrections. The total effect of the
power corrections is quite moderate due to the partial cancellation in the
total sum, see Eq.(139). Therefore the total numerical effect from the next-
to-leading power correction varies from $1.7\%$ for $\Xi$ up to $16\%$ for
nucleon.

Table 2:  Numerical value of the power corrections for different baryons, see Eq.(139). | $N$ | $\Lambda$ | $\Sigma$ | $\Xi$
---|---|---|---|---
$m_{B}^{2}/M_{\psi}^{2}\times 100\%$ | $9.2$ | $13.0$ | $14.8$ | $18.1$
$\delta A_{1}^{\text{nlo}}/\delta A_{1}^{\text{lo}}\times 100\%$ | $-25.4$ | $-16.4$ | $-23.1$ | $-19.8$

The numerical results for branching fractions and ratios $\gamma_{B}$ are
presented in Table 3. The obtained values of $\gamma_{B}$ include the error
bars, which describe the estimate of the uncertainty occurring from to the
interference with the electromagnetic FF’s as described by the coefficient
$R_{em}$.

Table 3:  The numerical results for branching ratios and for the ratio
$\gamma_{B}$ in comparison with experimental data. The values of the
electromagnetic FF’s $|G^{B}_{M}|$ and the branching rations are given in
units of $10^{-2}$ and $10^{-3}$, respectively. For more explanations see
discussion in the text.
xxx | $|G^{B}_{M}|$ | $\cos\phi_{M}$ | $Br[J/\psi\to B\bar{B}]$ |  $Br$[exp] | $\gamma_{B}$ |  $\gamma_{B}$[exp]
---|---|---|---|---|---|---
$p$ | $3.47$ | $0.70$ | $2.10$ | $2.12(3)$ | $0.82^{+0.10^{\phantom{X}}}_{-0.12_{\phantom{X}}}$ |  $0.83(2)$[9]
$n$ | $2.10$ | $0.80$ | $1.94$ | $2.09(2)$ | $0.82^{+0.10^{\phantom{X}}}_{-0.10_{\phantom{X}}}$ |  $0.95(6)$[9]
$\Lambda$ | $2.29$ | $0.80$ | $1.94$ | $1.89(9)$ | $0.85^{+0.15^{\phantom{X}}}_{-0.07_{\phantom{X}}}$ |  $0.83(4)$ [12]
$\Sigma^{0}$ | $2.00$ | $0.80$ | $1.10$ | $1.17(3)$ | $1.85^{+0.10}_{-0.24}$ |  $2.11(5)$ [12]
$\Sigma^{+}$ | $1.50(3)$ |  $2.27(5)$ [13]
$\Xi^{+}$ | $1.60$ | $0.50$ | $1.38$ | $0.97(8)$ | $0.65^{+0.10}_{-0.03}$ |  $0.61(5)$ [10]
$\Xi^{0}$ | $1.16(44)$ |  $0.66(3)$ [11]

The given choice of the different parameters for baryon DAs allows one to
describe the data for the $\gamma_{B}$ within the accuracy $20\%$, which is
sufficiently reasonable taking into account various higher order effects. The
largest differences are obtained for neutron and for the $\Sigma$ decay
channels. In considered description, the neutron channel differs from the
proton one only in the contribution of electromagnetic form factors. However
this interference can not provide a sufficient numerical impact. On the other
hand the neutron data in Ref.[9] have very large systematic error, which gives
$\gamma_{n}=0.95(6)\pm 0.27_{syst}$. Therefore the obtained value of
$\gamma_{n}$ can be accepted as quite reasonable result. For $\Sigma$ decay
channel, it is possible to obtain the larger value of $\gamma_{\Sigma}$ but
corresponding parameter set gives a worse description of the branching ratio
because of relatively large and negative power correction.

Comparing the current analysis with that in Ref.[20] we conclude that the NLP
corrections improve the phenomenological description and make it possible to
better constrain the unknown parameters of baryon DAs. The negative effects of
the power contribution reduces the absolute value of the amplitude
$\mathcal{G}_{M}$ and increases the value of the ratio $\gamma_{B}$ that
improves the phenomenological description. However, this effect cannot be very
large, since this can significantly worsen the description of the branching
fractions.

##  5 Conclusions

In conclusion, in this work we continue to develop the analysis proposed in
Ref.[20]. The decay amplitudes $J/\psi$ into baryon-antibaryon pairs are
considered within the NRQCD and collinear factorisation approach, power
suppressed corrections and electromagnetic contributions are taken into
account. This allows us to better appreciate some of the theoretical
uncertainties.

The next-to-leading power correction, which is described by twist-4 and
twist-5 three-quark baryon DAs is calculated and discussed for the first time.
It is shown that this contribution is well defined and this allows one to
study their effects in a more systematic way. From the phenomenological point
of view the NLP correction allows to improve the description of the data and
provide more rigorous constrains on the baryons DAs parameters.

In the given analysis we also take into account the interference with the
electromagnetic amplitude, which depends on the time-like baryon form factors.
Unfortunately, these values are little known due to the lack of experimental
information. Therefore, an evaluation of their numerical impact can only be
made under certain assumptions.

The existing experimental data for decays $J/\psi\to B\bar{B}$ provide an
interesting information about two observables: branching fractions and angular
distribution coefficients $\alpha_{B}$. The value of $\alpha_{B}$ is
determined by the ratio of the two independent decay amplitudes, which is
especially interesting because it does not contain many theoretical
uncertainties that are significant for the branching fractions. This also
provides an attractive opportunity to study the baryon distribution
amplitudes.

One of the main goal of the given phenomenological analysis is to verify that
the current data can be described within the factorisation framework. In order
to constrain the non-perturbative baryon DAs we use results obtained in the
framework of QCD sum rules [25, 5] and from the lattice calculations [28].
Some of unknown moments have been estimated using flavour $SU(3)$-symmetry and
from data fit. It is found that under reasonable assumptions about the non-
perturbative parameters, the data can be described with an accuracy of
$5\%-20\%$.

It is shown that the power corrections arising from higher twist baryonic DAs
are moderate. This conclusion is, of course, sensitive to the models of
baryonic DAs. The maximal numerical effect is about $16\%$ is obtained for
nucleon only. It is interesting to note that the model ABOI proposed in
Ref.[25] for describing nucleon electromagnetic form factors also works well
for describing charmonium decay, if we use the value for the normalisation
couplings $f_{N}$ obtained from the QCD sum rules [5]. The values of the
branching fractions are very sensitive to the normalisation scale, a reliable
description is obtained only for the relatively low scale $\mu^{2}=1.5\ $GeV2.

Note also that a large value of $\gamma\Sigma\sim 2$ leads to a large
numerical effect, which is crucially important for the correct estimate of the
decay width. Unlike other baryons, in the case of $\Sigma$ the contribution of
the power suppressed term in the expression for the width, see the equation
(160 ), is of the order of one, because it is enhanced by a factor of four due
to large $\gamma_{\Sigma}$.

The angular coefficient $\gamma_{B}$ is quite sensitive to the interference
with the electromagnetic decay amplitudes. Under the described assumptions,
the corresponding numerical effect is estimated in the range of $10\%-25\%$.
This definitely reduces the accuracy in the constraining of DA parameters.
Nevertheless, the obtained results allow us to conclude that the calculated
hadronic amplitudes make a dominant contribution and therefore the QCD
factorisation provides a suitable basis for understanding dynamics of these
decays.

##  6 Acknowledgements

This work is supported by the Deutsche Forschungsgemeinschaft (DFG, German
Research Foundation) Project-ID 445769443.

## Appendix

##  Appendix A  Analytical expressions for the models of baryon DAs

In this Appendix we provide the explicit expressions for the models of the
basic light-cone baryon DAs, which define the QCD nonperturbative input in the
compact form. These DAs are used in order to compute the auxiliary DAs, which
are convenient for the calculation of the decay amplitudes. Recall, that in
this work we follow the definitions from Ref.[25].

We start from the description of the twist-3 DAs. For the nucleon DA
introduced in Eq.(30) we use the following expression

$\displaystyle\varphi^{N}_{3}(x_{i})=120x_{1}x_{2}x_{3}\left(1+\phi^{N}_{10}\mathcal{P}_{10}(x_{i})+\phi^{N}_{11}\mathcal{P}_{11}(x_{i})+\phi^{N}_{20}\mathcal{P}_{20}(x_{i})+\phi^{N}_{21}\mathcal{P}_{21}(x_{i})+\phi^{N}_{22}\mathcal{P}_{22}(x_{i})\right).$
(A.1)

For other baryons $B\neq N$ the DAs are defined in Eqs.(19) and (20).
Corresponding models read

$\displaystyle\varphi_{3+}^{B}(x_{i})$ $\displaystyle=$ $\displaystyle
120x_{1}x_{2}x_{3}~{}\left(1+\phi_{11}^{B}\mathcal{P}_{11}(x_{i})\right),\,\varphi_{3-}^{B}(x_{i})=120x_{1}x_{2}x_{3}\,\phi_{10}^{B}\mathcal{P}_{10}(x_{i}).$
(A.2) $\displaystyle\Pi_{3}^{B}(x_{i})$ $\displaystyle=$ $\displaystyle
120x_{1}x_{2}x_{3}\left((1-\delta_{\Lambda
B})\\{1+\pi_{11}^{B}P_{11}(x_{i})\\}+\delta_{\Lambda
B}\,\pi^{B}_{10}\mathcal{P}_{10}(x_{i})\right).$ (A.3)

The orthogonal polynomials $\mathcal{P}_{21}(x_{i})$ read

$\mathcal{P}_{10}(x_{i})=21(x_{1}-x_{3}),~{}~{}\mathcal{P}_{11}(x_{i})=7(x_{1}-2x_{2}+x_{3}),$
(A.4)
$\displaystyle\mathcal{P}_{20}(x_{i})=\frac{63}{10}\left[3(x_{1}-x_{3})^{2}-3x_{2}(x_{1}+x_{3})+2x_{2}^{2}\right],~{}$
(A.5)
$\displaystyle\mathcal{P}_{21}(x_{i})=\frac{63}{2}(x_{1}-3x_{2}+x_{3})(x_{1}-x_{3}),$
(A.6)
$\displaystyle\mathcal{P}_{22}(x_{i})=\frac{9}{5}\left[x_{1}^{2}+9x_{2}(x_{1}+x_{3})-12x_{1}x_{3}-6x_{2}^{2}+x_{3}^{2}\right].$
(A.7)

The moments $\phi_{ik}^{B}$ are multiplicatively renormalisable and the
corresponding anomalous dimensions can be found in Ref.[25].

The set of the models for the twist-4 DAs are defined as following. The
nucleon DAs from Eqs.(31) and(32) have the following twist decomposition

$\Phi_{4}(x_{123})=f_{N}\Phi_{4}^{(3)}(x_{123})+\lambda_{1}\bar{\Phi}_{4}(x_{123}),$
(A.8)
$\Psi_{4}(x_{123})=f_{N}\Psi_{4}^{(3)}(x_{123})-\lambda_{1}\bar{\Psi}_{4}(x_{123}),$
(A.9)

The kinematical twist-3 contributions read

$\displaystyle\Phi_{4}^{(3)}(x_{123})=$ $\displaystyle
40x_{1}x_{2}\left(1-2x_{3}\right)-20x_{1}x_{2}\sum_{k=0,1}\phi^{N}_{1k}\left(3-\frac{\partial}{\partial
x_{3}}\right)x_{3}\mathcal{P}_{1k}(x_{123})$
$\displaystyle-12x_{1}x_{2}\sum_{k=0,1,2}\phi^{N}_{2k}\left(4-\frac{\partial}{\partial
x_{3}}\right)x_{3}\mathcal{P}_{2k}(x_{123}).$ (A.10)
$\Psi_{4}^{(3)}(x_{123})=40x_{1}x_{3}\left(1-2x_{2}\right)-20x_{1}x_{3}\sum_{k=0,1}\phi^{N}_{1k}\left(3-\frac{\partial}{\partial
x_{2}}\right)x_{2}\mathcal{P}_{1k}(x_{213})$ (A.11)
$-12x_{1}x_{3}\sum_{k=0,1,2}\phi^{N}_{2k}\left(4-\frac{\partial}{\partial
x_{2}}\right)x_{2}\mathcal{P}_{2k}(x_{213}).$ (A.12)

Notice that the differentiations must be computed with the unmodified
expressions of the polynomials $\mathcal{P}_{nk}(x_{i})$ in Eqs.(A.4)-(A.7)
and only after that one can apply the condition $x_{1}+x_{2}+x_{3}=1$.

The genuine twist-4 DAs in Eqs.(A.8), (A.9) and (33) reads

$\bar{\Phi}_{4}(x_{123})=24x_{1}x_{2}\left(1+\eta^{N}_{10}\mathcal{R}_{10}(x_{312})-\eta^{N}_{11}\mathcal{R}_{11}(x_{312})\right),$
(A.13)
$\bar{\Psi}_{4}(x_{123})=24x_{1}x_{3}\left(1+\eta^{N}_{10}\mathcal{R}_{10}(x_{231})+\eta^{N}_{11}\mathcal{R}_{11}(x_{231})\right),$
(A.14)
${\Xi}_{4}(x_{i})=24x_{2}x_{3}\left(1+\frac{9}{4}\xi^{N}_{10}\mathcal{R}_{10}(x_{132})\right).\
$ (A.15)

where the orthogonal polynomials $\mathcal{R}_{1i}$ read

$\mathcal{R}_{10}(x_{1},x_{2},x_{3})=4\left(x_{1}+x_{2}-3/2x_{3}\right),~{}\mathcal{R}_{11}(x_{1},x_{2},x_{3})=\frac{20}{3}\left(x_{1}-x_{2}+x_{3}/2\right).$
(A.16)

For other baryons we use definitions from Eqs.(37)-(40). Recall, that DA
$\Upsilon_{4}$ from Eq.(41) is not independent and can be obtained from Eq.(70
). The kinematical twist-3 contributions in Eqs.(37)-(39) read

$\Phi_{4+}^{B(3)}(x_{123})=40x_{1}x_{2}\left(1-2x_{3}\right)-20x_{1}x_{2}\phi^{B}_{11}\left(3-\frac{\partial}{\partial
x_{3}}\right)x_{3}\mathcal{P}_{11}(x_{123}),$ (A.17)
$\Phi_{4-}^{B(3)}(x_{123})=-20x_{1}x_{2}\
\phi^{B}_{10}\left(3-\frac{\partial}{\partial
x_{3}}\right)x_{3}\mathcal{P}_{10}(x_{123}).$ (A.18)
$\Pi_{4}^{B(3)}(x_{123})=40x_{1}x_{2}\left(1-2x_{3}\right)-20x_{1}x_{2}\pi_{11}^{B}\left(3-\frac{\partial}{\partial
x_{3}}\right)x_{3}\mathcal{P}_{11}(x_{123}).$ (A.19)
$\Pi_{4}^{\Lambda(3)}(x_{123})=-20x_{1}x_{2}\
\pi_{10}^{\Lambda}\left(3-\frac{\partial}{\partial
x_{3}}\right)x_{3}\mathcal{P}_{10}(x_{123}).$ (A.20)

The genuine twist-4 contributions read

$\displaystyle\bar{\Phi}^{B}_{4+}(x_{123})$ $\displaystyle=$ $\displaystyle
24x_{1}x_{2}\left(-\eta^{B}_{11}\right)\mathcal{R}_{11}(x_{312}),$ (A.21)
$\displaystyle\bar{\Phi}^{B}_{4-}(x_{123})$ $\displaystyle=$ $\displaystyle
24x_{1}x_{2}\left(1+\eta^{B}_{10}\mathcal{R}_{10}(x_{312})\right).$ (A.22)
$\displaystyle\bar{\Xi}^{B}_{4+}(x_{123})$ $\displaystyle=$ $\displaystyle
24x_{2}x_{3}\left(1+\frac{9}{8}\zeta^{B}_{10}(\mathcal{R}_{10}(x_{132})+\mathcal{R}_{10}(x_{123}))\right),$
(A.23) $\displaystyle\bar{\Xi}^{B}_{4-}(x_{123})$ $\displaystyle=$
$\displaystyle 27x_{2}x_{3}(4^{-\delta_{\Lambda
B}})\zeta^{B}_{10}(\mathcal{R}_{10}(x_{132})-\mathcal{R}_{10}(x_{123})).$
(A.24) $\displaystyle\bar{\Pi}^{B}_{4}(x_{123})$ $\displaystyle=$
$\displaystyle
24x_{1}x_{2}\left(-\zeta^{B}_{11}\right)\mathcal{R}_{11}(x_{312}),$ (A.25)
$\displaystyle\bar{\Pi}^{\Lambda}_{4}(x_{123})$ $\displaystyle=$
$\displaystyle
24x_{1}x_{2}\left(1+\zeta^{\Lambda}_{10}\mathcal{R}_{10}(x_{312})\right).$
(A.26)

The basic twist-5 DAs are used in order to obtain the auxiliary DAs, which are
discussed in Sec.2.2. To obtain them, only the asymptotic expressions for DAs
of twist-3 and twist-4 are used, i.e. all higher moments except
$f_{B},\,f^{B}_{\perp}$, $\lambda^{B}_{1,2}$ and $\lambda^{B}_{\bot}$ are
neglected in the above expressions. For the twist-5 basic DAs we used the
following input

$\displaystyle\Phi^{B(3)}_{5+}(x_{i})=40x_{3}(1-x_{2}-2x_{1}-2x_{1}x_{2}),\quad\Phi^{B(3)}_{5-}(x_{i})=0.$
(A.27)
$\displaystyle\Phi^{B(4)}_{5+}(x_{i})=0,\quad\Phi^{B(4)}_{5-}(x_{i})=8x_{3}(1-x_{2}).$
(A.28) $\displaystyle\Xi^{B(4)}_{5+}(x_{i})$ $\displaystyle=$ $\displaystyle
4x_{1}(1+x_{1}),\quad\Xi^{B(4)}_{5-}(x_{i})=-12x_{1}(x_{2}-x_{3}),$ (A.29)
$\displaystyle\Upsilon_{5}^{B(4)}(x_{i})$ $\displaystyle=$
$\displaystyle(1-\delta_{B\Lambda})4x_{1}\left(1+x_{1}\right)-\delta_{B\Lambda}12x_{1}(x_{2}-x_{3}).$
(A.30) $\displaystyle\Pi_{5}^{B(3)}(x_{i})$ $\displaystyle=$
$\displaystyle(1-\delta_{B\Lambda})40x_{3}(1-2x_{1}-x_{2}-2x_{1}x_{2}),$
(A.31) $\displaystyle\Pi_{5}^{B(4)}(x_{i})$ $\displaystyle=$
$\displaystyle\delta_{B\Lambda}8x_{3}(1-x_{2}).$ (A.32)

Numerical values of all moments are given in the Tab.1. All these moments are
multiplicatively renomalisable and the corresponding anomalous dimensions can
be found in Refs.[23, 25].

##  Appendix B  Relations between the twist-5 DAs

Here we discuss the technical details that clarify the derivation of
relationships between the various DAs introduced in Sec.2.

Let us start from discussion of the chiral even DAs. There are two ways to get
additional equations, which relate the auxiliary and the basic DAs. One way is
to consider the light-cone operator with one derivative and to repeat the same
analysis as in Sec.2. The second method is to use the off light-cone
correlation function (CF). Below we use this one. For simplicity we consider
the nucleon case. We also do not write explicitly the colour indices and the
Wilson lines.

Consider the following CF ($x^{2}\neq 0,\ y^{2}\neq 0$)

$\displaystyle\left\langle
0\right|u(x)C\gamma^{\sigma}u(y)d_{3}(0)\left|B(k)\right\rangle=-k^{\sigma}\left[\gamma_{5}N\right]_{3}\text{FT}\left[V^{B}_{1}\right]\phantom{\left\langle
0\right|u(x)C\gamma^{\sigma}u(y)d_{3}(0)\left|k\right\rangle=}$
$\displaystyle-\frac{m_{N}}{2}\left[\gamma^{\sigma}\gamma_{5}N\right]_{3}\text{FT}\left[V^{B}_{3}\right]-m_{N}~{}ik^{\sigma}\left[\gamma_{\alpha}\gamma_{5}N\right]_{3}\left(x^{\alpha}\text{FT}\left[\mathcal{V}^{B}_{1}\right]+y^{\alpha}\text{FT}\left[\mathcal{V}^{B}_{2}\right]\right)$
$\displaystyle-\frac{m_{N}^{2}}{4}k^{\sigma}~{}\left[\gamma_{5}N\right]_{3}\left\\{x^{2}\text{FT}\left[V^{B}_{1xx}\right]+y^{2}\text{FT}\left[V^{B}_{1yy}\right]+2(xy)\text{FT}\left[V^{B}_{1xy}\right]\right\\}$
$\displaystyle-\frac{m_{N}^{2}}{2}k^{\sigma}x_{\alpha}y_{\beta}\left[i\sigma^{\alpha\beta}\gamma_{5}N\right]_{3}\
\text{FT}\left[T^{B}_{xy}\right]-\frac{m_{N}^{2}}{2}~{}\left[\gamma_{5}N\right]_{3}\left\\{x^{\sigma}\text{FT}\left[V^{B}_{x4}\right]+y^{\sigma}\text{FT}\left[V^{B}_{y4}\right]\right\\}$
$\displaystyle-\frac{m_{N}^{2}}{4}~{}\left[i\sigma^{\alpha\sigma}\gamma_{5}N\right]_{3}~{}\left\\{x_{\alpha}\text{FT}\left[V^{B}_{x5}\right]+y_{\alpha}\text{FT}\left[V^{B}_{y5}\right]\right\\}.$
(B.1)

Performing the expansion of the operator in lhs of Eq.(B.1) around the light-
cone direction up to twist-5 terms and expanding the rhs  of Eq.(B.1) one
obtains the set of relations.

The twist-5 operators resulting from the expansion can be divided into three
groups. The first group includes the operators with two transverse
derivatives, which arise from the expansion of the arguments of the fields

$\displaystyle u(x)$
$\displaystyle\rightarrow\frac{1}{2}x_{\bot}x_{\bot}\partial_{\bot}\partial_{\bot}u(x_{-}),\text{
}u(y)\rightarrow\frac{1}{2}y_{\bot}y_{\bot}\partial_{\bot}\partial_{\bot}u(y_{-}),\
$ (B.2) $\displaystyle\ u(x)u(y)$ $\displaystyle\rightarrow
x_{\bot}\left[\partial_{\bot}u(x_{-})\right]y_{\bot}\left[\partial_{\bot}u(y_{-})\right].$
(B.3)

There second group includes the terms with one transverse derivative and one
small field $\eta$

$u(x)\rightarrow x_{\bot}\partial_{\bot}\eta(x_{-}),\ \ u(x)u(y)\rightarrow
x_{\bot}\left[\partial_{\bot}u(x_{-})\right]\eta(y_{-}),...\ .$ (B.4)

The third group includes the operators with two small fields $\eta$

$u(x)u(y)\rightarrow\eta(x_{-})\eta(y_{-}),...\ .$ (B.5)

All operators with fields $\eta$ can be converted into the operators with two
transverse derivatives using Eq.(88) and corresponding matrix elements can be
obtained using Eqs. (75)-(2.2).

The matrix elements of the operators from the first group simply reproduce the
definitions in Eqs. (75), (76) and (2.2).

The required relations can be obtained from the matrix elements of the
operators from the second group (B.4). Consider the terms with $x_{\bot}$.
There are three twist-5 light-cone operators, which appear in this case

$\left\\{O_{21},O_{22},O_{23}\right\\}=\left\\{\left[\partial_{\bot}^{\alpha}\xi\right]C\gamma_{\bot}^{\sigma}\eta\xi_{3},\,\left[\partial_{\bot}^{\alpha}\xi\right]C\gamma_{\bot}^{\sigma}\eta\xi_{3},\,\left[\partial_{\bot}^{\alpha}\xi\right]C\gamma^{\sigma}\xi\eta_{3}\right\\},$
(B.6)

where we assume $qqq_{3}\equiv q(x_{-})q(y_{-})q_{3}(0)$. The matrix elements
of these operators can be obtaned using QCD EOM (88) and definitions Eqs.
(75), (76) and (2.2). On the other side the same matrix elements can be
obtained from the expansion of rhs in Eq.(B.1). Comparing these results one
obtains

$\left\langle
0\right|O_{23}^{\sigma}\left|B(k)\right\rangle:-2i(py)\text{FT}\left[T_{xy}\right]=\text{FT}\left[\frac{1}{x_{3}}\left(V^{B}_{1xx}+V^{B}_{1xy}+T^{B}_{xy}\right)-2\mathcal{V}^{B}_{1}-V^{B}_{x5}\right],$
(B.7) $\left\langle
0\right|O_{21}^{\sigma}+O_{22}^{\sigma}\left|B(k)\right\rangle:V^{B}_{x5}=\frac{1}{x_{2}}\left(A^{B}_{1xy}+T^{B}_{xy}\right)-\frac{1}{x_{1}}A^{B}_{1xx}.$
(B.8)

Excluding $V^{B}_{x5}$ one obtains

$-2i(ky)\text{FT}\left[T^{B}_{xy}\right]=\text{FT}\left[\frac{1}{x_{3}}\left(V^{B}_{1xx}+V^{B}_{1xy}+T^{B}_{xy}\right)+\frac{1}{x_{1}}A^{B}_{1xx}-\frac{1}{x_{2}}\left(A^{B}_{1xy}+T^{B}_{xy}\right)-2\mathcal{V^{B}}_{1}\right],$
(B.9)

where we assume $2(ky)=k_{+}y_{-}$.

The similar analysis for the operators with $y_{\bot}$

$\left\\{O_{21},O_{22},O_{23}\right\\}=\left\\{\xi
C\gamma_{\bot}^{\sigma}\left[\partial_{\bot}^{\alpha}\eta\right]\xi_{3},\,\xi
C\gamma_{\bot}^{\sigma}\left[\partial_{\bot}^{\alpha}\eta\right]\xi_{3},\,\xi
C\gamma^{\sigma}\left[\partial_{\bot}^{\alpha}\xi\right]\eta_{3}\right\\},$
(B.10)

allows one to obtain one more relation

$2i(kx)\text{FT}\left[T^{B}_{xy}\right]=\text{FT}\left[\frac{1}{x_{3}}\left(V^{B}_{1yy}+V^{B}_{1xy}-T^{B}_{xy}\right)-\frac{1}{x_{2}}A^{B}_{1yy}+\frac{1}{x_{1}}\left(A^{B}_{1xy}+T^{B}_{xy}\right)-2\mathcal{V}^{B}_{2}\right].$
(B.11)

Consider now the axial CF

$\displaystyle\left\langle
0\right|u(x)C\gamma^{\sigma}\gamma_{5}u(y)d_{\gamma}(0)\left|B(k)\right\rangle=-k^{\sigma}\left[N\right]_{3}\text{FT}\left[A^{B}_{1}\right]\phantom{\left\langle
0\right|u(x)C\gamma^{\sigma}\gamma_{5}u(y)d_{\gamma}(0)\left|B(k)\right\rangle=}$
$\displaystyle-\frac{m_{N}}{2}\left[\gamma^{\sigma}N\right]_{3}\text{FT}\left[A^{B}_{3}\right]-ik^{\sigma}m_{N}~{}\left[\gamma^{\alpha}N\right]_{3}\left(x_{\alpha}\text{FT}\left[\mathcal{A}^{B}_{1}\right]+y_{\alpha}\text{FT}\left[\mathcal{A}^{B}_{2}\right]\right)$
$\displaystyle-\frac{m_{N}^{2}}{4}k^{\sigma}\left[N\right]_{3}\left\\{x^{2}\text{FT}\left[A^{B}_{1xx}\right]+y^{2}\text{FT}\left[A^{B}_{1yy}\right]+2(xy)\text{FT}\left[A^{B}_{1xy}\right]\right\\}$
$\displaystyle-\frac{m_{N}^{2}}{2}k^{\sigma}\left[i\sigma^{\alpha\beta}N\right]_{3}x^{\alpha}y^{\beta}\text{FT}\left[G^{B}_{xy}\right]-i\frac{m_{N}^{2}}{2}~{}\left[N\right]_{3}\left\\{x^{\sigma}\text{FT}\left[A^{B}_{4x}\right]+y^{\sigma}\text{FT}\left[A^{B}_{4y}\right]\right\\}$
$\displaystyle-i\frac{m_{N}^{2}}{4}~{}\left[i\sigma^{\alpha\sigma}N\right]_{3}~{}\left\\{x^{\alpha}\text{FT}\left[A^{B}_{5x}\right]+y^{\alpha}\text{FT}\left[A^{B}_{5y}\right]\right\\}.$
(B.12)

The analogous consideration of terms with $x_{\bot}$ and $y_{\bot}$ yields the
two relations

$\displaystyle
2i(ky)\text{FT}\left[G^{B}_{xy}\right]=\text{FT}\left[-\frac{1}{x_{3}}\left(A^{B}_{1xx}+A^{B}_{1xy}+G^{B}_{xy}\right)-2\mathcal{A}^{B}_{1}+\frac{1}{x_{2}}\left(V^{B}_{1xy}+G^{B}_{xy}\right)-\frac{1}{x_{1}}V^{B}_{1xx}\right],$
(B.13) $\displaystyle
2i(kx)\text{FT}\left[G^{B}_{xy}\right]=\text{FT}\left[\frac{1}{x_{3}}\left(A^{B}_{1yy}+A^{B}_{1xy}-G^{B}_{xy}\right)+2\mathcal{A}^{B}_{2}+\frac{1}{x_{1}}\left(V^{B}_{1xy}+G^{B}_{xy}\right)-\frac{1}{x_{2}}V^{B}_{1yy}\right].$
(B.14)

It is convenient to combine these two relations with ones in Eqs.(B.9) and
(B.11) and simplify the obtained expressions using the relations (102)–(105).
This gives the relations given in Eqs.(110) and (111).

Consider now the auxiliary chiral-odd DAs. The set of these functions is
defined in Eqs.(81)-(83) and reads

$T^{B}_{0xx},T^{B}_{0yy},T^{B}_{2xx},T^{B}_{2yy},T^{B}_{0xy},T^{B}_{2xy},\tilde{T}^{B}_{2xy}.$
(B.15)

The algebraic relations can be obtained in the same way as described for the
chiral-even DAs in Sec.2.2. Therefore we skip the details and only provide a
schematical description. The consideration of the scalar operators $\eta
C\xi\eta_{3}$ and $\xi C\eta\eta_{3}$ yields

$T^{B}_{0xx}+2T^{B}_{2xx}+T^{B}_{0xy}+2T^{B}_{2xy}-\tilde{T}^{B}_{2xy}=2x_{1}x_{3}\left(S^{B}_{2}-T^{B}_{4}\right),$
(B.16)
$T^{B}_{0yy}+2T^{B}_{2yy}+T^{B}_{0xy}+2T^{B}_{2xy}+\tilde{T}^{B}_{2xy}=2x_{2}x_{3}\left(-S^{B}_{2}-T^{B}_{4}\right),$
(B.17)

respectively. From the analysis of the pseuvdoscalar operators $\eta
C\gamma_{5}\xi\eta_{3}$ and $\xi C\gamma_{5}\eta\eta_{3}$ one finds

$T^{B}_{0xx}-2T^{B}_{2xx}+T^{B}_{0xy}-2T^{B}_{2xy}-\tilde{T}^{B}_{2xy}=-2x_{1}x_{3}\left(P^{B}_{2}+T^{B}_{8}\right),$
(B.18)
$T^{B}_{0yy}-2T^{B}_{2yy}+T^{B}_{0xy}-2T^{B}_{2xy}+\tilde{T}_{2xy}=2x_{2}x_{3}\left(P^{B}_{2}-T^{B}_{8}\right).$
(B.19)

One more relation follows from the consideration of the tensor operator $\eta
C\bar{n}\gamma_{\bot}\eta\xi_{3}$:

$T^{B}_{2xy}=-x_{1}x_{2}T^{B}_{5}.$ (B.20)

Using the same consideration as in Eq.(107) we consider the matrix element of
the operator $(\bar{n}\partial)\left[\xi
Cn\gamma_{\bot}^{\nu}\xi\xi_{3}\right]$. This yields

$-2x_{1}x_{2}x_{3}T^{B}_{1}=x_{2}\bar{x}_{2}T^{B}_{0xx}+x_{1}\bar{x}_{1}T^{B}_{0yy}+2x_{1}x_{2}T^{B}_{0xy}.$
(B.21)

All other operators give the linear combinations of these six relations.

In order to get one more relation we use the method of CF, which is described
above. Consider the following CF

$\displaystyle\left\langle
0\right|u(x)C\sigma^{\mu\nu}u(y)d_{3}(0)\left|B(k)\right\rangle=-~{}i(k^{\nu}g^{\alpha\mu}-k^{\mu}g^{\alpha\nu})\left[\gamma^{\alpha}\gamma_{5}N\right]_{3}\text{FT}\left[T^{B}_{1}\right]\phantom{\langle
0|u(x)C\sigma^{\mu\nu}u(y)}$
$\displaystyle-~{}2m_{N}\left[\sigma^{\mu\nu}\gamma_{5}N\right]_{3}\text{FT}\left[\mathcal{T}^{B}_{3}\right]-m_{N}\left(k^{\nu}g^{\alpha\mu}-k^{\mu}g^{\alpha\nu}\right)\left[\gamma_{5}N\right]_{3}\left(x^{\alpha}\text{FT}\left[\mathcal{T}^{B}_{2x}\right]+y^{\alpha}\text{FT}\left[\mathcal{T}^{B}_{2y}\right]\right)$
$\displaystyle~{}-m_{N}\left(k^{\nu}g^{\alpha\mu}-k^{\mu}g^{\alpha\nu}\right)\left[i\sigma^{\alpha\beta}\gamma_{5}N\right]_{3}\left(x^{\beta}\text{FT}\left[\mathcal{T}^{B}_{4x}\right]+y^{\beta}\text{FT}\left[\mathcal{T}^{B}_{4y}\right]\right)$
$\displaystyle~{}-m_{N}^{2}\left[\gamma^{\alpha}\gamma_{5}N\right]_{3}i(p^{\nu}g^{\alpha\mu}-p^{\mu}g^{\alpha\nu})\left(x^{2}\text{FT}\left[T^{B}_{1xx}\right]+y^{2}\text{FT}\left[T^{B}_{1yy}\right]+2(xy)\text{FT}\left[T^{B}_{1xy}\right]\right)$
$\displaystyle-
m_{N}^{2}\left[\gamma^{\alpha}\gamma_{5}N\right]_{3}i(p^{\nu}g^{\beta\mu}-p^{\mu}g^{\beta\nu})\left(x^{\alpha}x^{\beta}\text{FT}\left[T^{B}_{2xx}\right]+y^{\alpha}y^{\beta}\text{FT}\left[T^{B}_{2yy}\right]\right)$
$\displaystyle-
m_{N}^{2}\left[\gamma^{\alpha}\gamma_{5}N\right]_{3}i(p^{\nu}g^{\beta\mu}-p^{\mu}g^{\nu\beta})\left\\{x^{\alpha}y^{\beta}+x^{\beta}y^{\alpha}\right\\}\text{FT}\left[T^{B}_{2xy}\right]$
$\displaystyle-
m_{N}^{2}\left[\gamma^{\alpha}\gamma_{5}N\right]_{3}i(p^{\nu}g^{\beta\mu}-p^{\mu}g^{\nu\beta})\left\\{x^{\alpha}y^{\beta}-x^{\beta}y^{\alpha}\right\\}\text{FT}\left[\tilde{T}^{B}_{2xy}\right]$
$\displaystyle-
m_{N}^{2}\left[\gamma^{\alpha}\gamma_{5}N\right]_{3}\left(g^{\alpha\nu}g^{\mu\beta}-g^{\alpha\mu}g^{\beta\nu}\right)\left(x^{\beta}\text{FT}\left[T^{B}_{5x}\right]+y^{\beta}\text{FT}\left[T^{B}_{5y}\right]\right)~{}$
$\displaystyle-
m_{N}^{2}\left[\gamma^{\alpha}\sigma^{\mu\nu}\gamma_{5}N\right]_{3}\left(x^{\alpha}\text{FT}\left[T^{B}_{8x}\right]+y^{\alpha}\text{FT}\left[T^{B}_{8y}\right]\right),$
(B.22)

where the rhs includes all terms up to twist-5 accuracy. Expanding the
operator in the lhs we consider the following linear in $x_{\bot}$ and
$y_{\bot}$ twist-5 terms

$\displaystyle\left.u(x)C\sigma^{\mu\nu}u(y)d_{3}(0)\right|_{\text{tw5}}=x_{\bot\alpha}\mathcal{O}_{x\bot}^{\alpha\mu\nu}+y_{\bot\alpha}\mathcal{O}_{y\bot}^{\alpha\mu\nu}\phantom{\left.u(x)C\sigma^{\mu\nu}u(y)d_{3}(0)\right|_{\text{tw5}}}$
$\displaystyle+\frac{i}{2}\left(\bar{n}^{\mu}g_{\bot}^{\nu\sigma}-\bar{n}^{\nu}g_{\bot}^{\mu\sigma}\right)\left(x_{\bot\alpha}\mathcal{O}_{x}^{\alpha\sigma}+y_{\bot\alpha}\mathcal{O}_{y}^{\alpha\sigma}\right)+\ldots,$
(B.23)

where

$\mathcal{O}_{x\bot}^{\alpha\mu\nu}=\left[\partial_{\bot}^{\alpha}\eta(x_{-})\right]C\sigma_{\bot\bot}^{\mu\nu}\xi(y_{-})\xi_{3}(0)+\left[\partial_{\bot}^{\alpha}\xi(x_{-})\right]C\sigma_{\bot\bot}^{\mu\nu}\eta(y_{-})\xi_{3}(0),$
(B.24)
$\mathcal{O}_{y\bot}^{\alpha\mu\nu}=\eta(x_{-})C\sigma_{\bot\bot}^{\mu\nu}\left[\partial_{\bot}^{\alpha}\xi(y_{-})\right]\xi_{3}(0)+\xi(x_{-})C\sigma_{\bot\bot}^{\mu\nu}\left[\partial_{\bot}^{\alpha}\eta(y_{-})\right]\xi_{3}(0),$
(B.25)
$\mathcal{O}_{x}^{\alpha\sigma}=\left[\partial_{\bot}^{\alpha}\xi(x_{-})\right]C\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\gamma_{\bot}^{\sigma}\xi(y_{-})\eta_{3}(0),$
(B.26) $\mathcal{O}_{y}^{\alpha\sigma}=\xi(x_{-})C\hbox to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\gamma_{\bot}^{\sigma}\left[\partial_{\bot}^{\alpha}\xi(y_{-})\right]\eta_{3}(0).$
(B.27)

The consideration of the corresponding matrix elements yields

$\left\langle
0\right|\mathcal{O}_{x\bot}^{\alpha\mu\nu}\left|B(k)\right\rangle:T^{B}_{5x}+T^{B}_{8x}=\frac{1}{8}\frac{1}{x_{1}}\left(T^{B}_{0xx}-2T^{B}_{2xx}\right)+\frac{1}{8}\frac{1}{x_{2}}\left(T^{B}_{0xy}-2T^{B}_{2xy}+\tilde{T}^{B}_{2xy}\right),$
(B.28) $\left\langle
0\right|\mathcal{O}_{y\bot}^{\alpha\mu\nu}\left|B(k)\right\rangle:T^{B}_{5y}+T^{B}_{8y}=\frac{1}{8}\frac{1}{x_{2}}\left(T^{B}_{0yy}-2T^{B}_{2yy}\right)+\frac{1}{8}\frac{1}{x_{1}}\left(T^{B}_{0xy}-2T^{B}_{2xy}-\tilde{T}^{B}_{2xy}\right).$
(B.29)

The sum of these two relations gives

$\displaystyle T^{B}_{5y}+T^{B}_{8y}+T^{B}_{5x}+T^{B}_{8x}$
$\displaystyle=\frac{1}{8}\frac{1}{x_{2}}\left(T^{B}_{0yy}-2T^{B}_{2yy}+T^{B}_{0xy}-2T^{B}_{2xy}+\tilde{T}^{B}_{2xy}\right)$
$\displaystyle+\frac{1}{8}\frac{1}{x_{1}}\left(T^{B}_{0xx}-2T^{B}_{2xx}+T^{B}_{0xy}-2T^{B}_{2xy}-\tilde{T}^{B}_{2xy}\right).$
(B.30)

Using Eqs.(B.18) and (B.19) one finds

$T^{B}_{5y}+T^{B}_{8y}+T^{B}_{5x}+T^{B}_{8x}=-\frac{1}{2}x_{3}T^{B}_{8}.$
(B.31)

The consideration of the matrix elements $\ \left\langle
0\right|\mathcal{O}_{x}^{\alpha\sigma}\left|B(k)\right\rangle$ and $\
\left\langle 0\right|\mathcal{O}_{y}^{\alpha\sigma}\left|B(k)\right\rangle$
yields

$\displaystyle
i(kx)\text{FT}\left[T^{B}_{2xx}\right]+i(ky)\text{FT}\left[T^{B}_{2xy}-\tilde{T}^{B}_{2xy}\right]=\text{FT}\left[-x_{1}\left(S^{B}_{2}-T^{B}_{4}\right)-4\left(T^{B}_{5x}+T^{B}_{8x}\right)+2\mathcal{T}^{B}_{2x}\right],$
(B.32) $\displaystyle
i(ky)\text{FT}\left[T^{B}_{2yy}\right]+i(kx)\text{FT}\left[T^{B}_{2xy}+\tilde{T}^{B}_{2xy}\right]=\text{FT}\left[x_{2}\left(S^{B}_{2}+T^{B}_{4}\right)-4\left(T^{B}_{5y}+T^{B}_{8y}\right)+2\mathcal{T}^{B}_{2y}\right],$
(B.33)

respectively. The sum of these two equations gives

$\displaystyle ip(x-y)\text{FT}\left[\tilde{T}^{B}_{2xy}\right]$
$\displaystyle=-i(kx)\text{FT}\left[T^{B}_{2xx}+T^{B}_{2xy}\right]-i(ky)\text{FT}\left[T^{B}_{2yy}+T^{B}_{2xy}\right]$
$\displaystyle\text{FT}\left[2\left(\mathcal{T}^{B}_{2x}+\mathcal{T}^{B}_{2y}\right)+2x_{3}T^{B}_{8}\right]+\text{FT}\left[\left(x_{1}+x_{2}\right)T^{B}_{4}+(x_{2}-x_{1})S^{B}_{2}\right],$
(B.34)

where we used Eq.(B.31). The rhs of this equation can be rewritten in terms of
basic DAs. Combaining Eqs.(B.16) and (B.18) one finds

$T^{B}_{2xx}+T^{B}_{2xy}=\frac{1}{2}x_{1}x_{3}\left(S^{B}_{2}-T^{B}_{4}+P^{B}_{2}+T^{B}_{8}\right).$
(B.35)

Similar combination Eqs.(B.17) and (B.19) gives

$T^{B}_{2yy}+T^{B}_{2xy}=\frac{1}{2}x_{2}x_{3}\left(-S^{B}_{2}-T^{B}_{4}-P^{B}_{2}+T^{B}_{8}\right).$
(B.36)

Substituting this into Eq.(B.34) one obtains the relation, which allows one to
find the DA $\tilde{T}_{2xy}$. The expressions for other DAs in the list
(B.15) can be derived from the linear relations (B.16)-(B.21). The
corresponding results are given in Eqs.(125)-(131).

Let us also notice that consistency of Eq.(B.34) implies (setting $x=y=0$)

$\int_{0}^{1}dx_{1}\int_{0}^{1-x_{1}}dx_{2}\left\\{2\left(\mathcal{T}^{B}_{2x}+\mathcal{T}^{B}_{2y}\right)+2x_{3}T^{B}_{8}+\left(x_{1}+x_{2}\right)T^{B}_{4}+(x_{2}-x_{1})S^{B}_{2}\right\\}(x_{123})=0,$
(B.37)

where we assume $x_{3}=1-x_{1}-x_{2}$. This provides a good check of the
obtained expressions.

Finally, let us also mention that the additional relation for the twist-4
chiral-odd DAs given in Eq.(69), can be obtained from the consideration of the
matrix element of the operator
$\xi(x_{-})C\sigma^{\mu\nu}\xi(y_{-})\eta_{3}(0)$.

##  Appendix C  On the determination of quark and quark-gluon contributions
for the next-to-leading power correction

In this section, we discuss a simple example that demonstrates a difficulty in
the determination of the pure quark contribution for the subleading power
corrections. We assume the same kinematical notations as for the baryon with
momentum $k$. Let us consider the light-cone matrix element for the vector
meson state ($J^{PC}=1^{--}$), which can be parametrised as

$\left\langle
0\right|\bar{\psi}(x_{-})\gamma^{\sigma}\psi(0)\left|V(k)\right\rangle=f_{V}m_{V}\
e_{+}\frac{\bar{n}^{\sigma}}{2}\text{FT}\left[\phi\right]+f_{V}m_{V}e_{\sigma}^{\bot}\text{FT}\left[g_{v}\right]-f_{V}m_{V}^{2}\frac{n^{\sigma}}{2}\frac{e_{+}}{k_{+}^{2}}\text{FT}\left[g_{3}\right],$
(C.1)

where $e_{\mu}\equiv e^{\lambda}_{\mu}(k)$ denotes the polarisation vector, FT
denotes the appropriate Fourier transformation

$\text{FT}\left[\phi\right]=\int_{0}^{1}due^{-i(1-u)(x_{-}k_{+})/2}\phi(u).$
(C.2)

All defined DAs are normalised to one

$\int_{0}^{1}du~{}\phi(u)\equiv\left\langle\phi\right\rangle=\left\langle
g_{v,3}\right\rangle=1.$ (C.3)

The twist-3 DA $g_{v}$ and the twist-4 DA $g_{3}$ in the rhs have been studied
in Refs.[41, 42]. The obtained results read

$g_{v}(u)=\frac{1}{2}\int_{0}^{u}dv\frac{\phi(v)}{1-v}+\frac{1}{2}\int_{u}^{1}dv\frac{\phi(v)}{v}+\bar{g}_{v}(u),$
(C.4) $g_{3}(x)=\phi(x)+\bar{g}_{3}(x),$ (C.5)

where $\bar{g}_{v}$ and $\bar{g}_{3}$ denote the contributions of quark-gluon
matrix elements.

The local limit in Eq.(C.1) gives

$\left\langle
0\right|\bar{\psi}(0)\gamma^{\sigma}\psi(0)\left|V(k)\right\rangle=f_{V}m_{V}e_{+}\frac{1}{2}\bar{n}^{\sigma}\left\langle\phi\right\rangle+f_{V}~{}m_{V}e_{\sigma}^{\bot}\left\langle
g_{v}\right\rangle-
f_{V}m_{V}\frac{1}{2}n^{\sigma}m_{V}^{2}\frac{e_{+}}{k_{+}^{2}}\left\langle
g_{3}\right\rangle$ (C.6)
$=f_{V}m_{V}\left\\{e_{+}\frac{1}{2}\bar{n}^{\sigma}+e_{\sigma}^{\bot}-\frac{1}{2}n^{\sigma}m_{V}^{2}\frac{e_{+}}{k_{+}^{2}}\right\\}=f_{V}m_{V}e^{\sigma},$
(C.7)

where it is used that

$e\cdot k=0\Rightarrow
e_{-}=-e_{+}k_{-}/k_{+}=-e_{+}k_{-}/k_{+}=-e_{+}\frac{m_{V}^{2}}{k_{+}^{2}}.$
(C.8)

The off light-cone CF can be written as [41]

$\displaystyle\left\langle
0\right|\bar{\psi}(x)\gamma^{\sigma}\psi(0)\left|V(k)\right\rangle=f_{V}m_{V}\frac{(ex)}{(kx)}k^{\sigma}\left(\text{FT}\left[\phi\right]+\frac{m^{2}_{V}x^{2}}{16}\text{FT}\left[\mathbf{A}\right]\right)\phantom{\left\langle
0\right|\bar{\psi}(x)\gamma^{\sigma}\psi(0)\left|k\right\rangle}$
$\displaystyle+f_{V}m_{V}\left(e_{\sigma}-k^{\sigma}\frac{(ex)}{(kx)}\right)\text{FT}\left[g_{v}\right]-f_{V}m_{V}\frac{1}{2}x^{\sigma}m_{V}^{2}\frac{(ex)}{(kx)^{2}}\text{FT}\left[g_{3}+\phi-2g_{v}\right].$
(C.9)

The derivation of the expression for the DA $\mathbf{A}$ is given in Ref.[41],
the result reads

$\mathbf{A}(u)=32\int_{0}^{u}dv\int_{0}^{v}d\omega(g_{v}-\phi)(\omega)+\bar{\mathbf{A}}(u),$
(C.10)

where $\bar{\mathbf{A}}$ again denotes the quark-gluon contributions. Such
notation for the contributions of quark-gluon matrix elements is also accepted
below. All results in Eqs.(C.4),(C.5) and (C.10) are obtained with the help of
differentiation of the CF, see details in Refs.[41, 42].

Let us now consider the twist-4 DAs with help of the technique discussed in
this work. The twist-4 operators, which appear in the expansion of the off
light-cone operator $\bar{\psi}(x)\gamma^{\sigma}\psi(0)$ read

$\displaystyle\left[\bar{\psi}(x)\gamma^{\sigma}\psi(0)\right]_{\text{tw4}}=\frac{1}{4}x_{\bot}^{2}\left[\partial_{\bot}^{2}\bar{\xi}(x_{-})\right]\gamma^{\sigma}\xi(0)+\frac{1}{2}x_{+}\left[(\bar{n}\partial)\bar{\xi}(x_{-})\right]\gamma^{\sigma}\xi(0)$
$\displaystyle~{}\ \ \
+x_{\bot}^{\alpha}\left[\partial_{\bot}^{\alpha}\bar{\xi}(x_{-})\right]\gamma_{\bot}^{\sigma}\eta(0)+x_{\bot}^{\alpha}\left[\partial_{\bot}^{\alpha}\bar{\eta}(x_{-})\right]\gamma_{\bot}^{\sigma}\xi(0)+\bar{\eta}(x_{-})\gamma^{\sigma}\eta(0)+\mathcal{O}_{4}^{\sigma}[\mathcal{A}],$
(C.11)

where, remind, the collinear fields $\xi$ and $\eta$ are defined in Eq.(87).
The operator $\mathcal{O}_{4}^{\sigma}[\mathcal{A}]$ denotes the contribution
of the quark-gluon operators.

The matrix element of the twist-4 operator with two derivatives can be
obtained from the expansion of the CF

$\left\langle 0\right|\left[\partial_{\bot}^{2}\bar{\xi}(x_{-})\right]\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\xi(0)\left|V(k)\right\rangle=f_{V}m_{V}^{3}\frac{e_{+}~{}}{4}~{}\text{FT}\left[\mathbf{A}+\bar{G}_{2}\right].$
(C.12)

Consider now the operators with the longitudinal derivative
$i(\bar{n}\partial)$ in (C.11). Using QCD EOM

$\left\langle 0\right|\left[i(\bar{n}\partial)\bar{\xi}(x_{-})\right]\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\xi(0)\left|V(k)\right\rangle=\left\langle
0\right|\left[(i\bar{n}\partial)^{-1}\partial_{\bot}^{2}\bar{\xi}(x_{-})\right]\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\xi(y_{-})\left|V(k)\right\rangle+\ldots$ (C.13)
$=f_{V}m_{V}^{3}\frac{e_{+}}{k_{+}}~{}\frac{1}{4}\text{FT}\left[\frac{1}{\bar{u}}\mathbf{A}(u)\right]+\ldots\
,$ (C.14)

where dots denote the quark-gluon contributions. Similarly

$\left\langle 0\right|\left[\bar{\xi}(x_{-})\right]\hbox to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\left[i(\bar{n}\partial)\xi(0)\right]\left|V(k)\right\rangle=f_{V}~{}\frac{m_{V}^{3}}{4}\frac{e_{+}}{k_{+}}~{}\text{FT}\left[\frac{1}{u}\mathbf{A}(u)\right]+\ldots\
.$ (C.15)

The sum of these terms gives the operator m.e. of the operator
$i(\bar{n}\partial)\left[\bar{\xi}(x_{-})\hbox to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\xi(0)\right]$ with the total derivative, which
can be computed from Eq.(C.9). This gives the following relation

$\mathbf{A}(u)=4u\bar{u}\phi(u)+\bar{G}_{1}(u),$ (C.16)

where $\bar{G}_{1}(u)$ again schematically denotes the contributions of quark-
gluon DAs.

Consider now the m.e.

$\displaystyle\left\langle 0\right|\bar{\eta}(x_{-})\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\eta(0)\left|V(k)\right\rangle=\left\langle
0\right|\left[(in\partial)^{-1}\partial_{\bot}^{\alpha}\bar{\xi}(x_{-})\right]\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\left[(in\partial)_{\bot}^{-1}\partial_{\bot}^{\alpha}\xi(y_{-})\right]\left|k\right\rangle+\dots$
$\displaystyle=-\left\langle
0\right|\left[(in\partial)^{-1}\partial_{\bot}^{2}\bar{\xi}(x_{-})\right]\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\left[(in\partial)^{-1}\xi(y_{-})\right]\left|k\right\rangle+\dots$
$\displaystyle=-f_{V}m_{V}\frac{m_{V}^{2}}{4}\frac{e_{+}}{k_{+}^{2}}~{}\text{FT}\left[\frac{1}{u\bar{u}}\mathbf{A}(u)\right]+\dots,$
(C.17)

where we show the two-particle quark-antiquark operator only. On the other
hand

$\left\langle 0\right|\bar{\eta}(x_{-})\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\eta(y_{-})\left|V(k)\right\rangle=\left\langle
0\right|\bar{\psi}(x_{-})\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\psi(y_{-})\left|k\right\rangle=-f_{V}m_{V}^{3}\frac{e_{+}}{k_{+}^{2}}\text{FT}\left[g_{3}\right].$
(C.18)

Therefore this gives

$g_{3}=\frac{1}{4}\frac{\mathbf{A}(u)}{u\bar{u}}+\bar{G}_{3}(u),$ (C.19)

where $\bar{G}_{3}$ denotes a contribution from quark-gluon DAs. Using (C.16)
one obtains

$g_{3}(u)=\phi(u)+\frac{\bar{G}_{1}(u)}{4u\bar{u}}+\bar{G}_{3}(u).$ (C.20)

Neglecting the quark-gluon terms $\bar{G}_{1}\sim\bar{G}_{3}\rightarrow 0$ one
finds the result

$g_{3}(u)\simeq\phi(u),$ (C.21)

that agrees with Eq.(C.5) assumming $\bar{g}_{3}\rightarrow 0.$ At the same
time the results for $\mathbf{A}$ in Eqs.(C.10) and (C.16) do not agree if one
simply neglects the quark-gluon contributions

$4u\bar{u}\phi(u)\neq\left.32\int_{0}^{u}dv\int_{0}^{v}d\omega(g_{v}-\phi)(\omega)\right|_{\bar{g}_{3}=\bar{g}_{v}=0}.$
(C.22)

This result can be explained by the fact that the neglected quark-gluon
contributions still implicitly depend on the twist-2 moments. For the local
matrix elements this can be seen from the QCD equation of motions, see e.g.
Ref.[41]. Therefore the simple recipe prescribing just to neglect the
quark–gluon contributions, as for the twist-3 operators, is not applicable for
the twist-4 calculations. Notice however that Eq.(C.22) is valid for one
particular case: for the asymptotic shape $\phi^{as}(u)=6u\bar{u}$. The
corresponding local operator can not be related with any quark-gluon one.

This observation shows that consistent twist-4 calculation, which only account
the two-particle quark operator requires some prescription, which makes it
possible to unambiguously and systematically distinguish between the twist-4
quark and the quark-gluon contribution.

In fact, the Lorentz symmetry of the CF (C.9) implies that DA $\mathbf{A}$
must satisfy to a certain condition. Let us briefly discuss this point.
Consider the following matrix element

$\left\langle 0\right|\bar{\psi}(0)\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\psi(0)\left|V(k)\right\rangle=-~{}f_{V}m_{V}^{3}~{}\frac{e_{+}}{k_{+}^{2}},$
(C.23)

where we used Eq.(C.1). Now the operator in lhs can be rewritten as

$\displaystyle\bar{\psi}(0)\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\psi(0)$
$\displaystyle=\bar{\eta}(0)\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\eta(0)={i^{2}}\bar{\xi}\overleftarrow{\hbox
to0.0pt{\hbox
to8.55695pt{\hfil/\hfil}\hss}D}_{\bot}(i\overleftarrow{D}n)^{-1}\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n(in\overrightarrow{D})^{-1}\overrightarrow{\hbox
to0.0pt{\hbox to8.55695pt{\hfil/\hfil}\hss}D}_{\bot}\xi$ (C.24)
$\displaystyle\simeq\left[(in\partial)^{-1}i{\partial}_{\alpha\bot}\bar{\xi}\,\right]\gamma_{\bot}^{\alpha}\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\gamma_{\bot}^{\beta}\left[(in\partial)^{-1}i{\partial}_{\bot\beta}\xi\right],$
(C.25)

where we again skip the quark-gluon operators. The matrix element of this
operator gives

$\left\langle
0\right|\left[(in\partial)^{-1}i{\partial}_{\alpha\bot}\bar{\xi}\,\right]\gamma_{\bot}^{\alpha}\hbox
to0.0pt{\hbox
to6.00235pt{\hfil/\hfil}\hss}n\gamma_{\bot}^{\beta}\left[(in\partial)^{-1}i{\partial}_{\bot\beta}\xi\right]\left|V(k)\right\rangle\simeq-\frac{e_{+}}{k_{+}^{2}}~{}f_{V}m_{V}^{3}\int_{0}^{1}du~{}\frac{\mathbf{A}_{\bar{q}q}(u)}{4u\bar{u}},$
(C.26)

where $\mathbf{A}_{\bar{q}q}$ denotes the pure quark part of the DA
$\mathbf{A}$

$\mathbf{A}=\mathbf{A}_{\bar{q}q}+\bar{\mathbf{A}}.$ (C.27)

Therefore one finds

$\left\langle 0\right|\bar{\psi}(0)\hbox to0.0pt{\hbox
to5.00002pt{\hfil$\bar{n}$\hfil}\hss}/\psi(0)\left|V(k)\right\rangle=-\frac{e_{+}}{k_{+}^{2}}~{}f_{V}m_{V}^{3}\int_{0}^{1}du~{}\frac{\mathbf{A}_{\bar{q}q}(u)}{4u\bar{u}}+\bar{G}_{4},$
(C.28)

where $\bar{G}_{4}$ schematically denotes the quark-gluon contributions.
Comparing with Eq.(C.23) one finds the following relation

$\int_{0}^{1}du~{}\frac{\mathbf{A}_{\bar{q}q}(u)}{4u\bar{u}}+G_{4}=1,$ (C.29)

One can assume that consistently defined quark contribution
$\mathbf{A}_{\bar{q}q}$ must satisfy

$\int_{0}^{1}du~{}\frac{\mathbf{A}_{\bar{q}q}(u)}{4u\bar{u}}=1,$ (C.30)

that ensures the correct Lorentz symmetry properties in the local limit, even
if the quark-gluon contributions are not taken into account, see Eq(C.7). For
instance, neglecting the quark-gluon contributions in Eq.(C.16) one finds the
expression

$\mathbf{A}_{\bar{q}q}(u)=4u\bar{u}\phi(u),$ (C.31)

which satisfy to Eq.(C.30).

In conclusion, the given consideration demonstrates that in general case the
definition of the pure quark contribution for the twist-4 DAs is not universal
due to the implicit contributions from the matrix element of the quark-gluon
operators. In the given example the universal part can only be defined in
terms of the asymptotic part of the twist-2 DA $\phi$ because the
corresponding local operator do not appear from the quark-gluon ones.

The similar arguments also applicable for the three-quark operators. Therefore
in our calculations we only consider the asymptotic contributions of the
twist-3 and twist-4 3-quark matrix elements in order to define the
corresponding universal three-quark components of twist-5 DAs.

##  Appendix D  Isospin relations between the auxiliary twist-5 DAs

The standard analysis of the light-cone matrix elements with two transverse
derivatives allows one to obtain many relations between the auxiliary DAs
defined in the Sec.2.2. Here we consider such relations for a nucleon state.
For simplicity, in this section we omit the baryon subscript “B=N” assuming
for all DAs $F^{N}\equiv F$. Let us define the operator basis as following

$\displaystyle\epsilon^{ijk}\left[\partial_{\bot}^{\alpha}u^{i}_{1}(x_{-})\right]\left[\partial_{\bot}^{\beta}u^{j}_{2}(y_{-})\right]\
d^{k}_{3}(z_{-}),$ (D.1)
$\displaystyle\epsilon^{ijk}\left[\partial_{\bot}^{\alpha}\partial_{\bot}^{\beta}u^{i}_{1}(x_{-})\right]u^{j}_{2}(y_{-})\
d^{k}_{3}(z_{-}),$ (D.2)
$\displaystyle\epsilon^{ijk}u^{i}_{1}(x_{-})\left[\partial_{\bot}^{\alpha}\partial_{\bot}^{\beta}u^{j}_{2}(y_{-})\right]\
d^{k}_{3}(z_{-}).$ (D.3)

In these formulas we assume that each quark field $q={u,d}$ is projected onto
the large component $\xi$ as defined in Eq.(87). Recall, that we use the short
notation for the Dirac indices $q_{\alpha_{i}}\equiv q_{i}$. From the
technical point of view the derivation of such relations is the same as in
Ref.[21] . This method gives the following three relations for the light-cone
matrix elements of defined operators

$\displaystyle\left\langle
0\right|\left[\partial_{\bot}^{\alpha}u_{1}(x_{-})\right]\left[\partial_{\bot}^{\beta}u_{2}(y_{-})\right]\
d_{3}(z_{-})\left|k\right\rangle+\left\langle
0\right|\left[\partial_{\bot}^{\alpha}u_{1}(x_{-})\right]\ u_{3}(z_{-})\
\left[\partial_{\bot}^{\beta}d_{2}(y_{-})\right]\left|k\right\rangle$
$\displaystyle+\left\langle 0\right|u_{3}(z_{-})\
\left[\partial_{\bot}^{\beta}u_{2}(y_{-})\right]\
\left[\partial_{\bot}^{\alpha}d_{1}(x_{-})\right]\left|k\right\rangle=0.$
(D.4) $\displaystyle\left\langle
0\right|\left[\partial_{\bot}^{\alpha}\partial_{\bot}^{\beta}u_{1}(x_{-})\right]u_{2}(y_{-})\
d_{3}(z_{-})\left|k\right\rangle+\left\langle
0\right|\left[\partial_{\bot}^{\alpha}\partial_{\bot}^{\beta}u_{1}(x_{-})\right]\
u_{3}(z_{-})\ \left[d_{2}(y_{-})\right]\left|k\right\rangle$
$\displaystyle+\left\langle 0\right|u_{3}(z_{-})\ u_{2}(y_{-})\
\left[\partial_{\bot}^{\alpha}\partial_{\bot}^{\beta}d_{1}(x_{-})\right]\left|k\right\rangle=0,$
(D.5) $\displaystyle\left\langle
0\right|u_{1}(y_{-})\left[\partial_{\bot}^{\alpha}\partial_{\bot}^{\beta}u_{2}(x_{-})\right]\
d_{3}(z_{-})\left|k\right\rangle+\left\langle
0\right|u_{1}(x_{-})u_{3}(z_{-})\
\left[\partial_{\bot}^{\alpha}\partial_{\bot}^{\beta}d_{2}(y_{-})\right]\left|k\right\rangle$
$\displaystyle+\left\langle 0\right|u_{3}(z_{-})\
\left[\partial_{\bot}^{\alpha}\partial_{\bot}^{\beta}u_{2}(y_{-})\right]d_{1}(x_{-})\left|k\right\rangle=0.$
(D.6)

Substituting in these equations the definitions of the auxiliary matrix
elements from Sec.2.2 and using the Fierz identities one derives various
relations between the auxiliary DAs. It is convenient to divide the operators
$O_{i}^{\alpha\beta}$ in Eqs.(D.1)-(D.3) into three independent groups:
symmetric ($O_{i}^{\alpha\beta}=g_{\bot}^{\alpha\beta}O_{i}$), symmetric
traceless $g_{\bot}^{\alpha\beta}O_{i}^{\alpha\beta}$ and and antisymmetric
$O_{1}^{\beta\alpha}=-O_{1}^{\alpha\beta}$. Each of these groups can be
considered independently.

For simplicity, in what follow we use the following convenient notations

$\displaystyle F(x_{1},x_{2},x_{3})\equiv
F,\,F(x_{2},x_{1},x_{3})=\hat{P}_{12}F,\,F(x_{3},x_{2},x_{1})=\hat{P}_{13}F\,\text{
and so on}.$ (D.7)

For the DAs, which describe the symmetric traceless matrix elements the set of
relations read

$\displaystyle T_{2xy}\
-\hat{P}_{23}\left(T_{2xy}+T_{2xx}\right)-\hat{P}_{13}\left(T_{2xy}+T_{2yy}\right)\
=0.$ (D.8) $\displaystyle
T_{2xx}+\hat{P}_{23}T_{2xx}+\hat{P}_{13}\left(T_{2xx}+T_{2yy}+2T_{2xy}\right)\
=0.$ (D.9) $\displaystyle
T_{2yy}+\hat{P}_{13}T_{2yy}+\hat{P}_{23}\left(T_{2xx}+T_{2yy}+2T_{2xy}\right)=0.$
(D.10)

For the DAs, which describe the antisymmetric matrix elements the set of
relations read

$\displaystyle
T_{xy}+G_{xy}-\hat{P}_{23}\left(T_{xy}+G_{xy}\right)+\hat{P}_{13}\,\tilde{T}_{2xy}=0,$
(D.11) $\displaystyle
T_{xy}-G_{xy}+\hat{P}_{23}\,\tilde{T}_{2xy}-\hat{P}_{13}\left(T_{xy}-G_{xy}\right)=0,$
(D.12)
$\displaystyle\tilde{T}_{2xy}+\hat{P}_{23}\left(T_{xy}-G_{xy}\right)+\hat{P}_{13}\left(T_{xy}+G_{xy}\right)=0.$
(D.13)

The relations for the symmetric matrix elements are described by the matrix
elements of the operators $\left[\partial_{\bot}^{2}u\right]ud$,
$u\left[\partial_{\bot}^{2}u\right]d$ and
$\left[\partial_{\bot}^{\alpha}u\right]\left[\partial_{\bot}^{\alpha}u\right]d$.
The analysis of the operator $\left[\partial_{\bot}^{2}u\right]ud$ gives

$\displaystyle
V_{1xx}+A_{1xx}+\hat{P}_{23}\left(V_{1xx}+A_{1xx}\right)-\hat{P}_{13}\left(T_{0xx}+T_{0yy}+2T_{0xy}\right)=0,$
(D.14) $\displaystyle V_{1xx}-A_{1xx}-\hat{P}_{23}T_{0xx}+\hat{P}_{13}\left(\
V_{1xx}+V_{1yy}+2V_{1xy}-A_{1xx}-A_{1yy}-2A_{1xy}\right)=0,$ (D.15)
$\displaystyle\hat{P}_{23}\left(V_{1xx}-A_{1xx}\right)+\hat{P}_{13}\left(V_{1xx}+V_{1yy}+2V_{1xy}+A_{1xx}+A_{1yy}+2A_{1xy}\right)-T_{0xx}=0.$
(D.16)

The analysis of the operator $u\left[\partial_{\bot}^{2}u\right]d$ gives

$\displaystyle
V_{1yy}+A_{1yy}+\hat{P}_{23}\left(V_{1xx}+V_{1yy}+2V_{1xy}+A_{1xx}+A_{1yy}+2A_{1xy}\right)-\hat{P}_{13}T_{0yy}=0.$
(D.17) $\displaystyle
V_{1yy}-A_{1yy}+\hat{P}_{13}\left(V_{1yy}-A_{1yy}\right)-\hat{P}_{23}\left(T_{0xx}+T_{0yy}+2T_{0xy}\right)=0.$
(D.18)
$\displaystyle\hat{P}_{13}\left(V_{1yy}+A_{1yy}\right)+\hat{P}_{23}\left(V_{1xx}+V_{1yy}+2V_{1xy}-A_{1xx}-A_{1yy}-2A_{1xy}\right)-T_{0yy}=0.$
(D.19)

The analysis of the operator
$\left[\partial_{\bot}u\right]\left[\partial_{\bot}u\right]d$ gives

$\displaystyle
V_{1xy}+A_{1xy}-\hat{P}_{23}\left(V_{1xy}+V_{1xx}+A_{1xy}+A_{1xx}\right)+\hat{P}_{13}\left(T_{0xy}+T_{0yy}\right)=0,$
(D.20) $\displaystyle
V_{1xy}-A_{1xy}+\hat{P}_{23}\left(T_{0xy}+T_{0xx}\right)-\hat{P}_{13}\left(V_{1xy}+V_{1yy}-A_{1xy}-A_{1yy}\right)=0,$
(D.21) $\displaystyle