Dark Matter in Split SUSY with Intermediate Higgses

The searches for heavy Higgs bosons and supersymmetric (SUSY) particles at the LHC have left the minimal supersymmetric standard model (MSSM) with an unusual spectrum of SUSY particles, namely, all squarks are beyond a few TeV while the Higgs bosons other than the one observed at 125 GeV could be relatively light. In light of this, we study a scenario characterized by two scales: the SUSY breaking scale or the squark-mass scale $(M_S)$ and the heavy Higgs-boson mass scale $(M_A)$.We perform a survey of the MSSM parameter space with $M_S<10^{10}$ GeV and $M_A<10^4$ GeV such that the lightest Higgs boson mass is within the range of the observed Higgs boson as well as satisfying a number of constraints. The set of constraints include the invisible decay width of the $Z$ boson and that of the Higgs boson, the chargino-mass limit, dark matter relic abundance from Planck, the spin-independent cross section of direct detection by LUX, and gamma-ray flux from dwarf spheroidal galaxies and gamma-ray line constraints measured by Fermi LAT. Survived regions of parameter space feature the dark matter with correct relic abundance, which is achieved through either coannihilation with charginos, $A/H$ funnels, or both. We show that future measurements, e.g., XENON1T and LZ, of spin-independent cross sections can further squeeze the parameter space.


Introduction
Supersymmetry (SUSY) is one of the most elegant solutions, if not the best, to the gauge hierarchy problem. SUSY provides an efficient mechanism to break the electroweak symmetry dynamically with a large top Yukawa coupling. Another virtue is that the lightest SUSY particle (LSP) is automatically a dark matter (DM) candidate to satisfy the relic DM abundance assuming the R-parity conservation. The fine-tuning argument in the gauge hierarchy problem requires SUSY particles at work at the TeV scale to stabilize the gap between the electroweak scale and the grand unified theory (GUT) scale or the Planck scale. With this scale the gauge coupling unification is also naturally achieved in renormalization group equation (RGE) running.
Although SUSY has quite a number of merits at least theoretically, the biggest drawback of SUSY is that so far we have not observed any sign of SUSY. Nevertheless, we have observed a light standard model (SM) like Higgs boson, which is often a natural prediction of SUSY. The null results for all the searches of SUSY particles have pushed the mass scale of squarks beyond a few TeV [1]. Such a high-scale SUSY scenario also draws more and more attention on CP problems [2], cosmological problems [3], and DM search [4]. On the other hand, the searches for the SUSY Higgs bosons provide the less stringent mass limits and it still seems possible to find them in the range of a few hundred GeV [5]. Consequently, we are left with an unusual spectrum of SUSY particles and Higgs bosons: (i) all squarks are heavy beyond a few TeV [1], (ii) the gluino is heavier than about 1 TeV [6], (iii) neutralinos and charginos can be of order O(100 − 1000) GeV, (iv) heavy Higgs bosons can be of order O(100 − 1000) GeV [5], and (v) a light Higgs boson with a mass 125 GeV [7]. The spectrum is somewhat similar to the proposal of split SUSY [8], except that the heavy Higgs bosons need not be as heavy as those of split SUSY. We name the scenario the "modified split SUSY" framework, with two distinct scales: the SUSY breaking scale M S and the heavy Higgs-boson mass scale M A . In the following, for simplicity we call this "modified split SUSY" as scenario A in which M S and M A are independent, while the original split SUSY as scenario B in which M A and M S are set to be equal.
We perform a survey of the parameter space of the minimal supersymmetric standard model (MSSM) characterized by two scales: (i) the SUSY breaking scale M S with M S < ∼ 10 10 GeV, and (ii) the heavy Higgs-boson mass scale (M A ) with M A < ∼ 10 4 GeV, such that the lightest Higgs boson mass with large radiative corrections from heavy squarks is within the range of the mass of the observed Higgs boson. We choose M A smaller than or at most equal to M S . Specifically, we assume the MSSM above the SUSY breaking scale M S . Then we do the matching at the scale M S while we decouple all the sfermions. We evolve from M S down to M A with a set of RGEs comprising of two-Higgs doublet model (2HDM), gauge couplings, and gaugino couplings. For this purpose, we derive the RGEs governing the range between M S and M A and present them in Appendix A. Then we do the matching at the scale M A while we decouple all the heavy Higgs bosons. We evolve from M A down to the electroweak scale with a set of RGE comprising of the SM and the gauginos. The matching is then done at the electroweak scale. Once we obtain all the relevant parameters at the electroweak scale, we calculate all the observables and compare to experimental data.
In this work the LSP of the MSSM is the DM candidate, which is the lightest neutralino in the current scenario. Since we are strongly interested in DM, we include a number of other existing constraints on SUSY particles and DM: 1. the invisible decay width of the Z boson and that of the Higgs boson, 2. the chargino-mass limit, 3. dark matter relic abundance from Planck, 4. the spin-independent cross section of direct detection by LUX, and 5. gamma-ray flux from dwarf spheroidal galaxies (dSphs) and gamma-ray line constraints measured by Fermi LAT.
Due to multidimensional model parameters involved in this work, it will be advantageous to adopt a Monte Carlo sampling technique to perform a global scan. In order to assess the robustness of our Monte Carlo results, we investigate both Bayesian maps in terms of marginal posterior (MP) and frequentist ones in terms of the profile likelihood (PL) technique. However, the likelihood functions of experimental constraints are the same for both approaches.
The organization is as follows. In the next section, we describe the theoretical framework of the modified split SUSY, including the matching conditions at the scales of M S and M A , and the corresponding interactions of the particles involved. In Sec. 3, we list the set of constraints from collider and dark matter experiments that we use in this analysis. In Sec. 4, we present the results of our analysis using the methods of PL and MP. We discuss and conclude in Sec. 5.

