QCD Effects on Direct Detection of Wino Dark Matter

We complete the calculation of the wino-nucleon scattering cross section up to the next-to-leading order in $\alpha_s$. We assume that the other sparticles are decoupled and wino interacts with the Standard Model particles via the weak interaction. As a result, the uncertainties coming from the perturbative QCD are significantly reduced to be smaller than those from the nucleon matrix elements. The resultant scattering cross section is found to be larger than the leading-order one by about 70%, which is well above the neutrino background. In the limit of large wino mass the spin-independent scattering cross section with proton turns out $\sigma_{\text{SI}}^p = 2.3~{}^{+0.2}_{-0.3}~{}^{+0.5}_{-0.4}\times 10^{-47}~ \text{cm}^2$ (errors come from perturbative calculation and input parameters, respectively). The computation for a generic SU(2)$_L$ multiplet dark matter is also presented.


Introduction
Weakly-interacting massive particles (WIMPs) are promising candidates for dark matter (DM) in the Universe. Many theoretical models beyond the Standard Model predict WIMPs and it is known that the thermal WIMP scenario can explain the present energy density of dark matter in those models. The early stage of the experiments at the Large Hadron Collider, however, has found no evidence for new physics near the electroweak scale so far. In particular, the experiments give severe bounds on new colored particles, such as gluino and squarks in the supersymmetric (SUSY) models [1]. This situation may imply that most of the new particles in the high energy theory have masses much larger than the electroweak scale and only a WIMP, probably accompanied with some non-colored particles, is accessible in the TeV-scale experiments.
The current experimental consequences would fit in with a simple SUSY breaking scenario. If the SUSY breaking is induced by non-singlet chiral superfields (as is often the case with the dynamical SUSY breaking [2]), gaugino masses are induced by the anomaly mediation mechanism [3,4] and thus suppressed by a loop factor compared with the gravitino mass. A generic Kähler potential gives masses of the order of the gravitino mass to scalar particles and higgsino. In this framework, the neutral wino is found to be the lightest SUSY particle and thus becomes a candidate for dark matter in the Universe. Actually, its thermal relic abundance explains the observed energy density of DM if the wino DM has a mass of 2.7 − 3.1 TeV [5]. For relatively light wino DM, on the other hand, the non-thermal production via the late time decay of gravitino could be invoked to provide the correct abundance of DM [6,7]. As this scenario [4,8] requires the SUSY breaking scale to be much higher than the electroweak scale, a relatively heavy mass for the Higgs boson is predicted [9], which is consistent with the observed value m h 125 GeV [10]. Such a high SUSY-breaking scale is phenomenologically desirable since it relaxes the SUSY flavor and CP problems [11], the dimension-five proton decay problem [12], and some cosmological problems [13]. Gauge coupling unification is found to be still preserved with good accuracy [14]. For these reasons, the wino DM scenario attracts a lot of attention, and its phenomenology has been studied widely.
A lot of efforts have been dedicated to searching for the wino DM. A robust constraint is provided by the Large Hadron Collider experiment; charged winos with a mass of 270 GeV or less have been excluded at 95% C.L. [15]. For prospects of the wino search in future collider experiments, see Ref. [16]. On the other hand, signal of the wino DM may be detected in cosmic ray observations. Since the wino DM has large annihilation cross section [17,18], cosmic rays from annihilating winos are promising tools to detect the wino DM indirectly. The mass of wino DM M is constrained as 320 GeV ≤ M ≤ 2.25 TeV and 2.43 TeV ≤ M ≤ 2.9 TeV at 95% C.L. [19] by using gamma ray data from dwarf spheroidal galaxies provided by Fermi-LAT collaboration [20]. Gamma rays from the Galactic center provided by the H.E.S.S. [21] may give a strong limit on the wino DM, though the consequences are quite dependent on the DM density profile used in the analysis [22]. Developments in both theory [23] and observation enable us to probe a wide range of mass region of the wino DM in future indirect detection experiments.
Direct detection of dark matter is another important experiment to study the nature of dark matter. Currently the most stringent limits are provided by the LUX experiment [24]; it sets an upper limit on the spin-independent (SI) WIMP-nucleon elastic scattering cross section as σ SI < 7.6 × 10 −46 cm 2 at a WIMP mass of 33 GeV. Moreover, various future projects with ton-scale detectors are now ongoing and expected to have significantly improved sensitivities. To test the wino DM scenario in the direct detection experiments, one needs to evaluate the wino-nucleon scattering cross section precisely, with the theoretical uncertainties being sufficiently controlled. This scattering is induced by loop diagrams if the higgsino and squarks are much heavier than wino [25]. At present, the leading order (LO) calculation for the scattering cross section is given in the literature [26][27][28][29]; in these works, the SI scattering cross section with a nucleon is evaluated as σ SI ∼ 10 −47 cm 2 . For other relevant works, see Refs. [30,31]. Since the predicted scattering rate of the wino DM is larger than those of the neutrino backgrounds [32], one expects that the future direct detection experiments may eventually catch a signal of the wino DM. However, it was pointed out by the authors of Ref. [31] that the present calculation may suffer from large uncertainties. They further found that these uncertainties mainly come from the neglect of the higher order contribution in perturbation theory, not from the error of the nucleon matrix elements, which may alter the SI cross section by a factor. To reduce the uncertainties, therefore, we need to go beyond the LO calculation.
In this paper, we complete this calculation up to the next-to-leading order (NLO) in the strong coupling constant α s . For this purpose, we first reformulate the computation based on the effective theoretical approach. The relevant interactions are expressed in terms of the effective operators, whose Wilson coefficients are given up to the NLO with respect to α s . The coefficients are evolved down to the scale at which the nucleon matrix elements of the effective operators are evaluated, by means of the renormalization group equations (RGEs). This procedure allows us to include the NLO QCD effects systematically.
The rest of the paper is organized as follows. In the next section, we describe the formulation mentioned above. All of the matching conditions as well as the RGEs are presented here. Then, in Sec. 3, we show our results for the SI scattering cross section and discuss the uncertainties of the calculation. In Sec. 4, we also present the results for a generic SU(2) L multiplet DM. Those who are interested in a quick reference may find our results in these two sections. Section 5 is devoted to conclusion and discussion.

