# Dispersion relation analysis of the radiative corrections to $g_{A}$ in the
neutron $\beta$-decay

Mikhail Gorchtein  Institut für Kernphysik, PRISMA Cluster of Excellence
Johannes Gutenberg-Universität, Mainz, Germany<EMAIL_ADDRESS>Chien-
Yeah Seng  Helmholtz-Institut für Strahlen- und Kernphysik and Bethe Center
for Theoretical Physics,
Universität Bonn, 53115 Bonn, Germany

###### Abstract

We present the first and complete dispersion relation analysis of the inner
radiative corrections to the axial coupling constant $g_{A}$ in the neutron
$\beta$-decay. Using experimental inputs from the elastic form factors and the
spin-dependent structure function $g_{1}$, we determine the contribution from
the $\gamma W$-box diagram to a precision better than $10^{-4}$. Our
calculation indicates that the inner radiative corrections to the Fermi and
the Gamow-Teller matrix element in the neutron $\beta$-decay are almost
identical, i.e. the ratio $\lambda=g_{A}/g_{V}$ is almost unrenormalized. With
this result, we predict the bare axial coupling constant to be
$\mathring{g}_{A}=-1.2754(13)_{\mathrm{exp}}(2)_{\mathrm{RC}}$ based on the
PDG average $\lambda=-1.2756(13)$.

##  I Introduction

The recent emergence of an apparent deficit in the top-row Cabibbo-Kobayashi-
Maskawa (CKM) matrix unitarity Zyla et al. (2020),

$|V_{ud}|^{2}+|V_{us}|^{2}+|V_{ub}|^{2}=0.9985(5)$ (1)

has triggered a renewed interest in precise experimental studies of various
$\beta$-decay processes giving access to $|V_{ud}|$. Superallowed
$0^{+}-0^{+}$ nuclear decays have long been regarded as the best avenue for
such purpose. Recent works, however, pointed out that current theory
uncertainties in the nuclear-structure corrections may have been significantly
underestimated Seng et al. (2018, 2019); Gorchtein (2019); Hardy and Towner
(2020). Reducing these uncertainties requires novel ab-initio nuclear theory
calculations that are not yet available. As a consequence, the role of
alternative channels such as the $\beta$-decays of the free neutron, mirror
nuclei and pion becomes increasingly important. With the future improvement in
the experimental precision, these other $\beta$ decay processes will offer
competitive determinations of $|V_{ud}|$ and complementary sensitivity to
possible beyond standard model (BSM) signals.

Free neutron $\beta$-decay currently provides the second best determination of
$|V_{ud}|$ through the following master formula Czarnecki et al. (2004); Seng
et al. (2020a):

$|V_{ud}|^{2}=\frac{4903.1(1.1)~{}\mathrm{s}}{\tau_{n}(1+3\lambda^{2})}~{},$
(2)

where the uncertainty in the numerator arises from the Standard Model (SM)
theory input. The two required experimental inputs are the neutron lifetime
$\tau_{n}$, and the decay parameter $\lambda\equiv g_{A}/g_{V}$ which is the
ratio between the neutron axial and vector coupling constant. This parameter
is renormalized by electroweak radiative corrections (RCs), and these latter
are the primary focus of this article.

The parameter $\lambda$ can be measured either via the P-even correlation
$\vec{p}_{e}\cdot\vec{p}_{\nu}$ (the $a$ coefficient), or the P-odd
$\hat{e}_{s}\cdot\vec{p}_{e}$ (the $A$ coefficient) and
$\hat{e}_{s}\cdot\vec{p}_{\nu}$ (the $B$ coefficient) ones, with $\hat{e}_{s}$
the unit 3-vector along the neutron polarization. The current best measurement
reported by the PERKEO III collaboration
$\lambda=-1.27641(45)_{\mathrm{stat}}(33)_{\mathrm{sys}}$, with a 0.04%
precision Märkisch et al. (2019) (in this paper we pick the sign convention
$\lambda<0$, also adopted in the Particle Data Group (PDG) review). However,
the current PDG average reads $\lambda=-1.2756(13)$ Zyla et al. (2020), where
the much larger uncertainty is due to a scale factor of 2.6 that accounts for
the large discrepancy between the results before Bopp et al. (1986);
Erozolimsky et al. (1997); Liaud et al. (1997); Mostovoi et al. (2001) and
after 2002 Schumann et al. (2008); Mund et al. (2013); Darius et al. (2017);
Brown et al. (2018); Märkisch et al. (2019) (see Ref.Czarnecki et al. (2018)
for more discussions). Future improvements are expected from the Nab Fry et
al. (2019) and PERC Dubbers et al. (2008); Wang et al. (2019) collaborations,
both aiming at an accuracy level of $10^{-4}$.

The exact value of $\lambda$ not only serves for extracting $V_{ud}$, but is
also interesting in itself. The “bare” (i.e. without electroweak corrections)
neutron axial coupling $\mathring{g}_{A}$ is one of the simplest hadronic
matrix elements and has received much attention. Unlike its vector counterpart
$\mathring{g}_{V}$ which remains non-renormalized due to the conserved vector
current (CVC), the bare axial coupling is not protected and must be
calculated, e.g. using lattice Quantum Chromodynamics (QCD) Khan et al.
(2006); Lin et al. (2008); Capitani et al. (2012); Horsley et al. (2014); Bali
et al. (2015); Abdel-Rehim et al. (2015); Alexandrou et al. (2017); Capitani
et al. (2019); Edwards et al. (2006); Yamazaki et al. (2008, 2009); Bratt et
al. (2010); Green et al. (2014); Yamanaka et al. (2018); Liang et al. (2018);
Ishikawa et al. (2018); Ottnad et al. (2018); Bhattacharya et al. (2016);
Berkowitz et al. (2017); Chang et al. (2018); Gupta et al. (2018); Walker-Loud
et al. (2020). The most recent FLAG average Aoki et al. (2019) reads:

$\displaystyle N_{f}=2+1+1:$ $\displaystyle\mathring{g}_{A}=-1.251(33)$
$\displaystyle N_{f}=2+1:$ $\displaystyle\mathring{g}_{A}=-1.254(16)(30)$
$\displaystyle N_{f}=2:$ $\displaystyle\mathring{g}_{A}=-1.278(86)~{},$ (3)

but individual calculations have achieved higher precision. For instance,
Ref.Chang et al. (2018) reported a percent-level determination of
$\mathring{g}_{A}=-1.271(10)(7)$ using an unconventional method inspired by
the Feynman-Hellmann theorem, and follow-up works are aiming for sub-percent
precision Walker-Loud et al. (2020). Such a rapid development makes $\lambda$
a powerful tool for searching for redBSM physics. By comparing first-
principles calculations of $\mathring{g}_{A}$ to the experimental results for
$\lambda$ one thereby constraints the strength of possible BSM contributions
that could modify $g_{A}$, in particular the right-handed currents González-
Alonso and Martin Camalich (2016); Alioli et al. (2017); Gonzalez-Alonso et
al. (2019); Falkowski et al. (2020).

When the lattice precision reaches $10^{-3}$, a valid comparison between
$\mathring{g}_{A}$ and $\lambda$ will require precise
$\mathcal{O}(\alpha_{em}/\pi)$ RCs that bring
$\mathring{g}_{A}/\mathring{g}_{V}$ to $g_{A}/g_{V}$. In particular, we need
to deal with sizable hadronic uncertainties originating from the $\gamma
W$-box diagram. The latter can be written as a $Q^{2}$-integral, and
performing operator product expansion (OPE) it is well-known that the
corrections to $g_{V}$ and $g_{A}$ coming from large $Q^{2}$ (which carries a
large electroweak logarithm) are the same. Therefore, it was believed that the
difference between $\mathring{g}_{A}$ and $\lambda$ is numerically small
Sirlin (1967); Garcia and Queijeiro (1983); Kurylov et al. (2002, 2003).
However, for a long time there was no serious attempt to understand the RC to
$g_{A}$ from the low-$Q^{2}$ part of the integral which is also of the order
$10^{-3}$, comparable to the high-$Q^{2}$ contribution. It includes the
elastic contributions that are fixed by the nucleon form factors, as well as
the inelastic contributions that are governed by non-perturbative QCD. The
first attempt for a complete analysis was performed recently in Refs.Hayen
(2021a, b). Anticipating the results of this work, we found that those Refs.
originally contained algebraic mistakes in the computation of the elastic
contribution, invalidating their numerical results. These mistakes were later
corrected in the published version of Ref.Hayen (2021a). Additionally, the
inelastic contribution residing at low $Q^{2}$ was obtained based on a
holographic QCD model, following Ref.Czarnecki et al. (2019) where the RC to
$g_{V}$ was addressed. The model-dependent nature of this approach makes a
rigorous estimation of the theoretical uncertainty complicated.

In this article we improve on both points. We perform a novel analysis of the
RC to $g_{A}$ based on the dispersion relation (DR) approach. It is a powerful
tool which has proved successful in the treatment of the RCs to the Fermi
amplitude in the free neutron and superallowed $\beta$-decays Seng et al.
(2018, 2019); Gorchtein (2019); Seng et al. (2020a). In this formalism, the
$\gamma W$-box diagram is expressed as a dispersion integral over structure
functions that are directly or indirectly related to experimental data. This
ultimately allows for a fully data-driven analysis of this RC. For $g_{A}$,
the required input relies on the spin-dependent structure functions $g_{1}$
and $g_{2}$, well-studied quantities in deep inelastic scattering (DIS)
experiments. We utilize high-precision world data on $g_{1,2}$ to evaluate the
dispersion integral, and fix the forward $\gamma W$-box diagram correction to
$g_{A}$ to an unprecedented precision better than $10^{-4}$. We observe that
the RCs to $g_{V}$ and $g_{A}$ are numerically very close and largely cancel
in the ratio, which practically removes any distinction between
$\mathring{g}_{A}$ and $\lambda$ down to $2\times 10^{-4}$.

The contents in this paper are arranged as follows. In Section II we define
our notation and introduce the starting point for the discussion of the RC. We
introduce the $\gamma W$-box diagram in Section III, and derive its dispersive
representation in Section IV. The elastic (Born) and inelastic contributions
to the box diagram are computed in Section V and VI respectively. The final
results and discussions are presented in Section VII.

##  II General framework

We start by defining the hadronic currents relevant to the $\beta$-decay of
the free neutron:

$\displaystyle J_{em}^{\mu}$ $\displaystyle=$
$\displaystyle\frac{2}{3}\bar{u}\gamma^{\mu}u-\frac{1}{3}\bar{d}\gamma^{\mu}d$
$\displaystyle J_{W}^{\mu}$ $\displaystyle=$
$\displaystyle\bar{u}\gamma^{\mu}(1-\gamma_{5})d~{}.$ (4)

Their single-nucleon matrix elements are given by:

$\displaystyle\langle N(p_{f},s_{f})|J_{em}^{\mu}|N(p_{i},s_{i})\rangle$
$\displaystyle=$
$\displaystyle\bar{u}_{s_{f}}(p_{f})\left[F_{1}^{N}\gamma^{\mu}+\frac{i}{2M}F_{2}^{N}\sigma^{\mu\nu}(p_{f}-p_{i})_{\nu}\right]u_{s_{i}}(p_{i})$
$\displaystyle\langle p(p_{f},s_{f})|J_{W}^{\mu}|n(p_{i},s_{i})\rangle$
$\displaystyle=$
$\displaystyle\bar{u}_{s_{f}}(p_{f})\left[F_{1}^{W}\gamma^{\mu}+\frac{i}{2M}F_{2}^{W}\sigma^{\mu\nu}(p_{f}-p_{i})_{\nu}+G_{A}\gamma^{\mu}\gamma_{5}-\frac{G_{P}}{2M}\gamma_{5}(p_{f}-p_{i})^{\mu}\right]u_{s_{i}}(p_{i})~{},$
(5)

where $N=p,n$ and $M=(M_{p}+M_{n})/2\approx 939$ MeV. All the form factors
above are functions of $-(p_{i}-p_{f})^{2}$. We may also define the isospin
combinations $F_{1,2}^{S}=F_{1,2}^{p}+F_{1,2}^{n}$ and
$F_{1,2}^{V}=F_{1,2}^{p}-F_{1,2}^{n}$. The values of the vector and the axial
charged weak form factors at zero momentum transfer define the “bare” vector
and axial coupling constants: $F_{1}^{W}(0)=\mathring{g}_{V}$,
$G_{A}(0)=\mathring{g}_{A}$, which represent the Fermi and Gamow-Teller matrix
element in neutron $\beta$-decay respectively. In particular,
$\mathring{g}_{V}=1$ from isospin symmetry, and the correction due to the
strong isospin-breaking effects is negligible due to the Behrends-Sirlin-
Ademollo-Gatto theorem Behrends and Sirlin (1960); Ademollo and Gatto (1964).
On the other hand, $\mathring{g}_{A}$ is not protected by any exact symmetry.
We do not include the isospin-breaking correction to $\mathring{g}_{A}$
separately because it is already included in the respective first-principles
calculations.