Theoretical Framework
In the case under consideration, we have the two characteristic scales: the high SUSY scale M S and the Higgs mass scale M A . The relevant phenomenology may be described by the effective Lagrangians depending on scale Q as follows:

Interactions for M A < Q < M S
At the scale M S all the sfermions decouple when we assume that they are heavier than or equal to the scale M S . We are left with the spectrum of the Higgs sector of the 2HDM, gauginos, and higgsinos.
In this work, we take the general 2HDM potential as follows: with the parameterization and v d = v cos β = vc β , v u = v sin β = vs β , and v ≃ 245 GeV. Then we have The wino(bino)-Higgsino-Higgs interactions are given by where σ a are the Pauli matrices. We note

Matching at M S
The couplings of the interactions when M A < Q < M S are determined by the matching conditions at M S and the RGE evolution from M S to Q. Assuming that all the sfermions are degenerate at M S , the quartic couplings at the scale M S are given by * with We note that the quartic couplings at M S consist of its tree level values and the threshold corrections induced by the A and µ terms. We further observe λ 5,6,7 vanish without including the threshold corrections.
On the other hand, for the wino(bino)-Higgsino-Higgs couplings at the scale M S , we haveg We note the relation g ′ = 3/5 g 1 .
The threshold corrections to the gauge and Yukawa couplings at M S also vanish in the framework under consideration or when all the sfermions are degenerate at M S .

Interactions for Q < M A
When the scale drops below M A , all the heavy Higgs bosons decouple. We are left with the SM particles, a light Higgs boson, gauginos, and higgsinos.
The SM Higgs potential is given by where G ±,0 denotes the would-be Goldstone bosons and h the physical neutral Higgs state. We note m 2 h = 2λv 2 . The wino(bino)-Higgsino-Higgs interactions are then given by

Matching at M A
The couplings of the interactions when Q < M A are determined by the matching conditions at M A and the RGE evolution from M A to Q.

Matching at the electroweak scale
Matching at the electroweak scale is exactly the same as in the original split SUSY framework. We closely follow Ref. [9] to include the threshold corrections to the gauge couplings at the electroweak scale and to calculate the pole masses for the Higgs boson and the top quark.

Experimental Constraints and Likelihoods
In this section, we describe how to construct the likelihood functions involved with experimental constraints which are used in both MP and PL approaches. For the experimental constraints considered in this work, we assume either half-Gaussian or Gaussian distribution when the central values µ, experimental errors σ, and theoretical errors τ are available. Otherwise, we take Poisson distributions.
In Table 1, in the second last column, we show the likelihoods of each experimental constraint. Here "hard cut" means we apply the 95% upper limits instead of constructing its likelihood. For the details of our statistical treatment, we refer to Appendix B. In the following subsections, we give more details of the constraint and likelihood of each measurement.

