Spin-dependent ${\mu \to e}$ Conversion on Light Nuclei

The experimental sensitivity to $\mu \to e$ conversion will improve by four or more orders of magnitude in coming years, making it interesting to consider the"spin-dependent"(SD) contribution to the rate. This process does not benefit from the atomic-number-squared enhancement of the spin-independent (SI) contribution, but probes different operators. We give details of our recent estimate of the spin dependent rate, expressed as a function of operator coefficients at the experimental scale, and explore the prospects for distinguishing coefficients by using different targets. For this purpose, a geometric representation of different targets as vectors in coefficient space is introduced. It is found that comparing the rate on isotopes with and without spin could allow to detect spin dependent coefficients that are at least a factor of few larger than the spin independent ones. Distinguishing among the axial, tensor and pseudoscalar operators that induce the SD rate would require calculating the nuclear matrix elements for the second two. Comparing the SD rate on nuclei with an odd proton vs odd neutron could allow to distinguish operators involving $u$ quarks from those involving $d$ quarks; this is interesting because the distinction is difficult to make for SI operators.


Introduction
Charged Lepton Flavour Violation (CLFV) is New Physics that must exist; only the rates are unknown. In this paper, we consider µ ↔ e flavour change, and assume that it can be parametrised by contact interactions involving Standard Model particles. Flavour change µ ↔ e can be probed in the decays µ → eγ [1] and µ → eēe [2], in µ → e conversion [3] and in various meson decays such as K →μe [4]. In µ → e conversion, a beam of µ − impinges on a target, where the µ is captured by a nucleus, and can convert to an electron while in orbit. The COMET [5] and Mu2e [6] experiments, currently under construction, plan to improve the sensitivity by four orders of magnitude, reaching a branching ratio ∼ 10 −16 . The PRISM/PRIME proposal [7] aims to probe ∼ 10 −18 . These exceptional improvements in experimental sensitivity motivate our interest in subdominant contributions to µ → e conversion.
Initial analytic estimates of the µ → e conversion rate were performed by Feinberg and Weinberg [8], for promising operators and nuclei. A wider range of nuclei were studied numerically by Shanker [9], and estimates for many operators and nuclei can be found in the review [10]. Relativistic effects relevant in heavier nuclei were included in [11]. The matching of CLFV operators constructed with quarks and gluons, onto operators constructed with nucleons, was performed in [13]. The current state of the art is the detailed numerical calculations of Kitano, Koike and Okada (KKO) [12], who studied all the CLFV nucleon operators that contribute coherently to µ → e conversion, for nuclei from Helium to Uranium. In such processes, the amplitude for µ → e conversion on each nucleon is coherently summed over the whole nucleus. Like "spin-independent"(SI) dark matter scattering, the final rate therefore is enhanced by a factor ∼ A 2 , where A is the atomic number of the nucleus. However, other conversion processes are possible. For instance, incoherent µ → e conversion, where the final-state nucleus is in an excited state, has been discussed by various people [14,9],and is expected to be subdominant with respect to the coherent process.
In a previous letter [15], some of us noted that "spin-dependent"(SD) µ → e conversion can also occur, if the target nuclei have spin(as is the case for Aluminium, the target of the upcoming COMET and Mu2e experiments). Although this process does not benefit from the ∼ A 2 enhancement associated to SI rates, it has the interest of being mediated by different CLFV operators from the coherent process.
The aim of this manuscript is to give details of our calculation, and explore whether the SD process could help distinguish models or operators, should µ → e conversion be observed. The operators which could induce SD µ → e conversion are listed in section 2. The conversion rate in Aluminium is estimated in section 3, and the extrapolation to other nuclei is discussed in subsection 3.2. The theoretical uncertainties in our estimates are briefly discussed in section 4. Section 5 explores the consequences of including the SD contribution to the µ → e conversion rate, both in the perspective of obtaining constraints on operator coefficients from an upper bound on the branching ratio, and for discriminating models when µ → e conversion is observed. This section comes in three parts: we study three leptoquark models which induce SD and SI conversion, then consider the same operators but with arbitrary coefficients, and calculate a covariance matrix. Finally, we allow all possible operators with arbitrary coefficients. We summarise in section 6. In our previous letter [15], we showed that the SI and SD operator coefficients mix under Renormalisation Group(RG) evolution between the experimental and weak scales. The effects of this mixing are significant: the largest contribution to the µ → e conversion rate from an "SD" coefficient at the weak scale, would be via the RG mixing to an SI coefficient (for example, a tensor coefficient at the weak scale induces a SI contribution to the rate which is ∼ A 2 larger than the SD contribution). In this paper, we focus on operator coefficients at the experimental scale, only including the RG evolution in the leptoquark models of section 5.1. The RG evolution of the operator coefficients is summarised in Appendix C.

