Connecting Dark Matter UV Complete Models to Direct Detection Rates via Effective Field Theory

Direct searches for WIMPs are sensitive to physics well below the weak scale. In the absence of light mediators, it is fruitful to apply an Effective Field Theory (EFT) approach accounting only for dark matter (DM) interactions with Standard Model (SM) fields. We consider a singlet fermion WIMP and effective operators up to dimension 6 which are generated at the mass scale of particles mediating DM interactions with the SM. We perform a one-loop Renormalization Group Evolution (RGE) analysis, evolving these effective operators from the mediators mass scale to the nuclear scales probed by direct searches. We apply our results to models with DM velocity-suppressed interactions, DM couplings only to heavy quarks, leptophilic DM and Higgs portal, which without our analysis would not get constrained from direct detection bounds. Remarkably, a large parameter space region for these models is found to be excluded as a consequence of spin-independent couplings induced by SM loops. In addition to these examples, we stress that more general renormalizable models for singlet fermion WIMP can be matched onto our EFT framework, and the subsequent model-independent RGE can be used to compute direct detection rates. Our results allow us to properly connect the different energy scales involved in constraining WIMP models, and to combine information from direct detection with other complementary searches, such as collider and indirect detection.


Introduction
The nature of dark matter (DM) is one of the main open questions in particle physics. Among many candidates [1][2][3], a Weakly Interacting Massive Particle (WIMP) with relic abundance obtained through thermal freeze-out [4][5][6][7] is quite appealing. Motivated frameworks for physics beyond the Standard Model (SM) naturally have a WIMP candidate [8][9][10][11][12][13], and it is suggestive that the same theory addressing the hierarchy problem also provides us with a DM candidate. Another exciting feature of the WIMPs is the fact that their typical couplings to SM particles are in the correct ballpark to give signals at colliders, direct and indirect detection experiments. Each of these searches is more sensitive to a certain parameter space region, so the WIMP paradigm can be tested with multiple and complementary methods.
In this work we develop a formalism to connect DM models to nuclear scales probed by direct detection. As shown in Refs. [57][58][59][60][61][62][63][64][65][66], there are examples where this large separation of scales has remarkable implications when a comparison with experiments is attempted. We focus on models for fermion DM with no SM gauge charge, and we assume that all the non-SM particles (with the possible exception of the DM itself) are above the weak scale. Fermion singlets cannot communicate with the SM at a renormalizable level, thus DM interactions must be necessarily mediated by these heavy particles. Once they are integrated out, the subsequent study to make a connection with nuclear scales does not depend on the specific model we started from. This general analysis is the goal of our paper, where we connect physics at the mediator mass scales with direct detection observables. A generic model of singlet fermion DM, Dirac or Majorana, can be matched onto our EFT framework at the mass scale of the mediator particles.
Our setup is sketched in Fig. 1. At high energy scales we imagine the DM field χ embedded in an ultraviolet (UV) complete model. We neglect the UV details and consider the low-energy EFT with only χ and the SM fields in the spectrum, which is suitable for our model-independent analysis. As emphasized in Ref. [65], a systematic study allows us to identify mixing among operators and bound interactions that are poorly constrained otherwise. For this reason we start from the most general basis of operators up to dimension 6, defined at the EFT cutoff Λ, which corresponds to the mediators mass scale. The operators are evolved via a proper Renormalization Group (RG) analysis down to the ElectroWeak Symmetry Breaking (EWSB) scale, where the W and Z gauge bosons, the Higgs boson and the top quark are integrated out. This procedure defines a different EFT, with only strong and electromagnetic gauge interactions and 5 quark flavors. We perform the RG analysis in this EFT as well, taking into account threshold corrections to the Wilson coefficients from integrating out the b and c quarks and the τ lepton, and evolving down to the nuclear scale µ ∼ 1 − 2 GeV at which hadronic matrix elements are evaluated (see e.g. [21,67,68]).
Why is this study relevant? After all, RG corrections are of the order of log(Λ/1 GeV) multiplied by a loop factor, and pushing Λ to 10 TeV barely changes the order of magnitude for the rate. One may question the usefulness of evaluating the precise numerical value of the cross section in a pre-discovery era. We are certainly not after such a precision, and our focus is rather on models where loop effects are the dominant contribution. The only DM interactions at the nuclear scale relevant for direct detection involve the u, d, s quarks, gluons and photons. However, many motivated models have mediator fields coupling the DM particle to heavy SM states, and in this case the main contribution to direct detection rates comes from loop effects. Furthermore, different light quarks couplings yield direct detection cross sections which could differ by orders of magnitude, as Goodman and Witten showed in their seminal paper [69]. If the mediator fields induce suppressed couplings to light quarks (e.g. DM velocity-suppressed and/or spin-dependent interactions), loop-induced couplings to non-suppressed operators are again the dominant contribution. The best current experimental limits come from XENON100 [70] and LUX [71], and will be significantly improved soon by SCDMS [72] and XENON1T [73]. They rule out electroweak processes with Z boson exchange by orders of magnitude, and are therefore powerful enough to put constraints even on loop-induced processes.
The paper is structured as follows. The bases of independent operators for both the EFTs in Fig. 1 as well as matching conditions at the EWSB scale are discussed in Sec. 2. The RGE equations in both EFTs are presented in Sec. 3, with details on loop calculations contained in App. B. The reader only interested in our results, not in their derivation, can safely jump from Sec. 2 to Sec. 4, where we present the applications of our results to spin-independent searches. Consistently with the spirit of this work, we focus on examples where the DM has either suppressed couplings to light quarks or couplings only to heavy SM states. In these cases our loop effects are the main contribution to spin-independent direct detection rates. In App. D we give a straightforward recipe that allows one to apply our results and constrain UV complete fermion WIMP models that give rise to dimension 6 effective operators. Sec. 5 contains our conclusions.