Invisible decay widths
The invisible decay width of the Z boson was accurately measured by taking the difference between the total width and the visible width, and is well explained by the three light active neutrino species of the SM. Any additional invisible decays of the Z boson are strongly constrained by this data. In the current framework, the additional invisible width comes from Z → χ 0 1 χ 0 1 . With the invisible width given in the PDG [10], Γ Z inv = 499 ± 1.5 MeV, we can constrain Z → χ 0 1 χ 0 1 . If the neutralino mass is below m h /2, the Higgs boson can decay into a pair of neutralinos, thus contributing to an invisible width of the Higgs boson. From a global fit using the Higgs-boson data at the 7 and 8 TeV runs of the LHC, the invisible width of the Higgs boson is constrained to be Γ h inv < 0.6 MeV [11] at 1-σ level if all other parameters are fixed at their SM values. If other parameters are allowed to vary, the Γ h inv would have a more relaxed limit, which is about the same as the bound from the direct search on the invisible mode of the Higgs boson, which has a branching ratio about 50% [19]. Nevertheless, we use Γ h inv < 0.6 MeV in this work, as shown in Table 1.

Chargino mass
The mass limits on charginos come either from direct search or indirectly from the constraint set by the non-observation of χ 0 2 states on the gaugino and higgsino MSSM parameters M 2 and µ. For generic values of the MSSM parameters, limits from high-energy e + e − collisions coincide with the highest value of the mass allowed by phase space, namely m χ ± < ∼ √ s/2.
The combination of the results of the four LEP collaborations of LEP2 running at √ s up to 209 GeV yields a lower mass limit of m χ ± 1 ≥ 103.5 GeV, which is valid for general MSSM models. However, it could be weakened in certain regions of the MSSM parameter space where the detection efficiencies or production cross sections are suppressed, e.g., when the mass difference m χ ± 1 − m χ 0 1 becomes too small. Regardlessly, we simply employ the mass limit of m χ ± 1 ≥ 103.5 in this work. We do not use the LHC constraint since it is more model dependent and does not give any bounds when m χ 0 1 > ∼ 70 GeV [20]. Furthermore, for m χ 0 1 < 70 GeV region, the H/Z resonance region (see next subsection) is not sensitive to this search [21]. Note that the χχ ± coannihilation is strongly forbidden by this limit especially when m χ 0 To deal with the chargino mass limit without detector simulations, we adopt the half Gaussian distribution when m χ ± 1 < 103.5 GeV to describe the tail of the chargino mass likelihood function. For the likelihood, we assume ∼ 1% theoretical uncertainty. When m χ ± 1 ≥ 103.5 GeV, we always assume the maximum likelihood.

Relic abundance
The half-Gaussian distribution for relic abundance likelihood in Table 1 suits the wellmotivated moduli decay scenario [22,23,24,25,26,27,28]. In this scenario, the relic abundance can be reproduced by moduli decay after the freeze-out, which is different from the usual multi-component DM scenario, in which the total relic abundance is shared among a few DM candidates, such as the axion. In the moduli decay scenario, all the DM is still assumed to be the neutralino, and the DM local density need not be rescaled with respect to the neutralino fraction as implemented in the multi-component DM scenario, so that the DM direct and indirect detection constraints will be stronger.
Very often, the neutralino DM in most of the MSSM parameter space over-produces the relic abundance, because the annihilation in the early Universe is too inefficient. Generally speaking, by opening the W + W − final state the wino-like neutralino can very efficiently reduce relic abundance for wino mass up to 3 − 4 TeV, e.g. see Ref. [27,29,30]. However, it requires some specific mechanisms for bino-like, Higgsino-like, or mixed neutralinos to fulfill correct relic abundance. Sometimes more than one mechanisms are needed. In most cases the (non-wino) regions both of correct relic abundance and still allowed by the current LHC direct searches in our modified split SUSY parameter space are: • The Z/h resonance region, where the neutralinos annihilate through the resonance with the Z boson at m χ 0 1 ∼ 45 GeV and Higgs boson at m χ 0 1 ∼ 62.5 GeV. In this region, neutralinos are governed mainly by the bino fraction but with a small mixing with the higgsino fraction.
• The chargino-neutralino coannihilation region, where the µ parameter is usually closed to gaugino parameters M 1 or M 2 so that the χ 0 1 , χ ± 1 , and χ 0 2 are almost degenerate. If the masses between χ 0 1 and χ ± 1 or χ 0 1 and χ 0 2 are very close to each other, the number densities of the next-to-lightest supersymmetric particle(s) (NLSP(s)) have only slight Boltzmann suppression with respect to the LSP number density. Therefore, all the interactions among the LSP and NLSP(s), such as χ 0 1 − χ ± 1 , χ 0 1 − χ 0 2 and χ ± 1 − χ 0 2 , play important roles to reduce the relic abundance. Note that χ 0 1 in this region shall have nonnegligible fractions of wino or higgsino in order to coannihilate with χ ± 1 and χ 0 2 .
• The A/H funnel region, where neutralinos annihilate through the resonance of the pseudoscalar Higgs boson A or the heavy scalar Higgs boson H. In the original split SUSY framework with M A = M S , because of the large mass of A/H as well as their large decay width, this mechanism becomes irrelevant. On the other hand, in our modified split SUSY scenario with light M A , this A/H funnel can still play a significant role in reducing the relic abundance. Nevertheless, we shall see later that the A/H-funnel for m χ 0 1 > 1 TeV is not efficient enough to reduce the relic abundance because of the larger A/H decay width.
In split SUSY scenario, because of the very heavy sfermion masses, all thef − χ coannihilation channels have been closed. On the other hand, the chargino annihilation is still allowed but the chargino mass must be above the LEP limit, m χ ± > 103.5 GeV. We found in our viable parameter space the majority of bino-like neutralino and chargino is always close to each other (χχ ± coannihilation on). Besides, χχ annihilation can have a few other choices. Lowering M A to less than 1 TeV, the A/H-funnel region can be important, especially for higgsino and mixed neutralino. For m χ 0 1 < 100 GeV, Z-and h-resonances can also significantly reduce relic abundance. Finally, the wino-like neutralinos can easily annihilate into the W + W − final state, which can sufficiently reduce relic abundance as well.