Operators
We are interested in contact interactions that can mediate µ → e conversion on nuclei, at a scale µ N ∼ 2 GeV. The focus of this manuscript is the subset of "spin-dependent" interactions, but for completeness, all QED×QCD invariant operators that mediate µ → e conversion on nuclei are included. The relevant operators in the quark-level Lagrangian are [12,13]: where the two-lepton operators are and O ′ ∈ {V, A, S, P, T } labels 2-lepton 2-quark operators in a basis where only the lepton currents are chiral: with σ αβ = i 2 [γ α , γ β ] and P L = (1 − γ 5 )/2. This choice of non-chiral quark currents is convenient for matching onto nucleons. However, often an operator basis with chiral quark currents is added to the Lagrangian as δL = [10,16,17], where for instance, O qq V,Y X ≡ (eγ α P Y µ)(qγ α P X q). In this case, the coefficients are related as (recall that O qq T,LR vanishes-see appendix C of [16]) In eqn (1), the coefficients and operators are evaluated close to the experimental scale, at µ N ≃ 2 GeV. The scale is relevant, because Renormalisation Group running mixes the tensor and axial vector operators (that induce SD µ → e conversion) into the scalar and vector operators (who mediate the SI process) [15] 1 . This is reviewed in Appendix C. Throughout the paper, coefficients without an explicit scale are assumed to be at µ N .
To compute the rate for µ → e conversion, the operators containing quarks should be matched at the scale µ N onto CLFV operators involving nucleons and mesons. The relevant nucleon operators are the four-fermion operators of eqn (3) with q → N and N ∈ {n, p}. As discussed below, rather than include mesons in the Lagrangian, we approximate their effects by form factors for some nucleon operators and two additional operators given in eqn (10). So the nucleon-level Lagrangian will be where O ′′ ∈ {V, A, S, P, T, Der}. At zero momentum transfer ( P f − P i → 0), we match onto operators with nucleon currents, by replacinḡ obtained at zero-recoil are given in appendix A, and we will assume that they are an acceptable approximation at the momentum-transfer of µ → e conversion, which is | P f − P i | 2 = m 2 µ . Various mesons are present in the low energy theory at µ N , so in principle the quark operators of eqn (1) should be also matched onto meson operators. χPT [26] involving nucleons (see e.g. the review [19]) would be the appropriate formalism for this calculation, and has been used to calculate WIMP scattering on nuclei [20,21,22,23], neutrinolessdouble-beta-decay [24], and SI µ → e conversion [25]. However, to avoid more notation, here we just give results for the simple diagrams of interest. We only consider the CLFV decays of pions, because the effects of heavier mesons would be suppressed by their masses, and diagrams where a pion is exchanged between two nucleons are suppressed by more propagators, and would require two nucleons in the initial and final states 2 . Pion decay can contribute to µ → e conversion via the second diagram of figure 1, in the presence of a pseudoscalar or axial vector quark current. We follow the notation of [26,19] in matching the axial vector and pseudoscalar quark currents onto the pion, at in order to obtain the usual expectation values 0|ū( MeV. Later in the manuscript, the matrix element for µ → e conversion on a nucleon, M(µ + N (P i ) → e(k) + N (P f )) will be required. In the case of vector, scalar or tensor interactions, it is is straightforward because conversion proceeds via a 2-nucleon-2-lepton contact interaction. In the case of axial vector and pseudoscalar interactions, there is a pion exchange contribution, as illustrated in figure 1, so we give the matrix elements here. The pion-nucleon interaction term in the Lagrangian is taken as ig πN N N γ 5 τ · πN , and the Goldberger-Treiman relation gives g πpp ≃ (G pu A − G pd A )m p /f π . In the following two equations, u N = (u p , u n ) represents a vector of spinors in isospin space. The matrix element M(µ + N (P i ) → e X (k) + N (P f )) mediated by the axial up quark current, can be written [27,20] where q = (0, − q) = P f − P i , the first term is written in terms of iso-scalar and iso-vector contributions (a 0 + a 1 )/2 = C uu A,X G p,u A , (a 0 − a 1 )/2 = C uu A,X G n,u A , whereas the pion contribution is only isovector. In the case of the pseudoscalar operator O uu P,Y , the pion exchange diagram is non-vanishing at | q| 2 = 0, so at finite momentum transfer, only the additional contribution ∝ 1/(| q| 2 + m 2 π ) − 1/m 2 π should be included. This gives: In summary, the axial vector and and pseudoscalar quark operators could equivalently have been matched at µ N to an EFT without pions, but with a q 2 -dependent "form factor" for the pseudoscalar nucleon operator, and an additional dimension seven derivative operator In this extended basis, the nucleon coefficients are where C N N Der,Y was evaluated at q 2 = −m 2 µ , and the scalar nucleon coefficients, to which contribute also gluon operator of eqn (2), are given in [13].
To obtain the µ → e conversion rate, the expectation values of the nucleon operators in the nucleus are required. This is discussed in the next section. We were unable to find nuclear expectation values of the tensor and pseudoscalar operator, so O N N P,Y will be neglected, and the tensor included in the scalar and axial operators, as described in eqn (19).

Estimating the SD and SI rate in light nuclei
In our previous paper [15], we gave analytic estimates of the SI and SD conversion rates on Aluminium. The aim of section 3.1 is to give details of the calculation in the notation of relativistic, second-quantised Field Theory. The results can then be matched onto the nuclear physics calculations of [12] (for SI conversion), and SD WIMP scattering [27,28,29,20] (for SD conversion). In subsection 3.1.3, the estimates are mapped onto the numerical results of KKO [12], and SD conversion in heavier targets is discussed in section 3.2.

Estimating the SD and SI rate in Aluminium
We define the bound state of momentum P i composed of an Aluminium nucleus and a muon in the 1S orbital as ≡ |Alµ(P i ) . We are interested in the S-matrix element for Alµ(P i ) → Al(P f ) + e − X (q) induced either by the dipole operator (which we discuss later), or by a four-fermion operator (e X Γ l µ)(N Γ n N ). To be concrete, we consider the S-matrix element where the nucleon N is a proton: where s is the spin of the electron selected by the chiral projector P X , field operators wear hats, and Γ n ∈ {I, γ 5 , γ α , γ β γ 5 , σ αβ }, Γ l ∈ {I, γ α , σ αβ }.

four-fermion operators
1. A first step is to write the motionless bound state |Alµ(0) as where w is the spin of the muon, the square-root prefactor accounts for one vs two-body normalisation of states in Lorentz-covariant field theory conventions where states are normalised ∝ √ 2E [30], and ψ µ ( l) = d 3 ze −i l· z ψ µ ( z) is the fourier transform of the Schrodinger wavefunction ψ µ ( z) for a muon in a central potential of charge Z.
For Zα ≪ 1, the unit-normalised wavefunction, for either spin state, can be approximated [31,32,33] as We approximate the outgoing electron as a free particle (plane wave), which should be acceptable for an Aluminium target. For heavy nuclei, the Dirac equation for the electrons outgoing in the field of the nucleus should be solved [11], allowing to express the electron as a superposition of free states. This approach was followed in KKO [12].
2. In the same non-relativistic bound state formalism (see e.g., Appendix B of [28] for more details), the Aluminium nucleus, of spin J A , can be written as a bound state composed of a proton of spin t, with another state M 1 of mass M 1 and spin J M containing Z − 1 protons and A − Z neutrons: 3. The fermion operators can be expanded as [30] and act on states asμ(y)| k, w = u w k e −ik·y |0 , where the spinors are normalised as u † k u k = 2k 0 . The S-matrix element of eqn (12) can then be evaluated as where the spinors subscripts are particle names rather than momenta, and To obtain this approximation, the states were taken to be non-relativistic, the wavefunctions expressed in position space, the proton wavefunction was assumed independent of the proton spin, and the dependence of spinors on three-momenta was neglected in many integrals. Notice the M Al /m p enhancement factor that arises automatically for both spin-dependent and spin-independent interactions, and that the usual (2π) 4 δ 4 (P i − P f − q), which accounts for four-momentum conservation, appears despite that there is a spatial integral over the nucleus. In the following, we drop the spin indices in the nucleon distribution in the nucleus |f N | 2 .
4. The leptonic spinor contraction is independent of x and can be factored out of the spatial integral in eqn (17).
In light nuclei such as Aluminium, the muon wavefunction can also be factored out [8], because the muon wavefunction decreases on the scale ∼ 1/(Zαm µ ), which is larger than the radius of the Aluminium nucleus, given in [34] as < ∼ 6 fm. On the other hand, the first zero of the electron plane wave (the e −i q· x of eqn (17)) would occur at r ∼ π/(m µ ) ∼ 6 fm.
5. The nucleon spinor contractions, in the non-relativistic limit, can be written (see eqn (47) of [35]): where the spin vector of the nucleon is defined as 2 S N = u † N Σu N /2E N , and the rotation generator The momentum transfer q = P i − P f has been neglected, except in the case of the pseudoscalar, where the leading term is O( q · S N ), and in the case of the tensor, where the there is a "spinindependent" contribution ∝ q.
These spinor identities allow the tensor interaction involving nucleons to be absorbed into the scalar and axial vector coefficients. Following [15], we define where in both cases the 2 arises from the two antisymmetric contributions of the tensor, the unprimed Cs are defined in eqn (11), X, Y ∈ {L, R}, and X = Y because only operators with electrons of the same chirality can interfere. Notice that there is an error in [15], where is written It remains to evaluate the expectation value of the nucleon currents in the nucleus.
(a) In the case of the scalar or vector operators, the matrix element of eqn (17) becomes where the sum over protons in the nucleus will give a factor Z, we drop the spin indices because the sum and average give one, and assume a spherically symmetric nucleon distribution |f p (r)| 2 in the nucleus, which allows to replace 3 e −i q· x → sin(qr) qr . The "form factors" are defined in eqns (29) and (30) of [12]: F p (m µ ) ∼ .53 for Al, and ∼ .35 for Ti.
(b) The expectation value of the axial current in Aluminium (A = 27, Z = 13, J Al = 5/2) was calculated by Engel et.al [29] and Klos et.al [20] using the shell model. In the zero-momentum transfer limit, where the spin expectation values S A N are defined by: they obtain S Al n = 0.0296, S Al p = 0.3430. (J k A is a quantum mechanical operator, to be evaluated in the ground state of the nucleus A). At finite momentum transfer, references [29,20] include the nucleon axial vector operators O N N A,X and the pion exchange operator O N N Der,X , in the combination induced by axial vector quark operators. The various terms in the matrix-element-squared have different spin sums, so the finite momentum transfer correction depends on C pp ′ A,X and C nn ′ A,X , and is quoted as a multiplicative factor S A (m µ )/S A (0) in the rate (see eqn (26)). Neglecting S Al n ≪ S Al p , the results of Engel et. al for Aluminium give [29] S Al (k) ∝ (0.31500480 − 1.857857y + 4.86816y 2 − 5.422770y 3 ) where y = (m µ b/2) 2 and b =1.73 fm. This gives S Al (m µ )/S Al (0) = 0.29. 3 Recall that a plane wave can be expanded on spherical harmonics as (c) At zero momentum transfer, the nuclear expectation value of tensor operators O N N T,X is proportional to that of axial vector operators, as accounted for in eqn (19). However, at finite momentum transfer, there is no pion exchange contribution for the tensor operator (while pion exchange induces O N N Der,X in the presence of the axial vector quark operators), so the redefinition of eqn (19) is not valid. Indeed, the tensor and axial vector operators are distinct at finite momentum transfer. However, we did not find nuclear calculations of SD scattering on Aluminium mediated by the tensor operator. We can try to estimate the error from using the axial results for the tensor: at q 2 = −m 2 µ , the pion exchange contribution to the matrix element in eqn (8) is comparable to the four-fermion contact interaction. Also, the finite-momentum-transfer suppressions of the axial and scalar rates on Aluminium are comparable (S Al (m µ )/S Al (0) ≃ |F N (m µ )| 2 ), despite that one might expect the oscillations of the electron wavefunction to suppress the SD rate more than the SI rate, because spin-carrying nucleons are likely to be at large radii. So we interpret that axial matrix element is amplified by a factor ∼ 2 at q 2 = −m 2 µ (due to the pion), and suppressed by an extra factor ∼ 1/2 (as compared to the scalar matrix element) due to the oscillations of the electron wavefunction, and estimate that the identification of eqn (19) could overestimate the tensor contribution to the branching ratio by a factor ∼ 2 → 4 (depending on whether the pseudoscalar and axial matrix elements interfere). (d) The pseudoscalar operator O N N P,X is proportional to the nucleon spin, is only present at finite momentum transfer, and at q 2 = −m 2 µ , is enhanced by a pion exchange contribution of comparable magnitude. Since the magnitude of the pseudoscalar spinor contraction in eqn (18) is suppressed with respect to the axial vector by ∼ m µ /2m N , its contribution to the SD branching ratio could be ∼ m 2 µ /4m 2 N × the axial vector contribution. However, the identification C P,X does not work, because the spin sums suppress the axial-pseudoscalar interference term. A dedicated nuclear calculation would seem required for both the pseudoscalar and tensor operators. 7. To obtain the matrix-element-squared, the lepton spinor part can be evaluated by Dirac traces. Then to perform the nuclear spin sums in the SD case, the identity can be used.
8. Finally, the conversion rate is obtained as where |M| 2 is averaged over the incident spins, and dΠ gives the integration over the final state phase space of the nucleus and electron.
These steps give an analytic estimate for the four-fermion contributions to the SI conversion rate on a nucleus of atomic number A and charge Z: where the F N are defined in eqn (21) and related to the overlap integrals of KKO in (34), the contribution of the dipole operator (estimated in subsection 3.1.2) was also included, and where Γ capt is the rate for the Standard Model process of muon capture [12,36]. Similarly, the SD conversion rate on a nucleus of atomic number A, charge Z and spin J A is where the spin expectation values S A N and the finite momentum tranfer correction S A (k) are given for Aluminium at eqn (23).

the dipole
In the case of the dipole operator of eqn (2), the S-matrix element can be written where the 2 under the integral is to account for E i = F 0i = F i0 , and the magnitude of the radial electric field induced by the nucleus is [12] To estimate the dipole matrix element analytically, we suppose that the electric field only contributes at radii within the first zero of the electron wavefunction r e , because outside the rapid oscillation of the electron wavefunction gives an approximate cancellation in M. The muon wavefunction is approximately constant at such radii. The radius of the Aluminium nucleus is comparable to r e , but if we nonetheless approximate the nucleon distribution |f p (r)| 2 as a constant for r < r e , we obtain wherer is a radial unit vector. The "matrix element" M neglects recoil of the nucleus, so the final state phase space in the rate is only one-body, and we reproduce the analytic estimate of [12] for light nuclei (D ∼ 8eS p given above eqn (29) of [12]): This estimate uses r 3 dr/3 ≃ r 2 dr, and applies in the absence of other contributions; the dipole coefficient sums with the scalar and vector coefficients in the amplitude, as given in eqn (25).

Comparing to KKO
This section compares our estimates to the more exact formulae of [12] (KKO). Our estimates use a solution of the Schrodinger equation for the muon, a plane wave for the electron, and chiral γ-matrices. KKO solve the Dirac equation in the potential of the nucleus, both for the electron and muon, use Bjorken and Drell γ-matrix conventions, and give the branching ratio as: where Γ cap is the rate for the muon to transform to a neutrino by capture on the nucleus (see [12,36]), and the nucleus-and nucleon-dependent "overlap integrals" V (N ) X , D (N ) correspond to the integral over the nucleus of the lepton wavefunctions and the appropriate nucleon density (vector, scalar, electric field for the dipole operator; the definitions and numerical values are given in KKO [12]). The numerical coefficient in eqn (33) differs from the result given in KKO, because 4 C| here =g| KKO .
Our unit-normalised nuclear density |f N (r)| 2 can be identified with the similarly normalised density ρ N (r) of KKO [12]. Our Schrodinger approximation for the muon wavefunction can be identified to the upper component (in Bjorken and Drell γ conventions) of the Dirac wavefunction obtained by [12]. Then the normalisation conventions of eqn (5) and (7) of [12] identify In the limit of massless electron, the upper (g e ) and lower components (if e ) of the electron wave function of [12] are comparable. The electron normalisation condition given in eqn (8) of [12] then implies that we can identify our electron plane wave as In the approximation where the muon wavefunction is constant in the nucleus, the overlap integrals of [12] can be identified to our approximations as as given in eqns (29) -(31) of KKO.

Spin-dependent conversion in other light nuclei
In this section we consider how the estimates of the previous section could be applied to other nuclei. Recall that light nuclei are interesting for SD detection, because the SD rate is relatively suppressed by 1/A 2 compared to the SI rate: the ratio Γ SD /Γ SI is largest for light nuclei.
The matrix element given in eqn (17) for SD µ → e conversion contains the integral of the axial current over the nucleus, weighted by the lepton wavefunctions. In the case of light nuclei (Z < ∼ 20), as discussed in the previous section, the muon wavefunction can be taken constant in the nucleus, and the electron can be treated as a plane wave. This allows to use the results of nuclear calculations [27] of matrix elements for spin-dependent WIMP scattering at finite-momentum-transfer. The zero-momentum-transfer matrix elements (spin expectation values; see eqn (22)), have been calculated for a wide variety of nuclei [37], and finite momentum transfer results also been obtained for some nuclei [38]. For µ → e conversion in heavier nuclei, a dedicated nuclear calculation would be required to obtain the expectation values of the SD operators weighted by the lepton wavefunctions.
An interesting light nucleus for SD µ → e conversion could be Titanium (Z=22) 4 , because it has isotopes with and without spin, so targets of different isotopic abundances could allow to distinguish SD from SI operators. Titanium has a spin-zero isotope with A = 48 and 74% natural abundance [39], an isotope with A = 47, J = 5/2, 7.5% abundance, and another isotope with A = 49, J = 7/2, 5.4 % abundance. These natural abundances of more than 5 % are large enough to make sufficiently-enriched sample targets.
In the Odd Group Model, Engel and Vogel [40]  However, we observe that in Aluminium, the SI and SD form factors are comparable: 0.28 = |F p (m µ )| 2 ≈ S Al (m µ )/S Al (0) = 0.29. A similar relation appears to hold [12,38] 36. This suggests that for light nuclei, the spin-expectation-squared at | q| 2 = 0 (that is, S A (m µ )), is similar to the square of the spin-expectation-value at zero momentum transfer, multiplied by the square of the SI | q| 2 = 0 form-factor. Or taking the square root: So we apply this approximation to Titanium, and estimate S T i (m µ )/S T i (0) ≈ 0.12.

Parametric expansions and uncertainties
Once µ → e conversion is observed, the aim will be to determine (or constrain) as many operator coefficients as possible.
This would require at least as many "independent" observations as operators, where observations are independent if, in spite of uncertainties, they depend on a different combination of coefficients. So the purpose of this section, is to estimate the uncertainties in relating the conversion rate to operator coefficients. The inputs for this calculation, (equivalently, the theoretical parameters to be extracted from data) are the coefficients of either the quark operators (see eqn 1), or of the nucleon operators (see eqns 11,19), in both cases at the experimental scale µ N . So uncertainties associated to the Renormalisation Group evolution from the New Physics scale to the experimental scale are not considered. In the remainder of this paper, we will sometimes use the quark operator basis, and sometimes the nucleon basis. As discussed in point 1 below, there are significant uncertainties in some of the {G N,q O }, which are required to extract the coefficients of the quark operators, but can be avoided by using the nucleon operators.
1. There are uncertainties in some of the matching coefficients that relate quark to hadron operators (see eqn (6) and appendix A). The G N,q V are from charge conservation, so should be exact. For the axial and scalar coefficients, the determinations from data (see eqn (62)) and from the lattice(63,65) are quoted with smaller uncertainties than their differences (this is especially flagrant for the G N,q S , whose lattice and data values differ by 30-50%, and are both quoted with < ∼ 10% uncertainties). First, it can be hoped that these discrepancies will be reduced in the future. Secondly, in some models (or equivalently, for some choices of coefficients), these factors can be cancelled by taking ratios. Finally, if we are only interested in discriminating SD from SI contributions to the rate, this distinction exists at the nucleon level, so the matching to quark operators is not required.
2. The lepton interactions with nucleons are calculated at Leading Order (LO) in χP T . At NLO, arise pion loops as well as processes with two nucleons in the initial and final states which exchange a pion that interacts with the leptons. For the case of WIMP scattering, such NLO contributions for the scalar quark operator have been discussed [22,21,41] and reference [21] estimates them to be a higher order effect ( < ∼ 10%), provided there are no cancellations among the LO contributions. The two-nucleon contributions were also calculated to be unexpectedly small for WIMP scattering on few-nucleon nuclei [42]. However, after this manuscript was completed, appeared a study of the µ → e conversion rate mediated by the scalar and vector interactions [25], where the authors estimate that the NLO effects associated to pion exchange between two nucleons can reduce the scalar matrix element by 20→ 30%(NLO corrections vanish for the vector). We will account for these nucleon/χPT uncertainties by including them in the uncertainties in the overlap integrals.
3. The µ → e conversion matrix element, expressed as a function of nucleon operator coefficients, relies on many perturbative expansions, among which an expansion in the finite-momentum-transfer | q| 2 = m 2 µ . Naively such corrections are O(m 2 µ /m 2 N ) (so negligible), however in practise there are various effects which are not so suppressed. First, the finite momentum transfer gives a significant suppression of the matrix element. In our analytic approximations, where the muon is at rest and the electron momentum k = q, this is encoded in the form factors F N (see eqn (21)), which are ∼ .2 → 0.7. KKO include this effect more accurately, by solving the Dirac equation for the leptons. Secondly, finite momentum transfer effects can change the nucleon and lepton spinor algebra. This is discussed for Dark Matter in [28,35], and gives the O(m µ /m N ) contribution of the tensor to the scalar coefficient given in eqn (19). We include this correction, because the tensor operator at zero momentum transfer contributes to the SD matrix element (suppressed by 1/A), whereas this (m µ /m N )-suppressed contribution gains a relative factor A because it contributes to the SI matrix element. The ratio of these contributions to the conversion rate is estimated in appendix B. Finally, pion exchange becomes relevant at | q| 2 = m 2 µ for the axial vector and pseudoscalar operators (see eqns (8,9)), and is included in the nuclear matrix elements of [29] that we use for the axial vector in Aluminium. Pion contributions at | q| = 0 to the SI rate are discussed above. We hope that these are the dominant finite-momentum-transfer corrections, such that any other effects are negligible (< 10%) corrections. 4. In our calculation of the SD matrix element, the velocity of the initial muon was neglected. This may seem doubtful, by analogy with the extended basis of WIMP scattering operators constructed in [28], because these authors expand in both the momentum transfer between the WIMP and nucleon, and the incoming velocity difference. However, in our case, the muon velocity is parametrically smaller: writing the binding energy of the 1s state as πZαm µ ∼ m v 2 , gives | v| ∼ √ Zα. We neglect any effects related to this velocity.
5. There could be nuclear uncertainties in the SI overlap integrals S N , V N , D, in addition to the effects discussed in point 2 above. These were estimated by [12] to be ∼ a few % in most cases, < ∼ 10% in the case of some heavier nuclei.
Consider first the uncertainty on the SI rate, because, when µ → e conversion is observed in a nucleus with spin, the SD conversion rate can only be observed, if it is larger than the uncertainty in the ubiquitous A 2 -enhanced SI rate. The uncertainty in Γ SI , written as a function of the quark operator coefficients C qq O,X , would arise from the G N,q O , from the overlap integrals S N , V N , D of [12], and from NLO contributions in χP T : sums on N ∈ {n, p} and q ∈ {u, d, s} are implicit, the gluon contribution to the scalar [13] was neglected, for simplicity a common uncertainty δIA IA was assigned to the overlap integrals in nucleus A, except for the effect of neglecting pion exchange between two nucleons [21,25] (discussed in point 2 above), which is parametrised as an uncertainty [δS N ] N LO in the scalar overlap integrals. Expressed this way, the uncertainty depends on the quark coefficients present: for C qq S,X ≫ C qq V,X , C D,X , the current discrepancies in the determination of the G N,q S and [δS N ] N LO give an O(1) uncertainty on the conversion rate, whereas if only the C qq V,X and C D,X were present, the rate uncertainty would come from the overlap integrals. The G N,q S uncertainties can be avoided by expressing the rate in terms of the coefficients of the nucleon Lagrangian; if in addition, [δS N ] N LO /S N < 10%, then the uncertainty in the SI rate comes from the overlap integrals. From the KKO discussion, 2 δIA IA < ∼ 10% in most cases, < 20% in all cases. In order to be concrete, we assume in the remainder of this paper, that the uncertainty on the SI rate expressed in terms of coefficients on nucleons, is < ∼ 10%. This suggests that the SD rate would need to be > ∼ 10 − 20% of the SI rate, to be observed.
A better sensitivity to the SD rate could be obtained by using isotopes with and without spin as targets: consider for instance, 48 Titanium (without spin), and 47 Titanium (with spin), whose SI matrix elements differ by one neutron. Using the analytic approximation of eqn (25), the ratio of the SI conversion rates, for real coefficients and left-handed electrons, is where the second term 5 is O(1/A). The theoretical uncertainty in this ratio will arise from the overlap integrals (equivalently, form factors F N ), so should be of order 1 A δIT i IT i < ∼ 0.002. This greatly improves the sensitivity to the SD rate, although it is unlikely to allow as good a sensitivity to SD as SI coefficients, because the SD rate is parametrically suppressed as 1/A 2 which is < ∼ 1 A δIT i IT i . Some prospects for distinguishing among SI operators by using different targets will be discussed in section 5.3. For this, the various targets need to probe different combinations of operator coefficients, and this difference needs to be larger than the theoretical uncertainty. In section 5.3, targets are parametrised as vectors in coefficient space, whose components are the overlap integrals (see eqn (54)), and targets are assumed to probe different combinations of operator coefficients if the angle between the vectors is > ∼ 10% > ∼ δIA IA . This estimate can be obtained in the 2-dimensional plane of the vectors, where the uncertainty on the angle θ of a point (I 1 ± δI 1 , I 2 ± δI 2 ) is

Implications of including the SD rate
The aim of this section is to explore the implications of including the SD contribution to µ → e conversion. At first sight, it appears to be of limited interest: the ratio of SD to SI rates is so for a SI operator coefficient C SI comparable to C SD , the SD contribution to the branching ratio is much smaller than the ∼ 10% theory uncertainty of the SI contribution, estimated in the previous section. Furthermore, as discussed in [15], renormalisation group running between the New Physics scale and low energy mixes the tensor and axial vector ("SD") operators to the scalars and vectors, so even in a New Physics model that only induces SD operators, their dominant contribution to µ → e conversion is via the SI operators that arise due to RG running. This perspective that SD conversion can be ignored is illustrated in section 5.1, where we consider three leptoquark models. They give negligeable SD branching ratios, but we explore the prospects of distinguishing them by comparing the SI rates in various nuclei. The SD conversion rate is nonetheless interesting, because it is an independent observable, that can be observed by comparing targets with and without spin. As in the case of dark matter, it is sensitive to different operator coefficients (evaluated at the experimental scale) from the SI process, so provides complementary information. In section 5.2 we allow C SD ≫ C SI such that the SD rate can be observable, and discuss the constraints that could be obtained from upper bounds on µ → e conversion. Finally in section 5.3, we allow arbitrary coefficients to all the operators of the nucleon-level Lagrangian, and explore the prospects for identifying coefficients should µ → e conversion be observed.

Leptoquarks
We consider three possible leptoquark scenarios, each containing an SU(2) singlet leptoquark, whose mass M > ∼ few TeV respects direct search constraints [43], and which has only one coupling to electrons and one to muons. The scenarios are represented by adding to the Standard Model the following Lagrangians where D µ is the appropriate covariant derivative of QCD and QED. At the leptoquark mass scale, we match these scenarios onto the SM extended by QED*QCD invariant operators, which mediate µ → e conversion. The coefficients and operators are given in table 1.  In each scenario, we translate the coefficients down to the experimental scale µ N =2 GeV via an approximate analytic solution to the one-loop RGEs of QED and QCD [17,16]: where λ = αs(M) αs(µN ) ≃ 1/3 for M = TeV, and I, J represent the super-and subscripts which label operator coefficients. The a I describe the QCD running and are only non-zero for scalars and tensors. We suppose five quark flavours for the running, which gives a I = Γ s II 2β0 = {− 12 23 , 4 23 } for I = S, T . Γ e is the one-loop QED anomalous dimension matrix, Γ e is this matrix with an additional factor multiplying the T S and ST entries [44,45] in order to account for the QCD running: We neglect the RG mixing out of our operator basis, because it is small: tensor mixing to the dipole is suppressed by light quark masses, and the mixing via the penguin diagram to vector operators O f f V,X is a few %, and does not generate operators interesting to us here. The RG evolution is described in more detail in appendix C.
This formalism allows to predict the ratio of SD to SI contributions to the branching ratio for µ → e conversion. In Aluminium, we find for the three scenarios, taking M = 1 TeV: so we see that the SD rate is smaller than the current ∼10% uncertainties on the SI rate, so cannot be observed in these models. This is as expected, because the leptoquark model imposes that the tensor/axial coefficients are comparable to the scalar/vector coefficients, so the SD rate is relatively suppressed with respect to the SI rate by 1/(AG N,q S ) 2 for tensor coefficients, and 1/A 2 for axial vector coefficients.
It is interesting to explore whether the three leptoquark scenarios could be distinguished by comparing the SI rates in various nuclei. We imagine that µ → e conversion has been observed in Aluminium(Z Al =13, A Al = 27), the target of the upcoming COMET and Mu2e experiments. We wish to identify alternative target materials, which could allow to distinguish our leptoquark scenarios.
A simple distinction between the leptoquarks S and S, is that the former couples to u quarks, and the latter to d quarks. To identify an appropriate target (A, Z), where the µ → e conversion rates induced by S and S would be significantly different (subject to the constraint that both reproduced the Aluminium observations), we consider the double ratio: where the operator coefficients cancel because we compare two models that each induce a single SI operator. This ratio should differ from 1 by > ∼ 20%, in order to unambiguously distinguish the S from S, given the ∼ 10% uncertainties on the theory calculation. The first approximate equality in eqn (45), applies for light nuclei, where the conversion rate can be written as eqn (25). The second equality uses the KKO conversion rate given in eqn (33) in terms of the overlap integrals V (N ) , and applies for all nuclei. The continuous green line (with stars) of figure 2 is the ratio of µ → e conversion rates mediated by S and S, assuming equal operator coefficients. It corresponds to the second fraction in the products appearing in eqn (45), so the double ratio of eqn (45) is simply obtained by dividing by the ratio for Aluminium. The stars are the light nucleus approximation, the green continous line is the ratio of overlap integrals. This shows that the approximation is very similar to the numerical results of KKO, and that a target with Z > ∼ 40 could allow to distinguish the first and third leptoquark scenarios. In the following, we take Niobium (Nb,Z=41,A=93) as a Z > ∼ 40 target.
It is also interesting to explore the prospects of distinguishing scalar operators involving u vs d quarks. So we also plot in figure 2, as a dashed red line, the ratio of µ → e conversion rates mediated upstairs by O uu S,X and downstairs by O dd S,X : For the G N,q S values given in appendix A, the scalar ratio is close to one (because G p,q S ≃ G n,q S ), suggesting that changing the target in µ → e conversion does not help distinguish O uu S,X from O dd S,X . The first and second leptoquark scenarios respectively induce scalar and vector operators. As discussed in [12,13], these can be distinguished by comparing the conversion rates in light and heavy targets. This is illustrated in figure  2, by the blue dotted line, which gives the double ratio normalised to Niobium  We see that measuring the µ → e conversion rate on Aluminium, some intermediate target around Z ∼ 40 and on a heavy nucleus like lead or gold (Z = 79), could distinguish the three leptoquark scenarios, that is a vector operator involving ds, vs vector operator involving us, vs a scalar operator involving us.

Bounds on arbitrary coefficients of four operators
This section considers the operators induced by the second and third leptoquark models (see equations (40), (41)) which are added simultaneously to the Lagrangian with arbitrary coefficients: This is clearly an incomplete basis (the complete basis of dimension six operators at µ N is given in eqns (1,3)); however, it is sufficient for our purpose 6 , which is to explore which constraints can be obtained on quark-level operators from the non-observation of µ → e conversion in targets with and without spin. We suppose that µ → e conversion has not been observed on Aluminium, Titanium (enriched in isotopes with spin) and Lead targets. These targets are chosen because heavy and light targets have different sensitivities to vector and scalar coefficients, and because the spin of Titanium and Aluminium is respectively associated to an odd neutron and an odd proton. In order to check that upper bounds on these branching ratios can constrain all the operator coefficients which we consider, we set the branching ratios to zero, and check that this forces the coefficients to vanish.
Setting the SD conversion rates in Titanium and Aluminium to zero gives two equations: so even allowing for a 10% theory uncertainty on the coefficients, it is clear that the only solution is for both coefficients to vanish. This is because the spin of Titanium isotopes arises from the odd number of neutrons, whereas in Aluminium the spin is from the odd proton, so these two conversion rates probe the SD coefficients C ′ N N A,X for both neutrons and protons. Then, since the matching coefficients G N u A,X and G N d A,X (equivalently G N u T,X and G N d T,X ) are of opposite sign and different magnitude, C uu A,X + 2C uu T,X and C dd A,X + 2C dd T,X can be distinguished. It is straightforward to check that setting the SI rates on Al, Ti and Pb to zero, forces C dd V,R , C uu S,L → 0. A more informative way to present the constraints on coefficients arising from the experimental bounds is to give the covariance matrix. We suppose an upper bound of BR (for instance, 10 −14 ) on the SI branching ratios on Pb and Al, and on the SD branching ratios on Al and Ti. The tensor operator gives comparable contributions to both SI and SD processes (see Appendix B), so the 4 × 4 covariance matrix does not split into 2 × 2 subblocks. Nonetheless, it is interesting to give the covariance matrices for different cases, in order to see the variation of the bounds, when different theoretical information is included.
First, the tensor contribution to the SI rates is neglected, in which case the covariance matrices for (C dd V,R , C uu S,L ) and (C uu T,L , C dd A,R ) are: So, for instance, |C uu S,L | is excluded above √ 0.0007 × BR, and |C dd A,R | < √ 73.6 × BR. If now the SD rates are neglected, but the tensor contribution to SI is included, then the covariance matrix for which shows that the exclusions become weaker due to potential cancellations between a large C uu T,L and the vector or scalar coefficients. Finally the full covariance matrix arising from imposing BR SI (µP b → eP b) ≤ 10 −14 , BR SI (µT i → eT i) ≤ 10 −14 , BR SD (µ 47 T i → e 47 T i) ≤ 10 −14 , BR SI (µAl → eAl) ≤ 10 −14 , and BR SD (µAl → eAl) ≤ 10 −14 , for the coefficients (C dd V,R , C uu S,L , C uu T,L , C dd A,R ), is Comparing to the bounds of eqn (51), shows that the tensor contribution to the SI rate is of little importance, provided the SD bounds are included. However, if the SD bounds are neglected, including the tensor in the SI rate significantly weakens the constraints, as can be seen in eqn (52). We also checked that including BR SI (µAu → eAu) ≤ 10 −14 only changes a few entries by about 25%, as expected, because Al, Ti and Pb were chosen as targets for their discriminating power.

Reconstructing nucleon coefficients
We now suppose that µ → e conversion is observed in Aluminium, where there can be SI and SD contributions to the rate, and that the New Physics is described by the nucleon-level Lagrangian of eqn (5) with arbitrary operator coefficients. It is interesting to consider which subsequent targets, in what order, would be required to distinguish the SD and SI contributions, and then to discriminate among the SI operators?
We first introduce a geometric representation of models and targets, which allows to visualise the ability of various targets to discriminate among models. A New Physics scenario can be represented as a two 5-dimensional vectors, each composed of SI coefficients which interfere C X ≡ (C D,X , C ′ pp S,X , C pp V,Y , C ′ nn S,X , C nn V,Y ), and two two-component vectors of SD coefficients ( C ′ nn A,X , C ′ pp A,X ). For simplicity, we focus on X = L, and drop this electron chirality subscript. Then we focus on discriminating among SI operators, because the spin of target nuclei is usually associated to either an unpaired n or p, giving an order of magnitude better sensitivity to the coefficient corresponding to the unpaired nucleon (see, e.g. the spin expectation values given after eqn (22)). This means that discriminating C ′ nn A,X vs C ′ pp A,X should be a straightforward matter of using targets with an unpaired p and n.
For the spin independent process, a target nucleus (Z, A) can be envisaged as a vector in the five-dimensional coefficient space, whose components are the appropriate overlap integrals. (In the following, the vectors and components are indiscriminately labelled by A or Z because we use the overlap integrals of KKO, obtained for a single abundant isotope.) The matrix element for µ → e conversion on target A, mediated by a combination of coefficients C, is proportional to C · v A , and target nucleus A allows to probe coefficients in the direction v A . If we define the unit-normalisedê A = v A /| v A |, then target A probes the same combination of coefficients as Aluminium if e A is parrallel toê Al , and the difference gives an invariant measure of whether the target A has sensitivity to an orthogonal direction in coefficient space. In eqn (55), θ is the angle betweenê A andê Al . Figure 3 givesê A ·ê Al as a function of Z. From eqn (38), the uncertainty in the direction ofê A is < ∼ 10%, so target A is indistinguishable from Aluminium forê A ·ê Al > ∼ 0.995, or Z < ∼ 25 − 30.
Perhaps a more transparent measure of the change of direction ofê A in coefficient space, is given in figure 4 by where with comparable coefficients could be difficult to distinguish from O pp V .) This normalised ratio of overlap integrals is interesting, because the normalisation "factors out" the growth with Z shared by all the overlap integrals, so this ratio parametrises the difference in direction in coefficient space, which allows different targets to discriminate amoung coefficients. This ratio also indicates that targets of Z < ∼ 25 cannot distinguish operators, if one admits a theory uncertainty of ∼10% in the calculation of the components e O A . Assisted by the measures of discriminating power given in eqns (55) and (56), we now speculate on a possible series of targets. A light nucleus without spin could be an interesting second target, because it would allow to distinguish whether the rate in Aluminium was dominantly SD or SI. In particular, the SI rate in Aluminium could be approximately predicted from the the rate observed in another spinless light nucleus. This is because the SI rate in all targets with Z < ∼ 20 is sensitive to a similar linear combination of operator coefficients, as illustrated in figures 3 and 4.
An interesting choice for the second target could be Titanium (Z=22, A = 48). As illustrated in figures 3 and 4, it of sufficiently low Z that the SI rate probes the same combination of operator coefficients as the SI rate in Aluminium. So measuring the SI rate in Titanium-48 would allow to determine whether there was a significant SD contribution to the µ → e conversion rate observed on Aluminium.
If there is indication for an SD contribution in Aluminium, then it could be interesting to measure the rate on a Titanium target enriched with the spin-carrying isotopes 47 and 49. This would give complementary information on the quark flavour of the tensor and/or axial vector operators, because the spin of Aluminium is largely due to the odd proton, whereas for Titanium, there is an odd neutron. So the SD rate in Aluminium is mostly sensitive to C ′ pp A,X , whereas the SD rate in Titanium depends on C ′ nn A,X . Finally, if there is no evidence of an SD rate in Aluminium, a heavy target such as lead could be interesting to discriminate the scalar vs vector coefficients in the SI rate.  (55), of the misalignment in coefficient space of the target with respect to Aluminium.