The nucleon mass difference $\Delta=M_{n}-M_{p}\approx 1.3$ MeV and the
electron mass $m_{e}\approx 0.511$ MeV are much smaller than $M$. Therefore,
the tree-level amplitude of the decay process $n(p_{n})\rightarrow
p(p_{p})e(p_{e})\bar{\nu}_{e}(p_{\nu})$ is given by:

$\mathcal{M}_{\mathrm{tree}}=-\frac{G_{F}}{\sqrt{2}}L_{\lambda}\bar{u}_{s^{\prime}}(p)\gamma^{\mu}(\mathring{g}_{V}+\mathring{g}_{A}\gamma_{5})u_{s}(p)+\mathcal{O}(\Delta^{2})~{},$
(6)

where $G_{F}=1.1663787(6)\times 10^{-5}$ GeV-2 is the Fermi constant measured
from the muon decay Zyla et al. (2020), $p=(p_{n}+p_{p})/2$ is the average
nucleon momentum, and
$L_{\lambda}=\bar{u}(p_{e})\gamma_{\lambda}(1-\gamma_{5})v(p_{\nu})$ is the
lepton piece. The recoil corrections scale as $\Delta/M\sim 10^{-3}$, which
are small but important in precision physics. They were studied in detail with
both conventional methods and effective field theory (EFT) Holstein (1974);
Wilkinson (1982); Ando et al. (2004); Gudkov (2008); Ivanov et al. (2013,
2021), and will not be discussed here.

RCs of the order $\mathcal{O}(\alpha_{em}/\pi)$ must be included for a precise
extraction of the weak coupling parameters. In the usual nomenclature, they
are divided into the “outer” and “inner” corrections The former is a function
of $\\{E_{e},m_{e}\\}$ calculable within Quantum Electrodynamics (QED) and
independent of details of strong interaction. The latter is instead a constant
in $E_{e}$ but depends on details of the hadronic structure. The squared
amplitude for the decay of a polarized neutron (to unpolarized final states)
after the inclusion of the $\mathcal{O}(\alpha_{em}/\pi)$ RCs reads:

$\displaystyle|\mathcal{M}|^{2}$ $\displaystyle=$ $\displaystyle
16G_{F}^{2}|V_{ud}|^{2}M_{n}M_{p}E_{e}(E_{m}-E_{e})g_{V}^{2}(1+3\lambda^{2})F(\beta)\left(1+\frac{\alpha_{em}}{2\pi}\delta^{(1)}\right)\left\\{1+\left(1+\frac{\alpha_{em}}{2\pi}\delta^{(2)}\right)a_{0}\frac{\vec{p}_{e}\cdot\vec{p}_{\nu}}{E_{e}E_{\nu}}\right.$
(7)
$\displaystyle\left.+\hat{e}_{s}\cdot\left[\left(1+\frac{\alpha_{em}}{2\pi}\delta^{(2)}\right)A_{0}\frac{\vec{p}_{e}}{E_{e}}+B_{0}\frac{\vec{p}_{\nu}}{E_{\nu}}\right]\right\\}+\mathcal{O}(\Delta^{3})~{},$

with

$a_{0}=\frac{1-\lambda^{2}}{1+3\lambda^{2}}~{},\quad
A_{0}=-\frac{2\lambda(\lambda+1)}{1+3\lambda^{2}}~{},\quad
B_{0}=\frac{2\lambda(\lambda-1)}{1+3\lambda^{2}}~{}.$ (8)

Here, $E_{e}=p_{e}\cdot p/M$ is the electron energy,
$E_{m}=(M_{n}^{2}-M_{p}^{2}+m_{e}^{2})/(2M_{n})\approx M_{n}-M_{p}$ is the
electron end-point energy, and $\beta=\sqrt{1-m_{e}^{2}M^{2}/(p_{e}\cdot
p)^{2}}$ is the electron speed, all in the nucleon’s rest frame. With these
notations, the functions $\delta^{(1,2)}$ describe the outer corrections
Sirlin (1967); Garcia and Maya (1981):

$\displaystyle\delta^{(1)}$ $\displaystyle=$ $\displaystyle
3\ln\frac{M_{p}}{m_{e}}-\frac{3}{4}+4\left(\frac{1}{\beta}\tanh^{-1}\beta-1\right)\left({\ln\frac{2(E_{\mathrm{m}}-E_{e})}{m_{e}}}+\frac{E_{m}-E_{e}}{3E_{e}}-\frac{3}{2}\right)$
$\displaystyle-\frac{4}{\beta}\mathrm{Li}_{2}\left(\frac{2\beta}{1+\beta}\right)+\frac{1}{\beta}\tanh^{-1}\beta\left(2+2\beta^{2}+\frac{(E_{m}-E_{e})^{2}}{6E_{e}^{2}}-4\tanh^{-1}\beta\right)$
$\displaystyle\delta^{(2)}$ $\displaystyle=$ $\displaystyle
2\left(\frac{1-\beta^{2}}{\beta}\right)\tanh^{-1}\beta+\frac{4(E_{m}-E_{e})(1-\beta^{2})}{3\beta^{2}E_{e}}\left(\frac{1}{\beta}\tanh^{-1}\beta-1\right)$
(9)
$\displaystyle+\frac{(E_{m}-E_{e})^{2}}{6\beta^{2}E_{e}^{2}}\left(\frac{1-\beta^{2}}{\beta}\tanh^{-1}\beta-1\right)~{},$

whereas $F(\beta)\approx 1+\alpha_{em}\pi/\beta$ is the Fermi’s function that
incorporates the Coulomb interaction between the final-state proton and the
electron Fermi (1934). The function $\delta^{(1)}$ is also known as Sirlin’s
function $g(E_{e},E_{m})$.

The axial to vector coupling constants’ ratio, parameter
$\lambda=g_{A}/g_{V}$, is understood as fully renormalized by the inner RCs.
In near-degenerate semileptonic $\beta$-decay processes, the inner RCs are
most conveniently studied in Sirlin’s representation Sirlin (1978) (see also
Refs.Seng et al. (2020b); Feng et al. (2020); Seng et al. (2020c) for a
detailed account). In this formalism, most of the
$\mathcal{O}(\alpha_{em}/\pi)$ electroweak RCs are either exactly known from
current algebra, or give rise to the outer corrections in Eq.(9) and the
Fermi’s function. As a result, the renormalized vector and axial coupling
constants read:

$\displaystyle g_{V}$ $\displaystyle=$
$\displaystyle\mathring{g}_{V}\left\\{1+\frac{\alpha_{em}}{4\pi}\left[3\ln\frac{M_{Z}}{M_{p}}+\ln\frac{M_{Z}}{M_{W}}+\tilde{a}_{g}\right]+\frac{1}{2}\delta_{\mathrm{HO}}^{\mathrm{QED}}+\Box_{\gamma
W}^{V}\right\\}$ $\displaystyle g_{A}$ $\displaystyle=$
$\displaystyle\mathring{g}_{A}\left\\{1+\frac{\alpha_{em}}{4\pi}\left[3\ln\frac{M_{Z}}{M_{p}}+\ln\frac{M_{Z}}{M_{W}}+\tilde{a}_{g}\right]+\frac{1}{2}\delta_{\mathrm{HO}}^{\mathrm{QED}}+\Box_{\gamma
W}^{A}\right\\}~{},$ (10)

where $\tilde{a}_{g}$ is a pQCD correction factor and
$\delta_{\mathrm{HO}}^{\mathrm{QED}}$ summarizes the leading-log higher-order
QED effects Marciano and Sirlin (1993); Erler (2004). One observes that the
fractional corrections to $\mathring{g}_{V}$ and $\mathring{g}_{A}$ are mostly
identical and cancel in the ratio $g_{A}/g_{V}$. The only exceptions are the
constants $\Box_{\gamma W}^{V}$ and $\Box_{\gamma W}^{A}$ that describe the
inner RCs originated from the $\gamma W$-box diagrams (see Fig.1), which are
the focus of this paper111In the existing literature, e.g. Refs.Seng et al.
(2018, 2019); Feng et al. (2020), the quantity $\Box_{\gamma W}^{V}$ was
written as $\Box_{\gamma W}^{VA}$, where the superscript indicates that it
involves the product of a vector current and an axial current. In this paper,
the superscript carries a different meaning, namely which weak coupling
constant (vector or axial) they are correcting.. With the above, we obtain:

$\lambda=\frac{\mathring{g}_{A}}{\mathring{g}_{V}}\left[1+\Box_{\gamma
W}^{A}-\Box_{\gamma W}^{V}\right]~{}.$ (11)

##  III $\gamma W$-box diagram

Figure 1: The forward $\gamma W$-box diagrams in the free neutron
$\beta$-decay.

The $\gamma W$-box correction is of the natural size $\alpha_{em}/\pi\sim
10^{-3}$. Taking into account recoil corrections $\sim\Delta/M,\,m_{e}/M$ on
top of the overall $\alpha_{em}/\pi$ factor would bring us to accounting for
effects in the $10^{-6}$ range that exceed the precision goal by two or three
orders of magnitude. This defines the level of the detalization that is needed
in our analysis. We will consistently set $\Delta=m_{e}=0$ throughout the
calculation below, as well as the proton recoil. This approximation also leads
to the neglect of the pion pole due to the partially-conserved axial current
(PCAC) hypothesis: the pion pole contribution, when contracted with the lepton
tensor, results in lepton mass terms which, as stated are neglected. This
precision level is supported by the fact that the lowest hadronic state is
separated by the pion mass $\sim 140$ MeV which is about hundred times larger
than $\Delta$. Notice however that this approximation may not be as safe for
nuclear $\beta$-decay where the available $Q$-values may be as large as 15-20
MeV which are comparable to the energy level of nuclear excitations.

The part of the $\gamma W$-box diagram amplitude that contributes to the inner
correction must involve an antisymmetric tensor that stems from the lepton
spinor structure. It reads:

$\mathcal{M}_{\gamma
W}^{\mathrm{inner}}=-\frac{G_{F}}{\sqrt{2}}L_{\lambda}I_{\gamma
W}^{\lambda}~{},$ (12)

where

$I_{\gamma
W}^{\lambda}=\bar{u}_{s}(p)\gamma^{\lambda}\left[\mathring{g}_{V}\Box_{\gamma
W}^{V}+\mathring{g}_{A}\Box_{\gamma
W}^{A}\gamma_{5}\right]u_{s}(p)=ie^{2}\int\frac{d^{4}q}{(2\pi)^{4}}\frac{M_{W}^{2}}{M_{W}^{2}-q^{2}}\frac{\epsilon^{\mu\nu\alpha\lambda}q_{\alpha}}{(q^{2})^{2}}T_{\mu\nu}^{\gamma
W}$ (13)

with $\epsilon^{0123}=-1$ in our convention. The forward generalized Compton
tensor describing the $W^{+}n\rightarrow\gamma p$ process, is defined as:

$\displaystyle T^{\mu\nu}_{\gamma W}=\int dxe^{iqx}\langle
p|T[J^{\mu}_{em}(x)J^{\nu}_{W}(0)]|n\rangle.$ (14)

To extract $\Box_{\gamma W}^{V}$ and $\Box_{\gamma W}^{A}$, we use following
identities:

$\frac{1}{2}\bar{u}_{s}(p)\gamma^{\mu}u_{s}(p)=p^{\mu}~{},\>\>\frac{1}{2}\bar{u}_{s}(p)\gamma^{\mu}\gamma_{5}u_{s}(p)=S^{\mu}~{},$
(15)

where the spin vector $S^{\mu}$ is fixed by $S^{2}=-M^{2}$ and $S\cdot p=0$,
we obtain:

$\displaystyle\Box_{\gamma W}^{V}$ $\displaystyle=$
$\displaystyle\frac{ie^{2}}{2M^{2}\mathring{g}_{V}}\int\frac{d^{4}q}{(2\pi)^{4}}\frac{M_{W}^{2}}{M_{W}^{2}-q^{2}}\frac{\epsilon^{\mu\nu\alpha\lambda}q_{\alpha}p_{\lambda}}{(q^{2})^{2}}T_{\mu\nu}^{\gamma
W}$ $\displaystyle\Box_{\gamma W}^{A}$ $\displaystyle=$
$\displaystyle-\frac{ie^{2}}{2M^{2}\mathring{g}_{A}}\int\frac{d^{4}q}{(2\pi)^{4}}\frac{M_{W}^{2}}{M_{W}^{2}-q^{2}}\frac{\epsilon^{\mu\nu\alpha\lambda}q_{\alpha}S_{\lambda}}{(q^{2})^{2}}T_{\mu\nu}^{\gamma
W}~{}.$ (16)