Formalism
In this section, we give a formalism to evaluate the SI scattering cross section of the wino DM with a nucleon. We will carry out the calculation up to the NLO in the strong coupling constant α s . The formalism given here is based on the method of effective field theories, which consists of the following three steps. Firstly, we obtain the effective operators at the electroweak scale µ W m Z (m Z is the mass of the Z boson) by integrating out heavy particles whose masses are not less than the electroweak scale. This step is carried out in terms of the operator product expansions (OPEs). Secondly, we evolve the Wilson coefficients of the effective operators using the RGEs down to the scale at which the nucleon matrix elements of the operators are evaluated. Finally, we express the SI effective coupling of a wino DM with a nucleon in terms of the Wilson coefficients and the nucleon matrix elements. From this effective coupling, one readily obtains the SI scattering cross section.

Effective Lagrangian
First let us formulate the effective Lagrangian which gives rise to the SI interactions of the wino DM with quarks and gluon. The effective Lagrangian comprises two types of the higher dimension operators-the scalar and the twist-2 type operators-as follows [33]: Here χ 0 , q, and G a µν denote the wino DM, quarks, and the field strength tensor of gluon field, respectively; m q are the masses of quarks; M is the mass of the wino DM; O q µν and O G µν are the twist-2 operators of quarks and gluon, respectively, which are defined by #1 with D µ the covariant derivative. These effective operators are renormalized at the electroweak scale µ W m Z with N f = 5 active quarks (q = u, d, s, c, b). The Wilson coefficients of the operators are to be determined below. Notice that we have included the strong coupling constant α s /π in the definition of the gluon scalar-type operator O G S [35]. We will discuss the validity in the next subsection.

Nucleon matrix elements
In order to compute the scattering cross section of the wino DM with a nucleon, we need the nucleon matrix elements of the scalar and twist-2 type quark and gluon operators #1 We have changed the definition of O G µν by a factor of −1 from those in Refs. [26][27][28]33]. We follows the conversion in Ref. [34].  [36,37].

Proton
Neutron 0.009 (22) presented above. Since these two types of the operators do not mix with each other under the renormalization group (RG) flow, it is possible to consider these two types separately. For the scalar-type quark operators, we use the results from the QCD lattice simulations. The values of the mass fractions of a nucleon N (= p, n), which are defined by f (N ) are shown in Table 1. Here m N is the nucleon mass. They are taken from Ref. [29], which are computed with the recent results of the lattice QCD simulations [36,37]. The nucleon matrix element of O G S , on the other hand, is evaluated by means of the trace anomaly of the energy-momentum tensor in QCD [38]: Here the beta-function β(α s ) and the anomalous dimension γ m are defined by the following equations: whose explicit forms will be given in Eqs. (2.37) and (2.38), respectively. By putting the operator (2.5) between the nucleon states at rest, we obtain This formula is obtained with N f = 3 quark flavors. Notice that the relation (2.5) is an operator equation and thus scale-invariant. This is because the energy-momentum tensor is corresponding to the current of the four momentums, which is a physical quantity and thus not renormalized. As a consequence, Eq. (2.7) should hold at any scales. We will evaluate the matrix element at the hadronic scale µ had 1 GeV in the following calculation.
Since β(α s ) = O(α 2 s ), the r.h.s. of Eq. (2.7) have a size of O(m N ). Namely, although we include a factor of α s /π in the definition of O G S , its nucleon matrix element is not suppressed by the factor. It should be also noted that the scalar-type quark operator Table 2: Second moments of the PDFs of proton evaluated at µ = m Z . We use the CJ12 next-to-leading order PDFs given by the CTEQ-Jefferson Lab collaboration [39].
m qq q is scale-invariant to all orders in perturbation theory (in a mass-independent renormalization scheme) and then the matrix element is independent of the scale at the LO in α s . This is another reason for our definition of O G S . Finally, the nucleon matrix elements of the twist-2 operators are given by the second moments of the parton distribution functions (PDFs): Here q (N ) (x, µ),q (N ) (x, µ) and g (N ) (x, µ) are the PDFs of quark, antiquark and gluon in nucleon at the scale µ, respectively. Contrary to the case of the scalar matrix elements, we have the values of the PDFs at various scales. In Table 2, for example, we present the second moments at the scale of µ = m Z . Here we use the CJ12 next-to-leading order PDFs given by the CTEQ-Jefferson Lab collaboration [39]. It turns out that with the definition of the gluon twist-2 tensor given in Eq.