The Effective Theories for Singlet Fermion Dark Matter
Our conceptual starting point is a renormalizable model for a fermion DM field χ that is a SM gauge singlet. Interactions between χ and the SM degrees of freedom ψ SM are due to the exchange of mediator fields Φ. The typical mass of the Φ's is assumed to be greater than the Fermi scale, and at such scales the full SM gauge symmetry SU (3) c × SU (2) L × U (1) Y is unbroken. The Lagrangian of the UV complete model schematically reads In the remaining part of this Section we give a basis of independent operators for both EFTs up to mass dimension 6, as well as a prescription for how to match SM χ EFT onto EMSM χ EFT at the EWSB scale.

SM χ Effective Theory
Right below the mediator scale Λ all the SM degrees of freedom are in the spectrum, and the SU (3) c × SU (2) L × U (1) Y SM gauge group is unbroken. Integrating out the mediators Φ in Eq. (2.1) generates an infinite tower of higher dimensional operators Our conventions for the SM Lagrangian L SM are summarized in App. A. In particular, since SM fermions are in a chiral representation of the gauge group, we use the matter fields The index i runs over the three different SM fermion generations, and the gauge quantum numbers are assigned as in Table 1. The index α runs over all gauge invariant operators of a given dimension d, with the dimensionless Wilson coefficients c (d) α encoding unresolved dynamics. These coefficients are renormalization-scale dependent, and we will quantify this dependence in the next Section.