In what follows, we use $\mathring{g}_{A}\approx\lambda=-1.2756$ as a
normalization in the second expression. Since $\Box_{\gamma W}^{A}\sim
10^{-3}$, the error induced by the ambiguity of $\mathring{g}_{A}$ is well
below our precision goal.

Only those components in $T_{\gamma W}^{\mu\nu}$ that contain an antisymmetric
tensor contribute to Eq.(16). These are

$\displaystyle T^{\mu\nu}_{\gamma W}$ $\displaystyle=$
$\displaystyle-\frac{i\epsilon^{\mu\nu\alpha\beta}q_{\alpha}p_{\beta}}{2(p\cdot
q)}T_{3}+\frac{i\epsilon^{\mu\nu\alpha\beta}q_{\alpha}}{(p\cdot
q)}\left[S_{\beta}S_{1}+\left(S_{\beta}-\frac{(S\cdot q)}{p\cdot
q}p_{\beta}\right)S_{2}\right]+...$ (17)

The spin-independent, parity-violating amplitude $T_{3}$ and spin-dependent,
parity-conserving amplitudes $S_{1,2}$ are functions of two invariants,
$\nu=(p\cdot q)/M$ and $Q^{2}=-q_{\mu}q^{\mu}$. Plugging Eq.(17) into Eq.(16)
gives:

$\displaystyle\Box_{\gamma W}^{V}$ $\displaystyle=$
$\displaystyle\frac{e^{2}}{2M\mathring{g}_{V}}\int\frac{d^{4}q}{(2\pi)^{4}}\frac{M_{W}^{2}}{M_{W}^{2}+Q^{2}}\frac{1}{(Q^{2})^{2}}\frac{\nu^{2}+Q^{2}}{\nu}T_{3}$
$\displaystyle\Box_{\gamma W}^{A}$ $\displaystyle=$
$\displaystyle\frac{e^{2}}{M\mathring{g}_{A}}\int\frac{d^{4}q}{(2\pi)^{4}}\frac{M_{W}^{2}}{M_{W}^{2}+Q^{2}}\frac{1}{(Q^{2})^{2}}\left\\{\frac{\nu^{2}-{2Q^{2}}}{3\nu}{S_{1}}-\frac{Q^{2}}{\nu}{S_{2}}\right\\}~{},$
(18)

where we have used the following identities:

$\displaystyle\int\frac{d^{4}q}{(2\pi)^{4}}q_{\alpha}F(\nu,Q^{2})$
$\displaystyle=$
$\displaystyle\int\frac{d^{4}q}{(2\pi)^{4}}\frac{\nu}{M}p_{\alpha}F(\nu,Q^{2})$
$\displaystyle\int\frac{d^{4}q}{(2\pi)^{4}}q_{\alpha}q_{\beta}F(\nu,Q^{2})$
$\displaystyle=$
$\displaystyle\int\frac{d^{4}q}{(2\pi)^{4}}\left[-\frac{\nu^{2}+Q^{2}}{3}g_{\alpha\beta}+\frac{4\nu^{2}+Q^{2}}{3}\frac{p_{\alpha}p_{\beta}}{M^{2}}\right]F(\nu,Q^{2})$
(19)

that hold for any Lorentz scalar function $F(\nu,Q^{2})$.

To evaluate the loop integrals, we need to discuss the symmetry properties of
the Compton amplitudes. We start by considering the isospin structure of
amplitudes $T_{i},\,S_{i}$. Electromagnetic interaction does not conserve
isospin and contains both isoscalar ($I=0$) and isovector ($I=1$) components.
Therefore, each amplitude $A=T_{i},\,S_{i}$ can be decomposed into components
contributed by the isoscalar and isovector electromagnetic current
respectively:

$\displaystyle A=A^{(0)}+A^{(1)}.$ (20)

The two isospin amplitudes have a different behavior under
$\nu\rightarrow-\nu$:

$\displaystyle A_{i}^{(I)}(-\nu,Q^{2})=\xi_{i}^{(I)}A_{i}^{(I)}(\nu,Q^{2}),$
(21)

with $\xi_{i}^{(I)}=\pm 1$ and $\xi_{i}^{(0)}=-\xi_{i}^{(1)}$. It is easy to
show that $\xi_{i}^{(0)}=-1$ for $T_{3}$ and $S_{1,2}$, so only the $I=0$
component of these amplitudes survives in the integrals in Eq.(18).

##  IV Dispersion representation of the forward Compton amplitudes

Figure 2: The singularities of the forward Compton amplitudes on the
complex-$\nu$ plane, and the contour chosen to derive the DRs in Eq.(24).

Forward Compton amplitudes have singularities along the real axis $\nu$: poles
due to a single nucleon intermediate state in the $s-$ and $u$-channels at
$\nu=\pm\nu_{B}=\pm Q^{2}/(2M)$, respectively, and unitarity cuts at
$\nu\geq\nu_{\pi}$ and $\nu\leq-\nu_{\pi}$ where
$\nu_{\pi}=(2Mm_{\pi}+m_{\pi}^{2}+Q^{2})/(2M)$ is the pion production
threshold (see Fig.2). The discontinuity of the forward Compton tensor in the
$s-$channel (i.e. $\nu\geq\nu_{B}$) is given by the generalization of the on-
shell hadronic tensor to the $\gamma W$-interference:

$\mathrm{Disc}T^{\gamma W}_{\mu\nu}(\nu)\equiv T^{\gamma
W}_{\mu\nu}(\nu+i\varepsilon)-T^{\gamma W}_{\mu\nu}(\nu-i\varepsilon)=4\pi
W_{\mu\nu}^{\gamma W}~{},$ (22)

where

$\displaystyle W_{\mu\nu}^{\gamma W}$ $\displaystyle=$
$\displaystyle\frac{1}{4\pi}\sum_{X}(2\pi)^{4}\delta^{(4)}(p+q-p_{X})\langle
p|J_{\mu}^{em}(0)|X\rangle\langle X|J_{\nu}^{W}(0)|n\rangle$ (23)
$\displaystyle=$
$\displaystyle{-\frac{i\epsilon_{\mu\nu\alpha\beta}q^{\alpha}p^{\beta}}{2(p\cdot
q)}F_{3}+\frac{i\epsilon_{\mu\nu\alpha\beta}q^{\alpha}}{(p\cdot
q)}\left[S^{\beta}g_{1}+\left(S^{\beta}-\frac{(S\cdot q)}{p\cdot
q}p^{\beta}\right)g_{2}\right]+...}$

The structure functions $F_{3}$ and $g_{1,2}$ can be decomposed similarly to
$I=0,1$ components just like Eq.(20).

According to the crossing behavior established earlier and noticing that they
cannot diverge faster than $\nu$ when $\nu\rightarrow\infty$, the amplitudes
entering Eq. (18) have the following dispersion representation:

$\displaystyle T_{3}^{(0)}(\nu,Q^{2})$ $\displaystyle=$
$\displaystyle-4i\nu\int_{0}^{\infty}d\nu^{\prime}\frac{F_{3}^{(0)}(\nu^{\prime},Q^{2})}{\nu^{\prime
2}-\nu^{2}}$ $\displaystyle S_{1}^{(0)}(\nu,Q^{2})$ $\displaystyle=$
$\displaystyle-4i\nu\int_{0}^{\infty}d\nu^{\prime}\frac{g_{1}^{(0)}(\nu^{\prime},Q^{2})}{\nu^{\prime
2}-\nu^{2}}$ $\displaystyle S_{2}^{(0)}(\nu,Q^{2})$ $\displaystyle=$
$\displaystyle-4i\nu\int_{0}^{\infty}d\nu^{\prime}\frac{g_{2}^{(0)}(\nu^{\prime},Q^{2})}{\nu^{\prime
2}-\nu^{2}}=-4i\nu^{3}\int_{0}^{\infty}d\nu^{\prime}\frac{g_{2}^{(0)}(\nu^{\prime},Q^{2})}{\nu^{\prime
2}(\nu^{\prime 2}-\nu^{2})}~{},$ (24)

where the $\nu$-integration is extended down to 0 to include the Born
contribution. Notice that in the last line we have slightly modified the DR of
$S_{2}^{(0)}$ using the Burkhardt-Cottingham (BC) sum rule Burkhardt and
Cottingham (1970):

$\int_{0}^{1}dxg_{2}(x,Q^{2})=0~{},$ (25)

where $x=Q^{2}/(2M\nu^{\prime})$ is the Bjorken variable. This sum rule is a
superconvergence relation and is expected to hold at all $Q^{2}$. The benefit
of this treatment will become apparent in the later section. Substituting
Eq.(24) into Eq.(18) and using the following Wick rotation formula Marciano
and Sirlin (1975),

$\int\frac{d^{4}q}{(2\pi)^{4}}F(\nu,Q^{2})=\frac{i}{8\pi^{3}}\int_{0}^{\infty}dQ^{2}Q^{2}\int_{-1}^{+1}du\sqrt{1-u^{2}}F(iQu,Q^{2})~{},$
(26)

we can integrate the variable $u$ analytically to obtain our final dispersive
representation of $\delta g_{V,A}^{\gamma W}$ as follows:

$\displaystyle\Box_{\gamma W}^{V}$ $\displaystyle=$
$\displaystyle\frac{\alpha_{em}}{\pi\mathring{g}_{V}}\int_{0}^{\infty}\frac{dQ^{2}}{Q^{2}}\frac{M_{W}^{2}}{M_{W}^{2}+Q^{2}}\int_{0}^{1}dx\frac{1+2r}{(1+r)^{2}}F_{3}^{(0)}(x,Q^{2})$
$\displaystyle\Box_{\gamma W}^{A}$ $\displaystyle=$
$\displaystyle-\frac{2\alpha_{em}}{\pi\mathring{g}_{A}}\int_{0}^{\infty}\frac{dQ^{2}}{Q^{2}}\frac{M_{W}^{2}}{M_{W}^{2}+Q^{2}}\int_{0}^{1}\frac{dx}{(1+r)^{2}}\left[\frac{5+4r}{3}g_{1}^{(0)}(x,Q^{2})-\frac{4M^{2}x^{2}}{Q^{2}}g_{2}^{(0)}(x,Q^{2})\right]~{},$
(27)

with $r=\sqrt{1+4M^{2}x^{2}/Q^{2}}$. As a useful crosscheck, the two-photon
exchange correction to the hyperfine splitting in ordinary and muonic atoms is
expressed through analogous two-fold integrals over electromagnetic spin
structure functions Carlson et al. (2008, 2011).

The quantity $\Box_{\gamma W}^{V}$, relevant for the extraction of $V_{ud}$
from superallowed $\beta$-decays is well-studied within the dispersive
approach Seng et al. (2018, 2019, 2020a); Shiells et al. (2020) and is not
addressed here. The main obstacle in those studies is the absence of direct
experimental data of the structure function $F_{3}^{(0)}$. This forces one to
either rely on data of $F_{3}$ from a different isospin channel (whose
relation to $F_{3}^{(0)}$ contains a residual model-dependence), or from
indirect lattice QCD data. On the other hand, despite having received much
less attention, a high-precision dispersive analysis of $\Box_{\gamma W}^{A}$
is in fact much more robust because it depends on the parity-conserving, spin-
dependent structure functions $g_{1,2}^{(0)}$. The isospin symmetry
unambiguously relates them to $g_{1,2}^{N}$ measured in ordinary DIS:

$g_{1,2}^{(0)}=\frac{1}{2}\left\\{g_{1,2}^{p}-g_{1,2}^{n}\right\\}~{},$ (28)

the latter are defined via

$\displaystyle W_{\mu\nu}^{\gamma\gamma,N}$ $\displaystyle=$
$\displaystyle\frac{1}{4\pi}\sum_{X}(2\pi)^{4}\delta^{(4)}(p+q-p_{X})\langle
N|J_{\mu}^{em}(0)|X\rangle\langle X|J_{\nu}^{em}(0)|N\rangle$ (29)
$\displaystyle=$
$\displaystyle{\frac{i\epsilon_{\mu\nu\alpha\beta}q^{\alpha}}{(p\cdot
q)}\left[S^{\beta}g_{1}^{N}+\left(S^{\beta}-\frac{(S\cdot q)}{p\cdot
q}p^{\beta}\right)g_{2}^{N}\right]+...}$

Therefore, it is possible to perform a fully data-driven analysis of
$\Box_{\gamma W}^{A}$ without introducing further model-dependence at low
$Q^{2}$. We will perform such an analysis in the sections below. Following
Sirlin’s notation Sirlin (1967), we express our result as

$\Box_{\gamma
W}^{A}=\frac{\alpha_{em}}{2\pi}d=\frac{\alpha_{em}}{2\pi}\left[d_{B}+d_{1}+d_{2}\right]~{},$
(30)

where $d_{B}$, $d_{1}$ and $d_{2}$ represent the elastic (Born) contribution,
the inelastic contributions from $g_{1}^{(0)}$ and the inelastic contributions
from $g_{2}^{(0)}$ respectively, which we will evaluate separately in the
following sections.