Wilson coefficients
Now we evaluate the Wilson coefficients of the effective operators at the electroweak scale µ W to the NLO in α s /π. We use the MS scheme in the following calculation. The scattering of a pure neutral wino χ 0 with a nucleon is induced via the weak interactions accompanied by the charged winos χ ± . The interaction Lagrangian is given by where g 2 and W µ are the SU(2) L gauge coupling constant and the W boson, respectively. Since the winos do not couple to the Higgs field directly and the mass difference ∆M between the neutral and charged winos is radiatively generated after the electroweak symmetry breaking, ∆M is much smaller than the DM mass itself or other masses which enter into our computation; according to the recent NLO computation given in Ref. [40], ∆M 165 MeV. Therefore, we safely neglect it in the following discussion. Before looking into the details of the calculation, we first summarize the procedure of the computation as well as the approximations we have used in the calculation. In Fig. 1, we show the diagrams which induce the couplings of wino DM with quarks and gluon, respectively [26][27][28][29]. These diagrams are classified into two types; one is the Higgs exchange type like the upper two diagrams and the other is the box diagrams corresponding to the lower two. We separately discuss each two type.
The Higgs contribution only induces the scalar-type operators. For the NLO-level calculation, we need to evaluate the two-and three-loop diagrams for the quark and gluon scalar-type operators, respectively.
For the box-type contribution, on the other hand, the NLO-level calculation requires us to determine the Wilson coefficients of the operators m qq q, αs π G a µν G aµν , and O i µν to O(α s /π). We first carry out the OPEs of the correlation function of the electroweak Table 3: Number of loops in diagrams relevant to the O(α s /π) calculation for each operator. We also show where we neglect the third generation contribution at the NLO. Here "−" means that there is no contribution or the contribution vanishes.

Operators
Higgs Box currents, as described in Refs. [27,28]. For the scalar operators, the NLO contribution to the OPEs of the correlation functions of vector and axial-vector currents is evaluated in Ref. [41] in the degenerate quark mass limit for each generation. The results are directly applicable to the contribution of the first two generations in our calculation since all of the quarks of the generations may be regarded as massless. Concerning the third generation contribution, the mass difference between top and bottom quarks is significant, and thus the mere use of the results in Ref. [41] is not justified. Their contribution is, however, found to be small compared with those of the first two generations. In our calculation, we neglect the NLO contribution of the third generation, and take into account the effects as a theoretical uncertainty. The Wilson coefficients of the twist-2 operators are evaluated in Ref. [42] to O(α s /π) in the massless limit. It is again not possible to use the results for the contribution of the third generation, and thus we will drop the contribution and estimate the effects as a theoretical uncertainty. By evaluating the W boson loop diagrams with this correlation function, we then obtain the Wilson coefficients of the operators in Eq. (2.2). As a result, C q S , C i T 1 , and C i T 2 are computed at the two-loop level, while C G S is evaluated at the three-loop level. In Table 3, we summarize the number of loops in diagrams relevant to the NLO calculation for each contribution. They complete the NLO matching condition for each Wilson coefficient at the electroweak scale µ W . In addition, we show in the table where we ignore the third generation contribution. As we will see below, the effect of dropping the NLO third-generation contribution is actually negligible.

Higgs exchange
The Higgs exchange processes are induced by the effective coupling of the wino DM with the Higgs boson. They only give the scalar-type interactions as we show in Table 3.
In the case that the wino DM is close to the electroweak eigenstate, the coupling is generated at one-loop level: where w ≡ m 2 W /M 2 with m W and M being the masses of W boson and wino, respectively, and c H (w) = Here g H (x) is a mass function presented in Ref. [26]. #2 By using the effective coupling we readily obtain the LO matching condition for the scalar-type quark operators as Here m h is the mass of the Higgs boson and α 2 ≡ g 2 2 /(4π). To evaluate the NLO matching condition, one needs to evaluate the QCD corrections in the full and effective theories at two-and one-loop levels, respectively. These corrections turn out to be equivalent, and thus the matching condition does not differ from the above equation, i.e., to the NLO in perturbation theory. For the scalar-type gluon operator, the one-loop long-distance contribution by the scalar-type quark operators is subtracted from the two-loop contribution in the full theory so that only the top-quark contribution is included in C G S . Then, we have [38] At the NLO, the above expression is modified to [43,44] Notice that it contains no logarithmic terms like those containing a factor of ln(m t /µ W ). This is because αs π G a µν G aµν is renormalization-group invariant up to this order in perturbation theory. #2 The mass functions used in text are collected in Appendix A.

