An updated view on the ATOMKI nuclear anomalies

In view of the latest experimental results recently released by the ATOMKI collaboration, we critically re-examine the possible theoretical interpretation of the observed anomalies in terms of a new BSM boson $X$ with mass $\sim17\;$MeV. To this end we employ a multipole expansion method and give an estimate for the range of values of the nucleon couplings to the new light state in order to match the experimental observations. Our conclusions identify the axial vector state as the most promising candidate, while other spin/parity assignments seems disfavored for a combined explanation. This results is however based on an order of magnitude estimate for the, currently unknown, axial nuclear matrix element of the $^{12}$C transition, that needs then to be evaluated before being able to draw a definite conclusion. Intriguingly, an axial vector state can also simultaneously accommodate other experimental anomalies, {\emph{i.e.}} the KTeV anomaly in $\pi^0 \to e^+ e^-$ decay while being compatible with the conflicting measurements of the anomalous magnetic moment of the electron $(g-2)_e$ and other constraints on the electron couplings of the $X$ boson. The PADME experiment will completely cover the relevant region of the parameter space, thus allowing for a strong test of the existence of the $X$ particle.


Introduction
The possibility that New Physics (NP) will manifest itself in the form of light and weakly coupled new states is nowadays raising more and more interest. This is mainly due to the null results from the LHC in searching for signs of TeV scale beyond the Standard Model (BSM) physics motivated by, e.g., the standard paradigms of compositeness or supersymmetry. While high−p T searches will continue to investigate these scenarios and, in the case of no discovery, further constrain them, it is of paramount importance to exploit any other experiment that can test the validity of the SM and ultimately falsify it. In this respect low energy and/or high intensity experiments provide an important probe.
Among the various processes that can be investigated, rare nuclear transition can provide a good handle to observe NP appearing at the MeV scale, since they can significantly be affected by BSM physics even if this is very weakly coupled. A nuclear transition occurs when an excited nucleus decays into a lower energy level of the same nucleus. Within the SM only electromagnetic (EM) interactions can mediated nuclear transition, which can mainly proceed through the following channels • γ−emission, where the nucleus decays emitting a real photon, • Internal Pair Creation (IPC), where the nucleus emits a virtual photon which then decays to an e + e − pair.
In recent years the ATOMKI collaborations has reported various anomalous measurements in the IPC decays of excited 8 Be [1,2], 4 He [3,4] and, more recently, 12 C [5] nuclei. These anomalies appear as bumps for both the invariant mass and the angular opening of the e + e − pairs and have a high statistical significance, well above 5σ. The ATOMKI collaboration has proposed to interpret them as due to the on-shell emission of a new boson X from the excited nuclei, subsequently decaying to an e + e − pair. The best fit mass for the hypothetical new particles is estimated to be ∼ 17 MeV. Although to this day no independent confirmation of these results has arrived, given the multitude of processes in which these anomalies have been observed the ATOMKI results have attracted a considerable attention from the particle physics community. Many theoretical interpretation of the X boson in terms of a new scalar or vector degree of freedom have been put forward, possibly unaccounted SM effects have been investigated and experimental searches have been proposed and/or are taking data with the goal of further investigating this anomaly.
In view of the latest experimental results recently released by the ATOMKI collaboration, we critically re-examine the possible theoretical interpretation of the anomaly in terms of a new BSM state. To this end we employ a multipole expansion method and give an estimate for the range of values of the nucleon couplings to the new light state in order to match the experimental observations. We will focus on the 8 Be and 4 He anomalies and comment on how the measurement of the anomalous signal in 12 C transitions impact our results. Our conclusions identify the axial vector state as the most promising candidate, while other spin/parity assignments seems disfavored for a combined explanation. However the axial nuclear matrix element of the 12 C transition is currently unknown and, as we will show, our findings regarding the compatibility of an axial vector candidate with the 12 C anomalous transition are based upon an order of magnitude estimate. Before being able to draw a definite and solid conclusion, the relevant matrix element must be evaluated. Intriguingly, an axial vector state can also simultaneously accommodate other experimental anomalies, i.e. the KTeV anomaly in π 0 → e + e − decay while being compatible with the conflicting measurements of the anomalous magnetic moment of the electron (g − 2) e and other constraints on the electron couplings of the X boson. The PADME experiment will completely cover the relevant region of the parameter space, thus allowing for a strong test of the existence of the X particle.
The paper is structured as the following. In Sec. 2 we review the anomalous measurement of the ATOMKI experiment and the theoretical interpretation proposed so far. In Sec. 3 we describe the multipole expansion formalism while in Sec. 4 we present the results of the decay rates of the 8 Be, 4 He and 12 C resonances for the various spin-parity assignment of the X boson. We present our results in Sec. 5 and we then conclude in Sec. 6. We also add some appendices with more technical details. In App. A we report some tables with the numerical values used throughout our analysis. In App. B we apply the multipole formalism to the electromagnetic case of real γ emission and IPC, while in App. C we derive the effective nuclear couplings from the X boson couplings to quarks. In App. D we present how our results change by considering also the 8 Be(17.64) excited state. In App. E we report useful formulae for the cross section of nuclear resonance production. Finally in App. F we collect the experimental bounds for the spin-1 case relevant for our study.