LUX: spin-independent cross section
At present the most stringent 90% C.L. limit on the spin-independent component of the elastic scattering cross section comes from LUX [15]. However, it did not take into account the systematic uncertainties from nuclear physics and astrophysics, otherwise the constraint becomes much less straightforward. The astrophysical uncertainties mainly come from our poor knowledge of the DM local density and velocity distribution. In order to account for the uncertainties of all the astrophysical parameters, we adopt the phase-space density factor and its associated error bars as computed in Ref. [31]. Nuclear physics uncertainties enter the systematic uncertainties through the nuclear matrix elements, mainly the pionnucleon sigma term Σ πN and the strange quark content of the nucleon f T s , which promote the spin-independent cross sections from quark level into nucleon level. In Table 2, we treat the Σ πN and f T s as nuisance parameters and distribute as Gaussian with central values and error bars obtained by recent lattice QCD calculations. Regarding the reconstruction of the LUX likelihood including the astrophysical and nuclear uncertainties, we refer to Ref. [14] for more detailed explanations.

Continuous gamma ray from dSphs
The most luminous gamma-ray source is the Galactic Center (GC) in the Milky Way, but it is also subject to higher astrophysical backgrounds. Better constraints were obtained from the diffuse gamma rays from the dSphs of the Milky Way. They are less luminous and dominated by DM, with little presence of gas or stars. Recently, the Fermi LAT Collaboration improved significantly the previous sensitivities to DM searches from dSphs [17].
Unlike the published limit from the Fermi LAT collaboration, we only include the eight classical dSphs in our analysis, because the DM halo distribution in the classical dSphs is measured with a higher accuracy from the velocity dispersion of the luminous matter [32]. We use the 273 weeks Fermi-LAT data and the Pass-7 photon selection criteria, as implemented in the FermiTools. The energy range of photons is chosen from 200 MeV to 500 GeV, and the region-of-interest is adopted to be a 14 • × 14 • box centered on each dSphs. The J-factors are taken from Table-I in Ref. [17].
In the likelihood analysis, the Fermi-LAT data are binned into 11 energy bins logarithmically spaced between 0.2 and 500 GeV, and we calculate the likelihood map of Fermi-LAT dSphs on the Ebin-flux plane following the method developed in Ref. [16].

Fermi photon line measured from GC
The experimental signature of monochromatic lines over the continuous spectrum is a clean signal of DM annihilation. In MSSM, the annihilation of χ 0 1 χ 0 1 into photons induced by loop diagrams also provides stringent constraints on parameter space, especially when χ 0 1 is wino-like and the annihilation cross section is enhanced. However, we do not reconstruct the likelihood for the Fermi-LAT photon line experiment but simply take the published limit at 5 GeV < m χ 0 1 < 300 GeV. In addition, we adopted the Isothermal profile since it is known to be more conservative than NFW or Einasto profile [18].

