Global constraints on vector-like WIMP effective interactions

The search for Dark Matter (DM) in the form of thermal relics represents one of the most active lines of research in astro-particle and particle physics. Indeed, there is an overwhelming number of dedicated experimental searches for DM, most of them concentrating on the so-called Weakly Interacting Massive Particles (WIMPs) paradigm Bergstrom (2000). These are classified into three different categories: (1) indirect detection searches, where the DM would annihilate or decay into SM particles which can be detected, (2) direct detection searches, where the DM would scatter against the protons and neutrons in a detector, producing an observable recoil, and (3) collider searches, where the DM would be produced in high-energy collisions, thus leading to events with missing momentum.¹¹1A complementary search strategy is the search for DM self-interactions which would impact structure formation as well as stellar evolution in particular scenarios (see e.g. Tulin et al. (2013); Hochberg et al. (2014); Vincent et al. (2015a, b)). Since the focus of this work is to probe the interactions between DM and visible particles these constraints will not be considered here.

In light of all these searches it is essential to look at the global picture of where we stand concerning the WIMP scenario. This is certainly a rather complicated task, given the enormous zoo of models available in the literature with an associated plethora of free parameters. Even for a single model, brute-force scans of the corresponding parameter space represent a significant computational challenge. Consequently, the usual practice is to either rely on simplified realisations of those models or to fix some of the free parameters, thus scanning only hyper-surfaces in the model parameter space. ²²2The literature in this respect is quite vast. Some examples relevant to our work are Beltran et al. (2009, 2010); Fox et al. (2011, 2012); Mambrini and Zaldivar (2011); Cheung et al. (2012a), in the context of effective field theory (see below).

On the other hand, Monte Carlo methods constitute an efficient alternative to scan over a multi-dimensional parameter space. However, for WIMPs, most of the Monte Carlo scans are performed in the context of supersymmetric models, where the neutralino is usually selected as the preferred DM candidate (e.g., Refs. Cabrera et al. (2013); Roszkowski et al. (2012); Strege et al. (2014), however see Cheung et al. (2012b) as an example of a non-SUSY search).

The main goal of this work is to explore the present status of our knowledge of DM couplings to the different SM constituents in as much generality as possible. In order to gain in model independence, it is thus interesting to be more agnostic about the DM interactions with the SM constituents and parametrize them through an effective field theory (EFT) approach. In particular, this approach can reveal how constrained the DM interactions with the different SM fields are in the light of currently available data, when no further assumptions are made. Since allowing all type of Lorentz structures for these effective operators would provide too much freedom with hundreds of operators that cannot be bounded, we rather choose to exemplify the different interactions between a Dirac fermion DM and the SM by operators of the form

with independent coefficients $c_{i,P}$ for all SM particle species $f_{i}$ with chiralities $P\equiv P_{L},P_{R}$. This type of interaction is well motivated by scenarios where the DM communicates with the SM via a vector portal. Moreover, scalar or tensor interactions typically require Higgs insertions, leading to higher dimensional and hence more suppressed operators. Additional dimension six operators with a vector-like structure involving the Higgs, such as $i\bar{\chi}\gamma^{\mu}\chi H^{\dagger}\overleftrightarrow{D}_{\mu}H$, could be included as well. For DM heavier than the Higgs, this operator would provide an additional annihilation channel relevant for indirect searches and relic abundance constraints, via $\bar{\chi}\chi\rightarrow hh$. However, as we will see, above 100 GeV the global upper bounds for any individual coupling correspond to the ones obtained from relic abundance alone (with no constraint for the top quark coupling, since that channel is not kinematically allowed). Therefore, we expect a similar behaviour for the Higgs, that is, no constraints for DM masses below the Higgs mass and recovering the relic abundance constraint for higher DM masses. Even with this restriction, if independent couplings for all SM constituents are allowed, 15 independent parameters remain to be constrained. Thus, given the large dimensionality of the parameter space, we make use of a nested sampling Monte Carlo algorithm to scan it efficiently.