Symbol
Operator Without the need of specifying the responsible symmetry, we make sure the DM field is stable by requiring that every operator contains at least two χ fields. As an example, if DM is stabilized by a Z 2 symmetry, only operators with an even number of χ fields are allowed. Furthermore, our focus is on DM elastic scattering off target nuclei, thus we only need to consider operators with two DM fields. A complete basis of DM bilinears O αχ is O αχ = χχ , χγ 5 χ , χγ µ χ , χγ µ γ 5 χ , χσ µν χ . ( αSM . (2.5) The part involving only SM fields O αSM has mass dimension d − 3, is a SM gauge singlet but not necessarily a Lorentz singlet. We derive a complete basis of operators for the SM χ EFT up to dimension 6 following this strategy: we first identify all possible gauge singlets O αSM by employing the same procedure described in Refs. [74,75], then we take all allowed Lorentz invariant contractions with DM bilinears in Eq. (2.4).
The first operators to look for are at dimension 5, which implies that we need gauge invariant SM operators O (2) αSM with mass dimension 2. The following options are available where H and B µν are the Higgs doublet and the hypercharge field strength, respectively. The Lorentz invariant combinations with DM bilinears are listed in Table 2. Since they are the lowest dimensional non-renormalizable operators, they do not mix onto other ones. As is well known, the resulting long-range interaction arising from the dipole operators severely constrain their Wilson coefficients [76][77][78]. The dimension 5 Higgs portal operator induces interactions with the gluon field strength once heavy quarks are integrated out [79], and current experiments are probing cross sections in its typical range [80][81][82][83][84]. DM interactions to the Higgs also yields mono-Higgs events at colliders [85][86][87], and for light enough DM (m χ < m h /2) they contribute to the invisible Higgs decay width [88][89][90][91]. No interesting mixing takes place in this dimension 5 sector [65], and for this reason our RG analysis will focus on dimension 6 operators, which we now identify. For dimension 6 operators, the relevant SM structures of dimension 3 are the currents  Table 3. Basis of dimension 6 operators for the SM χ EFT. The first two columns have three different replicas, corresponding to the SM generations. We consider a generic χ Γ µ χ, which can be either a vector (Γ µ = γ µ ) or an axial (Γ µ = γ µ γ 5 ) DM current or any linear combination of them.
where we do not assume any flavor violation The index i runs over the three different fermion generations, thus the above vector has 5 × 3 + 1 = 16 components. The doublearrow derivative entering the Higgs current reads Lorentz invariant operators can be obtained by contracting the currents in Eq. (2.7) with a DM current χ Γ µ χ, where both vector Γ µ = γ µ and axial Γ µ = γ µ γ 5 currents are possible. This gives a total of 16 × 2 = 32 independent operators. However, since χ is a singlet, the DM current χ Γ µ χ is invariant under RG evolution, thus we can study two 16-dimensional sectors separately. The basis for dimension 6 operators with a specific DM current χ Γ µ χ is shown in Table 3. For future convenience, we introduce a 16-dimensional vector of Wilson coefficients where c α is associated with the operator O α in Table 3. The solid double line divides DM interactions with the Higgs from the ones with SM fermions. The solid single lines divide different SM generations and within each generation quarks and leptons are divided by a dashed line.
We stress that the dimension 6 operator does not need to be included in our list since it can be expressed as a linear combination of the ones listed in Table 3 by using classical equation of motion [92] for the hypercharge field strength (see Eq. (A.10)). More specifically, the effect of this operator can be absorbed Table 4. SM matter fields and their gauge quantum numbers in the broken phase. In this case q i , u i and e i are Dirac fermions. The index i = 1, 2, 3 runs over the threes different generations. The top quark (i.e. u 3 ) is not included.
into the following shifts of the Wilson coefficients

EMSM χ Effective Theory
The construction of the operator basis for the EMSM χ EFT is analogous. The Lagrangian reads (2.17) Details and conventions for the renormalizable EMSM part can be found in App. A. Since matter fields fill vector-like representations of the gauge group in this phase, we employ the set of Dirac fermions The index i runs again over the three different SM generations, but this time without the top quark (i.e. without u (3) ). Gauge quantum numbers are assigned as in Table 4. The effective operators are still of the form The top quark is not in the spectrum, thus we count 6 × 3 − 2 = 16 independent currents. Also in this case they can be contracted with either a vector or and axial DM current, giving a total of 32 independent operators. Each 16 dimensional sector shown in Table 6 can be studied separately. In analogy to Eq. (2.9), we define the vector ΓAe .
(2.22) Here, the solid double line is used to divide DM couplings to a vector or an axial SM current, whereas single solid lines divide quarks from leptons.
The redundant dimension 6 operator in this case is Equations of motion for the electromagnetic field strength (see Eq. (A.12)) translates this operator into a linear combination of the ones listed in Table 6, which equivalently amounts to this shift of the Wilson coefficients for the operators with SM vector currents The operators with SM axial currents are not affected, since the photon only couples to vector currents.

Matching the two EFTs at the EWSB scale
We conclude this Section by giving matching conditions between the two theories, namely the relations between the Wilson coefficients in Eq. (2.22) and those in Eq. (2.9), both evaluated at the EWSB scale, which is smaller than Λ in this setup. As we will see shortly, the leading contribution arises already at tree level, therefore we do not need to consider the subleading one-loop contributions. When performing the tree-level matching, going from left-and right-handed currents to vector and axial currents is straightforward. But this is not the full story, since the operator coupling the DM to the Higgs current leads to the following contribution obtained by giving the Higgs doublets an EWSB VEV. The result is to induce an effective tree-level coupling between the DM and the Z boson The Z boson also couples to SM fermions where the neutral current J µ 0 is defined as follows Here, T 3 f is the third component of the weak isospin, s w the sine of the weak mixing angle and Q f the fermion electromagnetic charge. The coefficients for the SM fermions explicitly read (2.32) Integrating out the Z boson gives rise to the Fermi Lagrangian for SM neutral currents Analogously, tree-level Z exchange gives a finite threshold corrections to the Wilson coefficients of the EMSM χ EFT. The complete matching conditions read

Renormalization Group Evolution
We present the complete one-loop RG equations in both EFTs. Here, we only show Feynman diagrams and quote final results. Regularization and renormalization at one loop in both EFTs are detailedly discussed in App.B. As explained in the previous Section, no interesting loop effect takes place among the dimension 5 operators, besides the well known heavy quark threshold contribution from the Higgs portal [79]. Thus we focus on dimension 6 operators.

From the messenger scale to the EWSB scale
The evolution of the Wilson coefficients in Eq. (2.9) is described by the differential equation where µ is the renormalization scale and γ SMχ is the anomalous dimension matrix. Our goal here is to fill out the 16 × 16 = 256 entries of the matrix γ SMχ . We start our one-loop analysis in this theory by considering external legs corrections. Since the DM field is a gauge singlet, these contributions only involve SM fields and interactions. We perform the field renormalizations where ψ i is any SM fermion, and we do it in such a way to subtract the infinite part from the residue of each one-loop propagator. There are only two possible sources for this effect, which are gauge and Yukawa interactions. As is well know, the Higgs quartic coupling does not induce a one-loop contribution to the wave-function renormalization. The relevant Feynman diagrams are shown in Figs. 2 and 3 for fermion and Higgs fields, respectively. When considering vertex corrections, one still has to deal only with these two interactions. We organize the presentation by fixing the external legs of a specific amplitude, and then identifying all the possible one-loop contributions. In other words, we fix a given effective operator from the ones in Table 6 and then look for operators mixing into it. We start from the loop corrections to the Wilson coefficient c L , which can be induced by gauge interactions (diagonal renormalization) and by Yukawa interactions (off-diagonal renormalization). The associated Feynman diagrams are shown in Fig. 4. Loop effects for c u R and c d R are analogous, with the important difference that right-handed fermions have no SU (2) L interactions, and therefore there are no diagrams with W i µ in the loop. The analysis of loop corrections to c H involves many more Feynman diagrams. The associated operator describes the DM interaction with two Higgs bosons, and we expect by gauge invariance also diagrams with two Higgses and one electroweak gauge boson. We have computed all the possible one-loop diagrams, both the ones with only two Higgs fields on the external legs and the ones with an additional W i µ or a B µ gauge boson on the external legs, and checked that they combine in a gauge invariant way to give the Higgs current defined in Eq. (2.8). We show in Fig. 5 only the diagrams with two external Higgs bosons, which get contributions from gauge and Yukawa interactions. For the latter, we checked that the associated diagrams with an external B µ combine in a gauge invariant way with the others only after using the SM hypercharge values, as it should be.
Finally, loop diagrams in Fig. 6 radiatively induce a contribution to the Wilson coefficient c B of the redundant operator in Eq. (2.10). By consistently using the equation of motion, this translates into a shift for the independent operators as in Eqs. (2.11)-(2.16). 1 .
The explicit expression for the one-loop amplitudes, as well as the consequent derivation of the RG equations, are presented in App. B. As shown there, the gauge interactions contribution to wave-function renormalization exactly cancels with the associated vertex corrections. This is consistent with the Ward identities in abelian gauge theories, and with the fact that non-abelian gauge interactions renormalize only axial currents starting at two loops. The final anomalous dimension matrix has two main pieces where the contribution proportional to the hypercharge gauge couplings is a consequence of the diagrams in Fig. 6. The explicit expressions read

4)
1 Alternatively, instead of dealing with a redundant operator, we can restrict ourselves to a minimal basis. In this case one has to compute one-loop corrections to Wilson coefficients coming from one-particlereducible (penguin-type) diagrams. We explicitly checked that the two procedures lead to the same results, as it should be. for the Yukawa contribution, and