Numerical Analysis
In this section, after describing the input parameters over which we perform the scan of the MSSM, we present the results of our numerical study.
To compute the DM observables such as the relic abundance Ω χ h 2 , DM-proton elastic scattering cross section σ SI p , annihilation cross section σv at the present time, and branching ratios of DM annihilation, we calculate couplings and mass spectra at the neutralinomass scale M χ ≡ √ µ × M 2 , where µ and M 2 denote the values at the scale M χ . First, we solve the RGEs from M S to M A with those given in Appendix A. For the evolution from M A to M χ , which is required when M χ < M A , we employ the split SUSY RGE code † . Then we generate the SLHA output and feed it into DarkSUSY 5.1.1 [34] to compute the DM observables. Finally, we use the DM annihilation information from DarkSUSY 5.1.1 to compute the likelihoods for direct and indirect detections by following the method developed in Ref. [14].
We perform the MSSM parameter space scan, including nuisance parameters, by use of MultiNest v2.18 [35] taking 15, 000 living points with a stop tolerance factor of 0.01 and an enlargement factor of 0.8.

Input Parameters
In this subsection, we provide detailed description of our MSSM input parameters and the nuisance parameters. For the SM input parameters we take the PDG values [10].

Range
Prior distribution bino mass 10 −2 < |M 1 |/TeV < 5 Log wino mass 9 × 10 −2 < |M 2 |/TeV < 5 Log µ 9 × 10 −2 < µ/TeV < 5 Log gluino mass  In Table 2, the input parameters, their prior ranges and types of prior distributions are shown. We take |M 1,2 | , |µ| < 5 TeV because it is hard to satisfy the relic abundance constraint with the LSP heavier than 3 − 4 TeV. We apply the same maximum value for the gluino mass parameter, which does not affect our results much. The smallest values of |M 2 | and µ are chosen by taking into account the LEP limit on the chargino mass. We are taking |M 3 | > 1 TeV because of the LHC limit on the gluino mass. We cover the range of tan β up to 62 and fix the trilinear parameter A 0 = µ cot β assuming the no-mixing scenario in the stop sector. The MSSM input parameters M 1,2,3 , µ, and A 0 are given at the scale M S while tan β is the value at the scale M A .
Note that, in this work, we are using m h as an input nuisance parameter and, accordingly, the value of the high SUSY scale M S is an output. Numerically, we solve the RGEs to find the value of M S which gives the input value of m h . The Higgs boson mass measurements in the diphoton decay channel now give m h = 125.4 ± 0.4 GeV (ATLAS) [40] and m h = 124.70 ± 0.31 (stat) ± 0.15 (syst) GeV (CMS) [41]. On the other hand, the theoretical error of Higgs mass is estimated to be around 2 − 3 GeV [42] which is much larger than the experimental errors of ∼ 0.4 GeV. Therefore, in this work, we are taking m h = 125.1 GeV with a Gaussian experimental uncertainty of σ = 2 GeV.
Depending on the relative size of M A to M S , we are taking two scenarios: In the scenario A, we are taking the maximum value of 10 TeV for M A , because the A/H-funnel (M A ∼ 2 m χ 0 1 ) mechanism becomes ineffective for neutralino annihilation when M A is beyond 10 TeV. Smaller values of M A may help to obtain the correct Higgs-boson mass when M S is too large to give m h ∼ 125 GeV in the original split SUSY framework. On the other hand, the choice of M A in scenario B is the same as in the original split SUSY framework. We note that the scenario B is a part of scenario A if M S < 10 TeV.
We further need inputs for the pion-nucleon sigma term Σ πN and the strange quark content of the nucleon f T s . To account for the systematic uncertainties involved in the evaluation of the relevant nuclear matrix elements, we also treat them as nuisance parameters, as mentioned before. The central values and errors are obtained by recent lattice QCD calculations.