The ATOMKI anomaly
The ATOMKI experiment [6] consists in a proton beam colliding a target nucleus A at rest, with the aim of producing an excited nucleus N * and measure its IPC transition to a ground state N , i.e. p + A → N * → N + e + e − . (2.1) The list of nuclei used by ATOMKI and their main properties are reported in Tab. 5 and Tab. 6 in App. A. In 2015 the ATOMKI group studied the IPC decay channel from the 8 Be(18.15) and 8 Be(17.64) excited energy levels of Beryllium nuclei [1]. To populate the two states, a beam of protons was prepared in order to collider with target 7 Li nuclei at rest. By varying the energy of the incident proton beam, the collaboration was able to scan across the 8 Be resonances. As a results they observed an anomalous peak corresponding to an opening angle for the e + e − pairs of ∼ 140 • for the IPC correlation distribution, in striking contrast with the QED prediction of a rapidly falling one. This has been interpreted as due to the decay of a short-lived neutral particles decaying into an e + e − pair, which would produce the observed peak at large angles. The observed deviation had a significance of 6.8σ. In order to confirm the anomalous origin of the signal, the collaboration repeated the measurement varying the energy of the incident proton beam. They found that the anomaly disappeared off the resonance peak, leading to the conclusion that it was probably originated by the decay of the 8 Be excited energy level. The best fit mass for N N * S π I Γ(keV) Γ γ (eV) 8 Be 0 + 0 5.57 ± 0.25 8 Be(18.15) 1 + 0 * 138 ± 6 1.9 ± 0.4 8 Be ( Table 1. Spin-parity J π and isospin I quantum numbers, total decay widths Γ and γ-decay widths Γ γ = Γ(N * → N γ) for the nuclei used in the ATOMKI experiment: 8 Be [7], 4 He [8,9] and 12 C [10,11] nuclei. Asterisks on isospin assignments indicate states with significant isospin mixing.
the hypothetical X neutral boson has been estimated to be m X = 16.70 ± 0.35 ± 0.5 MeV, where the former uncertainty corresponds to the statistical error, while the latter to the systematic one. In 2018 the ATOMKI collaboration repeated the experiment with an improved setup [2], which confirmed both the presence of the anomaly and the compatibility with the previous measurement. Later in 2019 the ATOMKI group replicated the experiment using 4 He nuclei [3], with the aim of searching for the anomalous signal in a difference source. The excitation energy was chosen to lie between two different resonances: the 4 He(20.21) state, with J π = 0 + , and the 4 He(21.01) state, with J π = 0 − . In this case the decay widths of the two excited states are large enough so that they can substantially overlap, so in the experiment both the excited states were populated, although off the resonance peak. In this case the IPC process was only possible for the 4 He(20.21) states, while it's forbidden for the 4 He(21.01), because of parity conservation in electromagnetic interactions. The group again found a rather sharp bump in the e + e − angular opening analogously to what has been observed in the 8 Be case, with a significance equal to 7.2σ. Interestingly, the peak was found to be located at an angle of ∼ 115 • which is compatible with the kinematics arising from the decay of the hypothetical X boson with a best fit mass of m X = 16.98 ± 0.16 ± 0.2 MeV. Also the anomaly in the 4 He channel was later confirmed by a second measurement, at different energies of the proton beam [4]. More recently, the group has released a new analysis, where the experiment has been replicated using now 12 C nuclei [5]. Also in this case an anomalous signal has been observed, with a peak at a larger value of the e + e − opening angle ∼ 150 • − 160 • , again compatible with the kinematic of the X particle, with a best fit mass of m X = 17.03 ± 0.11 ± 0.20 MeV. In Tab. 1 we list all the ground and excited states considered in the ATOMKI analyses, together with their main properties: spin-parity assignment S π , isospin I, total decay width Γ and γ-decay transition width Γ γ .
Clearly, new and independent measurements are needed in order to confirm, or disproof, the results of the ATOMKI collaboration and test the consistency of the X particle hypothesis. The MEG II experiment [12] at PSI has the possibility to repeat the ATOMKI measurement on 8 Be nucleus. At present time, the first dedicated data taking has recently been completed and data analysis is ongoing [13]. A similar experiment is also being set up at the Montreal Tandem accelerator [14] with data taking that should take place in early 2023, and at the Van-de-Graaff laboratory [15]. Finally, also the PADME experiment [16,17] in Frascati is in its data taking phase, and dedicated analyses will completely test the available parameter space for the coupling of the X boson to electrons relevant for the explanation of the ATOMKI anomalies [18,19].
These results and the interpretation given by the ATOMKI collaboration in terms of a new BSM particle comprehensibly attracted the attention of the theory community. However there also exist the possibility that the anomalous signal is due to unknown and/or underestimated SM effects. In this respect, after the publication of the 8 Be measurements, an attempt has been made in order to explain the anomaly with effects arising from known nuclear physics. In [20] it has been proposed an improved nuclear physics model of the experiment, inspired by the so-called Halo Effective Field Theory (EFT) framework [21], showing that the nuclear form factor needed to explain the anomaly suggests an unrealistic large length scale on the order of 10 fm for the 8 Be nucleus. Other Authors investigated the possibility of new exotic bound states [22][23][24] as possible SM explanation of the ATOMKI anomaly while, on less exotic lines, it has been claimed in [25] that the experimental results can be reproduced within the SM by carefully considering the full set of next-toleading-order corrections and the interference terms to the Born-level decay amplitudes, also proposing experimental improvements in order to test this hypothesis. All together it is fair to say that no firm explanation of the ATOMKI measurements in terms of SM effects has been established.
As regarding possible BSM interpretation, the observation of the 8 Be and 4 He transitions restrict the X boson to be either a vector, an axial vector or a pseudoscalar state, under the assumption of definite parity. These options have been all investigated in recent literature both from a model independent and/or effective parameterization and from a more ultraviolet (UV) completed perspective. Among all the possibilities the one of a spin-1 boson stands out as an appealing one, since it could be related to a new symmetry of Nature. This scenario has been deeply analyzed in [26][27][28], where, by working within an EFT framework, it has been found that a combined explanation of the 8 Be and 4 He anomalies in terms of a new vector states is possible, with the main constraint on this explanation coming from the search for a dark photon γ D in π 0 → γ γ D decay by the NA48/2 collaboration [29], whose non observation requires the X boson to be protophobic. Subsequently in [30] it has been pointed out that the contribution from a protophobic vector boson with mass around 17 MeV to direct proton capture processes, i.e. to processes which do not proceed through an intermediate resonance, would be dominant with respect to the contribution from the resonant 8 Be(18.15) state, in sharp contradiction with the experimental observation that the anomaly disappears off the nuclear resonance [1]. However a new experimental result from ATOMKI [31] claims to have observed the anomaly at different energies of the proton beams, opposite to their previous results. The group explains that this difference is due to a wrong estimate of the background in the previous analyses and then reopens the window for a protophobic vector scenario. The case of a pure axial vector has been investigated in [32] where the authors applied a multipole expansion method to the anomalous nuclear decay rates and evaluated the related nuclear matrix elements by ab-initio calculation using realistic nuclear forces. Their estimation concludes that the X axial coupling to quark should be order O(10 −4 − 10 −5 ) to explain the anomalous signal in the 8 Be transition. Vector with mixed parity have also attracted attention, especially in the case of more specific BSM UV construction, see e.g. [33][34][35]. Lastly, the possibility of a light pseudoscalar has been considered in [36], where the authors made a rough estimation of the range of the values of the Yukawa couplings, assuming a nuclear shell model for the 8 Be nucleus. The strongest constraints they reported come from flavor changing neutral current interactions as K → πX, which however can be satisfied simultaneously explaining the ATOMKI 8 Be results. The interesting possibility of the QCD axion being responsible for the ATOMKI anomalies has been entertained in [37,38], see also [39]. Here the Authors focus on a axion candidate with dominant coupling to the first generation of SM fermions and piophobic, e.g. with suppressed isovector coupling. It has subsequently however been pointed out in [40] that for such a scenario a large pion decay rate for π → 3X → 3e + 3e − of O(10 −3 ) is expected, exceeding the SM double-Dalitz decay by a factor of thirty. There has been no direct measurement of this process so far, although it is reasonable to assume that such a large decay rate would have been noticed.

Process kinematics
The ATOMKI anomalies show simple but well defined features, which are: • the excesses are resonant bumps located at the same e + e − invariant mass for all the 8 Be and 4 He transitions, • the e + e − opening angles of the anomalous peaks are around 140 • , 115 • and 155 • − 160 • , respectively, for the 8 Be, 4 He and 12 C, • the anomalous signal in the 8 Be transition have been observed only inside the kinematic region given by |y| < 0.5, where y is the energy asymmetry of the lepton pair, i.e. the ratio between the difference and the sum of their energies.
As we review below, closely following earlier results appeared in [27,28], these features are naturally explained by the hypothesis of resonant production of a new particle.
In the experimental setup the target nucleus A is at rest in the laboratory frame, while the proton beam energy E b is of the order of MeV, so that the colliding protons are mostly non relativistic. The Center of Mass (CM) energy E CM is then given by where m pA = (m −1 A + m −1 p ) −1 is the reduce mass of the proton-target system. The experiment calibrates the beam energy in order to populate the N * state, which is produced almost at rest in the CM frame, and then measure its IPC transition to the ground state N . Regardless on whether the resonance is (fully) populated, the collision between the proton and the target leads to the production of the N nucleus via the emission of a boson of mass m, which could be a real or virtual photon or an hypothetical BSM particle. In the CM frame the boson energy ω is given by where the threshold energy E th = m p + m A − m N is the energy gap between the N nucleus and the proton-target system. Note that the boson energy roughly only depends on the beam energy. Once produced the (real or virtual) boson decays into an e + e − pair, whose angular correlation is the main observable measured by the ATOMKI experiment. In the CM frame the total energy of the leptons is given by is the positron (electron) energy, while we label the opening angle between the leptons with θ ± . The decay of the boson into the e + e − pair is controlled by the energy asymmetry while the opening angle is given by where δ = 2m e /ω (0 < δ < 1) and is the boson velocity. For an hypothetical X boson with mass m X = 17 MeV, the maximum value of the energy asymmetry and the minimal value of the opening angle are respectively In Fig. 1 we plot the opening angle as a function of the energy asymmetry for different values of the X boson mass around 17 MeV for all the transitions studied by ATOMKI. Note that the signal region for the 8 Be case is all contained in |y| ≤ 0.5 in agreement with the ATOMKI experiment. We then show in Fig. 2 the normalized distribution of θ ± , integrated over the asymmetry y, for the three nuclei considered by ATOMKI in the spin-0 boson hypotheses case, where the distribution only depends on phase space quantities. In the spin-1 case there is a dynamical dependence due to the polarization state of the X boson which can modify this distribution, which is however expected to be qualitatively similar to the spin-0 case. As we see the opening angle distributions peak at the lowest end with peak values compatible with the ones reported by the experiments. We further report in Tab. 7 in App. A the predicted values for the opening angle at the signal peak for the various transitions investigated by ATOMKI. Lastly, from Eq. (2.5) the invariant where here cos θ ± is now a free parameter independent on the energy asymmetry.

Signal computation: overview
In this section we describe the multipole expansion formalism used in order to estimate the anomalous nuclear decay rates relevant for the ATOMKI experiment, see also [32].

Nuclear states and processes
We describe the interaction of a spin s X boson to nuclear matter through the Hamiltonian where the nuclear scalar density S and the nuclear current J µ = (J 0 , ⃗ J ) are quantum operators containing all the information of the matter fields. In the case of electromagnetic interaction the nuclear current J µ is replaced by the electromagnetic current J γ µ , which allows for γ−emission and IPC processes. At the lowest order in the interaction picture, the nuclear matrix element for the X emission is given by where ⃗ k is the boson momentum, the index a = 0, ±1 labels the polarization of the vector boson and |i⟩ and |f ⟩ indicate the nuclear matter initial and final states which are |f ⟩ = |N ; J f M f ⟩ and |i⟩ = |p + A; J p M p ; J A M A ; ⃗ p CM ⟩, where ⃗ p CM is the proton momentum in the CM frame. We employ the narrow width approximation and factorize the excited resonance production from its decay, thereby assuming that the initial state is described by |i * ⟩ = |N * ; J * M * ⟩. In the following we want to compute the decay widths of the excited N * states for real X emission in order to compare to the experimental results on this quantity reported by the ATOMKI collaboration. To make the calculation, it will turn out to be useful to expand the nuclear matrix elements in terms of spherical tensor operators, which will allow to also compute the electromagnetic real γ-emission and IPC processes, which we report in App. B. The diagrams for these three processes are shown in Fig. 3. To perform the calculation we will expand the nuclear matrix elements in terms of spherical tensor operators through a multipole expansion.

Multipole expansion
Spherical operators O JM are irreducible tensor operators which satisfy the Wigner-Eckart theorem [41] ⟨J The reduced matrix element ⟨J f ||O J ||J i ⟩ contains all the physical information of the operator while its behavior under rotation is completely set by the Clebsh-Gordan coefficient where r = |⃗ r| and expand the nuclear matrix elements as a sum of reduced matrix elements of spherical operators. In the case of interest of the emission of the X boson in the process where here the indices J and M denote the total angular momentum of the emitted boson (sum of its spin and relative angular momentum with the N nucleus) and its projection. The rotation D matrices play the role of the wave function 2 , whose moduli squared give the probability for the X boson to be emitted in the (ϕ, θ) direction with β defining a rotation along this direction 3 . An explicit calculation for the unpolarized decays gives (3.12) 1 YJM , jJ (x) and Y JℓM (r) are respectively the spherical harmonics, the vector Bessel functions and the vector spherical harmonics, see [41]. In the case where the vector boson is coupled to a conserved current, i.e. ∂ µ J µ = 0, a simplification occurs. By assuming the nuclear initial and final state to be eigenstates of the nuclear Hamiltonian the continuity equation 13) and the partial width for the emission of a vector X boson then reduces to (3.14) The above results are equally useful for the electromagnetic processes once we substitute the electromagnetic current in the spherical operators and put m γ = 0, i.e. k = ω.