When constraining the DM EFT parameter space, we consider bounds from all types of experiments where a WIMP signal is being actively sought for, i.e., direct detection (namely, LUX Akerib et al. (2014) and EDELWEISS Ahmed et al. (2011)), indirect detection (AMS Aguilar et al. (2013) positron fraction data and Fermi-LAT data for dwarf galaxies Ackermann et al. (2015)), cosmic microwave background (CMB) and relic density constraints from Planck Ade et al. (2015), and monojet and monophoton searches in colliders (from LHC Aad et al. (2015) and LEP data Fox et al. (2011)).³³3Contrary to the cases of direct and indirect detection strategies, for collider searches the use of the EFT framework -mainly at the LHC- may not be optimal. We take this issue into account when recasting the limits coming from LHC. See Sec.III. Also, note that other searches different from monojets can be interesting in the context of specific models (see e.g. Ref. Khachatryan et al. (2015).

Additional information on some annihilation channels such as $\tau\tau$, $bb$, or $\nu\nu$ come from the non-observation of indirect neutrino signals from the Sun in SuperKamiokande Choi et al. (2015) and IceCube Aartsen et al. (2016). However, the bound is only relevant if apart from a large DM coupling to the final state particle, a sizeable coupling to first family quarks is also present. Since we are assuming all couplings to be independent, this bound will always be satisfied in the fit by choosing a set of parameters with one of the relevant couplings being very small. Thus, the inclusion of these experiments would not change our results.

We derive bounds on the coefficients accompanying each effective operator as a function of the DM mass, as well as bounds on the DM mass itself. These constitute the most general constraints which can be derived assuming that DM interacts with the SM as in Eq. (1).

The outline of the work is as follows. In Sec. II we describe the set of effective operators that will be jointly analyzed, introducing the parameters to be constrained. Sec. III lists the set of experimental constraints considered, and how these have been implemented. Some details regarding the numerical tools employed in the fit are explained in Sec. IV. Finally, Sec. V summarizes our results and we present our conclusions in Sec. VI.

In order to explore how much freedom the present global data on DM allows for its couplings to the SM, we will define a series of working models characterized by a set of independent effective operators, described by the following effective Lagrangian

where the effective current $j_{\mu}^{\rm eff}$ is given by

the index $i$ denotes the quark generations such that $u_{1}=u$, $u_{2}=c$, $u_{3}=t$, $d_{1}=d$, $d_{2}=s$, and $d_{3}=b$, while the sum over $f$ runs over all right-handed (RH) SM fermions. The coefficients $c_{X}$ are the couplings of the operators in the effective Lagrangian. As the effective operators are of dimension six, these coefficients will have a mass dimension of minus two. In expr.(2), $\chi$ represents a Dirac fermion DM.

This set of operators provide enough freedom so as to parameterize DM interactions with different strength to the various SM particles while, at the same time, keeping the set of operators at a viable level. Indeed, allowing extra operators would not only imply a more challenging numerical analysis, but would also be rather uninformative since mostly any value could be fitted and particular UV completions tend to have a much more limited set of free parameters. Furthermore, for simplicity and due to the generally stronger constraints, we will not consider flavour-changing operators between the SM fermion bilinears Kamenik and Smith (2012). Despite these restrictions, 15 different operators fall into this category, parametrizing the DM interactions with the 3 quark and lepton doublets as well as the 3 singlets of up-type and down-type quarks and charged leptons. Starting from this general setup, we will also define more restrictive working models that can exemplify other interesting scenarios such as leptophilic, leptophobic or flavour-blind setups, characterized by different correlations among the couplings. The models considered in this work are the following:

1. 1.

   General model: The first model makes no additional assumptions regarding the coefficients $c_{X}$, which are all allowed and free to vary independently. This represents the least restrictive choice possible in our given EFT framework, and includes a total of 15 coefficients in addition to the mass of the DM.

2. 2.

   Flavour-blind: In this model, all operators involving fermions with the same gauge quantum numbers, e.g., all left-handed (LH) quarks, are assumed to also have the same couplings to DM and therefore have the same coefficients in the effective lagrangian. We are then left with five different coupling parameters which are related to those of the general model as

   	$$\displaystyle c_{Le}=c_{L\mu}=c_{L\tau}\equiv c_{L},\quad c_{Q1}=c_{Q2}=c_{Q3}\equiv c_{Q},$$
   	$$\displaystyle c_{uR}=c_{cR}=c_{tR},\quad c_{dR}=c_{sR}=c_{bR},\quad c_{eR}=c_{\mu R}=c_{\tau R}.$$		(4)

   For the latter three groups of coefficients, we will generally use the notation for the coefficient of the first generation to refer to it within this model.

3. 3.

   Family-oriented: In this model, we make the assumption that the effective couplings with DM are equal for all particles belonging to the same generation. This can be considered a quite crude proxy to flavour theories were the successive SM families are characterized by hierarchical couplings or charges in order to explain the observed mass hierarchy between them, following the philosophy of the Froggatt-Nielsen mechanism⁴⁴4See Calibbi et al. (2015) for a realisation in the context of dark matter. Froggatt and Nielsen (1979); Leurer et al. (1993, 1994); Ibanez and Ross (1994); Binetruy and Ramond (1995). This leaves only three independent operator coefficients,

   	$$\displaystyle c_{Le}=c_{eR}=c_{Q1}=c_{uR}=c_{dR}\equiv c_{1},\quad c_{L\mu}=c_{\mu R}=c_{Q2}=c_{cR}=c_{sR}\equiv c_{2},$$
   	$$\displaystyle c_{L\tau}=c_{\tau R}=c_{Q3}=c_{tR}=c_{bR}\equiv c_{3}.$$		(5)

   This is the model with the fewest number of free parameters which we will consider and, as such, the correlations among the couplings will allow to obtain strong constraints, particularly for the coupling $c_{1}$.

4. 4.

   Leptophobic: This model assumes that DM does not have significant interactions with any of the leptons. As a result, all of the coefficients associated to operators involving leptons are set to zero, i.e.,

   $$c_{L\ell}=c_{\ell R}=0,$$

   (6)

   with $\ell=e,\mu,\tau$. At the same time, no restrictions are imposed on the operators involving quarks and the model therefore contains nine free parameters.

5. 5.

   Leptophilic: In direct contrast to the leptophobic model, this model instead considers the situation where the only relevant DM interactions with the visible sector are those involving leptons. In this model, we therefore set all of the coefficients for operators involving quark fields to zero, i.e.,

   $$c_{Qi}=c_{u_{i}R}=c_{d_{i}R}=0\,,$$

   (7)

   where $i$ indicates the generation. In this situation, we are left with six independent coupling parameters for the DM interactions with leptons.

We wish to stress that these models are only meant to be phenomenological tools to assess our present global knowledge of how well constrained the DM couplings to the SM fermions are when not assumed to be universal or related through any particular UV completion. We will thus allow all couplings to vary freely when fitting the present data without questioning the apparent naturalness (or lack thereof) of the preferred regions. In particular, as we will see in Sec. V, most models prefer $2c_{Q1}+c_{uR}+c_{dR}\simeq 0$ so as to avoid the very stringent limits coming from direct detection searches. Thus, even if these cancellations may seem unnatural, we will not avoid them by adding artificial “naturalness priors” to guide the Monte Carlo in any way, in the spirit of allowing the Monte Carlo to choose the points in parameter space which provide the best fit to the data. Symmetry arguments could perhaps be invoked in a particular UV completion when trying to reproduce the best fit found in the effective description in a natural way.

Moreover, these working models are not intended to be self-consistent low energy descriptions of any particular UV completion. In fact, due to the limited amount of operators considered, we are not allowed to take our effective theory description beyond tree level processes. Indeed, radiative corrections could induce other operators with interesting phenomenology.⁵⁵5For recent works in an EFT framework including one loop processes, see e.g. Refs. Crivellin et al. (2014); D’Eramo and Procura (2015) For example, in the leptophilic model the lepton legs could be closed in a loop and, through the emission of a virtual photon, induce a coupling to first generation quarks. Therefore, the strong bounds from direct searches would also apply to these models and put a constraint on the lepton couplings (see Ref. Kopp et al. (2009)). In a similar fashion, other signals at the LHC such as dijet or dileptons could be generated from loop-induced 4-SM-fermion operators. Furthermore, when moving from the flavour basis to the mass basis, the DM couplings to the SM doublets will induce flavour-changing operators that, when the DM legs are closed in a loop, could contribute to the oscillations of neutral kaons or other FCNC processes, for which stringent experimental constraints would apply. In full consistency, these operators should have been included from the start, since they are compatible with the particle content and symmetries of the theory and they are required to renormalize the divergences stemming from the loop contributions. However, when doing so, new unbounded coefficients are introduced and, in order to avoid these stringent limits, the fit will prefer the points in the parameter space where these new coefficients exactly cancel the loop-induced contributions. Thus we will simply not consider these loop-induced constraints in our list of observables, since in our approach they would not imply additional constraints on our parameter space unless particular relations between the couplings are invoked.

In our numerical analysis, the DM contributions to the different DM observables are computed with micrOmegas Belanger et al. (2014). These are then compared to the current experimental bounds. A chi-square function is then computed for each observable independently. The total $\chi^{2}$ is obtained as the sum of all separate contributions

where $i$ is an index which runs over the observables and $c$ is a vector containing the coupling coefficients of the model being tested.

The contribution to the chi-square from each experiment is computed assuming that the limit on the physical quantity being bounded (i.e., relic abundance, the thermally averaged cross section, etc) by the experimental collaborations behaves like a gaussian. The chi-square can be expressed as a function of the number of events as:

where $N^{th}$ and $N_{bg}$ correspond to the number of signal and background events, respectively, and $\Delta N_{sys}$ refers to the systematic uncertainty on the number of events. Notice that, given the unfortunate lack of a positive DM signal in any of the observables described below with the exception of the relic abundance constraint, in this expression it has been assumed that the observed number of events corresponds with the expected background. In some cases, small upward or downward fluctuations over that background will be present and thus small deviations with respect to the scaling derived below can take place. However, notice that we will normalize this scaling with the data presented by each experimental collaboration itself. Thus, our $\chi^{2}$ function will always reproduce correctly the result obtained by the collaboration at the confidence level at which it is reported (typically 90 or 95%) and the deviations should be small for the purpose of our exploration as long as CL close to these limits are investigated.

Experimental collaborations present their results as bounds on a given quantity at a certain confidence level (CL). This can be translated into a certain number of events $N^{exp}$ allowed at that CL. For a bound at 68% CL, for instance, we get⁶⁶6When the experimental result is quoted at a different CL, it can be simply rescaled to 68% CL within a gaussian approximation. :

By taking the ratio between Eqs. (9) and (10), we get the following expression for the chi-square:

The expression above can be further simplified in the two following limits:

• •

  Statistically-dominated: in this case, the number of background events and the systematic errors can be neglected, and we can therefore write the chi-square contribution as:

  $$\chi^{2}\simeq\frac{N^{th}}{N^{exp}}\,.$$

  (12)

• •

  Background/Systematics-dominated: in this case, the signal events can be neglected with respect to the background and/or the systematic error, and the chi-square simplifies to:

  $$\chi^{2}\simeq\left(\frac{N^{th}}{N^{exp}}\right)^{2}\,.$$

  (13)

By making use of Eqs. (12) or (13), all the finer details of the detector response cancel in the ratio. We will therefore use the experimental bounds on the different observables reported by the collaborations instead of the number of events (which are not always publicly available).

All the experiments considered in this work fall in the background/systematics-dominated case, and therefore we will apply Eq. (13). The only exception to this will be the case of dwarf galaxies studied by the Fermi-LAT collaboration Ackermann et al. (2015). These are dark-matter–dominated objects and constitute extremely clean probes to test for dark matter interactions with SM particles. We have explicitly checked the behaviour of the chi-square for a subset of the observed dwarf galaxies, using publicly available data from the collaboration from Ref. Ackermann et al. (2015), and we find that the behaviour is in between the two extreme cases identified above. We therefore opt for the most conservative approach in this case, which corresponds to Eq. (12), when using the limit provided by the collaboration (which has been obtained using the full data sample).

In the following, we describe the full set of observables sensitive to interactions between the DM and the SM fermions which have been included in our simulations.

All the phenomenological models described in Sec. II have been implemented in CalcHEP Belyaev et al. (2013) using FeynRules Alloul et al. (2014). The result is then fed into micrOmegas Belanger et al. (2014) for the computation of the DM observables.

With the chi-square function implemented as described above, the parameter spaces of the different models need to be efficiently explored in order to derive constraints. While a simple grid scan may be viable in the case of the family-oriented model (with only three parameters), it quickly becomes prohibitively expensive for models with a larger number of free parameters. In particular, the general model with 15 free parameters would be unfeasible to scan in this fashion. Moreover, we expect the parameter space to be affected by several degeneracies. These further decrease the effectiveness of a grid scan and require the grid spacing to be very small in order to find the global minimum. For these reasons, we have opted to perform our scans using the nested sampling Skilling (2004) Monte Carlo algorithm implemented in the MultiNest software Feroz et al. (2009). This method is particularly designed for handling parameter space degeneracies as well as for preferentially scanning the regions of parameter space where the likelihood function $\mathcal{L}=\exp(-\chi^{2}/2)$ is large. Although the MultiNest software was initially designed for computing Bayesian evidence, it produces a sample of the parameter space with the corresponding likelihood values as a by-product. In this paper, we take a purely frequentist approach and only consider the $\chi^{2}$ values obtained from the likelihood. Thus, MultiNest is used only for its capability of efficiently sampling the regions with relatively low $\chi^{2}$ values.

In our MultiNest simulations, the number of live points was set to 750, the tolerance to 0.2 and the efficiency to 0.9. Constant efficiency mode was not used. The range of parameters scanned was adapted to the values of the coefficients required to reproduce the correct relic abundance, for each of the dark matter masses considered. The scan was done in linear scale over both positive and negative values for the coefficients accompanying each of the operators. The number of distinct samples obtained ranges between $10^{5}$ and $5\times 10^{5}$, depending on the particular model and DM mass considered.

This section describes the results obtained from a global fit using the models described in Sec. II. All constraints considered in Sec. III have been included in our simulations.

Note that the dimensionful couplings $c_{j}$ will depend on the effective theory mass scale $\Lambda$ as $\Lambda^{-2}$. In this section, we will use units of TeV^-2 for the couplings $c_{j}$ and the rough limit $\Lambda^{2}\lesssim 1/c_{j}$ may be used to assess the range of validity of the EFT assumption.

Using the $\chi^{2}$ function built as described in Sec. III, we have derived two main types of constraints on the parameters of the model:

• i)

  For a fixed DM mass $m_{\chi}$, we scan the parameter space for all couplings in each model. This results in a sample of points in the parameter space, which is generally concentrated in the regions corresponding to the lowest values of the $\chi^{2}$ function. We use these to perform parameter estimation and determine the allowed regions for the different couplings within a given model, for given values of $m_{\chi}$.