Summary
This paper gives some details of the calculation of the Spin Dependent (SD) µ → e conversion rate in light nuclei, previously outlined in [15]. Section 2 reviews the operators involving quarks and gluons that contribute [12,13] at the experimental scale (µ N = 2 GeV), and matches them onto the nucleon operators which enter the nuclear physics calculation. Some attempt is made to include pion exchange in this matching (it is relevant because the momentumtransfer is m 2 µ ). Section 3 calculates as much as possible of the conversion rate in the notation of relativistic, secondquantised, QFT [30]; in the last steps, the results of nuclear calculations are included. The final rate is given in equation (26). This section is not original; its purpose is to make the result accessible to affictionados of QFT. We recall the SD µ → e conversion is incoherent, like SD WIMP scattering, so it is best searched for in light nuclei, where the 1/A 2 suppression with respect to the coherent Spin Independent (SI) rate (given in eqns (25,33)) is less significant.
Our SD rate estimate relies on nuclear physics calculations of the expectation value of nucleon axial currents in the nucleus. The results we use were obtained for SD WIMP scattering, which are often at zero momentum transfer. As discussed in point 6 of section 3.1, additional nuclear calculations seem required to include tensor and pseudoscalar operators at finite momentum transfer, in light targets such as Aluminium and Titanium. In this paper, we did not discuss SD conversion on heavy nuclei; however, one can speculate that the nuclear expectation values could be of interest, because heavy nuclei could be sensitive to a different combination of tensor and axial operators from light nuclei. This is because the anti-lepton wavefunction contributes with opposite sign to the tensor vs axial operators, and is more relevant in heavy nuclei (this sign difference allows to discriminate scalar and vector operators in SI conversion on light and heavy nuclei [12]). Of course, the SD rate might be unobservably small (due to the 1/A 2 suppression), but heavy nuclei could nonetheless give an independent constraint on the many operator coefficients.
Both the SD and SI conversion rates depend on the modulus-squared of a sum of coefficients, weighted by nucleusdependent numbers-see equations (25,26,33). This allows for cancellations, making it difficult to constrain individual coefficients, or identify the operators responsable for µ → e conversion when it is observed. In the SI case, Kitano Koike and Okada (KKO) [12] pointed out that scalar vs dipole vs vector operators could be distinguished by changing the nuclear target. Section 5 explores, from various approaches, the prospects of distinguishing a wider variety of operators, including SD vs SI, and u-vs d-quark operators.
The prospects for discriminating vector or scalar operators involving either u or d quarks are illustrated in figure  2: vector operators involving u or d quarks could be distinguished by comparing the µ → e conversion rate in light (Z < ∼ 20) and intermediate (Z ∼ 40) targets, but distinguishing scalar u versus d operators seems difficult. Curiously, the u vs d distinction is more transparent in the SD rates, as discussed after eqn (50). So if both SD and SI conversion are observed, possibly the quark flavour could be extracted from the SD rates 7 .
The SD and SI contributions to the conversion rate could be distinguished (if the SD rate is large enough) by comparing the conversion rate in nuclei with and without spin. Section 4 reviews the theoretical uncertainties in the calculation of the µ → e conversion rate, in order to estimate the sensitivity to the subdominant SD process.
Comparing µ → e conversion on isotopes with and without spin would cancel the leading theory uncertainties, giving a sensitivity (see the discussion after eqn 37) to Γ SD /Γ SI > ∼