Box type
Let us move on to the contribution of the box diagrams. They induce both scalar-type and twist-2 operators. To compute the effective operators, we first consider the OPEs of the correlation function of the charged currents: with P L ≡ (1 − γ 5 )/2. We evaluate the Wilson coefficients of the scalar and twist-2 operators in the OPEs up to the NLO in α s /π. We first consider the scalar part. It is convenient to decompose the correlator into the transverse and the longitudinal parts as Here we give only the transverse part since the longitudinal one does not contribute to C q S and C G S [27,28]. As for the contribution to the scalar-type quark operators of the first two generations, there is no O(α 0 s ) term since the charged current J W µ is pure chiral (we take small quark mass limit for q = u, d, s, c, b). Thus, only the one-loop diagrams are relevant in this case. It readily follows from the results given in Ref. [41], in which the correlation functions for vector and axial currents are evaluated with the OPEs, that for q = u, d, s, c. On the other hand, the tree-level contribution of the bottom quark to the scalar-type operator does not vanish because of the large top-quark mass. We have (2.24) Here, as mentioned above, we neglect the NLO contribution and take its effects into account as a theoretical uncertainty. The gluon contribution of the first two generations is also obtained straightforwardly from Ref. [41]. The contribution of the third generation quarks, however, is not evaluated reliably by means of the method used in Ref. [41] due to the large mass of top quark. Here again, we neglect the NLO effects and consider them as a theoretical uncertainty. As a result, we obtain (2.25) where the first and second terms in bracket correspond to the contribution of the first two generations and the third generation, respectively. Next, we consider the twist-2 part. For the contribution of q = u, d, s, c to the quark twist-2 operators, the relevant parts are written as The Wilson coefficients c q W,2 and c q W,L are evaluated in Refs. [42] as follows: For the third generation contribution, on the other hand, we take into account top mass in the LO part and neglect the NLO part as mentioned above. As a result, we have (2.29) Note that we have included the logarithmic part though it is induced at the NLO; otherwise, the Wilson coefficient shows wrong dependence on the factorization scale µ W . Finally let us derive the gluon twist-2 operator. It is always induced at O(α s /π). For the contribution of massless quarks, we use the results given in Refs. [42]. The result is with a factor of four counting the number of the first two generation quarks. As before, we neglect the NLO contribution of the third generation quarks but keep its logarithmic part in order to guarantee the appropriate scale dependence. This reads Then, the sum of the above contributions gives the total twist-2 contribution: Our remaining task is to obtain the Wilson coefficients of the effective operators in Eq. (2.2) by computing another loop with the electroweak current correlator Π W µν (q). For the scalar-type operators, we have where τ ≡ m 2 t /M 2 . The mass function g B1 (x) is given in Ref. [26], and g top (x, y) and g btm (x, y) are equivalent to g (1) B3 (x, y) and g (2) B3 (x, y) in Ref. [28], respectively. These functions are also presented in Appendix A.
For the twist-2 type operators, on the other hand, we have where the functions g T i (x, y), h T i (x) and g log T i (x, y; µ W ) are given in Appendix A. g T i (x, 0) agrees with g T i (x) in, e.g., Ref. [26]. The terms proportional to g log T i (x, y) come from the , and g B1 (w) in red solid, blue dotted, green dash-dotted, and gray dashed lines, respectively, to show smallness of third generation contributions.
logarithmic terms in the OPEs of the correlation function of the charged currents, while the terms with g T i (x, y) and h T i (x) are from the non-logarithmic terms in c and c q/G W,L , respectively. The NLO contribution to the gluon twist-2 operator is also given in Ref. [31].
Here we note that to obtain the proper dependence of the above coefficients on the scale µ W , we need to include all of the NLO corrections. Otherwise, the mismatch in the scale dependence between the matching conditions and the RGEs causes large uncertainties.
To that end, it is important to appropriately perform the order counting with respect to α s /π. Especially, the two-loop contribution to C G T i should be regarded as the NLO in α s /π, #3 not the LO, which is contrary to the case of the gluon scalar operator C G S ; in this case, the two-loop contribution is the LO in α s /π. Our convention for the definition of the gluon operators clarifies this order counting.
As we have already commented several times, we neglect the NLO contribution of the third generation quarks. Indeed, we expect that its significance is quite small, and thus we safely regard it as a theoretical uncertainty. In Fig. 2 we compare the mass functions corresponding to the LO third generation contributions with those of the LO massless quark contributions, which corresponds to g top (w, τ )/g B1 (w), g T 1 (w, τ )/g T 1 (w, 0), and g T 2 (w, τ )/g T 2 (w, 0). g btm (w, τ )/g B1 (w) is also shown as its contribution to C G S via integration of the bottom quark is given by −C b S /12. It is found that the LO third generation contributions are smaller than those of the first and second generations by almost an order of magnitude. Hence, we expect that the NLO contributions of the third generation are also considerably small compared with those of the other two generations. This allows us #3 One may easily check that the logarithmic parts in the NLO contribution to the twist-2 operators reproduce the one-loop RGEs presented in Sec. 2.4. This justifies the order counting discussed here.
to ignore the third-generation NLO contribution, and treat it as a theoretical uncertainty.