Numerical Results
We are taking both the PL and MP methods and make comparisons where it is informative. We note that, when we present our result based on the MP method, the systematic uncertainties of the input parameters are automatically included by utilizing a Gaussian prior distribution, see the nuisance parameters in Table 2. On the other hand, when we are using the PL method, the systematic uncertainties are added to the likelihood function.
In Fig. 1 we show the scatter plot on the (M A , M S ) plane by varying input parameters as in Table 2, while requiring m h to be in the 2-σ range: 121.1 GeV < m h < 129.1 GeV. Different colors represent different tan β ranges. We observe that a larger M S is required for small values of tan β and also as M A decreases. When tan β > ∼ 10, M S becomes almost independent of M A and it lies between ∼ 3 TeV and ∼ 15 TeV. When M A = M S is taken as in the scenario B, the value of M S is smaller in order to achieve m h ∼ 125 GeV. Therefore, in the split-SUSY framework with the intermediate Higgses lighter than ∼ 10 TeV, M S is generally predicted to be higher especially when tan β is small.
In Fig. 2 we present the probability density functions (PDFs) for marginalized posterior and profiled likelihood in the (|M 2 |/µ, |M 1 |/µ) plane. All the experimental constraints in Table 1 are applied and we make comparisons of the scenarios A (left) and B (right). We represent the bino-like, wino-like, higgsino-like and mixed neutralinos in red, blue, green and gray, respectively. Precisely, we identify the lightest neutralino χ 0 1 as bino-, wino-or higgsino-like when the corresponding fraction g b > 0.9, g W > 0.9 or g h > 0.9, respectively. ‡ Otherwise we identify it is the mixed lightest neutralino. Comparing the scenarios A and B, we can see that the difference lies in the bino region. This is because the bino-like χ 0 1 ‡ The parameters g b,W,h are defined as g b = Z 2 bino , g W = Z 2 wino , and g h = Z 2 Hu + Z 2 H d when χ 0 1 is decomposed into bino, wino, and higgsinos as follows  Table 2 while requiring m h to be in the 2-σ range: 121.1 < m h / GeV < 129.1. The color scheme are: 2 < tan β < 3 (red circle), 3 < tan β < 5 (blue square), 5 < tan β < 10 (green triangle), and tan β < 10 (gray cross). In the pink region, M A > M S which is out of our current consideration.
can satisfy the relic abundance constraint only through Z/h-resonance in the scenario B, where A/H-funnel does not work because M A = M S > ∼ 3 TeV. In fact, the mechanism of Z/h-resonance requires a small fraction of higgsino but it cannot be too large because of the constraint from the Fermi dSphs gamma ray measurement. In particular, we find that the higgsino composition is between 0.06 to 0.1 in the h resonance region which leads to the ratio |M 1 |/µ ∼ 0.4.
Furthermore, we find that the chargino-neutralino coannihilation working in reducing the relic abundance in both scenarios. Being different from the original split SUSY framework (scenario B), one can obtain the correct relic abundance in scenario A without resorting to the coannihilation mechanism thanks to the intermediate Higgses A and H. To  address this point, we show in Fig. 3 the points with δχ 2 < 5.99 on the (m χ 0 1 , m χ ± ) plane for the scenario A (left) and B (right). In addition to the Z/h-resonance regions around m χ 0 1 ∼ 50 , 60 GeV and the chargino-neutralino coannihilation region along the m χ 0 1 = m χ ± line, we observe there are more points appearing in the scenarios A (left panel) due to the A/H-funnel. We find that the A/H-funnel region disappears when m χ 0 1 > 1 TeV, because the decay widths of A and H become too large and the Breit-Wigner resonance effect is not strong enough to reduce the relic abundance when M H,A > ∼ 2 TeV. In Fig. 4, we show the marginalized 2D posterior 2-and 3-σ credible regions (CRs) for the scenario A (left) and B (right) in the (m χ 0 1 , σv ) plane. We also show the PL 2-σ region (scattered points) for the bino-like (red) and mixed (gray) χ 0 1 in the upper frames and the wino-like (blue) and higgsino-like (gray) χ 0 1 in the lower frames. Here σv denotes the annihilation cross section at the present time which is relevant to the DM indirect detections and through which one may easily identify different mechanisms for the relic abundance.
When m χ 0 1 < 100 GeV, via the Z/h resonances, the marginalized posterior CRs are located at the bino-like neutralino region with a small amount of higgsino component in both scenarios (see the upper frames). Although the Z/h-resonance channels have very good likelihoods, they only fall into the 3σ (99.73%) CR owing to the small prior volume effect. The similar effect happens for the bino-like χ 0 1 when m χ 0 1 > 100 GeV and the correct relic abundance is obtained by the χχ ± coannihilation. The fact that more parameter space survives in the scenario A (left) than scenario B (right) is due to the A/H-funnel. Nevertheless, most of the additional parameter space is a result of the mixture mechanism between A/H-funnel and coannihilation. In the lower frames, we observe that the 2σ CR has the wino-like branch (blue) with the higher σv than the higgsino-like one (green). For the wino-like branch, the relic abundance is mainly reduced by the wino-like DM annihilation into W + W − pairs. However when m χ 0 1 > ∼ 3 TeV, the wino DM cannot give the correct relic abundance as is well known. This mass limit can be slightly extended if coannihilation is taken into account. Since the wino DM have higher annihilation cross sections, the indirect detection constraint is stringent. Indeed, the lower bound for the wino-like neutralino mass is about 300 GeV from the Fermi dSphs gamma ray constraints. Incidentally, the lower bound for the higgsino-like neutralino mass is about 100 GeV, set by the LEP limit of m χ ± 1 > 103.5. We further see there is no particular lower bound for the bino-like or mixed neutralino, as seen from the upper frames.
Finally, in Fig. 5 we show the marginalized 2D posterior PDF 2σ and 3σ contours in the (m χ 0 1 , σ SI p ) plane. The red solid line denotes the recent LUX result, the black dashed line the XENON1T projected sensitivity, and the blue dash-dotted line the LZ projected sensitivity [43]. The orange dashed line represents the approximate line below which the DM signal becomes hardly distinguishable from the signals from the coherent scattering of the 8 B solar neutrinos, atmospheric neutrinos and diffuse supernova neutrinos with nuclei. We observe that a part of 2-σ CR is below the LZ projected sensitivity. We can see that, , σv ) plane for the scenario A (left) and B (right). The inner (outer) contours bounded the 2(3)-σ CR. All the scatter points superimposing on the contours agree with likelihood in the criteria δχ 2 < 5.99. The red dots, blue squares, green stars, and gray triangle are for the bino-like, wino-like, higgsino-like, and mixed neutralino, respectively.
in the 2-σ CRs, there is no significant difference between the scenarios A and B. The 3-σ CRs are slightly different in the lower σ SI p region. Moreover, in both scenarios, the future 7-tons experiments, LZ, can set a lower limit on the neutralino DM at m χ 0 1 > 100 GeV.