• ii)

  The calculation described in (i) can be repeated for several values of the DM mass. By minimizing the global $\chi^{2}$ over all couplings $c$ in the model, the resulting $\Delta\chi^{2}$,

  $$\Delta\chi^{2}(m_{\chi})=\min_{c}\chi^{2}(c,m_{\chi})-\min_{m_{\chi},c}\chi^{2}(c,m_{\chi}),$$

  (31)

  allows parameter estimation of the DM mass $m_{\chi}$. By deriving confidence intervals around the global minimum of a given model, this will indicate whether there are any values of the DM mass which are disfavoured by current data. Since the simulations are computationally rather expensive, they have been performed for a few values of the DM mass, namely $m_{\chi}=10$, 50, 100, 200, 500 and 1000 GeV. Our minimization in the $m_{\chi}$ variable is therefore performed in a discrete fashion, using only this set of values. However, from $m_{\chi}\geq 200$ GeV on, we find essentially no change in the $\chi^{2}$ value.

Our main results are summarized in Figs. 1 and 2. Figure 1 shows the value of the $\Delta\chi^{2}$ (see Eq. (31)) obtained for the different models under consideration, as indicated in the legend, see item (ii) above. The dashed horizontal line shows the value of the $\Delta\chi^{2}$ corresponding to $3\sigma$ limit for 1 d.o.f.. As can be seen from this figure, DM masses around 40 GeV and below are disfavored at $3\sigma$ with respect to higher masses for most models under consideration. The only exceptions are the general and leptophilic models, for which only masses below 20 GeV are disfavored.