##  V Elastic (Born) contribution

Figure 3: The pole (a,b) and seagull (c) diagrams that contribute to
$T^{\mu\nu}_{\gamma W}$.

Substituting $X=p$ into Eq.(23) and using the elastic form factors defined in
Eq.(5) give us the Born contribution to the spin structure functions needed
for the evaluation of the box correction:

$\displaystyle
g_{1}^{(0),B}=\frac{F_{1}^{W}G_{M}^{S}+F_{1}^{S}G_{M}^{W}}{8}\delta(1-x)~{},\quad
g_{2}^{(0),B}=-\tau\frac{F_{2}^{W}G_{M}^{S}+F_{2}^{S}G_{M}^{W}}{8}\delta(1-x)~{},$
(31)

where $G_{E}\equiv F_{1}-\tau F_{2}$ and $G_{M}\equiv F_{1}+F_{2}$ are the
usual electric and magnetic Sachs form factors defined for both the
electromagnetic and charged weak form factors, with $\tau=Q^{2}/(4M^{2})$. All
the form factors above are functions of $Q^{2}$ that drop at high $Q^{2}$,
therefore one can neglect the $Q^{2}$ dependence of the $W$-boxon propagator
$\sim Q^{2}/M_{W}^{2}$. With this the Born contribution reads,

$d_{B}=-\frac{1}{2\mathring{g}_{A}}\int_{0}^{\infty}\frac{dQ^{2}}{Q^{2}}\frac{1}{(1+r_{B})^{2}}\left\\{\frac{5+4r_{B}}{3}\left[F_{1}^{W}G_{M}^{S}+F_{1}^{S}G_{M}^{W}\right]+\left[F_{2}^{W}G_{M}^{S}+F_{2}^{S}G_{M}^{W}\right]\right\\}~{},$
(32)

where $r_{B}\equiv r|_{x=1}$. We notice that $F_{1,2}^{W}=F_{1,2}^{V}$ by
isospin symmetry, which means Eq.(32) is fully determined by the four nucleon
electromagnetic form factors: $\\{G_{E}^{p},G_{M}^{p},G_{E}^{n},G_{M}^{n}\\}$.
Different parameterizations of these form factors Drechsel et al. (2003);
Lorenz et al. (2012, 2015); Ye et al. (2018); Lin et al. (2021a, b) all give
consistent results within their respective error bars. In particular, the
parametrization of Ref. Ye et al. (2018) leads to

$d_{B}=1.216(6)_{G_{E}^{p}}(9)_{G_{M}^{p}}(1)_{G_{E}^{n}}(2)_{G_{M}^{n}}=1.22(1)_{\mathrm{FF}}~{}.$
(33)

In particular, the central values of the contribution from $g_{1}^{(0),B}$ and
$g_{2}^{(0),B}$ are 1.17 and 0.04 respectively. We observe that the latter is
much smaller, which turns out to also be the case for the inelastic
contributions.

We pause here to comment on the Born contribution before moving on to the
inelastic contributions. One may also try to derive it by calculating the
Compton amplitudes from the first two Feynman diagrams in Fig.3 (the “pole
diagrams”) using the form factors in Eq.(5) as effective vertex functions, and
then plugging them into Eq.(18). The pole diagrams give:

$\displaystyle S_{1}^{(0),B}(\nu,Q^{2})$ $\displaystyle=$
$\displaystyle-\frac{iM\nu
Q^{2}}{Q^{4}-4M^{2}\nu^{2}-i\varepsilon}\left[F_{1}^{W}G_{M}^{S}+F_{1}^{S}G_{M}^{W}\right]+\frac{i\nu}{2M}F_{2}^{W}F_{2}^{S}$
$\displaystyle S_{2}^{(0),B}(\nu,Q^{2})$ $\displaystyle=$
$\displaystyle\frac{i\nu
Q^{4}}{4M\left(Q^{4}-4M^{2}\nu^{2}-i\varepsilon\right)}\left[F_{2}^{W}G_{M}^{S}+G_{M}^{W}F_{2}^{S}\right]-\frac{i\nu}{4M}\left[F_{2}^{W}G_{M}^{S}+G_{M}^{W}F_{2}^{S}\right]~{}.$
(34)

We split each expression into two term, where the first term contains a
singularity at $\nu^{2}=\nu_{B}^{2}$ and vanishes as $1/\nu$ when
$\nu\rightarrow\infty$, while the second term is regular and diverges as $\nu$
when $\nu\rightarrow\infty$. It is easy to see that, retaining only the first
term leads again to Eq.(32), apart from a numerically small difference
originating from our accounting for the BC sum rule in the DRs, effectively
redefining the $g_{2}$ contribution to $d_{B}$. The regular terms in Eq.(34)
lead to an extra small deviation from Eq.(32). The origin of this deviation
lies in the Gerasimov-Drell-Hearn sum rule Gerasimov (1965); Drell and Hearn
(1966) and its extension to finite $Q^{2}$ Deur et al. (2018) which relate the
regular low-energy term to an integral over the inelastic part of $g_{1}$.

This discussion simply means that the definition of the “elastic” contribution
is not exactly the same in the diagrammatic and the dispersive representation.
Of course, if we were able to calculate the full (i.e. pole + seagull)
$T_{\gamma W}^{\mu\nu}$ exactly at all values of $\\{\nu,Q^{2}\\}$ with the
diagrammatic approach, then the outcome must be identical to the DR analysis.
But since this is impossible, the dispersive representation provides a much
better starting point. We want to also point out that Refs.Hayen (2021a, b)
attempted to calculate the Born contribution from the pole diagrams in Fig.3
(let us call it $d_{B}^{\prime}$). 222In the earlier versions of these Refs.,
the author made some algebraic mistakes when dealing with the symmetric loop
integral of the form $\int d^{4}qq_{\alpha}q_{\beta}F(\nu,Q^{2})$ (i.e.
Eq.(19)). As a consequence, an incorrect analytic formula which gave an
unexpectedly large value of $d_{B}^{\prime}=2.64(3)$ was obtained. The
published version of Ref.Hayen (2021a) corrected these mistakes, but retained
only the term proportional to $(5+4r_{B})/3$. Should all terms be retained,
they would have obtained the following result:

$\displaystyle d_{B}^{\prime}$ $\displaystyle=$
$\displaystyle-\frac{1}{2\mathring{g}_{A}}\int_{0}^{\infty}\frac{dQ^{2}}{Q^{2}}\frac{1}{(1+r_{B})^{2}}\biggl{\\{}\frac{5+4r_{B}}{3}\left[F_{1}^{W}G_{M}^{S}+F_{1}^{S}G_{M}^{W}\right]+\left[F_{2}^{W}G_{M}^{S}+F_{2}^{S}G_{M}^{W}\right]-\frac{3\tau}{2}(1+r_{B})^{2}F_{2}^{W}F_{2}^{S}\biggr{\\}}$
(35) $\displaystyle=$ $\displaystyle 1.23(1)~{}$

which only differs from the DR’s definition of $d_{B}$ in Eq.(32) by a
numerically small term. For the benefit of interested readers, we also provide
in Appendix A an alternative derivation of Eq.(35) without making use of the
invariant amplitudes.

##  VI Inelastic contributions

Extensive measurements of the structure function $g_{1}^{N}$ were carried out
in SLAC Anthony et al. (1993); Abe et al. (1995, 1997), CERN Adams et al.
(1994); Alexakhin et al. (2007); Alekseev et al. (2010); Aghasyan et al.
(2018), DESY Ackerstaff et al. (1997) and JLab Deur et al. (2004); Wesselmann
et al. (2007); Deur et al. (2008); Guler et al. (2015); Fersch et al. (2017).
In particular, we utilize the results from the EG1b experiment at JLab that
measured the $g_{1}^{p}$ Fersch et al. (2017) and $g_{1}^{n}$ Guler et al.
(2015) in a wide range of $\\{x,Q^{2}\\}$, from which the moments
$\Gamma_{i}^{N}$ ($i=1,3,5$) were computed in bins of $Q^{2}$, common for $p$
and $n$, from 0.05 GeV2 to 3.5 GeV2. Full results are available in the
supplementary material of each respective paper333There is a more recent
measurement of $g_{1}^{p}$ from JLab at low $Q^{2}$ Zheng et al. (2021), but
unfortunately it does not measure $g_{1}^{n}$ simultaneously.. Data on $g_{2}$
are generally scarce Anthony et al. (1999, 2003); Amarian et al. (2004);
Wesselmann et al. (2007); Kramer et al. (2005); Fersch et al. (2017) and
insufficient for a fully data-based analysis. Fortunately, its contribution is
generally expected to be small due to the BC sum rule that forces the first
moment of $g_{2}$ to vanish identically when accounting for elastic and
inelastic contributions.

Since the data do not extend to an arbitrarily large $Q^{2}$ needed to
evaluate the integrals, we make use of perturbative QCD results which are well
under control theoretically above a separation scale $Q_{0}^{2}$, while
directly using the experimental data that contains both perturbative and
nonperturbative physics below that scale. Following our earlier works on the
vector RC Seng et al. (2018, 2019) we take $Q_{0}^{2}=2$ GeV2. For the vector
RC case, not only does $Q_{0}^{2}=2$ GeV2 mark the onset of pQCD regime, it
also corresponds to the scale, below which the quality of data deteriorates
severely leading to some sensitivity to $Q_{0}^{2}$. In the case of the axial
RC, this scale lies well within the region covered by data, hence shifting it
to a slightly higher value (not lower because pQCD description starts to break
down) does not change the result.

###  VI.1 Contribution of $g_{1}^{(0)}$

As shown in Eq.(28), the polarized structure function $g_{1}^{(0)}$ is simply
related to $\\{g_{1}^{N}\\}$ that are measurable in DIS experiments. We define
their moments as:

$\Gamma_{i}^{N}(Q^{2})\equiv\int_{0}^{x_{\pi}}x^{i-1}g_{1}^{N}(x,Q^{2})~{},$
(36)

where $x_{\pi}=Q^{2}/\left[(M+m_{\pi})^{2}-M^{2}+Q^{2}\right]$ is the pion
production threshold. Notice that the definition above excludes the elastic
contribution at $x=1$, which is a general convention adopted by most of the
experimental papers. The first moment
$\Gamma_{1}^{p-n}\equiv\Gamma_{1}^{p}-\Gamma_{1}^{n}$ is of particular
interest because it satisfies the polarized Bjorken sum rule Bjorken (1966,
1970) at $Q^{2}\rightarrow\infty$. However, at large but finite $Q^{2}$ it
receives a number of corrections Ji and Unrau (1994):

$\Gamma_{1,\mathrm{th}}^{p-n}(Q^{2})=\frac{|\mathring{g}_{A}|}{6}C_{\mathrm{Bj}}(Q^{2})+\sum_{i=2}^{\infty}\frac{\mu_{2i}^{p-n}}{Q^{2i-2}}~{},\quad\text{large
$Q^{2}$}~{}$ (37)

here the subscript “th” denotes the theory prediction (at large $Q^{2}$). The
first term at the right hand side is the Bjorken sum rule with a pQCD
correction factor444Please be reminded that one should not consider again the
running effect of the QED coupling constant in $\Box_{\gamma W}^{V}$ and
$\Box_{\gamma W}^{A}$, because it is already contained in the factor
$\delta_{\mathrm{HO}}^{\mathrm{QED}}$ in Eq.(10). So, throughout this paper we
always take $\alpha_{em}=7.2973525693(11)\times 10^{-3}$ as a constant., while
the second term summarizes the higher-twist (HT) effects starting from twist-
four. The pQCD correction factor is written as:

$C_{\mathrm{Bj}}(Q^{2})=1-\sum_{n=1}^{\infty}\tilde{c}_{n}\left(\frac{\alpha_{s}}{\pi}\right)^{n}~{},$
(38)

where $\alpha_{s}$ is the running strong coupling constant in the
$\overline{\mathrm{MS}}$ scheme, while the coefficients $\\{\tilde{c}_{n}\\}$
are calculated at present to $n=4$ Baikov et al. (2010a, b):

$\displaystyle\tilde{c}_{1}$ $\displaystyle=$ $\displaystyle 1$
$\displaystyle\tilde{c}_{2}$ $\displaystyle=$ $\displaystyle 4.583-0.333n_{f}$
$\displaystyle\tilde{c}_{3}$ $\displaystyle=$ $\displaystyle
41.44-7.607n_{f}+0.177n_{f}^{2}$ $\displaystyle\tilde{c}_{4}$ $\displaystyle=$
$\displaystyle 479.4-123.4n_{f}+7.697n_{f}^{2}-0.1037n_{f}^{3}~{},$ (39)