Discussion
In this work, we have studied a "modified split SUSY" scenario, characterized by two separate scales -the SUSY-breaking scale M S and the heavy Higgs-boson mass scale M A . This is different from the split SUSY scenario, in which the scale M A is also set at M S . The current scenario is motivated by (i) the absence of direct SUSY signals from the searches of scalar quarks up to a few TeV, (ii) the observed Higgs boson is somewhat on the heavy side which needs a large radiative correction to the tree-level mass from heavy stops, and (iii) absence of signals from heavy Higgs bosons A/H and H ± which can be as light as a few hundred GeV. Therefore, the choice of M A need not be as large as M S . We have studied two scenarios: (i) M A ≤ min(M S , 10 Tev) and (ii) M A = M S (the same as split SUSY).
Because of two distinct scales M S and M A the running of the soft parameters and couplings are separated in two steps. We start with the set of RGEs given in appendix A to run from M S down to M A and perform the matching at the scale M A . Then run from M A down to the electro-weakino scale M χ ≡ √ µ × M 2 with the set of RGEs of split SUSY.
Because of this two-step RGEs the predictions for DM observables and the Higgs boson mass are more reliable than just a single-step RGE.
We have scanned the MSSM parameter space characterized by the two scales: M S and M A subjected to many existing experimental constraints: invisible widths of the Z boson and the Higgs boson, the chargino mass limit, relic abundance of the LSP, spinindependent cross sections from direct detection, and the gamma-ray data from indirect detection. We found interesting survival regions of parameter space with features of either chargino-neutralino coannihilation, the A/H funnel, or wino-like. These regions survive because of the large enough annihilation to reduce the relic abundance to the observed values, as well as give a large enough Higgs boson mass to fit to the observed value. Finally, the survived parameter space can be further scrutinized by near future direct detection experiments such as XENON1T and LZ.
We offer a few important comments as follows.
1. We used the Higgs boson mass in the range range 121.1 < m h < 129.1 GeV to search for suitable M S . Since m h is on the rather heavy side, it requires a large radiative correction to the tree-level mass. This can be achieved by a large stop mass and/or large mixing in the stop sector. Since the radiative correction is proportional to some powers of tan β, a smaller tan β requires then a larger M S in order to achieve a large enough m h . Typically, M S > ∼ 10 5−6 GeV for tan β < 3. For large enough tan β the values of M S is more or less independent of M A .
2. On the other hand, if we set M A = M S as we do in scenario A, the allowed M S is rather short from about 10 3 − 10 4 GeV with large tan β (see Fig. 1).
3. An interesting region that satisfies the relic abundance constraint is characterized by nearly degenerate mass among the first two neutralinos and the lightest chargino, indicated by M 2 /µ ≈ M 1 /µ ≈ 1. The increased effective annihilation cross section can help reducing the relic abundance.
4. Another interesting region is the Z/h resonance region (m χ 0 1 ∼ 50 − 60 GeV), though it is relatively fine-tuned region because of the narrow width of the Z boson and the Higgs boson.
5. Yet, another interesting survival region is the A/H funnel region. If m χ 0 1 falls around the vicinity of m A/H /2 the resonance effect is strong, provided that the width is not too large. This can be achieved for m χ 0 In scenario B, where M S = M A , large values of M S then cannot be accepted because the A/H funnel is not working efficiently. However, in scenario A, where M A < M S , the A/H funnel can be very effective in reducing the relic abundance, thus more parameter space is allowed.
6. Both wino-like and higgsino-like LSPs have large annihilation cross sections. The allowed mass for m χ 0 1 ranges from about 300 GeV to 3 TeV for wino-like LSP while from about 100 GeV to 2 TeV for higgsino-like LSP.
7. The current allowed parameter space has a large region below the current LUX limit σ SI p < ∼ 10 −9 pb. Although the future XENON1T can improve the limit by an order of magnitude, there is still a sizable region below the XENON1T sensitivity. Yet, there still exist some allowed regions even with the future 7-tons size direct detection experiment LZ. Therefore, this modified split SUSY scenario is hard to be excluded in the future. 8. We have used both the methods of profile likelihood and marginal posterior. Though these two statistical approaches have very different methodology, the resulting 2-and 3-σ regions are quite consistent, as shown in the figures.