The fact that very light DM masses are generally disfavored can be understood from the complementarity between different data sets. As can be seen from Eq. (16), the annihilation cross section is proportional to $m_{\chi}^{2}$, while the energy density scales as $m_{\chi}$. Thus, models with very light DM masses will require large couplings to the SM in order to satisfy relic density constraints. However, since there is no positive signal from DM in collider, direct or indirect detection data, a significant tension between the different data sets occurs, which eventually increases the minimum value of the $\chi^{2}$. The tension is stronger in models where DM either does not couple to leptons (leptophobic), or in models where the lepton couplings are related to others (i.e., flavour-blind or family-oriented). For the general and leptophilic models, the tension is relaxed since the bounds which mainly constrain the couplings to leptons (LEP, CMB and AMS data) are not as strong as those constraining the quark couplings (direct detection, LHC and FermiLAT bounds). Finally, one can observe a slight preference in the fit for DM masses around 100 GeV for all models except for the leptophobic and family-oriented. This can be understood from the mild preference of AMS data for a DM signal around this mass, since smaller masses are too strongly constrained, while heavier masses do not produce an appreciable signal. It should be stressed out, however, that this preference is far from being statistically significant. This feature is absent in the leptophobic and family-oriented models since the required fermion couplings for this signal to take place are either absent or constrained to be very small by other experimental bounds.