From the EWSB scale to the nuclear scale
The Wilson coefficients given in Eq. (2.22) for the EFT below the EWSB scale evolve according to We now discuss how to obtain the 16 × 16 anomalous dimension matrix γ EMSMχ . The external leg corrections only come from the gauge sector. For strong interactions they are identical to the ones in the SM χ EFT, and for electromagnetic interactions they can be easily obtained from the analogous hypercharge diagrams. Their effect is again to cancel out against the associated vertex corrections. Also in this case there are two classes of vertex corrections. The first ones are due to the SM four-fermion interactions, in the way we show in Fig. 7. This diagram in the EMSM χ EFT is the analogous of the correction to c H discussed in the SM χ EFT, but the Z boson is integrated out in this phase of the theory. Despite the fact that these diagrams are suppressed by the Fermi constant, we keep them to be consistent with the analysis above the EWSB scale, since their contribution is proportional to G F m 2 ψ ∝ λ 2 ψ . The second effect is the radiative correction to the Wilson coefficient c F of the redundant operator in Eq. (2.23). The diagrams are analogous to the fermion loop in Fig. 6, but this time with an external photon. We evaluate them, and use the equations of motion to induce the evolution of the independent operators. 2 The anomalous dimension is The running driven by fermion masses reads 2 If a minimal operator basis is used, photon penguin diagrams have to be evaluated. We checked that we obtain the same results with both procedures. 3 We correct an overall sign typo in the 5 × 5 quark block given in Ref. [65].

Applications to Spin-Independent Searches
In Sec. 3 we presented the full 16×16 anomalous dimension matrices describing the one-loop RG evolution for dimension 6 operators above and below the EWSB scale. As promised in the introduction, these details can be skipped by a reader only interested in our final results. For the benefit of such a reader, we now briefly summarize the RG procedure. The boundary conditions for the RG system are the SM χ EFT Wilson coefficients at the cutoff Λ. In a generic UV complete model, they are obtained by integrating out the mediator fields. Then we RG evolve them down to the EWSB scale, which we take equal to the Z boson mass. It is convenient to introduce this dimensionless variable related to the renormalization scale µ, (4.1) In this notation the matching is performed at t = 0, whereas the Wilson coefficient c Λ are specified at the cutoff scale, t Λ = ln [Λ/m Z ]. The RG evolution in the SM χ EFT is obtained by solving the system of differential equations with the Wilson coefficients vector C SMχ defined in Eq. (2.9), and the explicitly expression for the anomalous dimension matrix γ SMχ given in Sec. 3.1. Once at t = 0, we perform the matching between the two theories as described in Sec. 2.3. The subsequent RG evolution for the Wilson coefficients C EMSMχ defined in Eq. (2.22) is described by The Λ-dependent evolution matrix U Λ is derived in App. C and for a user-friendly recipe we refer to App. D. The rest of this Section is devoted to applying Eq. (4.5) to limits from direct detection experiments. We focus on spin-independent searches, since they have much stronger bounds, and this has two implications. First, we need to consider effective operators with DM vector currents χγ µ χ. For pure elastic scattering this operator is non vanishing only for Dirac fermions, but our results are also valid for inelastic scattering of two splitted Majorana states [93]. Second, matrix elements of SM fermion currents have only contributions from valence quarks in the target nuclei, therefore the direct detection cross section reads (4.6) Here, c V V u and c In what follows, we consider specific choices of Wilson coefficients c Λ at the cutoff scale and we evolve them down to the nuclear scale as in Eq. (4.5). We compare the predicted rate as in Eq. (4.6) to the experimental limits, and extract bounds on the Wilson coefficients. The constraints we find are model independent, in the sense that every UV complete model generating that specific set c Λ when matched on the SM χ is subject to our constraints.