0.1
A , assuming a 10% uncertainty on Γ SI . Among the SD operators, it is not currently possible to distinguish pseudoscalar, axial and tensor coefficients, because only the nuclear expectation value of the axial operator has been calculated. However, as mentioned in the previous paragraph, it could be possible to discriminate SD operators involving u vs d quarks, because they contribute differently in nuclei where the odd nucleon is a proton or neutron.
The upcoming COMET and Mu2e experiments will initially search for µ → e conversion on Aluminium, a target which has spin -so if they observe a signal, it could be mediated by the SD or SI operators. So in section 5.3, we considered what series of subsequent targets could give information about the dominant coefficients. To this purpose, we represent a target material as a vector in the space of nucleon-level operators, whose components are numbers which multiply the operator coefficient in the rate (overlap integrals, in the SI case). Different targets can discriminate between operators, if they point in different directions of operator space. We plot in figures 3 and 4 two different measures of the misalignment between target vectors.
If µ → e conversion is observed on Aluminium, the following sequence of targets could be interesting: as second target, a light nucleus without spin, such as Titanium-48, would discriminate whether the dominant contribution was from the SD rate, because the SI rate in Titanium is comparable to Aluminium (see figures 3 and 4). If there is an SD contribution to the rate in Aluminium, then Titanium isotopes with spin, could be an interesting target: the spin of Titanium is related to the odd neutron (whereas in Aluminium there is an odd proton), so this could discriminate whether the SD operators involved u or d quarks. Finally, a heavy target such as gold or lead could allow to discriminate scalar vs vector operators, as pointed out in [12].
When the quark Lagrangian of eqn (1) is matched onto the nucleon Lagrangian, the coefficients of the nucleon operators can be computed as for O ∈ T, A, V, P ; for the scalar operator there is an ad-ditional gluon contribution as described in [13]. We take the G N,q O , defined at zero-momentum-tranfer such that , G p,s T = G n,s T = .008 (9) .
where the parenthese gives the uncertainty in the last figure(s). The axial G A are the results inferred in Ref. [46] by using the HERMES measurements [47]. The scalar G S induced by light quarks are from a dispersive determination [48], and an average of lattice results [49] is used for the strange quark; in all cases, the MS quark masses at µ = 2 GeV are taken as m u = 2.2 MeV, m d = 4.7 MeV, and m s = 96 MeV [50]. The nucleon masses are m p = 938 MeV and m n = 939.6 MeV. The pseudoscalar results were calculated from data in the large-N c approximation at q 2 ≈ 0 [51], and here extrapolated to neutrons using isospin. The tensor results for the neutron are the lattice results of Cirigliano etal [52], which are here extrapolated to protons using isospin. For comparaison, the G A have been obtained on the lattice; a recent determination [53] is The scalar G N,q S have also recently been obtained on the lattice [54]: We observe that there is a 50% discrepancy with respect to the results of [48], obtained from pionic atoms and π − N scattering [55]. Results similar to [48] were earlier obtained in [56], also using an effective theory.