The second main result of this paper is shown in Fig. 2, where the allowed regions at 1, 2 and 3$\sigma$ (for 1 d.o.f.) are shown for all couplings associated to the effective operators and for the five models under consideration, see Sec. II. For each coupling, our results are shown for three different values of the DM mass, 10 GeV (upper/blue bands), 100 GeV (middle/green bands) and 500 GeV (lower/red bands). Each panel also shows the constraint imposed by relic density on the weighted sum of the squares of all couplings, i.e., $\Sigma_{C}$ as defined in Eq. (16).

The tension between different data sets can be understood from the results found in Fig. 2. As already explained, the constraint from relic density tends to bring the couplings to larger values as the DM mass is decreased. This is observed when comparing the upper/blue and lower/red bands, for all models under consideration, and in particular for the allowed values of $\Sigma_{C}$. One can clearly see how the tension between different data sets can lead to having only one coupling as the major contributor to the relic density for a given model. For instance, for the general and leptophilic models with $m_{\chi}=50$ GeV, the relic density constraint is satisfied with a sizable coupling to the second family lepton doublet, and the preferred region for this coupling at $1\sigma$ does not include zero. This can be understood as follows. First, the bounds on couplings to quarks of all generations are very strongly constrained by direct detection and Fermi-LAT experiments. Moreover, on the leptonic side, the bounds on RH leptons from both AMS (electrons) and Fermi-LAT (taus) are very strong. In addition, LH couplings would also imply DM annihilations through muon neutrinos, so that the indirect searches are relaxed while the relic abundance constraint is more easily fulfilled. As such, the weakest global bound appears for the second lepton generation, which can accommodate the required relic abundance without being in tension with AMS or Fermi-LAT.