D5 and D7 operators
The connection between different DM negative searches is often expressed in terms of limits on the coefficients for the effective operators introduced in Ref. [24]. For a vector current of a fermion WIMP, the relevant operators are We now connect this description to the notation used in this paper, and explore the consequences of connecting EFT scales.
Keeping the complementarity among different searches in mind (e.g. between collider and direct searches as in Ref. [24]), we take the operators in Eqs. (4.7) and (4.8) as defined at the EFT cutoff Λ. In other words these are operators in the SM χ EFT. Considering flavor universal coupling to SM quarks, D5 and D7 are reproduced by this set of Wilson coefficients where Our results are shown in the four panels of Fig. 8. In the top-left panel we consider the case where only D5 is switched on, and plot current and projected experimental limits in the (m χ , Λ) plane for c 5 = 1. As is well known, huge scales for the mediator masses are necessary to be consistent with experimental exclusion bounds. We gain valuable information from this plot: given the extremely strong constraints on this operator, we are still likely to get useful limits on the scale Λ in other cases where the dominant contribution to direct detection rates is via D5 generated by SM loop effects. We deal with these cases in the next subsections, but we first complete the discussion of the (D5, D7) set.
The top-right panel of Fig. 8 shows the analogous case where only D7 is switched on at the scale Λ, with c 7 = 1. The limits on Λ are weaker than the case of D5, but still in the multi-TeV region, as pointed out in Ref. [65]. For both D5 and D7 we also plot the line that gives a correct thermal relic density, obtained using the annihilation cross section in Ref. [94]. The exclusion limits push the scale Λ well above the value for which χ is a cold thermal relic.
Both upper panels are for either c 5 or c 7 equal to 1. To go beyond this approximation one cannot just rescale the vertical axis by Λ → Λc −1/2 i . While it is true that the rate in Eq. (4.6) scales as c 2 i Λ −4 , there is an additional Λ dependence in the running Wilson coefficient c i (Λ). For this reasons, when presenting our results we specify both c i and Λ.
We go beyond the c i = 1 approximation in the left and right bottom panels of Fig. 8. In each panel we fix the DM mass value and explore the allowed region in the (c 5 , c 7 ) plane for three different values of Λ. If we ignore the running, the allowed regions would be the faded vertical bands, namely no restriction at all on c D7 . However, the RG evolution mixes different operators, and the actual experimental limits are the oblique bands.
To summarize, if SM loop effects are included, it is not consistent to assume that c 5 = 0 or c 7 = 0 at all energy scales. This is our main point in this subsection. Furthermore, any sensible UV completion is likely to generate both c 5 and c 7 , at least at one loop [95]. The constraints that such a model has to satisfy are the ones in the bottom panel of Fig. 8.

DM interacting with heavy quarks
Let us now focus on a different class of models, where the DM vector current only couples to heavy SM quarks The results for this case are shown in Fig. 9. In the top-left panel we set c Q = c U = c D = 1, or in other words we couple the DM to heavy quark vector currents at the EFT cutoff. The limits are shown again in the (m χ , Λ) plane. In the top-right panel we choose the Wilson coefficients such that the DM couples to heavy quark axial currents. We observe a distinctive feature of couplings to heavy quarks: the limits for interactions to axial currents are much stronger than the ones for vector currents, unlike for the case of D5 and D7 (see top panels in Fig. 8). When considering D5 in the previous case, we had couplings to vector currents of light quarks already at tree level, which explains why the limits were much stronger than the loop induced couplings when starting from D7. Here, couplings to light quarks are induced via loop effects in both cases, therefore we have to look at the RG equations and see how this is achieved. The largest coupling driving the mixing onto light quark currents is the top Yukawa, and its effect is encoded in the anomalous dimension matrix in Eq. (3.4). As it turns out, for coupling to vector currents of SM fermions there is no contribution to the running from Yukawa interactions, and the mixing is driven by the sub-leading hypercharge contribution (see Eq. (3.5)). On the contrary, for couplings to SM axial currents the top Yukawa contribution is maximal, and a substantial mixing onto the Higgs operator O V H is radiatively induced, which in turns gives interactions to light quarks once the Z boson is integrated out at the EWSB scale.
We also consider cases beyond the |c i | = 1 limit. In the bottom-left panel of Fig. 9 we fix the DM mass to m χ = 100 GeV, impose the isospin conserving condition c U = c D , and show the allowed region in the (c Q , c U ) plane for three different values of Λ. Unsurprisingly, the allowed region lies close to the diagonal line c Q = c D , since limits are weaker for coupling to vector currents. In the bottom-right panel we consider isospin violation by fixing c D = 0.5 and identifying the allowed region in the (c Q , c U ) plane. The bands are still close to the diagonal line, since the effect is driven by the top Yukawa which coupled the SM fields q  where for simplicity we consider flavor universal coupling to leptons. In such models there are many sources of couplings to light quarks currents. The Yukawa coupling of the τ induces a mixing onto the Higgs current, which in turns leads to a coupling to light quarks when the Z is integrated out. Hypercharge (electromagnetic) interactions above (below) the EWSB scale also induce mixing onto light quark currents.

Leptophilic Dark Matter
Our results for coupling to vector (axial) lepton currents to the vector DM current are shown in the top-left(right) panel of Fig. 10. As in the heavy quarks case, limits are stronger for axial couplings, although in this case the difference is milder since λ τ λ t . The allowed regions in the (c l , c e ) plane are shown in the bottom panels for two different values of the DM mass. Here, we do not have an order one coupling as λ t , therefore the bands are wider than the heavy quarks case. We still have the characteristic orientation along the diagonal, for the same reason as before.