B The tensor contribution to the SD and SI rates
We consider tensor operators at the experimental scale µ N , which contribute at finite-momentum-transfer to the SI conversion process (see eqn (19)), and also to the SD processes: The ratio of these contributions, for a single operator, is where we assumed that the form factors are comparable SA(mµ) SA(0) ≃ |F p (m µ )| 2 as is the case in Aluminium. Recall that G n,u T ∼ − 1 2 G p,u T , so there is a partial cancellation in the SI amplitude, whereas the SD process arises mostly from an odd proton S A p ≫ S A n , or mostly from an odd neutron S A p ≪ S A n .
The estimates of eqn (69) assume that only one tensor coefficient is non-zero, so they neglect interferences, which can easily enhance the SI rate. For instance, RG running of the tensor operator from the New Physics scale to the experimental scale generically generates a scalar operator with comparable coefficient. The scalar-tensor interference contribution to the SI rate would be relatively enhanced, with respect to the tensor-squared, by G N,q S /G N,q T ∼ 10, which would suppress the ratio in eqn (69) by another factor 1/10.

C RG Evolution
In this appendix, we review the Renormalisation Group evolution of operator coefficients from the leptoquark mass scale M (∼ TeV) down to the experimental scale µ N (2 GeV), via the one-loop RGEs of QCD and QED [17,16]. We consider the QED× QCD invariant operator basis discussed in section 2. We neglect matching onto the SMEFT basis [57,58] and running with the full SM RGEs [59], on the assumption that QED is a reasonable approximation if M is not much larger than m W .
After including one-loop corrections in the M S scheme, the operator coefficients will run with scale µ according to [16] µ ∂ ∂µ where I, J represent the super-and subscripts which label operator coefficients, Γ e and Γ s are the QED and QCD anomalous dimension matrices and − → C is a row vector that contains the QCD ×QED invariant operators coefficients listed in section 2.
In this work, we use the approximate analytic solution [15] given in eqn (42): where the factors are given after eqn (42) and log M µN ∼ 6.21. Only QED loops contribute to operators mixing, while QCD loops only rescale scalar and tensor operators. In figure 5, we present the QED diagrams required to compute the anomalous dimension γ of the four-fermion operators, where f 1 ∈ {e, µ} and f 2 ∈ {u, d, s, c, b, e, µ, τ }. Figure 5: Examples of one-loop gauge vertex corrections to 4-fermion operators. The first two diagrams are the penguins. The last six diagrams contribute to operator mixing and running, but can only change the Lorentz or gauge structure of the operators, not the flavour structure. Missing are the wave-function renormalisation diagrams; for V ± A Lorentz structure in the grey blob, this cancels diagrams 3 and 4.
The operators coefficients below the scale M are organized in the vector − → C as following : In the basis of − → C , the QED anomalous dimension matrix can be written Vector and axial operators The first penguin diagram and the last four give the following matrices : Using these anomalous dimension matrices and the RGEs give : where q ∈ {u, d}. We see that axial operators mix to vector operators and vice versa, but there is no rescaling for axial and vector operators.

Scalar operators
Combining the third and fourth diagrams of figure 5 with the wavefunction diagrams renormalize the scalars while the last four diagrams mix the tensors to the scalars :

Tensor operators
Similarly, the last four diagrams mix the scalars to the tensors. Only the wavefunction diagrams renormalize the tensors, because for the third and fourth diagrams γ µ σγ µ = 0. We obtain the following matrices : Finally, the coefficients at the experimental scale µ N are obtain via the matching condition :