with $n_{f}$ the number of active quark flavors, and we refer the reader to
Refs. Seng et al. (2018, 2019) for full detail of the pQCD contribution and
relevant discussions and references. In the meantime, only the twist-four term
among all the HT corrections needs to be included for our precision goal.
There are several recent determinations of the coefficient $\mu_{4}^{p-n}$
Deur et al. (2014); Kotlorz et al. (2017); Ayala et al. (2018) that are
largely consistent with each other. In this work we quote the value
$\mu_{4}^{p-n}=(-0.047\pm 0.020)M^{2}$ in Ref.Kotlorz et al. (2017). We find
that at 2 GeV2, the inclusion of the twist-four correction reduces the size of
$\Gamma_{1,\mathrm{th}}^{p-n}$ by about 13%, but its total contribution to
$d_{1}$ through the integral at $Q^{2}>2$ GeV2 is only about 1%. Coming back
to our problem, we write

$d_{1}=-\frac{3}{2\mathring{g}_{A}}\int_{0}^{\infty}\frac{dQ^{2}}{Q^{2}}\frac{M_{W}^{2}}{M_{W}^{2}+Q^{2}}\bar{\Gamma}_{1}^{p-n}(Q^{2})~{},\quad\bar{\Gamma}_{1}^{p-n}(Q^{2})\equiv\int_{0}^{x_{\pi}}dx\frac{4(5+4r)}{9(1+r)^{2}}\left\\{g_{1}^{p}(x,Q^{2})-g_{1}^{n}(x,Q^{2})\right\\}~{}.$
(40)

When $Q^{2}\rightarrow\infty$ the function $\bar{\Gamma}_{1}^{p-n}$ reduces to
$\Gamma_{1}^{p-n}$, but at low $Q^{2}$ the two are not identical due to target
mass corrections $\sim M^{2}/Q^{2}$ contained in the factor
$f(x,Q^{2})=4(5+4r)/(9(1+r)^{2})$, where $r=\sqrt{1+4x^{2}M^{2}/Q^{2}}$

Figure 4: Left panel: Comparison between the function $4(5+4r)/(9(1+r)^{2})$
(blue solid curve), its global fit approximation (orange dashed curve) and the
Taylor-expansion approximation (red solid curve) at $Q^{2}=1$ GeV2,
$0<x<x_{\pi}$. Right panel: Data points of $\bar{\Gamma}_{1}^{p-n}$
reconstructed from the EG1b experiment Fersch et al. (2017); Guler et al.
(2015) versus the theory prediction at large $Q^{2}$ with (red band) and
without (brown curve) the twist-four correction.

We use the following strategy to reconstruct the full $\bar{\Gamma}_{1}^{p-n}$
at low $Q^{2}$ from data. We fit $f(x,Q^{2})$ as function of $x$ at fixed
$Q^{2}$ as

$f(x,Q^{2})=a(Q^{2})+b(Q^{2})x^{2}+c(Q^{2})x^{4}+\dots,\quad 0<x<x_{\pi}~{},$
(41)

and obtain the three fitting parameters $a(Q^{2})$, $b(Q^{2})$ and $c(Q^{2})$
by first dividing $0<x<x_{\pi}$ into, say, 1000 equal intervals, evaluating
the 1000 respective discrete values for $f(x,Q^{2})$, and performing a three-
parameter fit with these discrete points using, e.g. Mathematica. We find that
this procedure allows for a very precise reproduction of the entire curve at
$0<x<x_{\pi}$, where the difference between the original and the fitted curve
is negligible for all practical purposes. Contrarily, a simple Taylor
expansion in powers of $x^{2}M^{2}/Q^{2}$,

$a(Q^{2})+b(Q^{2})x^{2}+c(Q^{2})x^{4}\ncong
1-\frac{10M^{2}x^{2}}{9Q^{2}}+\frac{7M^{4}x^{4}}{3(Q^{2})^{2}}$ (42)

is only applicable at high $Q^{2}$ and low $x$ but significantly deviates for
larger $x$. It breaks down completely at $Q^{2}=0$ where it becomes divergent,
whereas $f(x,0)$ remains finite. As an illustration, the two approximated
expressions evaluated at a representative value of $Q^{2}=1$ GeV2 are plotted
together with the analytic form in the left panel of Fig.4, and we clearly see
that Eq.(41) nicely reproduces the latter for the full range of $x$.

We thus reconstruct the full $\bar{\Gamma}_{1}^{p-n}$ in each bin of $Q^{2}$
in terms of the lowest Mellin moments,

$\bar{\Gamma}_{1}^{p-n}(Q^{2})=a(Q^{2})\Gamma_{1}^{p-n}(Q^{2})+b(Q^{2})\Gamma_{3}^{p-n}(Q^{2})+c(Q^{2})\Gamma_{5}^{p-n}(Q^{2})+\dots~{},$
(43)

where $\Gamma_{i}^{p-n}$ ($i=1,3,5$) are obtained from Refs.Fersch et al.
(2017); Guler et al. (2015). We find that the difference between
$\bar{\Gamma}_{1}^{p-n}$ and $\Gamma_{1}^{p-n}$ (i.e., the effect of higher
moments) does not exceed 3.5% for $Q^{2}=2~{}$GeV2, which implies a negligible
difference in the integral at $Q^{2}>2~{}$GeV2. Therefore we will not
distinguish between the two above 2 GeV2. In contrast, the higher-twist
correction due to $\mu_{4}^{p-n}$ reaches 13% and needs to be kept along.

The right panel in Fig.4 shows the reconstructed $\bar{\Gamma}_{1}^{p-n}$ data
points versus the large-$Q^{2}$ theory prediction using Eq.(37). We find that
the theory and experiment match well at $Q^{2}>2$ GeV2 (observe how the twist-
four correction is needed to reconcile the two), which justifies our choice of
$Q_{0}^{2}=2$ GeV2 as the separation scale between the perturbative and non-
perturbative regime. We therefore evaluate $d_{1}$ separately in these two
regions. At $Q^{2}<Q_{0}^{2}$, we fit three curves that correspond to the
upper bounds, central values and lower bounds of the discrete data points
respectively, and evaluate the $Q^{2}$-integral and its uncertainty by
integrating these three curves. Since the uncertainties of the data points are
mainly systematics, this prescription takes into account the possible positive
correlation effects. The resulting uncertainty is thus a conservative one; it
is likely that it can further be reduced, but this would require a dedicated
study of the systematic uncertainties of the data which lies beyond the scope
of the present work. Meanwhile, at $Q^{2}>Q_{0}^{2}$ we evaluate the integral
using the theory prediction in Eq.(37). The results are as follows:

$\displaystyle d_{1}^{<}$ $\displaystyle=$
$\displaystyle-\frac{3}{2\mathring{g}_{A}}\int_{0}^{Q_{0}^{2}}\frac{dQ^{2}}{Q^{2}}\frac{M_{W}^{2}}{M_{W}^{2}+Q^{2}}\left[\bar{\Gamma}_{1}^{p-n}(Q^{2})\right]_{\mathrm{data}}=0.30(4)_{\mathrm{data}}$
$\displaystyle d_{1}^{>}$ $\displaystyle=$
$\displaystyle-\frac{3}{2\mathring{g}_{A}}\int_{Q_{0}^{2}}^{\infty}\frac{dQ^{2}}{Q^{2}}\frac{M_{W}^{2}}{M_{W}^{2}+Q^{2}}\left[\frac{|\mathring{g}_{A}|}{6}C_{\mathrm{Bj}}(Q^{2})+\frac{\mu_{4}^{p-n}}{{Q^{2}}}\right]=1.83(1)_{\mathrm{HT}}~{},$
(44)

where we neglected the uncertainty associated with the leading twist
contribution compared to those coming from the data at low $Q^{2}$ and the HT
correction (i.e. the coefficient $\mu_{4}^{p-n}$) at high $Q^{2}$. The total
contribution of $g_{1}$ reads,

$d_{1}=d_{1}^{<}+d_{1}^{>}=2.14(4)_{\mathrm{data}}(1)_{\mathrm{HT}}~{}.$ (45)

We note that the integral below $Q^{2}\leq 0.05$ GeV2 not covered by the data
but making part of $d_{1}^{<}$ is controlled by the isovector GDH sum rule
Gerasimov (1965); Drell and Hearn (1966),
$d\Gamma_{1}^{p-n}(Q^{2}=0)/dQ^{2}=(\kappa_{n}^{2}-\kappa_{p}^{2})/8M^{2}$,
with $\kappa_{p,n}$ denoting the proton’s (neutron’s) anomalous magnetic
moment, respectively. Connecting the GDH-fixed value at $Q^{2}=0$ to the
lowest data point produces a negligible $d_{1}^{\mathrm{low}~{}Q^{2}}\lesssim
0.003$ contribution which is safely accommodated within the uncertainty.

###  VI.2 Contribution of $g_{2}^{(0)}$

In fact, we have already implemented this sum rule in the derivation of the DR
of $S_{2}^{(0)}$. We emphasize the importance of this procedure for a reliable
estimate of the contribution of $g_{2}^{(0)}$ to $\Box_{\gamma W}^{A}$: since
experimental data typically only cover the inelastic region, enforcing an
exact vanishing of the first moment of $g_{2}^{(0)}$ while operating with
phenomenological parametrizations of different pieces can be a delicate
matter. The explicit use of the BC sum rule thus precludes any numerically
significant mistake caused by an imperfection of these parametrizations. As a
result, the dispersive representation of $\Box_{\gamma W}^{A}$ only contains
higher moments of $g_{2}^{(0)}$, in which the non-perturbative physics at
small $x$ is suppressed. Additionally, since every extra power of $x^{2}$ is
accompanied by $1/Q^{2}$, the contribution of $g_{2}^{(0)}$ bears no large
logarithms. Due to the smallness of the $g_{2}$ contribution, we opt for an
approximate treatment and include this result in the estimate of the
systematic uncertainty.

The inelastic contribution coming from $g_{2}^{(0)}$ reads,

$\displaystyle d_{2}$ $\displaystyle=$
$\displaystyle\frac{2}{\mathring{g}_{A}}\int_{0}^{\infty}\frac{dQ^{2}}{Q^{2}}\frac{M_{W}^{2}}{M_{W}^{2}+Q^{2}}\int_{0}^{x_{\pi}}\frac{dx}{(1+r)^{2}}\frac{4M^{2}x^{2}}{Q^{2}}\left[g_{2}^{p}(x,Q^{2})-g_{2}^{n}(x,Q^{2})\right]~{},$
(46)

where the isospin relation in Eq.(28) is used. Rather than relying on data on
$g_{2}^{N}$, we decompose $g_{2}$ into twist-two and twist-three (and higher)
components Zyla et al. (2020),
$g_{2}^{N}=g_{2,\mathrm{tw}2}^{N}+g_{2,\mathrm{tw3}+}^{N}$, and use the
Wandzura-Wilczek relation Wandzura and Wilczek (1977) for the former,

$\displaystyle
g_{2,\mathrm{tw}2}^{N}(x,Q^{2})=-g_{1}^{N}(x,Q^{2})+\int\limits_{x}^{1}\frac{dy}{y}g_{1}^{N}(y,Q^{2})~{}.$
(47)

Notice that the relation above should be understood to not contain the elastic
contribution at $x=1$, because otherwise one could take $x_{\pi}<x<1$ at both
sides, and then the left hand side and the first term at the right hand side
would vanish but the second term at the right hand side would not, which is a
contradiction. With this in mind, we multiply both sides by $x^{i-1}$ and
integrate them at $0<x<x_{\pi}$ to obtain:

$\int_{0}^{x_{\pi}}x^{n-1}g_{2,\mathrm{tw}2}^{N}(x,Q^{2})=\frac{1-n}{n}\Gamma_{n}^{N}(Q^{2})~{},$
(48)

where $\Gamma_{n}^{N}$ is defined in Eq.(36). Therefore, we may use the
available information of $g_{1}$ to evaluate the twist-two contribution to
$d_{2}$.

Again, we discuss the integral at large and small $Q^{2}$ separately. For
$Q^{2}>Q_{0}^{2}$, we keep the leading term,

$d_{2,\mathrm{tw}2}^{>}\approx-\frac{4}{3\mathring{g}_{A}}\int_{Q_{0}^{2}}^{\infty}\frac{dQ^{2}}{Q^{2}}\frac{M_{W}^{2}}{M_{W}^{2}+Q^{2}}\frac{M^{2}}{Q^{2}}\Gamma_{3}^{p-n}(Q^{2})~{},$
(49)

where we have used Eq.(48) and set $r\rightarrow 1$. The $Q^{2}$ dependence of
the $W$ propagator can also be safely neglected as the integral converges.
There is no simple sum rule for $\Gamma_{3}^{p-n}$ at large $Q^{2}$, but we
may adopt a naïve valence quark picture that assumes each valence quark
carries 1/3 of the nucleon’s momentum. This gives:

$g_{1}^{p-n}(x)\approx\frac{|\mathring{g}_{A}|}{6}\delta\left(x-1/3\right)~{},\quad\text{large
$Q^{2}$}$ (50)