Nevertheless, a significant tension persists, as seen in Fig. 1. This is particularly the case for light DM masses, see, e.g., the bands corresponding to 50 GeV in Fig. 2. However, for 100 GeV mass (and heavier) the situation is different. For such large values, the indirect detection constraints are already weaker than the relic abundance one. Thus, since a coupling to LH fermions implies two annihilation channels, both contributing to decrease the DM density, we see that the LH couplings are more constrained than the corresponding RH couplings. Similarly, due to color multiplicity factors, couplings to quarks are more strongly constrained than those to leptons (in particular the LH ones).

Another interesting example is found in the panel for the flavour-blind model. Since the coupling to the muon is equal to that of the electron and tau (which are very strongly constrained) the relic abundance cannot be obtained out of annihilation to the second lepton doublet, as for the previous cases. For $m_{\chi}=50$ GeV, we see from the figure that the DM would prefer to annihilate to up-type RH quarks instead (the 1$\sigma$ region for this coupling is around the value of 1). However, this results in a strong tension with the Fermi-LAT data, and this mass is therefore disfavored at $\sim 3\sigma$ with respect to the best fit at higher mass (see Fig.1).

Finally, as was already mentioned in Sec. III, we have found two interesting degeneracies among different parameters in the models considered in this work. The first one is due to the constraints coming from direct detection experiments, which allows for a degeneracy between the couplings to quarks. This is shown in the left panel in Fig. 3, where the allowed confidence regions for the $c_{Q1}$, $c_{uR}$ and $c_{dR}$ couplings are depicted, showing the degeneracy explicitly. For $m_{\chi}=100$ GeV, the most stringent bound on the effective operators involving the first generation of quarks comes from direct detection experiments. Nevertheless, as can be seen from this figure, the allowed regions extend to rather large values of the couplings, as long as the relation $2.12c_{Q1}\simeq-c_{uR}-1.12c_{dR}$ is satisfied, given the values of $Z$ and $A$ of the Xenon material used by LUX (see Sec. III for details). While the inclusion of EDELWEISS-II, with a different target material, could in principle lift this degeneracy, the difference in the ratio of protons and neutrons is not large enough for this task. The final limits found in Fig. 3 for the degeneracy line are rather stemming from the relic abundance constraint.