Dimension 6 Higgs Portal
The last example we discuss is DM communicating with SM fields only via Higgs couplings both of these cases the thermal relic scenario is excluded. This conclusion persists in any Z portal model even if the Z couples to other SM fields besides the Higgs.

Discussion
In the last decade we witnessed an enormous progress by direct searches, which excluded a large parameter space fraction for WIMPs. In particular, models with tree-level vectorvector interactions with light quarks (i.e. D5 in Eq. (4.7)) are severely challenged, as shown in the left panels of our Fig. 8. This is not surprising, since spin-independent cross sections arising from tree-level Z boson exchange were ruled out years ago. This effort continues with the new LUX run and with future generation experiments like SCDMS and XENON1T, which will soon either discover DM or improve the exclusion limits by one or two orders of magnitude. At the same time, center of mass energies at the LHC will be increased up to a factor of two, and new regions of mass and couplings in DM simplified models will be probed. Last but not least, new generation ground-based and satellite experiments will look for products of DM reactions, constraining the Milky Way WIMPs annihilation rate. We will soon get access to many complementary results on WIMPs property.
When combining the negative outcome of different searches, it is crucial to properly handle the separation among the scales probed by the different experiments. This paper fills the gap between the mediator scales and the nuclear scales for singlet fermion DM models. As sketched in Fig. 1, we assumed that all non-SM particles (except possibly the DM) are heavier than the weak scale, and we considered the EFT obtained from integrating out these heavy degrees of freedom. It is certainly true that assuming contact interactions with SM fields is not always justified in a collider environment. However, given the typical recoil energy for DM-nucleus scattering, this is always the case for direct searches. Every DM model (possibly constrained by collider or indirect searches) can be matched onto an EFT with only DM interactions with the SM (the SM χ EFT in Fig. 1), and the subsequent model-independent RG evolution properly connects physics at the mediator scales to direct detection observables.
We used a systematic power counting in the suppression scale of the contact effective interactions, and we considered non-renormalizable operators up to dimension 6. Mixing among operators is a standard phenomenon in the dimension 6 sector, where SM interactions mix every operator into all the others. We identified the interactions responsible for these mixings, and computed the associated anomalous dimension matrices. Despite the size of the strong coupling constant α s , QCD interactions do not play any significant role since SM quark currents are not renormalized at one loop. QCD contributions in general cannot induce flavor-or parity-violating mixing, and unless there is a cancellation among subamplitudes [21,58,62,63] they can only modify the overall rate by a numerical factor of order one. The RGE in our case was instead driven by Yukawa and electroweak interactions. The former can mix SM fermion currents onto Higgs currents, which in turn generates DM interactions with SM fermion currents once the Z boson is integrated out. The latter mixes every SM current onto all the others. Given this rich mixing structure, for an arbitrary choice of DM interactions we are very likely to radiatively induce the strongly constrained D5 vector-vector operator defined in Eq. (4.7), and get significant bounds from experimental results.
In Sec. 4 we discussed applications of our results, focusing on examples where the operator mixing is the main source for direct detection rates. For dimensionless couplings of order one we found multi-TeV constraints on the mediator mass scale. These values are well above the mass scales ensuring a thermal relic, although one can always resort to a non-standard cosmological history [96,97]. For fixed DM and mediator masses, being consistent with the exclusion limits implies a peculiar alignment among the different dimensionless couplings. Besides the examples discussed here, further singlet fermion DM models generating dimension 6 operators at the mediator scale can be matched onto our framework, and exclusion bounds from direct searches can be derived.
Our work showed how to systematically account for loop effects in fermion singlet DM models. We envision several possible future directions. The next natural step for singlet fermion is to include dimension 7 operators, for which relevant mixing effects induced by electroweak interactions are already known [59,61,66]. Singlet scalar WIMPs belong to a separate chapter, since scalar and vector currents have different mass dimensions, and dimension 6 operators would involve also SM field strengths. Moreover, a singlet scalar WIMP can have renormalizable interactions with SM fields [98]. On a different route, one can apply the results of our paper to the simplified models recently proposed in Refs. [51,52]. As we eagerly await for new data, and hopefully for a DM discovery, it is important to stress that loop effects extend the notion of complementarity among different DM searches.

A Standard Model: Conventions and Results
In this Appendix we define our conventions and collect useful SM results. Our conventions for the space-time metric, the spinor and gauge fields are the same as in Ref. [99].
Above the EWSB scale the SM has an unbroken SU (3) c ×SU (2) L ×U (1) Y gauge symmetry, and the matter fields shown in Table 1. We divide the full Lagrangian into four contributions The Lagrangian for the gauge sector read where the indices A = 1, . . . , 8 and I = 1, 2, 3 run over the adjoint representations of SU (3) c and SU (2) L , respectively. We define a covariant derivative acting on matter fields as follows The Gell-Mann (λ A ) and Pauli (σ I ) matrices act on color and weak-isospin indices (if any), respectively, whereas the hypercharge Y is assigned as in Table 1. Then we have Finally, the Yukawa couplings between fermions and the Higgs doublet are For the purpose of this work we neglect flavor violating terms, therefore we take the Yukawa matrices to be diagonal. The SU (2) L × U (1) Y → U (1) em breaking Higgs VEV is chosen as which gives mass terms for the electroweak gauge bosons and the fermions.  Table 4. The Lagrangian has two contributions The Yang-Mills term is analogous to the one in Eq. (A.2), but this time with the photon field strength F µν in the electroweak sector. The fermion piece reads where we also include fermion masses. The fermion covariant derivative is defined in analogy to Eq. (A.3), but with only strong and electromagnetic interactions.