Renormalization group equations and matching conditions
The effective operators are scale dependent and their scale evolution is described by the RGEs. During the RG evolution, heavy quarks are integrated out around their mass scale. Thus we need to match the theories above and below the threshold. Here we summarize the RGEs and the matching conditions.
To begin with, we write down beta-function of α s and anomalous dimension of quark mass operator: is the number of colors, N f denotes the number of quark flavors in an effective theory and C F is the quadratic Casimir invariant defined by C F ≡ N 2 c −1 2Nc .) Here for the MS quark masses, we use the one-loop anomalous dimension since their effects first appear at the NLO level as we will see below soon. Now we give the RGEs for the Wilson coefficients of the above operators. First, we consider the RGEs for the scalar-type operators. To that end, notice that the quark mass operator is RG invariant in a mass-independent renormalization scheme like the MS scheme, i.e., To evaluate the evolution of the gluon scalar operator, we use the trace anomaly formula (2.5). Differentiating Eq. (2.5), we then obtain the differential equation for the gluonic scalar operator αs π G a µν G aµν . As a result, we have #4 where Γ S is a (N f + 1) × (N f + 1) matrix given by

41)
#4 In fact, we implicitly assume that the operators are to be evaluated between the on-shell states. As discussed in Refs. [43,45], during the RG flow, the scalar operators mix with other (gauge-variant) operators whose on-shell matrix elements vanish.
The solutions of the RGEs are given as follows: Next, we consider the RGEs for the twist-2 operators. The two-loop anomalous dimension matrix of the operators is evaluated as [46,47] (2.46) Finally we give the threshold corrections at the scale where heavy quarks are integrated out. For example, in the vicinity of the bottom-quark threshold µ b m b , we match the strong gauge coupling constant and the Wilson coefficients as and with q = u, d, s, c for the first and third equations. #5 In the following section, we estimate the uncertainties coming from the neglect of the higher order perturbation by varying the matching scale µ b around the µ b m b . We repeat a similar procedure for the charm-quark threshold around µ c m c .
Here we note that besides the above threshold corrections, the higher dimension operators suppressed by a power of the threshold quark mass are also generated in general. For instance, if the scalar-type quark operator is integrated out at a quark threshold m Q , then we will obtain the following dimension-nine operators at one-loop level [45,48]: where f abc is the SU(3) structure constant. In particular, those generated at the charmquark threshold give the largest effects. By using the naive dimensional analysis, we see that their contribution to the nucleon matrix element may give a correction by a factor of Λ 2 QCD /m 2 c = O(0.1), which could be additionally suppressed by the prefactors of these operators. Since we do not know precise values of the nucleon matrix elements of the operators in Eq. (2.49), we should also consider their effects as an uncertainty.

Results
Now we compute the wino-nucleon scattering cross section and evaluate the theoretical uncertainties. We first separately consider the scalar and twist-2 contributions to the wino-nucleon effective coupling in Sec. 3.1 and 3.2, respectively. Then, we show the result for the scattering cross section in the following subsection. In Table 4, we summarize the input parameters we use in our computation. For the mass of top quark, we use the pole mass as an input parameter, and convert it to the MS mass using the one-loop relation: The matching condition for C G Ti here differs from that given in Ref. [31].