Selection rules
The angular momentum conservation law, encoded in the Clebsh-Gordan coefficient of Eq. (3.3), states that the matrix element of the spherical operators vanishes unless the following conditions are satisfied Moreover, if the X boson has a definite parity π X , additional constraints on the matrix elements come from the requirement of parity conservation. By denoting the relative angular momentum between the boson and N as L, one has We report in Tab. 2 the relative angular momentum between the X boson and N in the various decay processes, based on the S π spin-parity assignments. One sees that a pure scalar solution to the 8 Be anomaly is excluded, while a pseudoscalar state can explain only the 8 Be and 4 He anomaly, if the latter is dominated by the 4 He(21.01) excited state transition, but not the 12 C one. On the other side a vector or axial-vector candidate can simultaneously explain all the three anomalies, but again only one of the two 4 He resonant states can contribute to the signal process.

Long wavelength approximation
The nuclear radius is approximately given by [42], which implies that in all the cases of interest the nucleus size is significantly smaller than the boson wavelength k −1 ∼ (10 MeV) −1 . We can thus expand the spherical Bessel function for small kr as  Table 2. Relative angular momentum between the X boson and N in the various decays, based on its possible parity-spin assignments. Note that parity conservation prohibits a pure scalar solution to the Beryllium anomaly.
with higher order corrections giving a contribution of order (kr) 2 ≃ 1% with respect to the leading one for the cases of interest, which can therefore be neglected.
An exception occurs for the monopole case L 00 , since it identically vanish at this order. The first contribution thus arises at the next order in the kr expansion and is given by The second order expansion is also needed for the M 00 monopole expression in the case of a conserved current. This is due to the fact that in this case the integral over space of J 0 (⃗ r) defines the generator Q of the symmetry associated with it. Then, with |i⟩ and |f ⟩ orthogonal eigenstates of the Hamiltonian, one as ⟨f |Q|i⟩ ∝ ⟨f |i⟩ = 0. It follows that the first contribution to the operator M 00 is given by Let's now consider the operators T el JM and T mag JM of Eq. (3.7) and Eq. (3.8). By using the identity