which automatically reproduces the free polarized Bjorken sum rule. This naïve
picture predicts $\Gamma_{3}^{p-n}\approx 0.11\Gamma_{1}^{p-n}$, which we may
check against the experimental data: At $Q^{2}\approx 3.4$ GeV2, Refs.Fersch
et al. (2017); Guler et al. (2015) give $\Gamma_{1}^{p-n}\approx 0.1558$,
$\Gamma_{3}^{p-n}\approx 0.0128$, i.e. $\Gamma_{3}^{p-n}\approx
0.08\Gamma_{1}^{p-n}$. So our naïve picture overestimates the size of
$\Gamma_{3}^{p-n}$ by some 30%, an acceptable uncertainty given our precision
goal. With the above, we obtain:

$d_{2,\mathrm{tw}2}^{>}\approx\frac{2}{81}\int_{Q_{0}^{2}}^{\infty}\frac{dQ^{2}}{Q^{2}}\frac{M^{2}}{Q^{2}}=\frac{2M^{2}}{81Q_{0}^{2}}=0.010(3)~{}.$
(51)

Next we turn to the small-$Q^{2}$ region where accounting for the leading
twist may not be sufficient. As before, the leading twist contribution

$\displaystyle d_{2,\mathrm{tw}2}^{<}$ $\displaystyle=$
$\displaystyle\frac{2}{\mathring{g}_{A}}\int_{0}^{Q_{0}^{2}}\frac{dQ^{2}}{Q^{2}}\int_{0}^{x_{\pi}}\frac{dx}{(1+r)^{2}}\frac{4M^{2}x^{2}}{Q^{2}}\left[g_{2,\mathrm{tw}2}^{p}(x,Q^{2})-g_{2,\mathrm{tw}2}^{n}(x,Q^{2})\right]~{},$
(52)

is reconstructed by using the WW relation and the $x$-integral related to
measured moments of $g_{1}$ by performing a two parameter fit of the
kinematical function $f^{\prime}(x,Q^{2})$ for fixed $Q^{2}$ (since the first
moment is removed by the BC sum rule, a two-parameter fit is already
sufficient for our precision goal).

$f^{\prime}(x,Q^{2})=\frac{1}{(1+r)^{2}}\frac{4M^{2}x^{2}}{Q^{2}}=b^{\prime}(Q^{2})x^{2}+c^{\prime}(Q^{2})x^{4}+\dots~{},\quad
0<x<x_{\pi}~{}.$ (53)

Using Eq.(48), we obtain

$\int_{0}^{x_{\pi}}\frac{dx}{(1+r)^{2}}\frac{4M^{2}x^{2}}{Q^{2}}\left[g_{2,\mathrm{tw}2}^{p}(x,Q^{2})-g_{2,\mathrm{tw}2}^{n}(x,Q^{2})\right]\approx-\frac{2}{3}b^{\prime}(Q^{2})\Gamma_{3}^{p-n}(Q^{2})-\frac{4}{5}c^{\prime}(Q^{2})\Gamma_{5}^{p-n}(Q^{2})~{}.$
(54)

With $\Gamma_{i}^{p-n}$ taken from Refs.Fersch et al. (2017); Guler et al.
(2015), this gives,

$d_{2,\mathrm{tw}2}^{<}\approx-\frac{2}{\mathring{g}_{A}}\int_{0}^{Q_{0}^{2}}\frac{dQ^{2}}{Q^{2}}\left[\frac{2}{3}b^{\prime}(Q^{2})\Gamma_{3}^{p-n}(Q^{2})+\frac{4}{5}c^{\prime}(Q^{2})\Gamma_{5}^{p-n}(Q^{2})\right]=0.01(1)~{},$
(55)

where we assigned a conservative 100% uncertainty to the entire contribution.

To quantify higher twist contributions we recall the definition of the “color
polarizability” Shuryak and Vainshtein (1982); Jaffe (1990)

${\mathbf{d}_{2}(Q^{2})=3\int_{0}^{1}dxx^{2}[g_{2}(x,Q^{2})-g_{2,\mathrm{tw2}}(x,Q^{2})],}$
(56)

of which we only consider the inelastic part $\bar{\mathbf{d}}_{2}(Q^{2})$
coming from the interval $0\leq x\leq x_{\pi}$ Alarcón et al. (2020) since the
elastic part is already taken into account. In terms of this polarizability
and neglecting higher moments, we obtain for the contribution of twist-three
and higher,

$\displaystyle d_{2,\mathrm{tw}3+}^{<}$ $\displaystyle=$
$\displaystyle{\frac{2M^{2}}{3\mathring{g}_{A}}\int_{0}^{Q_{0}^{2}}\frac{dQ^{2}}{Q^{4}}[\bar{\mathbf{d}}_{2}^{p}(Q^{2})-\bar{\mathbf{d}}_{2}^{n}(Q^{2})]}~{}.$
(57)

For numerical estimates, we rely on the recent analysis of generalized spin
polarizabilities of the nucleon in baryon chiral effective theory Alarcón et
al. (2020). We obtain, assigning a conservative 100% uncertainty,

$\displaystyle d_{2,\mathrm{tw}3+}^{<}$ $\displaystyle=$ $\displaystyle
0.03(3)~{}.$ (58)

Combining the various pieces we finally arrive at

$d_{2}=0.05(3)$ (59)

as our estimate of the total inelastic contribution from $g_{2}^{(0)}$. We
observe that it is two orders of magnitude smaller than $d_{1}$, following the
same hierarchy as in $d_{B}$.

The inelastic contribution in our DR analysis then reads
$d_{1}+d_{2}=2.19(4)_{\mathrm{data}}(1)_{\mathrm{HT}}(3)_{g_{2}}$. This is to
be compared with $2.31(9)$ from Refs. Hayen (2021a, b), and we see that the
two do not quite agree within error bars. While being identical in the
large-$Q^{2}$ treatment, their estimation of the low-$Q^{2}$ contribution is
largely model-based, raising questions about the reliability of the
uncertainty. In contrast, in our treatment the low-$Q^{2}$ contribution is
completely fixed by experimental data without any further assumption, apart
from the small $g_{2}$ correction for which only few assumption were made,
which will become testable as soon as new, higher-quality low-$Q^{2}$ data for
$g_{2}$ will become available.

##  VII Final discussions

Collecting all the results from Sec.V and VI gives:

$\Box_{\gamma
W}^{A}=3.96(1)_{\mathrm{FF}}(5)_{\mathrm{data}}(1)_{\mathrm{HT}}(3)_{g_{2}}\times
10^{-3}=3.96(6)\times 10^{-3}~{},$ (60)

where the uncertainties come from the elastic form factors, the low-$Q^{2}$
$g_{1}$ data, the HT-correction to $g_{1}$ and $g_{2}$, respectively. We
compare this to our recent update of $\Box_{\gamma W}^{V}$ using indirect
lattice inputs: $\Box_{\gamma W}^{V}=3.83(11)\times 10^{-3}$ Seng et al.
(2020a). These two numbers are very close to each other, and in fact their
difference is consistent with zero:

$\Box_{\gamma W}^{A}-\Box_{\gamma W}^{V}=0.13(11)_{V}(6)_{A}\times
10^{-3}~{}.$ (61)

Using Eqs.(11), (61) and the PDG average $\lambda=-1.2756(13)$ Zyla et al.
(2020), we obtain:

$\mathring{g}_{A}=-1.2754(13)_{\mathrm{exp}}(2)_{\mathrm{RC}}~{},$ (62)

which is consistent with the result from the current best lattice QCD
determination. Our result indicates that there is no practical distinction
between $\lambda$ and $\mathring{g}_{A}$, unless the experimental precision of
the former and the lattice precision of the latter have reached $3\times
10^{-4}$ or better.

We wrap up with some discussions of the future prospects. Within the same DR
framework, a much better precision is achieved for $\Box_{\gamma W}^{A}$ than
for $\Box_{\gamma W}^{V}$ thanks to the existence of high-quality data of the
structure function $g_{1}$ at $Q^{2}<2$ GeV2. On the other hand, the precision
of $\Box_{\gamma W}^{V}$ is limited by the low-quality data of the structure
function $F_{3}$ from neutrino (antineutrino)-nucleus scattering experiments
in the 80s Bolognese et al. (1983); Allasia et al. (1985). Better-quality data
may come from the Deep Underground Neutrino Experiment (DUNE) in the next
decade Acciarri et al. (2016); Alvarez-Ruso et al. (2018).

It was pointed out that a direct lattice QCD calculation of $\Box_{\gamma
W}^{V}$ is a promising way to proceed at the present stage Seng and Meißner
(2019). Several exploratory calculations of mesonic $\gamma W$-box diagrams
have shown great success Feng et al. (2020); Ma et al. (2021) and the same
technology is directly applicable to nucleon. At present, no such direct
calculation on the nucleon is available yet. A more involved comparison of the
DR result for $\Box_{\gamma W}^{V}$ with the lattice computation of the
respective quantity on the pion, amended with further phenomenological
ingredients shows a nearly perfect agreement Seng et al. (2020a). Even with
this reassuring agreement, it is not unthinkable of that a direct lattice
calculation could still disagree with the phenomenological, DR-based
evaluation. Examples of such an unexpected disagreement are the pion-nucleon
sigma term $\sigma_{\pi N}$ (see Ref.Aoki et al. (2019) and references
therein) and, more recently, the hadronic vacuum polarization contribution to
$g_{\mu}-2$ Borsanyi et al. (2020); Aoyama et al. (2020). They show that even
carefully-performed first-principles calculations or fully data-driven
analysis may still contain unknown, previously unanticipated systematic
effects which may seriously affect the implications of the corresponding
precision experiments. Given these precedents, it is always useful to cross-
check the lattice calculations with alternative methods. Our new result of
$\Box_{\gamma W}^{A}$ is perfectly up to this task as it is a solid
phenomenological determination with the uncertainty very well under control.

Further effort from the DR side should be dedicated to RCs to the GT strength
in nuclear mirror decays where a recent study Hayen and Severijns (2019)
revealed inconsistencies in the previous analyses. Removing these
inconsistencies led to a better agreement for the $V_{ud}$ extracted across
mirror and superallowed nuclear decays, as well as neutron decay. However,
Ref. Hayen and Severijns (2019) only partially accounted for the $\gamma
W$-box contribution. The dispersion formulation of the $\Box_{\gamma W}^{A}$
correction developed in this work can be directly applied to mirror systems.
Following Refs. Seng et al. (2019); Gorchtein (2019), nuclear modifications of
the universal free-neutron $\gamma W$-box correction can be computed, and we
defer this task to future work.

## Acknowledgements

We are extremely thankful to Alexandre Deur and Xiaochao Zheng for their
detailed explanations of the JLab experiments. We also appreciate Leendert
Hayen for many inspiring discussions. We are furthermore grateful to Vadim
Lensky and Vladimir Pascalutsa for providing their code for computing the
generalized spin polarizabilities. This work is supported in part by EU
Horizon 2020 research and innovation programme, STRONG-2020 project under
grant agreement No 824093 and by the German-Mexican research collaboration
Grant No. 278017 (CONACyT) and No. SP 778/4-1 (DFG) (M.G), by the Deutsche
Forschungsgemeinschaft (DFG, German Research Foundation) and the NSFC through
the funds provided to the Sino-German Collaborative Research Center TRR110
“Symmetries and the Emergence of Structure in QCD” (DFG Project-ID 196253076 -
TRR 110, NSFC Grant No. 12070131001) (C.Y.S).

##  Appendix A Alternative derivation of Eq.(35)

In this Appendix we outline a derivation of Eq.(35), namely the elastic
contribution in the diagrammatic representation, directly from the tensor
$T^{\gamma W}_{\mu\nu}$ without going through its invariant amplitudes. It
involves some interesting tricks to deal with fermionic spinors and thus is
worthwhile to be displayed for pedagogical purposes.

We start by computing $T^{\gamma W}_{\mu\nu}$ from the direct (D) and crossed
(C) pole diagrams in Fig.3. Using the elastic form factors as effective
vertices, we obtain:

$T_{\mu\nu}^{\gamma
W,B}=\frac{\bar{u}_{s}(p)\Gamma_{\mu\nu}^{\mathrm{D}}u_{s}(p)}{2M\nu-Q^{2}}+\frac{\bar{u}_{s}(p)\Gamma_{\mu\nu}^{\mathrm{C}}u_{s}(p)}{-2M\nu-Q^{2}}$
(63)

where