Useful Equations of Motion
The only equations of motion needed to get rid of redundant operators are the ones for abelian gauge bosons. In the unbroken phase, the hypercharge field strength satisfies where we define the current (A.11) Analogously, for the electromagnetic field strength in the broken phase where the electromagnetic current reads

B Loops and RG Equations
In this Appendix we give results for one-loop diagrams, from which we derive the RG equations. We regularize UV divergences by computing loop integrals in d = 4 − 2 dimensions, and subtract infinities using a mass-independent subtraction scheme. As done throughout the paper, we divide the discussion in two cases, correspondent to the two different EFTs.

RGE in the SM χ EFT
Here we discuss the renormalizion of the EFT described in Sec. 2.1, with a specific focus on the dimension 6 operators listed in Table 3. The wave-function renormalization factors, generated by the one-loop diagrams in Figs. 2 and 3 have the following 1/ -pole structure Here C 2 (N ) denotes the Casimir for the SU (N ) fundamental representation, and the hypercharge values are assigned as follows We have also define the Yukawa coupling vectors The Higgs quartic coupling does not induce any field renormalization. We now move to vertex corrections shown in Figs. 4 and 5. The sum of diagrams with gauge and Yukawa vertices gives the following one-loop corrections to the Wilson coefficients A further shift comes from the diagrams in Fig. 6. They generate a one-loop contribution to the redundant Wilson coefficient which becomes a shift of the independent operators upon using the equations of motion.
With all the loop amplitudes in hand, we can derive the RG equations. First of all we observe that the external legs and vertex corrections coming from gauge interactions have opposite divergent pieces, and therefore do not contribute to the running. The only leftover contribution from gauge interactions comes from the diagrams in Fig. 6 inducing the redundant operator. Then we make sure that every amplitude in the theory, obtained by the sum of tree-level and one-loop contributions, is finite by renormalizing the coefficients where we divide the matrix Z SMχ into identity and interacting parts, with the latter in turn given by Yukawa and hypercharge contributions. The δZ matrices explicitly read The RG equations are derived by imposing that the bare Wilson coefficients are renormalization scale independent, d d ln µ C SMχ ren = γ SMχ C SMχ ren , The last thing is to identify the µ dependence in the δZ i matrices. The key point to observe here is that we work in a mass independent scheme, namely the matrices δZ i do not have an explicit dependence on the mass scale but they only depend on µ through the Yukawa and gauge couplings. By defining at one loop the anomalous dimension matrix reads

RGE in the EMSM χ EFT
Below the EWSB scale wave-function renormalization is only due to gauge interactions. The results can be derived from Eqs.(B.1 -B.5) for fermions above the EWSB scale, with identical QCD contribution and electromagnetic factor obtained by replacing the hypercharges with the electric charges They all cancel again with the associated vertex corrections, therefore we do not need to further consider them. As discussed in the paper there are two classes of vertex corrections also in this case. The first corrections are due to the SM four-fermion interactions, as shown in Fig. 7. These diagrams do not vanish only if we start from an interaction between the DM bilinear and the axial current of SM fermions. The explicit shifts read where we find it convenient to isolate the common factor accounting for all the possible SM fermions in the loop of Fig. 7, The couplings g V f and g Af of SM fermions to the Z boson are given in Eq. (2.32).
The second contribution is the radiative correction to the Wilson coefficient c F of the redundant operator in Eq. (2.23), which results in (B.32) Using the equations of motion again leads to a shift for the independent operators.
Analogously to what we have done in App. B, the RG equations are derived renormalizing the Wilson coefficients, where we divide again the matrix Z EMSMχ into different contributions. They explicitly read 0 0 0 0 0 0 0 0 3m 2 u g Au g V e 3m 2 d g Ad g V e 3m 2 c g Au g V e 3m 2 s g Ad g V e 3m 2 b g Ad g V e m 2 e g Ae g V e m 2 µ g Ae g V e m 2 τ g Ae g V e 0 0 0 0 0 0 0 0 3m 2 u g Au g V e 3m 2 d g Ad g V e 3m 2 c g Au g V e 3m 2 s g Ad g V e 3m 2 b g Ad g V e m 2 e g Ae g V e m 2 µ g Ae g V e m 2 τ g Ae g V e 0 0 0 0 0 0 0 0 3m 2 u g Au g V e 3m 2 d g Ad g V e 3m 2 c g Au g V e 3m 2 s g Ad g V e 3m 2 b g Ad g V e m 2 e g Ae g V e m 2 µ g Ae g V e m 2 τ g Ae g V e 0 0 0 0 0 0 0 0 3m 2 u g Au g Au 3m 2 d g Ad g Au 3m 2 c g Au g Au 3m 2 s g Ad g Au 3m 2 b g Ad g Au m 2 e g Ae g Au m 2 µ g Ae g Au m 2 τ g Ae g Au 0 0 0 0 0 0 0 0 3m 2 u g Au g Ad 3m 2 d g Ad g Ad 3m 2 c g Au g Ad 3m 2 s g Ad g Ad 3m 2 b g Ad g Ad m 2 e g Ae g Ad m 2 µ g Ae g Ad m 2 τ g Ae g Ad 0 0 0 0 0 0 0 0 3m 2 u g Au g Au 3m 2 d g Ad g Au 3m 2 c g Au g Au 3m 2 s g Ad g Au 3m 2 b g Ad g Au m 2 e g Ae g Au m 2 µ g Ae g Au m 2 τ g Ae g Au 0 0 0 0 0 0 0 0 3m 2 u g Au g Ad 3m 2 d g Ad g Ad 3m 2 c g Au g Ad 3m 2 s g Ad g Ad 3m 2 b g Ad g Ad m 2 e g Ae g Ad m 2 µ g Ae g Ad m 2 τ g Ae g Ad 0 0 0 0 0 0 0 0 3m 2 u g Au g Ad 3m 2 d g Ad g Ad 3m 2 c g Au g Ad 3m 2 s g Ad g Ad 3m 2 b g Ad g Ad m 2 e g Ae g Ad m 2 µ g Ae g Ad m 2 τ g Ae g Ad 0 0 0 0 0 0 0 0 3m 2 u g Au g Ae 3m 2 d g Ad g Ae 3m 2 c g Au g Ae 3m 2 s g Ad g Ae 3m 2 b g Ad g Ae m 2 e g Ae g Ae m 2 µ g Ae g Ae m 2 τ g Ae g Ae 0 0 0 0 0 0 0 0 3m 2 u g Au g Ae 3m 2 d g Ad g Ae 3m 2 c g Au g Ae 3m 2 s g Ad g Ae 3m 2 b g Ad g Ae m 2 e g Ae g Ae m 2 µ g Ae g Ae m 2 τ g Ae g Ae 0 0 0 0 0 0 0 0 3m 2 u g Au g Ae 3m 2 d g Ad g Ae 3m 2 c g Au g Ae 3m 2 s g Ad g Ae 3m 2 b g Ad g Ae m 2 e g Ae g Ae m 2 µ g Ae g Ae m 2 τ g Ae g Ae and The RG equations read where the anomalous dimension can be computed analogously to what done in Eq. (B.23).