The second degeneracy stems from the strong constraints from relic abundance on $\Sigma_{C}$, which generally imply that all couplings must lie on the surface of a hyperellipsoid. In particular, for the family oriented scenario, this reduces to the ellipse displayed in the right panel in Fig. 3, where the allowed regions at 1, 2 and 3$\sigma$ are shown in the $c_{2}-c_{3}$ plane for 2 d.o.f.. Given that the strongest constraints from colliders and direct detection experiments apply to the first generation, and in this model all particles in the first generation have the same coupling, the constraints on $c_{1}$ are very strong since the direct detection degeneracy relation cannot be fulfilled. This is shown in Fig. 2, where it can be seen that the coupling to the first family is restricted to $c_{1}<0.001/{\rm TeV}^{2}$ at $2\sigma$ for this model. This model only has two additional couplings, $c_{2}$ and $c_{3}$. Therefore, in order to satisfy relic density constraints, either the coupling to the second or the third family (or both) have to be different from zero. As a consequence, the allowed region is shaped as an ellipse in the $c_{2}-c_{3}$ plane.

A final remark regarding the validity of the EFT is in order. As can be seen from eq. (16) and from Fig. 2, the weighted sum of the coefficients required to obtain the correct thermal relic abundance goes like $\Sigma_{c}\propto 1/m_{\chi}$. On the other hand, for DM annihilation, the EFT validity condition reads $\Lambda\gtrsim 2m_{\chi}$. Since $\Sigma_{c}$ is a combination of couplings, we can generically write it as $\Sigma_{c}\sim 1/\Lambda^{2}$. Therefore, by combining the validity condition with the relic abundance requirement, a maximum DM mass is obtained for which the EFT ceases to be valid about $\mathcal{O}(2-5$ TeV) . Notice however that we do not expect the data to strongly constrain thermal DM in this regime.