Scalar part
The spin-independent effective coupling of the wino with nucleon is defined by The contribution of the scalar operators to the coupling is given by where we take the hadron scale µ had = 1 GeV with N f = 3 active quarks. Fig. 3 shows f p scalar with various types of errors. In Fig. 3 (a) f p scalar at the LO (blue dashed) and NLO (red solid) with corresponding bands showing the theoretical error due to the perturbative calculation are shown. In the plot the uncertainty coming from lack of the NLO contribution of the third generation, which is multiplied by a factor of five just for the purpose of presentation, is also shown (gray band). For the evaluation of the error from the ignorance of higher order contribution in perturbation, we vary each matching scale by a factor of two; i.e., The prescription is, however, less effective for the scalar-type operators since these operators are almost scale-invariant. For this reason, when evaluating the error resulting from the quark threshold matching for the NLO (LO) calculation, we use the NNLO (NLO) matching conditions to artificially generate the logarithmic dependence of the Wilson coefficients on the scale by using the mismatch between the matching conditions and RGEs. The NLO matching conditions are given in Eq. (2.48), while the NNLO ones are found in Ref. [53]. In addition, for the LO contribution, we evaluate the uncertainty caused by the electroweak-scale matching by merely multiplying the LO contribution by a factor of α s /π. Since the scalar-type operators are scale-invariant at the LO, it is impossible to estimate the LO uncertainty from the electroweak-scale matching by varying the scale µ W . At the NLO, on the other hand, we are to estimate the uncertainty with the scale variation since the NLO RGEs yield the scale dependence of the scalar operators.  The error from the LO perturbative calculation is more than 5%, which reduces to a few % level with the NLO calculation. The upper errors smaller than the lower errors in the LO and NLO perturbative calculations in the Fig. 3 (a). This comes from difference between α s (m Q /2) and α s (2m Q ) for Q = b, c. On the other hand, as for the uncertainty due to the lack of the third-generation NLO contribution, we estimate its effect by multiplying the LO contribution by a factor of α s /π. From the figure, we find that the ignorance of the third-generation NLO contribution only gives a negligible effect on the resultant value. The effect is much smaller than the uncertainty due to the perturbative calculation. Fig. 3 (b) shows comparison of the uncertainty in the NLO perturbative QCD calculation (pink) with that from the errors in the input parameters we have used in the calculation (gray). Among them, the uncertainty coming from the Higgs mass error is especially shown in the dark red band. We see that thanks to the NLO calculation the perturbative error now becomes smaller than the error from the input parameters, though they are still of the same order of the magnitude.
Finally we plot the theoretical uncertainty which could arise due to the higher dimension operators induced at each quark threshold in Fig. 3 (c). To evaluate the effects of the higher dimension operators, we vary the scalar gluon contribution induced at the charm-quark threshold by 2%, which is expected from the naive dimensional analysis as discussed in Sec. 2.4. #6 Since the higher dimension operators generated at the bottomquark threshold are suppressed by the bottom quark mass, their effects are negligible. As seen from the figure, this uncertainty may be as large as the NLO perturbative QCD error. To reduce the uncertainty, one of the most efficient ways is to use the nucleon matrix elements computed above the charm-quark threshold, say, at the scale of 2 GeV. In this case, we need to evaluate the charm-quark content in nucleon, f (N ) Tc = N |m cc c|N /m N , as well. Currently, the QCD lattice simulations are not able to compute it accurately [54]. If this quantity is evaluated with good precision in the future, then the uncertainty due to the higher dimension operators will be significantly reduced. We expect that the perturbative QCD error will also decrease, since we do not need the charm-quark threshold matching procedure any more. Thus, we strongly encourage the development in this field.

Twist-2 part
Contrary to the scalar-type operators, the twist-2 operators have the scale dependence at the leading order in α s . Therefore, it is necessary to determine the appropriate scale for the matching of the full theory onto the effective theory in order not to suffer from large logarithmic factors. To that end, we require that the logarithmic dependent parts g log T i in the Wilson coefficients presented in Eq. (2.36) should not be large, say, within O(1). Since the terms proportional to g log T i come from the logarithmic terms in the OPEs of the correlation function of the charged currents, this condition guarantees the validity of the perturbative QCD expansion. In Fig. 4, we show g log T i (w, 0; µ W ) (i = 1, 2) as function of the factorization scale µ W . Here M = 3 TeV (solid) and 300 GeV (dashed). The vertical gray line shows µ W = m Z . It turns out that the size of these functions is within O(1) #6 Since the first (second) operator in Eq. (2.49) receives additional suppression by a factor of five (sixty) compared with the contribution of the scalar gluon operator, −α s /(12π)GGχ 0 χ 0 , we estimate the significance of the former contribution as ∼ 2% of that of the scalar-type gluon operator, while the latter contribution is negligible. if one takes the scale µ W to be around the electroweak scale. This consequence rarely depends on the DM mass. The absolute values for these functions are minimum at a scale of O(10) GeV, which is much smaller than the DM mass. This observation reflects the fact that the typical scale of the loop momentum flowing in the loop diagrams in Fig. 1 is around the electroweak scale, as pointed out in Ref. [25]. In the following calculation, we take µ W = m Z , which assures that g log T i is within O(1) and thus the perturbative expansion is justified.
To calculate the contribution of the twist-2 operators, we also need to choose the scale at which the nucleon matrix elements of the twist-2 operators are evaluated. As mentioned above, contrary to the case of the scalar-type operators, the twist-2 matrix elements are obtained at various scales. Since the result does not depend on the choice of the scale within the uncertainty of the calculation, it is desirable to choose the scale so that the error in calculation is reduced. Thus, we take it to be the same as the factorization scale, i.e., µ = m Z . This choice allows us to decrease the error which would arise from the process where the operators are evolved down to the low-energy region; for instance, if one evaluates the matrix elements at a scale µ < m b , the result suffers from the uncertainty resulting from the bottom-quark mass threshold. See Ref. [35] for further discussion. Now we evaluate the contribution of the twist-2 operators to the SI effective coupling in Eq. (3.51), which is given by where q runs over the active quarks (q = u, d, s, c, b for our choice of the scale µ = m Z ). In Fig. 5, we show f p twist2 as function of the wino mass. We compare the LO and NLO results in the left panel, shown in the blue dashed and red solid lines, respectively, with the corresponding bands representing the uncertainties. The uncertainties are evaluated by varying the scale µ W between m Z /2 and 2m Z . Besides, it is found that to drop the NLO contribution of the third generation quarks causes only the negligible effects, so we do not show the error due to the contribution. The O(1)% error in the LO computation now reduces to ∼ 0.5% when going to the NLO level, though the central value shifts more than expected, i.e. about 5% change. This is due to a large NLO term in C q T i of Eq. (2.36). In the large DM mass limit, the contributions of quarks and gluon at the NLO are 0.90 and −0.047 in 10 −9 GeV −2 unit, respectively, while the quark contribution at the LO is 0.82 in 10 −9 GeV −2 unit. #7 In the right panel of Fig. 5, we also illustrate the uncertainty resulting from the input error, which turns out to be as large as the NLO uncertainty. The uncertainty mainly comes from those of the PDFs, which we estimate following the method given in Ref. [39] with the χ 2 tolerance T taken to be T = 10. After all, in the case of the twist-2 contribution, both the NLO and input uncertainties are less #7 To be concrete, in C q Ti the NLO term summed over i = 1, 2 gives (α 2 2 α s /4πm 3 W ) × (41π/12) in the large DM mass limit. Here logarithmic term g log Ti is neglected for simplicity. See also Eqs. (A.18)-(A.31) for the mass functions in the large DM mass limit. Figure 6: Wino-proton SI scattering cross section. Blue dashed and red solid lines represent LO and NLO results, respectively, with corresponding bands show perturbative uncertainties. Gray band shows uncertainty resulting from the input error. Yellow shaded area corresponds to the region in which neutrino background overcomes DM signal [32]. than 1%, and thus well controlled compared to the scalar contribution.