$\displaystyle\Gamma_{\mu\nu}^{\mathrm{D}}$ $\displaystyle=$
$\displaystyle\left[F_{1}^{p}\gamma_{\mu}-\frac{iF_{2}^{p}}{2M}\sigma_{\mu\alpha}q^{\alpha}\right]i\left[\not{p}+\not{q}+M\right]\left[F_{1}^{W}\gamma_{\nu}+\frac{iF_{2}^{W}}{2M}\sigma_{\nu\beta}q^{\beta}+G_{A}\gamma_{\nu}\gamma_{5}-\frac{G_{P}}{2M}\gamma_{5}q_{\nu}\right]$
$\displaystyle\Gamma_{\mu\nu}^{\mathrm{C}}$ $\displaystyle=$
$\displaystyle\left[F_{1}^{W}\gamma_{\nu}+\frac{iF_{2}^{W}}{2M}\sigma_{\nu\beta}q^{\beta}+G_{A}\gamma_{\nu}\gamma_{5}-\frac{G_{P}}{2M}\gamma_{5}q_{\nu}\right]i\left[\not{p}-\not{q}+M\right]\left[F_{1}^{n}\gamma_{\mu}-\frac{iF_{2}^{n}}{2M}\sigma_{\mu\alpha}q^{\alpha}\right]$
(64)

are matrices in the Dirac space.

We can get rid of the nucleon spinors $\bar{u}_{s}(p)$, $u_{s}(p)$ in the
expression of $T_{\mu\nu}^{\gamma W,B}$ using the following trick. First, we
recall that any Dirac structure $\Gamma$ can be decomposed in terms of
standard Dirac basis $1$, $\gamma_{5}$, $\gamma^{\alpha}$,
$\gamma^{\alpha}\gamma_{5}$, $\sigma^{\alpha\beta}$ using the following
identity:

$\Gamma=\frac{1}{4}\mathrm{Tr}[\Gamma]+\frac{1}{4}\mathrm{Tr}[\gamma_{5}\Gamma]\gamma_{5}+\frac{1}{4}\mathrm{Tr}[\gamma_{\alpha}\Gamma]\gamma^{\alpha}-\frac{1}{4}\mathrm{Tr}[\gamma_{\alpha}\gamma_{5}\Gamma]\gamma^{\alpha}\gamma_{5}+\frac{1}{8}\mathrm{Tr}[\sigma_{\alpha\beta}\Gamma]\sigma^{\alpha\beta}.$
(65)

Next, we have the following identities when a Dirac basis is sandwiched
between $\bar{u}_{s}(p)$ and $u_{s}(p)$:

$\displaystyle\bar{u}_{s}(p)u_{s}(p)$ $\displaystyle=$ $\displaystyle 2M$
$\displaystyle\bar{u}_{s}(p)\gamma_{5}u_{s}(p)$ $\displaystyle=$
$\displaystyle 0$ $\displaystyle\bar{u}_{s}(p)\gamma^{\alpha}u_{s}(p)$
$\displaystyle=$ $\displaystyle 2p^{\alpha}$
$\displaystyle\bar{u}_{s}(p)\gamma^{\alpha}\gamma_{5}u_{s}(p)$
$\displaystyle=$ $\displaystyle 2S^{\alpha}$
$\displaystyle\bar{u}_{s}(p)\sigma^{\alpha\beta}u_{s}(p)$ $\displaystyle=$
$\displaystyle\frac{2}{M}\epsilon^{\alpha\beta\rho\sigma}p_{\rho}S_{\sigma}~{}.$
(66)

Combining Eqs.(65) and (66), we obtain:

$\bar{u}_{s}(p)\Gamma_{\mu\nu}^{i}u_{s}(p)=\frac{M}{2}\mathrm{Tr}[\Gamma_{\mu\nu}^{i}]+\frac{p^{\alpha}}{2}\mathrm{Tr}[\gamma_{\alpha}\Gamma_{\mu\nu}^{i}]-\frac{S^{\alpha}}{2}\mathrm{Tr}[\gamma_{\alpha}\gamma_{5}\Gamma_{\mu\nu}^{i}]+\frac{1}{4M}\epsilon^{\alpha\beta\rho\sigma}p_{\rho}S_{\sigma}\mathrm{Tr}[\sigma_{\alpha\beta}\Gamma_{\mu\nu}^{i}],$
(67)

where $i=\mathrm{D,C}.$ The right hand side is free from the nucleon spinors.
The trace of the Dirac matrices can be performed using various Mathematica
packages, so we do not display the explicit results here.

We then plug our expression of $T_{\mu\nu}^{\gamma W,B}$, now free from
nucleon spinors, into Eq.(16) and use Eq.(19) to simplify the integrals.
Finally, we perform the Wick rotation using Eq.(26) and integrate out the
variable $u$ analytically. This brings us exactly to Eq.(35).