23)
Spin-1 case Table 3. Leading term of the non relativistic expansion for the relativistic vector current, the relativistic axial current and the anomalous magnetic moment terms (upper table) and for the scalar and pseudoscalar density (lower table). δ ⃗ r,⃗ rj = δ(⃗ r − ⃗ r j ).
they can be rewritten at the first order in the kr expansion as [42] T mag where the vector current has been split into an irrotational field ⃗ J irr and a solenoidal field ⃗ ∇ × ⃗ µ as ⃗ J = ⃗ J irr + ⃗ ∇ × ⃗ µ in virtue of the Helmholtz's theorem. For a conserved current, the matrix element expression of T el J can be simplified to (3.26) by again using the continuity equation.

Non relativistic expansion for nuclear operators
Through statistical considerations [43], the maximal kinetic energy E c per nucleon in the nucleus is estimated to be around 30 MeV, implying that a nucleus can be then modeled as a quantum mechanical system of non relativistic point-like nucleons [42]. One can then take the non relativistic limit of the nuclear operator and write it in first quantization formalism. The nuclear operator is given by with the single particle operatorÔ in the nucleon point-like approximation. We want to match the expression of Eq. (3.27) with its relativistic counterpart, where the nucleons are described in terms of quantum fields p(x) and n(x) and the nucleon operators are bilinears in p(x) and n(x). We report in Tab. 3 the leading terms of the non relativistic expansion for the vector current, axial current, anomalous magnetic moment, scalar and pseudoscalar densities. For operators which are even or odd under parity, as it is in our case, higher order terms in the expansions are of order of respect to the leading one, and can then be safely neglected. For the specific case of the pseudoscalar density one also has that where u are the spinors which enter the quantum field expression. Hence the non relativistic expansion of the pseudoscalar current only contains terms proportional to ⃗ k = ⃗ p − ⃗ p ′ . Since for the effective non relativistic operator it holds the substitution [42] i ⃗ ∇ → ⃗ k (3.29) the non relativistic expansion of the pseudoscalar density only contains operators given by a total divergence, meaning that the pseudoscalar density is a derivative coupling. The monopole operator G 00 then vanishes at first order in the long-wavelength expansion and the leading contribution is then given by The techniques introduced in this section will be implemented in the next to derive the theoretical decay rates of the nuclear transitions.

Signal computation: X dynamics
With the formalism described in the previous section we are now ready to describe the BSM dynamics of e + e − emission from the X boson, illustrated in Fig. 3. We refer to App. B for the details of the SM processes of real γ emission and IPC. We parametrize the interaction of the X boson with the scalar density S for the spin-0 case and the nuclear current J µ for the spin-1 cases in terms of effective couplings as 4,5 4 Effective nucleon operators are in principle also function of form factors F (q 2 ). In all practical cases however the transferred momentum is much smaller than the hadron scale ΛQCD, so that we approximate the form factors as constants. 5 We neglect electric dipole moment (EDM) operators since, by naive counting analysis, they would contribute at a higher order than the axial currentψγ µ γ 5 ψ. Moreover, they can only be generated by EDM effective quark operators, so that they will generally be suppressed by loop effects.
see App. C. Although the pure scalar hypothesis is not able to explain the anomaly observed in the 8 Be decay, we'll present for completeness explicit expressions also in this case, since it can anyway affect the 4 He and 12 C decays, see Tab. 2 and it can be relevant in the mixed parity hypothesis. The effective matching between these effective interactions and the UV interactions of the X boson with quark and gluons are reported in App. C. Once the X boson is produced from the nuclear collision, it decays to an e + e − pair with a branching ratio which depends on the size of the X coupling to electrons.

Spherical operators
In the non relativistic and long wavelength approximation the spherical operators with J = 0, 1 for the various S π assignment for the X boson are given by Scalar case S π = 0 + For the J = 0 term we go here beyond the leading order in the long wavelength and non relativistic expansion since leading order term proportional to the identity gives an identically vanishing contribution. One has 6 where z s = z p (z s = z p ) if the s-th nucleon is a proton (neutron). A similar notation is adopted for the rest of the section.
Pseudoscalar case S π = 0 − For the pseudoscalar case the relevant spherical operators are where for G 00 we have split the expression among the isoscalar and isovector contributions.
In this case conservation of the vector current implies a relation between the operator L JM and M JM , c.f.r. Eq. (3.13), and one has Axial vector case S π = 1 + Finally the spherical operators for the axial vector case read (4.17)

Decay rates
We can now express the decay rates for the various spin-parity assignment of the X boson in the case of the 8 Be, 4 He and 12 C transitions. We report them Tab. 4, expressed in function of nuclear matrix element of the relevant operators involved in the transition. Symmetry consideration allow to express these matrix elements in function of known ones. We list here the relevant relations 8 Be He(20.21) Table 4. Decay rates for the 8 Be, 4 He and 12 C nuclear processes for the various spin assignment of the X boson. In the case of the Helium transition for each spin-parity possibility we indicate the 4 He excited state involved. In the expressions ⟨O⟩ represents the matrix element between the ground state and the excited nucleus of the corresponding operator, e.g. for the 12 C transition in the 0 + case ⟨||d (γ) ||⟩ = ⟨ 12 C||d (γ) || 12 C(17.23)⟩. For the 4 He and 12 C cases we only report the non vanishing isoscalar and isovector contributions respectively.