C Full Solution of the RG System
The evolution operator U Λ appearing in Eq. (4.5) connects the mediator with the nuclear scale. In this Appendix we derive an expression for this operator. As described extensively in the paper, this connection requires three different steps, encoded in three different contributions The matrix U SMχ (t Λ ) evolves the Wilson coefficients from Λ down to m Z , whereas the matrix U EMSMχ describes the evolution from m Z to µ ∼ 1 − 2 GeV. The intermediate matching at the EWSB scale is taken care of by the matrix U match . We start by setting up a general method to solve the RG system, which can be applied both above and below the EWSB scale. In both cases we always have to deal with a system The only scale dependence in γ(t) comes from SM running couplings, namely Here, the index j runs over all SM running couplings g j (t), whereas the matrices γ j are constant. The running of the SM couplings can be found by solving the RG equations given in Refs. [100][101][102]. In both EFTs t = 0 is one boundary of our integration range, thus we rewrite the SM running couplings as follows where t = 0 at the Z pole. The SM couplings at the Z pole can be found in Refs. [103,104]. The anomalous dimension matrix has then the form This setup is identical to a time-dependent perturbation theory problem in quantum mechanics with time replaced by the dimensionless variable t = ln(µ/m Z ). In particular, the matrix γ 1 (t) can be treated as a time-dependent perturbation to the constant matrix γ 0 , since their difference is at most a logarithmic running over two or three orders of magnitude. The procedure to set up a perturbative series is well known. We first define the "interaction picture We apply this result to derive the full evolution operator as defined in Eq. (C.1). The running from the mediator to the EWSB scale is obtained by applying the linear operator (C.11) The analogous evolution from the EWSB scale down to the nuclear scale is described by dt γ EMSMχ1I (t) . (C.12) Finally, the matching at the intermediate scale is achieved by using (C.13)

D A Simple Recipe for an Approximate Analytical Solution
In this Appendix we describe a simple analytical way to perform the RG evolution and compute direct detection cross sections based on our results. This is meant to be a useful recipe for practitioners. What is given here is the 0-th order of the Dyson series in Eq. (C.9), which neglects the running of the SM couplings. This approximation is quite satisfactory for our purposes since it induces errors at the level of 10% in the scale Λ if the DM couples to the top quark, and at the level of 1% without this coupling. Let us consider the situation where, starting from a UV complete model of singlet fermion WIMP, integrating out the mediators at the scale Λ generates a subset of the dimension 6 operators in Table 6. Then the spin-independent rate for the scattering off a target nucleus with Z protons and A neutrons is obtained through the following straightforward steps.
(i) Write down the Wilson coefficients at the scale Λ organized as in Eq. (2.9), and define where t Λ ≡ ln (Λ/m Z ); (ii) Evolve the vector of Wilson coefficients c Λ down to the nuclear scale using and finally the matrix (D.5) -32 -(iii) After obtaining c N , its first two components directly appear in the spin-independent WIMP-nucleus cross section which can be computed according to where m N is the mass of the target nucleus.