Scattering cross section
Finally, we evaluate the wino-nucleon SI scattering cross section, which is given by We plot σ p SI as function of the wino mass in Fig. 6. Additionally we indicate the parameter region where the neutrino background dominates the the DM-nucleon scattering [32] and then it becomes hard to detect the DM signal in the DM direct detection experiments (yellow shaded). Here we estimate each error by varying the scalar and twist-2 contributions within their uncertainties evaluated above. The result shows that the large uncertainty in the LO computation is significantly reduced once the NLO QCD corrections are included, which is now smaller than that from the input error. In the large DM mass limit, the SI scattering cross section converges to a constant value, where the first and second terms represent the perturbative and input uncertainties, respectively. As seen from Fig. 6, σ p SI has little dependence on the DM mass; its variation is actually within the uncertainties of the calculation, for the wino mass larger than 270 GeV. Both the scalar and twist-2 contributions depend on the DM mass when the mass is smaller than ∼ 1 TeV as shown in Figs. 3 and 5. However, the dependence in the cross section is accidentally canceled. The NLO result is found to be larger than the LO result by almost 70%. After all, the resultant scattering cross section is well above that of the neutrino background [32], and therefore the future direct detection experiments are promising to test the wino DM scenario.

Electroweakly-interacting DM
Although we have focused on the wino DM in this paper, a similar formalism may be constructed for a more general class of the DM candidates; i.e., an SU(2) L multiplet with hypercharge Y that contains a neutral component for DM, and their thermal relic may explain the observed DM density with O(1) TeV masses. For previous works on such DM candidates, see Refs. [55][56][57][58][59][60][61][62]. Some theories beyond the Standard Model actually predict this kind of DM. For example, the higgsino and wino in the SUSY models are representative of the SU(2) L multiplet DM. Moreover, such a particle may show up in grand unified theories [63][64][65], whose stability is explained by a remnant discrete symmetry of extra U(1) symmetries in the theories [66][67][68][69][70]. Before concluding our discussion, we give the results of the NLO calculation for this class of DM candidates. If the DM particle is a fermion, its interactions with quarks and gluon are completely determined by the electroweak gauge interactions, #8 so we consider the fermionic DM candidates in the following discussion. If Y = 0, the DM is a Dirac fermion, while a Majorana fermion if Y = 0. Pure Dirac fermion DM is, however, severely constrained by the direct detection experiments already, since the vector interactions via the Z boson exchange yield too large scattering cross section with nucleon. The constraint may be evaded if there are some new physics effects that give rise to the mass difference between the neutral components to split them into two Majorana fermions. If the mass difference is larger than O(100) keV, the scatterings with nucleon are not induced by the tree-level Z boson exchange. In what follows, we assume the presence of the mass difference and regard the lighter neutral component χ 0 as a DM candidate. The mass difference is assumed to be small enough to be neglected in the following calculation. In this case, the interactions including the neutral components are given by Here n is the number of the components in the DM SU(2) L multiplet, g Z ≡ g 2 Y + g 2 2 with g Y the U(1) Y gauge coupling constant, and η 0 and Z µ for the heavier neutral component and the Z boson, respectively. The LO calculation of the scattering cross section with a nucleon for this type of DM candidates is given in Ref. [28]. As in the case of the wino DM, we find that there is a significant cancellation among the contributions to the scattering amplitude. Therefore, the NLO corrections are of importance to evaluate the scattering cross section precisely. We compute the NLO scattering cross section in a similar manner to above discussion. The only difference is the electroweak matching conditions, which we summarize in Appendix B. Below the electroweak scale, the procedure is completely the same as before.
In Fig. 7 we plot the SI scattering cross sections for several SU(2) L multiplet DM candidates. Here the red solid, green dashed, and blue dash-dotted lines represent the (n, Y ) = (3, 0), (2, 1/2), and (5, 0) cases, respectively. The triplet case corresponds to the wino DM, while the doublet one is regarded as the higgsino DM. The (n, Y ) = (5, 0) fermion DM is the so-called minimal DM [55], for which the gauge symmetry guarantees its stability. Again, the yellow shaded area indicates the region in which neutrino background overcomes the DM signal [32]. We find that all of the scattering cross sections are almost constant in the mass region we are interested in, as already seen in the case of wino DM. In the heavy DM mass limit, the DM-proton effective coupling f p ≡ f p scalar + f p twist2 at the NLO is given by from which one readily obtains the SI scattering cross section for a generic SU(2) L DM candidate. It is seen that the (n, Y ) = (3, 0) and (5, 0) cases offer the SI scattering cross sections well above the neutrino background, while that of the (n, Y ) = (2, 1/2) case falls far below the background. Compared to the previous results in Ref. [28], slightly larger SI scattering cross sections are obtained for DM candidates with Y = 0. As for the (n, Y ) = (2, 1/2) case, on the other hand, we obtain a smaller SI scattering cross section. #9