In this work we have explored how strongly the combination of present data from very different experimental probes can constrain the different couplings between Dark Matter (DM) and the Standard Model (SM) particle content. The focus of the work was a bottom-up approach to understand what data itself is able to tell us regarding DM interactions, minimizing when possible any theoretical input or bias. To this end, we make use of an Effective Field Theory (EFT) approach to attain the desired model-independence at the expense of assuming that DM–SM interactions are mediated by heavy particles that can be integrated out of the theory. Furthermore, we have allowed independent couplings of DM to all SM fields so as to minimize any theoretical bias, but we have limited these interactions to flavour-conserving ones in the SM sector, since generally stronger constraints apply to them. Finally, if all possible Lorentz structures were simultaneously allowed, the parameter space would become too large to explore and constrain with present data. Thus, we restrict our analysis to Dirac DM fermions which interact with the SM constituents through independent operators of the form

with independent coefficients $c_{i,P}$ for all SM fermions $f_{i}$ and chiralities $P\equiv P_{L},P_{R}$.

This working model, which is not intended to reproduce any particular ultraviolet completion, provides the necessary parametrization to assess through a joint analysis the constraints on the individual interaction strengths from present data. We constrain 15 independent coefficients through a combination of 8 experimental probes comprising relic abundance, direct and indirect DM searches as well as collider limits. In addition to this general setup, we also explore how these constraints are affected when correlations between the different coefficients are allowed, or when some of the operators are forbidden. In particular, we studied the completely leptophilic, leptophobic and flavour-blind cases, as well a family-oriented scenario in which all members of the same generation share the same coupling to the DM particle.

From our global fits we find that for DM masses $m_{\chi}<20$ GeV the tension with respect to the best fit between the couplings necessary to reproduce the observed DM relic abundance and the upper bounds from null results exceeds the $3\sigma$ level, for the general and leptophilic scenarios considered, while the same happens for $m_{\chi}\lesssim 40$ GeV for the leptophobic, flavour blind and family oriented cases. For these three cases, the $\chi^{2}$ rises very steeply as the DM mass decreases such that for $m_{\chi}<\mathcal{O}(10-20)$ GeV this tension reaches the $5\sigma$ level, within our EFT framework.

Furthermore we find that, due to the slightly weaker present constraints, couplings to the second-generation lepton doublet are preferred, whenever they are available within the model under study. If this coupling is not available or is instead further constrained by an assumed correlation, as is the case of the flavour-blind scenario, a more sizable coupling to the right-handed, up-type quark singlet is instead preferred for similar reasons. Thus, it would be interesting to improve our present constraints on these couplings. Finally, we also find that, since all couplings are assumed to be independent and free to vary in the fit, the very stringent constraints stemming from direct search experiments such as LUX imply that either the first generation quark couplings to DM are extremely small or that they are related to each other such that the different contributions to these processes cancel against each other, i.e., $2c_{Q1}+c_{uR}+c_{dR}\simeq 0$.

DM searches worldwide are now probing and constraining essentially all possible interaction channels between DM and the known matter constituents through extremely different and complementary and search techniques. In this context, and given our lack of a unique theoretical DM paradigm, it is important to test different DM models against present data with bottom-up approaches where theoretical biases are not imposed. Unfortunately, true and complete model independence cannot be achieved. Indeed, while EFT offers an ideal frame for this sort of studies, its adoption already enforces some assumptions, such as the decoupling nature of the mediating particle or the uniqueness of the DM candidate. Moreover, true model independence would imply the inclusion of hundreds of independent operators with different flavour and chiral structures, rendering the analysis too general to actually provide any useful information. In this work we have reduced the number of theoretical assumptions taking a first step which implies some unavoidable restrictions to the number and types of effective operators included. It would be very interesting to supplement our results with additional complementary analyses, including different Lorentz structures, different DM fields (Majorana fermions, or scalars), or with extensions beyond the EFT approach by allowing one (or several) generic light mediator(s).