## References

  * Zyla et al. (2020) P. Zyla et al. (Particle Data Group), PTEP 2020, 083C01 (2020).
  * Seng et al. (2018) C.-Y. Seng, M. Gorchtein, H. H. Patel, and M. J. Ramsey-Musolf, Phys. Rev. Lett. 121, 241804 (2018), eprint 1807.10197.
  * Seng et al. (2019) C. Y. Seng, M. Gorchtein, and M. J. Ramsey-Musolf, Phys. Rev. D100, 013001 (2019), eprint 1812.03352.
  * Gorchtein (2019) M. Gorchtein, Phys. Rev. Lett. 123, 042503 (2019), eprint 1812.04229.
  * Hardy and Towner (2020) J. C. Hardy and I. S. Towner, Phys. Rev. C 102, 045501 (2020).
  * Czarnecki et al. (2004) A. Czarnecki, W. J. Marciano, and A. Sirlin, Phys. Rev. D 70, 093006 (2004), eprint hep-ph/0406324.
  * Seng et al. (2020a) C.-Y. Seng, X. Feng, M. Gorchtein, and L.-C. Jin, Phys. Rev. D 101, 111301 (2020a), eprint 2003.11264.
  * Märkisch et al. (2019) B. Märkisch et al., Phys. Rev. Lett. 122, 242501 (2019), eprint 1812.04666.
  * Bopp et al. (1986) P. Bopp, D. Dubbers, L. Hornig, E. Klemt, J. Last, H. Schutze, S. J. Freedman, and O. Scharpf, Phys. Rev. Lett. 56, 919 (1986), [Erratum: Phys.Rev.Lett. 57, 1192 (1986)].
  * Erozolimsky et al. (1997) B. Erozolimsky, I. Kuznetsov, I. Stepanenko, and Y. A. Mostovoi (1997), [Erratum: Phys.Lett.B 412, 240–241 (1997)].
  * Liaud et al. (1997) P. Liaud, K. Schreckenbach, R. Kossakowski, H. Nastoll, A. Bussiere, J. P. Guillaud, and L. Beck, Nucl. Phys. A 612, 53 (1997).
  * Mostovoi et al. (2001) Y. A. Mostovoi et al., Phys. Atom. Nucl. 64, 1955 (2001).
  * Schumann et al. (2008) M. Schumann, M. Kreuz, M. Deissenroth, F. Gluck, J. Krempel, B. Markisch, D. Mund, A. Petoukhov, T. Soldner, and H. Abele, Phys. Rev. Lett. 100, 151801 (2008), eprint 0712.2442.
  * Mund et al. (2013) D. Mund, B. Maerkisch, M. Deissenroth, J. Krempel, M. Schumann, H. Abele, A. Petoukhov, and T. Soldner, Phys. Rev. Lett. 110, 172502 (2013), eprint 1204.0013.
  * Darius et al. (2017) G. Darius et al., Phys. Rev. Lett. 119, 042502 (2017).
  * Brown et al. (2018) M. A. P. Brown et al. (UCNA), Phys. Rev. C 97, 035505 (2018), eprint 1712.00884.
  * Czarnecki et al. (2018) A. Czarnecki, W. J. Marciano, and A. Sirlin, Phys. Rev. Lett. 120, 202002 (2018), eprint 1802.01804.
  * Fry et al. (2019) J. Fry et al., EPJ Web Conf. 219, 04002 (2019), eprint 1811.10047.
  * Dubbers et al. (2008) D. Dubbers, H. Abele, S. Baessler, B. Maerkisch, M. Schumann, T. Soldner, and O. Zimmer, Nucl. Instrum. Meth. A 596, 238 (2008), eprint 0709.4440.
  * Wang et al. (2019) X. Wang et al. (PERC), EPJ Web Conf. 219, 04007 (2019), eprint 1905.10249.
  * Khan et al. (2006) A. A. Khan et al., Phys. Rev. D 74, 094508 (2006), eprint hep-lat/0603028.
  * Lin et al. (2008) H.-W. Lin, T. Blum, S. Ohta, S. Sasaki, and T. Yamazaki, Phys. Rev. D 78, 014505 (2008), eprint 0802.0863.
  * Capitani et al. (2012) S. Capitani, M. Della Morte, G. von Hippel, B. Jager, A. Juttner, B. Knippschild, H. B. Meyer, and H. Wittig, Phys. Rev. D 86, 074502 (2012), eprint 1205.0180.
  * Horsley et al. (2014) R. Horsley, Y. Nakamura, A. Nobile, P. E. L. Rakow, G. Schierholz, and J. M. Zanotti, Phys. Lett. B 732, 41 (2014), eprint 1302.2233.
  * Bali et al. (2015) G. S. Bali, S. Collins, B. Glässle, M. Göckeler, J. Najjar, R. H. Rödl, A. Schäfer, R. W. Schiel, W. Söldner, and A. Sternbeck, Phys. Rev. D 91, 054501 (2015), eprint 1412.7336.
  * Abdel-Rehim et al. (2015) A. Abdel-Rehim et al., Phys. Rev. D 92, 114513 (2015), [Erratum: Phys.Rev.D 93, 039904 (2016)], eprint 1507.04936.
  * Alexandrou et al. (2017) C. Alexandrou, M. Constantinou, K. Hadjiyiannakou, K. Jansen, C. Kallidonis, G. Koutsou, and A. Vaquero Aviles-Casco, Phys. Rev. D 96, 054507 (2017), eprint 1705.03399.
  * Capitani et al. (2019) S. Capitani, M. Della Morte, D. Djukanovic, G. M. von Hippel, J. Hua, B. Jäger, P. M. Junnarkar, H. B. Meyer, T. D. Rae, and H. Wittig, Int. J. Mod. Phys. A 34, 1950009 (2019), eprint 1705.06186.
  * Edwards et al. (2006) R. G. Edwards, G. T. Fleming, P. Hagler, J. W. Negele, K. Orginos, A. V. Pochinsky, D. B. Renner, D. G. Richards, and W. Schroers (LHPC), Phys. Rev. Lett. 96, 052001 (2006), eprint hep-lat/0510062.
  * Yamazaki et al. (2008) T. Yamazaki, Y. Aoki, T. Blum, H. W. Lin, M. F. Lin, S. Ohta, S. Sasaki, R. J. Tweedie, and J. M. Zanotti (RBC+UKQCD), Phys. Rev. Lett. 100, 171602 (2008), eprint 0801.4016.
  * Yamazaki et al. (2009) T. Yamazaki, Y. Aoki, T. Blum, H.-W. Lin, S. Ohta, S. Sasaki, R. Tweedie, and J. Zanotti, Phys. Rev. D 79, 114505 (2009), eprint 0904.2039.
  * Bratt et al. (2010) J. D. Bratt et al. (LHPC), Phys. Rev. D 82, 094502 (2010), eprint 1001.3620.
  * Green et al. (2014) J. R. Green, M. Engelhardt, S. Krieg, J. W. Negele, A. V. Pochinsky, and S. N. Syritsyn, Phys. Lett. B 734, 290 (2014), eprint 1209.1687.
  * Yamanaka et al. (2018) N. Yamanaka, S. Hashimoto, T. Kaneko, and H. Ohki (JLQCD), Phys. Rev. D 98, 054516 (2018), eprint 1805.10507.
  * Liang et al. (2018) J. Liang, Y.-B. Yang, T. Draper, M. Gong, and K.-F. Liu, Phys. Rev. D 98, 074505 (2018), eprint 1806.08366.
  * Ishikawa et al. (2018) K.-I. Ishikawa, Y. Kuramashi, S. Sasaki, N. Tsukamoto, A. Ukawa, and T. Yamazaki (PACS), Phys. Rev. D 98, 074510 (2018), eprint 1807.03974.
  * Ottnad et al. (2018) K. Ottnad, T. Harris, H. Meyer, G. von Hippel, J. Wilhelm, and H. Wittig, PoS LATTICE2018, 129 (2018), eprint 1809.10638.
  * Bhattacharya et al. (2016) T. Bhattacharya, V. Cirigliano, S. Cohen, R. Gupta, H.-W. Lin, and B. Yoon, Phys. Rev. D 94, 054508 (2016), eprint 1606.07049.
  * Berkowitz et al. (2017) E. Berkowitz et al. (2017), eprint 1704.01114.
  * Chang et al. (2018) C. C. Chang et al., Nature 558, 91 (2018), eprint 1805.12130.
  * Gupta et al. (2018) R. Gupta, Y.-C. Jang, B. Yoon, H.-W. Lin, V. Cirigliano, and T. Bhattacharya, Phys. Rev. D 98, 034503 (2018), eprint 1806.09006.
  * Walker-Loud et al. (2020) A. Walker-Loud et al., PoS CD2018, 020 (2020), eprint 1912.08321.
  * Aoki et al. (2019) S. Aoki et al. (Flavour Lattice Averaging Group) (2019), eprint 1902.08191.
  * González-Alonso and Martin Camalich (2016) M. González-Alonso and J. Martin Camalich, JHEP 12, 052 (2016), eprint 1605.07114.
  * Alioli et al. (2017) S. Alioli, V. Cirigliano, W. Dekens, J. de Vries, and E. Mereghetti, JHEP 05, 086 (2017), eprint 1703.04751.
  * Gonzalez-Alonso et al. (2019) M. Gonzalez-Alonso, O. Naviliat-Cuncic, and N. Severijns, Prog. Part. Nucl. Phys. 104, 165 (2019), eprint 1803.08732.
  * Falkowski et al. (2020) A. Falkowski, M. González-Alonso, and O. Naviliat-Cuncic (2020), eprint 2010.13797.
  * Sirlin (1967) A. Sirlin, Phys. Rev. 164, 1767 (1967).
  * Garcia and Queijeiro (1983) A. Garcia and A. Queijeiro, Phys. Rev. D 27, 2101 (1983).
  * Kurylov et al. (2002) A. Kurylov, M. J. Ramsey-Musolf, and P. Vogel, Phys. Rev. C 65, 055501 (2002), eprint nucl-th/0110051.
  * Kurylov et al. (2003) A. Kurylov, M. J. Ramsey-Musolf, and P. Vogel, Phys. Rev. C 67, 035502 (2003), eprint hep-ph/0211306.
  * Hayen (2021a) L. Hayen, Phys. Rev. D 103, 113001 (2021a), eprint 2010.07262.
  * Hayen (2021b) L. Hayen (2021b), eprint 2102.03458.
  * Czarnecki et al. (2019) A. Czarnecki, W. J. Marciano, and A. Sirlin (2019), eprint 1907.06737.
  * Behrends and Sirlin (1960) R. E. Behrends and A. Sirlin, Phys. Rev. Lett. 4, 186 (1960).
  * Ademollo and Gatto (1964) M. Ademollo and R. Gatto, Phys. Rev. Lett. 13, 264 (1964).
  * Holstein (1974) B. R. Holstein, Rev. Mod. Phys. 46, 789 (1974), [Erratum: Rev.Mod.Phys. 48, 673–673 (1976)].
  * Wilkinson (1982) D. H. Wilkinson, Nucl. Phys. A 377, 474 (1982).
  * Ando et al. (2004) S. Ando, H. W. Fearing, V. P. Gudkov, K. Kubodera, F. Myhrer, S. Nakamura, and T. Sato, Phys. Lett. B 595, 250 (2004), eprint nucl-th/0402100.
  * Gudkov (2008) V. P. Gudkov, Phys. Rev. C 77, 045502 (2008), eprint 0801.4896.
  * Ivanov et al. (2013) A. N. Ivanov, M. Pitschmann, and N. I. Troitskaya, Phys. Rev. D 88, 073002 (2013), eprint 1212.0332.
  * Ivanov et al. (2021) A. N. Ivanov, R. Höllwieser, N. I. Troitskaya, M. Wellenzohn, and Y. A. Berdnikov, Results Phys. 21, 103806 (2021), eprint 2010.14336.
  * Garcia and Maya (1981) A. Garcia and M. Maya, Phys. Rev. D 23, 2603 (1981).
  * Fermi (1934) E. Fermi, Z. Phys. 88, 161 (1934).
  * Sirlin (1978) A. Sirlin, Rev. Mod. Phys. 50, 573 (1978), [Erratum: Rev. Mod. Phys.50,905(1978)].
  * Seng et al. (2020b) C.-Y. Seng, D. Galviz, and U.-G. Meißner, JHEP 02, 069 (2020b), eprint 1910.13208.
  * Feng et al. (2020) X. Feng, M. Gorchtein, L.-C. Jin, P.-X. Ma, and C.-Y. Seng, Phys. Rev. Lett. 124, 192002 (2020), eprint 2003.09798.
  * Seng et al. (2020c) C.-Y. Seng, X. Feng, M. Gorchtein, L.-C. Jin, and U.-G. Meißner, JHEP 10, 179 (2020c), eprint 2009.00459.
  * Marciano and Sirlin (1993) W. J. Marciano and A. Sirlin, Phys. Rev. Lett. 71, 3629 (1993).
  * Erler (2004) J. Erler, Rev. Mex. Fis. 50, 200 (2004), eprint hep-ph/0211345.
  * Burkhardt and Cottingham (1970) H. Burkhardt and W. N. Cottingham, Annals Phys. 56, 453 (1970).
  * Marciano and Sirlin (1975) W. Marciano and A. Sirlin, Nucl. Phys. B 88, 86 (1975).
  * Carlson et al. (2008) C. E. Carlson, V. Nazaryan, and K. Griffioen, Phys. Rev. A 78, 022517 (2008), eprint 0805.2603.
  * Carlson et al. (2011) C. E. Carlson, V. Nazaryan, and K. Griffioen, Phys. Rev. A 83, 042509 (2011), eprint 1101.3239.
  * Shiells et al. (2020) K. Shiells, P. Blunden, and W. Melnitchouk (2020), eprint 2012.01580.
  * Drechsel et al. (2003) D. Drechsel, B. Pasquini, and M. Vanderhaeghen, Phys. Rept. 378, 99 (2003), eprint hep-ph/0212124.
  * Lorenz et al. (2012) I. T. Lorenz, H. W. Hammer, and U.-G. Meißner, Eur. Phys. J. A 48, 151 (2012), eprint 1205.6628.
  * Lorenz et al. (2015) I. T. Lorenz, U.-G. Meißner, H. W. Hammer, and Y. B. Dong, Phys. Rev. D 91, 014023 (2015), eprint 1411.1704.
  * Ye et al. (2018) Z. Ye, J. Arrington, R. J. Hill, and G. Lee, Phys. Lett. B 777, 8 (2018), eprint 1707.09063.
  * Lin et al. (2021a) Y.-H. Lin, H.-W. Hammer, and U.-G. Meißner, Phys. Lett. B 816, 136254 (2021a), eprint 2102.11642.
  * Lin et al. (2021b) Y.-H. Lin, H.-W. Hammer, and U.-G. Meißner (2021b), eprint 2106.06357.
  * Gerasimov (1965) S. B. Gerasimov, Yad. Fiz. 2, 598 (1965).
  * Drell and Hearn (1966) S. D. Drell and A. C. Hearn, Phys. Rev. Lett. 16, 908 (1966).
  * Deur et al. (2018) A. Deur, S. J. Brodsky, and G. F. De Téramond (2018), eprint 1807.05250.
  * Anthony et al. (1993) P. L. Anthony et al. (E142), Phys. Rev. Lett. 71, 959 (1993).
  * Abe et al. (1995) K. Abe et al. (E143), Phys. Rev. Lett. 74, 346 (1995).
  * Abe et al. (1997) K. Abe et al. (E154), Phys. Rev. Lett. 79, 26 (1997), eprint hep-ex/9705012.
  * Adams et al. (1994) D. Adams et al. (Spin Muon (SMC)), Phys. Lett. B 329, 399 (1994), [Erratum: Phys.Lett.B 339, 332–333 (1994)], eprint hep-ph/9404270.
  * Alexakhin et al. (2007) V. Y. Alexakhin et al. (COMPASS), Phys. Lett. B 647, 8 (2007), eprint hep-ex/0609038.
  * Alekseev et al. (2010) M. G. Alekseev et al. (COMPASS), Phys. Lett. B 690, 466 (2010), eprint 1001.4654.
  * Aghasyan et al. (2018) M. Aghasyan et al. (COMPASS), Phys. Lett. B 781, 464 (2018), eprint 1710.01014.
  * Ackerstaff et al. (1997) K. Ackerstaff et al. (HERMES), Phys. Lett. B 404, 383 (1997), eprint hep-ex/9703005.
  * Deur et al. (2004) A. Deur et al., Phys. Rev. Lett. 93, 212001 (2004), eprint hep-ex/0407007.
  * Wesselmann et al. (2007) F. R. Wesselmann et al. (RSS), Phys. Rev. Lett. 98, 132003 (2007), eprint nucl-ex/0608003.
  * Deur et al. (2008) A. Deur et al., Phys. Rev. D 78, 032001 (2008), eprint 0802.3198.
  * Guler et al. (2015) N. Guler et al. (CLAS), Phys. Rev. C 92, 055201 (2015), eprint 1505.07877.
  * Fersch et al. (2017) R. Fersch et al. (CLAS), Phys. Rev. C 96, 065208 (2017), eprint 1706.10289.
  * Zheng et al. (2021) X. Zheng et al. (CLAS), Nature Physics, https://doi.org/10.1038/s41567-021-01198-z, (2021), eprint 2102.02658.
  * Anthony et al. (1999) P. L. Anthony et al. (E155), Phys. Lett. B 458, 529 (1999), eprint hep-ex/9901006.
  * Anthony et al. (2003) P. L. Anthony et al. (E155), Phys. Lett. B 553, 18 (2003), eprint hep-ex/0204028.
  * Amarian et al. (2004) M. Amarian et al. (Jefferson Lab E94-010), Phys. Rev. Lett. 92, 022301 (2004), eprint hep-ex/0310003.
  * Kramer et al. (2005) K. Kramer et al., Phys. Rev. Lett. 95, 142002 (2005), eprint nucl-ex/0506005.
  * Bjorken (1966) J. D. Bjorken, Phys. Rev. 148, 1467 (1966).
  * Bjorken (1970) J. D. Bjorken, Phys. Rev. D 1, 1376 (1970).
  * Ji and Unrau (1994) X.-D. Ji and P. Unrau, Phys. Lett. B 333, 228 (1994), eprint hep-ph/9308263.
  * Baikov et al. (2010a) P. A. Baikov, K. G. Chetyrkin, and J. H. Kuhn, Nucl. Phys. B Proc. Suppl. 205-206, 237 (2010a), eprint 1007.0478.
  * Baikov et al. (2010b) P. Baikov, K. Chetyrkin, and J. Kuhn, Phys. Rev. Lett. 104, 132004 (2010b), eprint 1001.3606.
  * Deur et al. (2014) A. Deur, Y. Prok, V. Burkert, D. Crabb, F. X. Girod, K. A. Griffioen, N. Guler, S. E. Kuhn, and N. Kvaltine, Phys. Rev. D 90, 012009 (2014), eprint 1405.7854.
  * Kotlorz et al. (2017) D. Kotlorz, S. V. Mikhailov, O. V. Teryaev, and A. Kotlorz, Phys. Rev. D 96, 016015 (2017), eprint 1704.04253.
  * Ayala et al. (2018) C. Ayala, G. Cvetič, A. V. Kotikov, and B. G. Shaikhatdenov, Eur. Phys. J. C 78, 1002 (2018), eprint 1812.01030.
  * Wandzura and Wilczek (1977) S. Wandzura and F. Wilczek, Phys. Lett. B 72, 195 (1977).
  * Shuryak and Vainshtein (1982) E. V. Shuryak and A. I. Vainshtein, Nucl. Phys. B 201, 141 (1982).
  * Jaffe (1990) R. L. Jaffe, Comments Nucl. Part. Phys. 19, 239 (1990).
  * Alarcón et al. (2020) J. M. Alarcón, F. Hagelstein, V. Lensky, and V. Pascalutsa, Phys. Rev. D 102, 114026 (2020), eprint 2006.08626.
  * Bolognese et al. (1983) T. Bolognese, P. Fritze, J. Morfin, D. H. Perkins, K. Powell, and W. G. Scott (Aachen-Bonn-CERN-Democritos-London-Oxford-Saclay), Phys. Rev. Lett. 50, 224 (1983).
  * Allasia et al. (1985) D. Allasia et al., Z. Phys. C28, 321 (1985).
  * Acciarri et al. (2016) R. Acciarri et al. (DUNE) (2016), eprint 1601.05471.
  * Alvarez-Ruso et al. (2018) L. Alvarez-Ruso et al. (NuSTEC), Prog. Part. Nucl. Phys. 100, 1 (2018), eprint 1706.03621.
  * Seng and Meißner (2019) C.-Y. Seng and U.-G. Meißner, Phys. Rev. Lett. 122, 211802 (2019), eprint 1903.07969.
  * Ma et al. (2021) P.-X. Ma, X. Feng, M. Gorchtein, L.-C. Jin, and C.-Y. Seng (2021), eprint 2102.12048.
  * Borsanyi et al. (2020) S. Borsanyi et al. (2020), eprint 2002.12347.
  * Aoyama et al. (2020) T. Aoyama et al., Phys. Rept. 887, 1 (2020), eprint 2006.04822.
  * Hayen and Severijns (2019) L. Hayen and N. Severijns (2019), eprint 1906.09870.