B The statistical framework
To calculate the probability of MSSM parameter given the experimental data, one can employ Bayes's' theorem to compute the posterior probability density function, p(θ, φ|d) = L(d|θ, φ)π(θ, φ) Z(d) . (B.1) Here, we denote the MSSM parameters and DM direct detection nuisance parameters as θ and φ, respectively. The likelihood L(d|θ, φ) is the probability of obtaining experimental data for observables given the MSSM parameters. The prior knowledge of MSSM parameter space is presented as prior distribution π(θ, φ). Our MSSM prior ranges and distributions are tabulated in Table 2. Finally, the evidence of the model in the denominator can be merely a normalization factor, because we are not interested in model comparison.
The Bayesian approach allows us to simply get ride of the unwanted parameters by using marginalization. For example, if there would be n free model parameters, r i=1,...,n , but one is only interesting in the two-dimensional figure (r 1 , r 2 ), the marginalization can be written as p(r 1 , r 2 |d) = p(r 1 , ..., r n |d) An analogous procedure can be performed with the observables. One should keep in mind that a poor prior knowledge or likelihood function can raise a volume effect. In other words, some regions gain more weight from higher prior probability but fine-tuning regions such as resonance regions for relic abundance likelihood only have lower prior probability. Although this is the feature of Bayesian statistics, in order to manifest these fine-tuning regions, we still present both profile likelihood and marginal posterior method at the same time.
In Bayesian statistics, a credible region (CR) is the smallest region, R, in the best agreement with experiments bounded with the fraction ̺ of the total probabilities. For example at MSSM (M 1 , M 2 ) plane, the ̺ credible region can be written as where the normalization in the denominator is the total probability with R → ∞. In this paper, we have shown ̺ = 0.95 and ̺ = 0.9973 corresponding to 2σ and 3σ credible region. As the comparison, we also present the scatter points with selected criteria δχ 2 = −2 ln L/L max ≤ 5.99. This criteria is 2σ confidence region of Profile Likelihood method in 2 degrees of freedom. We can see from our result that most of 2σ confidence region of PL method is similar to the 3 σ credible region in MP method. We would like to note that the total profile likelihood here takes the likelihoods including the nuisance parameters distribution, which is the prior distribution in marginal posterior method.