Beryllium matrix elements
Assuming the static quark model κ p ≃ −κ n ≃ 2(C p − C n ), see App. C for the details.
From isospin symmetry then one has where the corresponding values are reported in App. B, while we take from [32] ⟨Be||σ (

Helium matrix elements
Isospin symmetry allows to relate where again the electromagnetic matrix element is reported in App. B while we take from [44] | ⟨ 4 He||d σ 0 || 4 He(21.01)⟩ | 2 ≃ 15.5 fm 2 ≃ 4 × 10 −4 MeV −2 (4. 22) but no uncertainty has been given. We will arbitrarily assume a 10% error on the matrix element in our calculation. Finally, to the best of our knowledge, the matrix element of the operatorK = s p 2 s /2m s has never been evaluated so far.

Carbon matrix elements
Isospin symmetry allows to relate whose values is reported in App. B, while the axial matrix elementD σ 3 has not been evaluated, to the best of our knowledge.

Results
We present in this section our main results, deriving the possible range of the nucleon couplings to the X particle that can explain both the 8 Be and 4 He anomalies, further commenting on the possibility of simultaneously explain the 12 C one. We analyze all the scenarios where the X boson has a definite parity, which implies that the pure scalar boson case is ruled out since it cannot explain the 8 Be anomaly. The best fit value for the anomalous decay rate for the 8 Be transition is [2] The ATOMKI collaboration observed no anomalous signal in the 8 Be(17.64) transition in the first experiment [1] but later they reported a non vanishing best fit for this anomalous decay rate in a contribution to the proceedings of International Symposium Advances in Dark Matter and Particle Physics 2016 [45]. In the following we will consider only the We now present our findings for the regions in the effective nucleon couplings parameter space for the various spin-parity assignments for the X boson. In presenting our results we assume, for simplicity, BR(X → e + e − ) = 1. For different BR assumptions the derived allowed space in the nucleon effective couplings will be rescaled according to Eq. (5.1), Eq. (5.4) and Eq. (5.5). We stress that our analysis relies on various assumptions, as the hypothesis of narrow width approximation for the nuclear production of the excited states. Other potential contribution, as for example direct capture processes, could potentially change the conclusions of our analysis.

Pseudoscalar and mixed parity scenario
We summarize the results for the pure pseudoscalar scenario in Fig. 4, where the shaded blue and orange areas represent the 1σ and 2σ compatibility regions with the ATOMKI 8 Be and 4 He anomalies respectively. We also overlay in red the region of parameter space satisfying the SINDRUM bound from π + → e + ν e X decay [37,46]. This is given in term of the pseudoscalar-pion mixing angle, linked to the isovector nucleon coupling as where g A ∼ 1.27 is axial nucleon factor and f π ∼ 93 MeV is the pion decay constant, and reads All together we see that the 8 Be and 4 He anomalies can be simultaneously satisfied for a range of effective nuclear coupling h n,p of O(10 −2 ). However, the recent observation of an anomalous signal in the 12 C transition would, if confirmed, exclude by itself the pure pseudoscalar scenario, see again Tab. 2. It's then interesting to entertain the possibility that the scalar X boson has both scalar and pseudoscalar couplings. As already mentioned, because of parity conservation the scalar contribution to the 8 Be transitions vanishes, so the latter processes only set a constraint on the range of the pseudoscalar couplings, which as we have shown in Fig. 4   so that the theoretical predictions for the 4 He transition is too large to match the ATOMKI results. We expect that this assertion holds even once ⟨ 4 He||K|| 4 He(20.21)⟩ is also included. Hence, one is forced to conclude that the scalar isoscalar coupling z p + z n is suppressed, at least respect the pseudoscalar one, leaving us with almost the same configurations as the pure pseudoscalar one of Fig. 4. This conclusion is in agreement with earlier results [28]. On the other side for the 12 C transition, the scalar isovector coupling would give in agreement with the order of magnitude of the ATOMKI fit (5.6) if z p − z n ≃ h p,n ≃ O(10 −2 ) and all the three anomalous measurements can be simultaneously satisfied.

Vector and axial vector scenarios
We summarize the results for the spin-1 cases in Fig. 5, with the same color code as Fig. 4 for the regions satisfying the ATOMKI anomalies. In the upper panels we show the results for the S π = 1 − assignment for the X boson. Here in the left and right plot we assume Regions of the C n,p effective nuclear couplings of a pure vector state where the 8 Be (blue), 4 He (orange) and 12 C (purple) anomalous ATOMKI transition can be explained at 1σ or 2σ. Inside the red and the gray region, respectively, the NA48 and the 208 Pb-n scattering bound are satisfied. In the left and right panel we assume ξ = 0 and ξ = 0.549 respectively, see App. B for details. Lower panels: Regions of the a n,p effective nuclear couplings of a pure axial vector state where the 8 Be (blue) and 4 He (orange) anomalous ATOMKI transition can be explained at 1σ or 2σ. In the green region the KTeV anomaly in π → e + e − decay can be satisfied, by assuming a positive (left panel) and negative (right panel) value for the C e A axial coupling of the X boson to electrons that can explain the anomalous (g − 2) e , see main text for more details. ξ = 0 and ξ = 0.549 respectively, where ξ represents non perturbative contribution to isospin breaking effects in the 8 Be case, see [27] and App. B for details. For the vector case the strongest bound comes from the non observation from the NA48 experiment of the π 0 → γX decay in dark photon searches [47]. This process receives two different contributions. One from the axial anomaly and a non anomalous one. The non anomalous contributions is proportional to the small quark masses [48,49] and can be neglected, while the anomalous one is proportional to the anomaly trace factor. One gets the bound |C p | × BR(X → e + e − ) ≲ 2.5 × 10 −4 , (5.12) which implies a protophobic nature for the X boson. We show in Fig. 5 in red the region of parameter space where the protophobia constraint is satisfied. Another relevant bound comes from observations of the angular dependence of neutron-lead scattering. The exchange of new, weakly-coupled boson produces a Yukawa potential 8 acting on the neutron, whose contribution has been constrained for the 208 Pb-n scattering as [50] |C n | 126 208 C n + 82 208 We show in Fig. 5 in gray the region of parameter space where the 208 Pb-n scattering constraint is satisfied. As it can be seen for the ξ = 0.549 assignment, a combined explanation of the 8 Be, blue region, and 4 He, orange region, anomalies at 1σ is in tension with the NA48 constraint, while it is possible at the 2σ level. Regarding the 12 C ATOMKI anomaly in the case of a S π = 1 − state the relevant matrix element is known. In this case one then gets Γ( 12 C(17.23) → 12 C + X) in agreement with the order of magnitude of the ATOMKI fit (5.6) if C p − C n ≃ O(10 −3 ). The 1σ and 2σ bands related to the 12 C transition are shown in purple in the upper panels of Fig. 5. Note that, if confirmed, the 12 C ATOMKI anomaly is in tension with a combined explanation of the 8 Be and 4 He anomalies and the protophobia constraint.
On the other side an axial vector S π = 1 + state can explain both the 8 Be and 4 He ATOMKI anomalies, as shown in the lower panels of Fig. 5, with axial couplings to the nucleon of O(10 −4 ). Within the green shaded area the KTeV anomaly in π 0 → e + e − decay can be explained for positive and negative values for the axial X coupling to electrons C e A , see Sec. 5.2.1 for details. As regarding the possibility of also explaining the 12 C ATOMKI anomaly the relevant nuclear matrix element, see Tab. 4, is currently unknown. While no definite claim can be made until it becomes available, we can make an order of magnitude estimate on the size of theD σ 3 and speculate on the possibility of a combined explanation of all the three ATOMKI anomalies with an axial vector state. We can expect that the isovector spin dipole would be of the order of the nuclear radius times the number of nucleons inside the nucleus. We can then estimate which is in order of magnitude accord with the ATOMKI result which predicts a value of 3.6(3) × 10 −6 for this rate [5]. We stress again that this conclusion strongly depends on our order of magnitude estimate of theD σ 3 matrix element, which seems to indicate that an axial vector state might be favored for a combined explanation. However to properly test its consistency with the 12 C anomalous transition, the relevant matrix element must be properly computed. Until then no definite conclusions can be drawn. In a general scenario where both vector and axial couplings to nucleons are present, the decay width for the real X emission is the direct sum of the two contributions. Assuming vector and axial couplings of the same order of magnitude, the axial contribution would typically dominate over the vector one.
Intriguingly, for the case a pure axial boson S π = 1 + , in the parameter space where the 4 He and 8 Be anomalies can be explained, other experimental anomalies can be simultaneously satisfied, while being compatible with current constraints on the electron couplings of the X boson. This is the case of the KTeV anomaly in π 0 → e + e − decay [51], inside the green region in Fig. 5, and the anomalous magnetic moment of the electron (g − 2) e , as we will explain in the following.
The discrepancy of 3.2σ 9 could be explained assuming a BSM contribution from a light axial boson [56]. The actual best fit from the data reads as [53] (a p − a n )C e where C e A is the axial coupling of the X boson to the electron, see Eq. (F.1). A light vector contributes to the anomalous magnetic moment of the electron. The SM prediction from the measurement of the fine structure constant α from Cs atoms [57] and the more recent prediction based on the measurement of α from Rb atoms [58] are in contradiction among themselves. By asking that the BSM contribution from the X boson given by [59] δa BSM Using the latest radiative corrections from [54,55], the full branching ratio extrapolated from the KTeV measurement is BR(π 0 → e + e − )exp = (6.85 ± 0.27 ± 0.23) × 10 −8 , thus reducing the discrepancy to 1.8σ. The Cs atoms SM prediction naturally suggests a pure axial boson and the discrepancy By fixing this value for C e A , we have shown in Fig. 5 the parameter space of the nucleon couplings which can explain the KTeV anomaly see Eq. (5.19), for the two distinct cases C e A > 0 and C e A < 0. Allowing instead for both a vector and axial contribution to the electron coupling, in the upper panel of Fig 6 we show the values of the isovector nucleon axial coupling a p − a n able to explain the KTeV anomaly at 1σ and assuming the discrepancy observed in the electron anomalous magnetic moment to be resolved for the Cs atoms (blue) and Rb atoms (orange) SM prediction for the two distinct cases C e A > 0 and C e A < 0. In the lower panel of the same figure we show instead the most relevant bounds on this scenario with generic X vector couplings to electrons, again for BR(X → e + e − ) = 1, which arise from the measurement of e + e − scattering from the KLOE experiment at DAΦNe collider [61], measurements on parity violation in Møller scattering at SLAC [62] and beam dump experiment at NA64 [63,64], see App. F for details. Interestingly, the PADME experiment will completely cover the region between the NA64 and the KLOE exclusions thus allowing for a strong test of the existence of the X boson [18,19].

Minimal SM extension with a new U (1) gauge symmetry
A common proposal in order to include a light vector X in a SM extension is by enlarging its gauge group In a minimal scenario, it is natural to assume that the SM Lagrangian is symmetric under U (1) X . It follows that the U (1) X charges are a linear combination of the hypercharge Y and all the accidental symmetries of the SM: the baryon number B and the three lepton family numbers L e,µ,τ . When the gauge group G SM × U (1) X is broken to SU (3) c × U (1) Q , the X field mixes with the other neutral gauge bosons in order to compose the physical states: the photon γ, the Z 0 boson and the light boson X. Moreover, once the symmetry is broken, the X charges acquire a contribution from the diagonal weak isospin T W 3 . In the quark sector, the X couplings are then a linear combination of three independent charge assignments: B, Y = Q − T W 3 and T W 3 . The baryon number B and the electric charge Q are vector symmetries so only the weak isospin T W 3 induces an axial coupling. Thus, the axial couplings of the light quarks satisfy while the vector coupling for up and down quarks are independent. The nucleon couplings are obtained from the quark ones by Eq. (C.3) and Eq. (C.11). However, the results from atomic parity violation experiments [68], see App. F, strongly constrain these couplings, and require for them (product of) values so small that the ATOMKI anomalies cannot be explained, see Eq. (F.9). A way to avoid this bound is to assume somehow a magical cancellation between the up and down vector couplings, which is Hence, we are left with two independent couplings, C p and C A . However, as we show in Fig. 7, there is no possible simultaneously explanation of the 8 Be and 4 He anomalies in the minimal BSM scenario considered here for both choices of the isospin breaking parameter ξ, which produce almost indistinguishable results.

Conclusions
Motivated by the latest experimental results recently released by the ATOMKI collaboration, we have critically re-examined the possible theoretical interpretation of the observed anomalies in 8 Be, 4 He and 12 C anomalies in terms of a BSM boson X with mass ∼ 17 MeV.
After having reviewed the current status of the ATOMKI results and the kinematic of the observed excesses we have employed a multipole expansion formalism to compute the anomalous decay rate for the decay of the excited nuclei into an e + e − pair via an intermediate on-shell BSM state. Our results identify an axial vector state as the most promising candidate to simultaneously explain all the three anomalous nuclear decay, while the other spin/parity assignments seems disfavored for a combined explanation. However, the axial nuclear matrix element of the 12 C transition is currently unknown and our conclusions are based on an order of magnitude estimate for its value. Before being able to make a definite claim regarding the possibility of combined explanation of the ATOMKI anomalies with an axial vector state, this matrix elements need to be evaluated. Intriguingly, the hypothesis of an axial vector state can also simultaneously accommodate other experimental anomalies, as the one observed by the KTeV experiment in π 0 → e + e − decay, while being compatible with the conflicting measurements of the anomalous magnetic moment of the electron (g − 2) e and other experimental constraints. The independent experimental tests that will be performed by the MEG II experiment [12] at PSI and the by the Montreal Tandem accelerator [14] will definitely answer the question regarding the nature of the ATOMKI results and tell us whether these anomalies are merely due to unaccounted SM and/or experimental effects or else are the first signs of the long sought new physics beyond the Standard Model.

Acknowledgments
The Authors thank Marco Nardecchia and Mauro Raggi for triggering and encouraging this work. Our gratitude to Luca Di Luzio, Gabriel Massoni Salla and Federico Ripani for useful comments and discussions. We are also grateful to Matheus Hostert and Maxim Pospelov who helped us to spot a mistake in the propagation of the errors for the best fit of the vector state which was present in the published version of our paper. The work of CT was supported in part by MIUR under contract PRIN 2017L5W2PT.

A Additional tables and input parameters
Here we report the numerical values for the nuclear transition of interest used throughout our analysis.  Table 5. Proton kinetic energy E b , mass of the target and of the excited nucleus, reduced protontarget mass and velocity of the excited nucleus produced by the collision in the lab frame at the resonance peak (E CM = m N * ) for all the measured N * states.  Table 6. N mass, threshold energy E th , boson energy ω in the CM frame the velocity v X of an hypothetical boson of mass m X ≃ 17 MeV emitted in the CM frame at the resonance peak (E CM = m N * ) for all the measured N * states.  Table 7. δ parameter, maximum energy asymmetry y max and minimal opening angle θ min ± for an hypothetical boson of mass m X = 17 MeV emitted in the CM frame.

B Electromagnetic dynamics
With the formalism described in Sec. 3 we can describe the dynamics for the SM processes pictured in Fig. 3 of real γ emission and IPC. The nuclear electromagnetic current, including the anomalous magnetic moments, is given by where κ p = +1.792847351(28), κ n = −1.9130427(5) [42] and Q p,n indicates the electric charge of the nucleon in units of the absolute electron charge. The magnetic momenta of the nucleons are µ p,n = (Q p,n + κ γ p,n )µ N , where µ N = e/2m N is the nuclear magneton. The conservation of the electromagnetic current implies that only three independent spherical operators have to be considered, c.f.r. Eq. (3.13). In the non relativistic and long wavelength approximation the spherical operators with J = 0, 1 are given by where we have defined the electromagnetic monopole ρ (γ) , the electric dipole d

Real γ emission
The rate for the process with a real γ emission can be readily computed from Eq. (3.14) by fixing m γ = 0. In this case, due to the transversality of the photon, only processes with J = 0 are allowed which are • electric type transitions EJ from the contribution of T el J with parity π(EJ) = (−1) J , • magnetic type transitions M J from the contribution of T mag J with parity π(M J) = (−1) J+1 .
The E1 and M 1 decay rates are equal to

Internal pair creation
At lowest order the IPC process involves the emission of a virtual photon decaying into an e + e − pair. The differential decay rate with respect to the energy asymmetry y and the opening angle θ ± is given by where the f M,T (y, cos θ ± , δ) functions are (B.14) Note that IPC processes acquire a contribution also from the longitudinal modes, absent in the real γ emission case. Hence, processes of pair production are of three types: • longitudinal type transition LJ from the contribution of M J with parity π(LJ) = (−1) J , • electric transition EJ from the contribution of T el J with parity π(EJ) = (−1) J , • magnetic transition M J from the contribution of T mag J with parity π(M J) = (−1) J+1 , with differential decay rates 11 We show in Fig. 8 the theoretical distributions of the angular correlation, obtained after integrating the above expressions over the asymmetry y. The reduced invariant mass of the lepton pair is a function of the energy asymmetry and opening angle. Given the (normalized) distribution function f (y, cos θ ± , δ) of the pairs over the plane defined by y and cos θ ± , the (normalized) distribution function g(s, δ) of the reduced invariant mass is given by where s = 2m 2 ee ω 2 (2δ 2 ≤ s ≤ 2) is the reduced invariant mass. For the J = 0, 1 multipoles we calculated, the integration over the energy asymmetry is easily performed since the integrand turns out to be polynomial in y. Hence, one finds We show in Fig. 9 their theoretical distributions.

Isospin mixing and electromagnetic nuclear matrix elements
The 8 Be(18.15) and 8 Be(17.64) states, close in energy and with same spin-parity assignment, presents a significant isospin mixing. In general, given a doublet of nuclear energy levels of spin J with mixed isospin, the physical states (denoted with a and b) are given by a linear combination of states with fixed isospin where a labels the lowest energy level between them. The coefficients α J and β J are real and satisfy α 2 J + β 2 J = 1. For the 8 Be nucleus, the values of the mixing coefficients have been evaluated through Quantum Monte Carlo simulation [69]. The result for the J = 1 doublet is According to this, we define the isospin magnetic strength M 1 γ I=0,1 by (B.27) whose values has been estimated to be [69] M 1 γ I=0 = 0.014(1), M 1 γ I=1 = 0.767 (9). (B.28) At this level, a direct comparison with the experimental values of the decay lengths shows a significant discrepancies with the theoretical prediction. Following [27], we'll consider the deficiency as due to isospin breaking effects we neglected in the first attempt. The inclusion of them is obtained trough the introduction of a ∆I = 1 spurion, whose effective result is to shift the nuclear matrix elements as The parameter ξ characterizes the strength of the spurion and its size is controlled by non-perturbative effects. One finds ξ = 0.549 by requiring that the resulting decay width Γ( 8 Be(17.64) → 8 Be+γ) reproduces its experimental value. For ξ = 0 the isospin breaking effects are simply neglected. Parity conservation prohibits electromagnetic interaction in the 4 He(21.01) transition to the ground state, thus the nuclear matrix element relative to this decay is equal to zero. Due to the massless nature of the photon, the γ-emission decay width of the 4 He(20.21) transition to the ground state also vanishes but the IPC process is still possible. Pair production is mediated by the monopole operator ρ (γ) , whose matrix element has been measured to be [9] √ 4π ⟨ 4 He||ρ The γ-emission decay length of the 12 C(17.23) transition to the ground state has been measured to be Γ γ = 44 eV [11]. Hence, from (B.10), one finds

C Nucleon effective couplings
We report in this section the matching between the effective interaction of the X boson with the nuclear matter and its interactions with the fundamental SM degrees of freedom, quark and gluons. Since nucleons are spin 1/2 particles, for an operator O composed by quark fields one has where Γ is a matrix with spinor indices and u p is the solution of the free Dirac equation 12 . Lorentz invariance as well C, P and CP T symmetries impose further constraints on this matrix element. As mentioned in the main text, since the transferred momentum in the considered processes is generally much smaller than Λ QCD we approximated the form factors that are in general present in these expressions as constants.

Vector interaction
In the UV an X vector boson interacts with a quark current of the form C q Vq γ µ q, which gives an effective interaction Conservation of the vector current implies g N = 0, while symmetry considerations fix Symmetry considerations don't allow to simplify the expression for the magnetic moments of the nucleons µ (X) where µ N is the Bohr magneton, since even the sea quarks can give a contribution to these quantities. However by using the static quark model one can make an estimation [70].
Working under the assumption that the valence quarks of the nucleons have mass equal to m u ≃ m d ≃ m N /3 13 , at lowest order the magnetic moment of each quark is given only by its charge and effective mass, i.e.
In the static quark model one thus finds For the electromagnetic couplings, C p = +e and C n = 0, one obtains values close to the experimental ones, µ em p | exp ≃ +2.792 µ N and µ em n | exp ≃ −1.913 µ N , within a 10%. The 12 This leaves the matrix element unchanged upon the substitution Γ(p ′ , p) → / p ′ +m 2m Γ(p ′ , p) / p+m 2m . 13 Here we are considering the effective mass of the quarks when they are bounded together by gluons, not their real mass.
interaction between quarks and the X particle might also come from an effective magnetic moment interaction L = κq Λ ∂ ν (qσ µν q)X µ generated, e.g., at loop level by integrating out some heavy particle in the low energy limit. This effective operator contributes to the magnetic moment of the nucleons through the substitution for each q. By lattice computation [71,72], it has been estimated for the light quarks contributions.

Axial interaction
In the UV an X axial vector boson interacts with a quark current of the form C q Vq γ µ γ 5 q which brings to an effective nucleon current CP conservation in QCD interactions forces d N = 0, while the term proportional to b N doesn't contribute to the considered processes when one has on-shell X 14 . The nucleon axial couplings a N are given by the sum of quark coupling a q weight by the fraction of the spin of the nucleon ∆ These fractions are given by integrals of helicity-dependent parton distributions and can be measured in lepton nucleon scattering. Their values are equal to [32,73] (27) , ∆ (p) s = ∆ (n) s = −0.026(4) , (C.12) while the contributions from heavy quark are small and can be neglected. As for the vector case, it's possible that the interaction between quarks and the X boson comes from an effective interaction like L = dq Λ ∂ ν (iqσ µν γ 5 q)X µ . This effective operator generates an electric dipole for the nucleons through the substitution where N = p, n, which again shifts the d N value as before. Unfortunately the values of δ are difficult to be measured and are poor known. Only recently [74] it has been measured the light quark contribution to the proton at Q 2 = 0.

Scalar interaction
The scalar interaction between quarks and a spin 0 particle is given by the scalar density operatorqq. The matching with the nucleon effective coupling is linked to the generation of nucleon masses [75]. From trace anomaly, the mass of the nucleons is given by where the β function at lowest order is β = −α 2 s /2π(11 − 2n f /3) and α s is the strong coupling constant. The heavy quark fields Q = c, b, t can be integrated out trough the expansion [76] so that We can now define the fractions of nucleon mass as We consider a scalar interaction term with the X particle defined by where v = 246 GeV is the Higgs vacuum expectation value. The last one is an effective interaction term that can be generated at loop level by massive particles in the low energy limit. The nucleon effective interaction reads are the effective scalar couplings of the nucleons. The values of the fractions of nucleon mass are given by [77] f (p) T u = 0.020 ± 0.004, f

Pseudoscalar interaction
The pseudoscalar density iqγ 5 q is proportional to the divergence of the axial currentqγ µ γ 5 q.
The matching with the nucleon effective operator is then done with the same ∆ (N ) q parameters already used. For the light quark contribution, we have while for the heavy quark fields it's enough to expand them as [76] The nucleon matrix element for the pseudoscalar gluon operator is given by [78] ⟨N | α s 8π By considering an interaction term for quarks and gluons given by which is equivalent to with C gg P = C g P + q C q P , the effective nucleon interaction then reads

Pseudoscalar scenario
We summarize the results for the pure pseudoscalar scenario in Fig. 10, where the shaded blue and orange areas represent the 1σ and 2σ compatibility regions with the ATOMKI 8 Be and 4 He anomalies respectively, where the former is a combination arising from both the 8 Be energy levels. We also overlay in red the region of parameter space satisfying the SINDRUM bound from π + → e + ν e X decay [37,46]. Note that a combined explanation of the 8 Be and 4 He anomalies is not anymore possible, once we include the constraint from the 8 Be(17.64) transition.

Vector and axial scenarios
We summarize the results for the spin-1 cases in Fig. 11, with the same color code as Fig. 10 for the regions satisfying the ATOMKI anomalies. In the upper panels we show the results for the S π = 1 − assignment for the X boson. As it can be seen for both ξ assignments, a combined explanation of the 8 Be, blue region, and 4 He, orange region, anomalies at 1σ is in tension with the NA48 constraint, while it is possible at the 2σ level. The 1σ and 2σ bands related to the 12 C transition are shown in purple in the upper panels of Fig. 5. Note that, if confirmed, the 12 C ATOMKI anomaly is in tension with a combined explanation of the 8 Be and 4 He anomalies and the protophobia constraint. On the other side an axial vector S π = 1 + state can explain both the 8 Be and 4 He ATOMKI anomalies at 2σ, as shown in the lower panels of Fig. 11, with axial couplings to the nucleon of O(10 −4 ). Thus the inclusion of the 8 Be(17.64) transition does not change drastically the conclusion for the spin-1 cases. Figure 11. Upper panels: Regions of the C n,p effective nuclear couplings of a pure vector state where the 8 Be (blue), 4 He (orange) and 12 C (purple) anomalous ATOMKI transition can be explained at 1σ or 2σ. Inside the red and the gray region, respectively, the NA48 and the 208 Pb-n scattering bound are satisfied. In the left and right panel we assume ξ = 0 and ξ = 0.549 respectively, see App. B for details. Lower panels: Regions of the a n,p effective nuclear couplings of a pure axial vector state where the 8 Be (blue) and 4 He (orange) anomalous ATOMKI transition can be explained at 1σ or 2σ. In the green region the KTeV anomaly in π → e + e − decay can be satisfied, by assuming a positive (left panel) and negative (right panel) value for the C e A axial coupling of the X boson to electrons that can explain the anomalous (g − 2) e . In all figures both the 8 Be(18.15) and the the 8 Be(17.64) transitions are considered, see main text for more details.

E Cross section for resonance production
For resonance production p + A → N * , the unpolarized cross section expression is given by where E 2 CM = (p A + p p ) 2 . Here we have treated the excited state N * as a bound states of p and A with mass m * and spin J * . The reverse process defines the decay width where λ(x, y, z) = x 2 + y 2 + z 2 − 2xy − 2xz − 2yz. The angular integration is trivial because the unpolarized squared matrix element of three momentum is just a function of the masses. Since the nuclear force, which mediate the interaction, is invariant under time reversal, it follows that M(N * → p + A) = M(p + A → N * ) .

F Experimental constraints on a spin-1 boson
A light boson coupled to first generation quarks and leptons is subject to a large variety of experimental constraints. In this section we recap the most relevant ones that affects a possible explanation of the ATOMKI anomaly through a BSM degree of freedom with mass ∼ 17 MeV for the spin-1 case. We parametrize the effective UV interactions of a spin-1 state X µ with leptons and quarks as where C f V,A = 1 2 (C f R ±C f L ) and we assume diagonal couplings in flavor space. The connection between the quark and nucleon couplings can be found in App. C. We list in the following the most relevant constraints for the spin-1 case coupling to the first generation of quarks and leptons.
e + e − → γX scattering By neglecting kinematic differences with respect to the pure dark photon case, we can recast the bound from the KLOE experiment at the DAΦNE collider [61] from which we obtain a bound (C e V ) 2 + (C e A ) 2 ≲ 6.1 × 10 −4 BR(X → e + e − ) . (F.2)

Parity violation
Parity violation in Møller scattering constraints the product of the vector and axial couplings. The most sensitive measurement arises from SLAC E158 [62] at Q 2 = 160 MeV 2 . The measurement in [62] has been recast in [79] and the obtained bound reads

Beam dump experiments
Beam dump experiments look for X production via bremsstrahlung from electrons scattering off target nuclei. For the X particle not to be seen in these experiments there are two possibilities: either the particle is not produced at all, or its decay products are caught in the dump, thus setting both an upper and lower limit for the couplings of the X boson. In the first case the stronger limit comes from the E137 experiment [80], see also [81], which is independent on the X decay rate and gives (C e V ) 2 + (C e A ) 2 ≲ 1.1 × 10 −8 , (F. 4) while in the second case the stronger limit comes from the NA64 experiment [63,64] for which we have (C e V ) 2 + (C e A ) 2 ≳ 3.6 × 10 −5 × BR(X → e + e − ) . (F.5)

Prompt decay in ATOMKI detector
The requirement of a prompt decay into the ATOMKI detector imposes now the constraint (C e V ) 2 + (C e A ) 2 ≳ 3 × 10 −7 × BR(X → e + e − ) , (F. 6) which is weaker than the bound from NA64.

Atomic parity violation
In atomic system, parity violation can be observed in the case, e.g, of and electric dipole transition between two atomic states with the same parity. The X µ gives additional contributions to these transitions due to the interaction between atomic electrons and the nucleus. In the effective operator where only the V × A part have been kept, only the A e × V u,d interaction give a relevant effect for parity violation observables. This is due to the fact this part of the interaction between the electron and the nucleus is coherent, and thus proportional to the total weak charge of the nucleus itself, while the A q × V e interaction adds incoherently. This effect is thus suppressed for heavy enough nuclei [82]. The BSM contribution to A e × V u,d can be expressed as a modification to the weak nuclear charge Q W [83] δQ . (F.8) The most accurate prediction comes from transition of 133 55 Cs [68] which, combined with the SM theoretical prediction [84], yields [85]