Conclusion and discussion
In this paper we have completed the calculation of the wino-nucleon scattering cross section up to the NLO in α s /π. It turns out that the inclusion of the NLO corrections allows us to reduce the theoretical uncertainty significantly, which is now O(10)% level. The NLO scattering cross section is larger than the LO one by about 70%. The resultant cross section is well above the neutrino background, and thus the DM direct detection experiment is a promising tool for examining the wino DM scenario. In addition, we give the NLO results for the cases with a generic SU(2) L multiplet DM, some of which may also be probed in future experiments. At present, the uncertainties from the input parameters, especially those of the scalar matrix elements, dominate the theoretical error. If future lattice simulations determine the charm-quark content in nucleon with good accuracy, the uncertainties are to be reduced considerably. We strongly anticipate the developments in the field.

A Mass functions
Here we list the mass functions used in text: are equal to g T 1 (x) and g T 2 (x) in Ref. [26], respectively. On the other hand, g log T i (x, y; µ W ) are given by the following integrals: We compute these integrals numerically. For the generic SU(2) L DM case, we further introduce the following functions: Again, these integrals are evaluated numerically. The functions f anl V (x, y) and f anl A (x, y) are given by functions in Ref. [28] as f anl V (x, y) = G t1 (x, y)/4 and f anl A (x, y) = G t2 (x, y)/4. In the large DM mass limit, i.e., x, y → 0 with the ratio y/x fixed, the above analytic functions are reduced to as follows: Here we have set the values for the masses of Z boson and top quark in z and τ , respectively.

B Results for the electroweak-interacting DM
In this Appendix, we summarize the electroweak matching conditions for generic SU(2) L multiplet DM.

B.1 Current correlator
To begin with, we consider the OPEs of the electroweak current correlators as in Sec. 2.3.2. The correlation function of the charged currents has been already discussed there. Here we give the OPEs of the neutral current correlator, for it is necessary to evaluate the Z boson contribution. The correlation function of the weak neutral current is defined by Let us first evaluate the Wilson coefficients of the scalar operators. For the scalar operators, the correlator is decomposed to the transverse and longitudinal parts as Π Z µν (q)| scalar = −g µν + q µ q ν q 2 Π Z T (q 2 ) + q µ q ν q 2 Π Z L (q 2 ) .

(B.35)
Again, only the transverse part is relevant to the calculation. The OPE coefficients are defined by Π Z T (q 2 ) = q c q Z,S (q 2 ; µ W )m qq q + c G Z,S (q 2 ; µ W ) α s π G a µν G aµν , (B.36) are then evaluated as follows [41]: with q = u, d, s, c, b. Here we drop the NLO contribution of top quark for simplicity. This contribution is also readily obtained from the results in Ref. [41]. The LO terms in the above equations agree with the results given in Ref. [28]. Next, we consider the twist-2 operators. Their contribution to the correlation function is written as [42] Π Z µν (q)| twist2 = g 2 Z i=q,G − g µρ g νσ q 2 − g µρ q ν q σ − q µ q ρ g νσ + g µν q ρ q σ (q 2 ) 2 c i Z,2 Here again, we have neglected the top-quark contribution to the NLO gluon coefficients.

B.2 Wilson coefficients
Now we calculate the electroweak-scale matching conditions. For the scalar-type quark operators, we have for q = u, d, s, c, and where θ W is the weak mixing angle and z ≡ m 2 Z /M 2 . #10 The Wilson coefficient of the scalar-type gluon operator is, on the other hand, computed as where the functions f V (x, y) and f A (x, y) are given in Appendix A. The twist-2 contribution is given by g log T i (w, 0; µ W ) + 9 4 g T i (w, 0) + 16 9 h T i (w) #10 Note that g S (x) in Ref. [28] is equal to 6g B1 (x).