Interpreting IceCube 6-year HESE data as an evidence for hundred TeV decaying Dark Matter

The assumption of a single astrophysical power-law flux to explain the IceCube 6-year HESE extraterrestrial events yields a large spectral index that is in tension with gamma-ray observations and the 6-year up-going muon neutrinos data. Adopting a spectral index belonging to the range $\left[2.0,2.2\right]$, which is compatible with the one deduced by the analysis performed on the 6-year up-going muon neutrinos data and with $p$-$p$ astrophysical sources, the latest IceCube data show an up to $2.6\,\sigma$ excess in the number of events in the energy range 40--200 TeV. We interpret such an excess as a decaying Dark Matter signal and we perform a likelihood-ratio statistical test to compare the two-component scenario with respect to the single-component one.

The assumption of a single astrophysical power-law flux to explain the IceCube 6-year HESE extraterrestrial events yields a large spectral index that is in tension with gamma-ray observations and the 6-year up-going muon neutrinos data. Adopting a spectral index belonging to the range [2.0, 2.2], which is compatible with the one deduced by the analysis performed on the 6-year upgoing muon neutrinos data and with p-p astrophysical sources, the latest IceCube data show an up to 2.6 σ excess in the number of events in the energy range  TeV. We interpret such an excess as a decaying Dark Matter signal and we perform a likelihood-ratio statistical test to compare the two-component scenario with respect to the single-component one.
In the last few years IceCube (IC) Neutrino Telescope has been collecting neutrino events in the TeV -PeV range [1][2][3][4][5][6] that are providing striking evidences for extraterrestrial neutrinos even though their origin still remains unknown. According to the IC veto implementation, the IC data are mainly divided in two samples: "Medium Energy Starting Events" (MESE) with an energy threshold of 1 TeV and "High Energy Starting Events" (HESE) for neutrino energies larger than 20 TeV. In both data sets, the neutrino interaction vertex is located inside the detector. The analysis of such events has triggered an interesting discussion in the scientific community about their possible astrophysical origin, as well as an intriguing connection with Dark Matter (DM). In particular, IC events have been related to the decay or annihilation of very heavy DM particles . A pure astrophysical explanation of neutrino events has been also extensively discussed in the literature where different sources have been proposed. In particular, among such sources one can quote the stellar remnants in starforming galaxies [38,39], extragalactic Supernovae and Hypernovae remnants [40], active galactic nuclei [41][42][43], and gamma-ray bursts [44].
At the ICRC 2017 Conference the IceCube Collaboration has released the 6-year HESE data [45], which consist of a sample of 82 events (track + shower). The analysis performed on the neutrinos with a deposited energy above 60 TeV provides a best-fit power-law E −γ with spectral index 1 equal to γ 6yr = 2.92 +0. 29 −0.33 . Note that such a best-fit spectral index is larger than the one obtained by the 4-year HESE data (γ 4yr = 2.58 ± 0.25) [5]. This is due to the fact that, during the last two years, IC experiment has not observed neutrinos with deposited energy larger than about 200 TeV.
According to multi-messenger analyses 2 , one would expect a spectral index γ < ∼ 2.2 in case of hadronic (pp) astrophysical sources [38,39,48,49]. On the other hand, photo-hadronic (p-γ) sources, that would produce a softer power-law neutrino flux, are strongly constrained by searches of spatial and temporal correlations with gamma-ray observations [50,51]. In this framework, if one takes a conservative approach, the ambiguity (p-p vs p-γ) that concerns the nature of dominant neutrino sources does not allow to fix a reliable bound on the spectral index.
Nevertheless, the analysis restricted to the 6-year upgoing muon neutrinos data requires a spectral index γ = 2.13 ± 0.13 as best-fit [6]. Therefore, assuming isotropy (extragalactic sources) and an equal neutrino flavor ratio at the Earth, one would reasonably expect such a hard power-law for the whole diffuse TeV-PeV neutrino flux.
In Figure 1 we show the residual in number of neutrino events with respect to the sum of the conventional atmospheric background and a power-law with spectral index 2.0 (upper panel) and 2.2 (lower panel). The power-law normalization has been obtained by a maximum likelihood procedure in both cases. We highlight the presence of an excess at low energy, in the energy bins between 40 and 200 TeV. In case of spectral index 2.0, the maximum local statistical significance of such an excess is 2.6 σ (6-year HESE), hence slightly larger than the one of ∼ 2 σ obtained in the previous analysis concerning 4-year HESE [26]. The maximum local statistical significance decreases to 2.1 σ in case of a power-law with spectral index 2.2. We note that in case of 2-year MESE data sample [4] the discrepancy is 2.3 σ (1.9 σ) for γ = 2.0 (γ = 2.2), as shown in Ref. [32]. Hence, independent of the spectral index adopted, the statistical evidence for an excess slightly improves once the 6-year HESE data are used instead of 4-year HESE or 2-year MESE.
In this Letter, we focus on the interpretation of the low-energy excess in terms of two components, an Astrophysical (Astro) power-law and a decaying Dark Matter (DM) signal. 3 Therefore, in addition to the conventional atmospheric background [61], the whole extraterrestrial differential neutrino flux is given by (1) The quantity φ Astro 0 is the normalization of the astrophysical neutrino flux, whereas the quantity φ DM 0 is the normalization of the DM neutrino flux, i.e. the inverse life-time 1/τ DM in case of decaying DM with mass m DM . For the sake of brevity we omit all the details of the above definitions. An extensive discussion of such a parametrization can be found in Ref. [32]. Note that the prompt atmospheric background (neutrinos produced by the decays of charmed mesons) [62] is here considered negligible, in agreement with the IceCube results contained in Ref. [2,4,6].
As shown in Eq. (1) the whole extraterrestrial neutrino flux depends on four free parameters, namely φ Astro 0 , γ, φ DM 0 and m DM . Hence, a statistical analysis to study the relevance of the additional DM component should consider the whole set of parameters. In this letter we follow the same approach of Ref. [32] that is strongly based on the prior 2.0 ≤ γ ≤ 2.2 previously discussed. Clearly, relaxing such a constraint the evidence for an additional component becomes almost statistically irrelevant. In particular, one would obtain the best-fit γ 6yr = 2.92 +0.29 −0.33 that is however in tension with the 6-year up-going muon neutrinos data. In the following analysis we consider two fixed values γ = 2.0, 2.2 that corresponding to the extreme cases give an idea of the dependence of the results on the spectral index. To further simplify the analysis, the remaining free parameter for the astrophysical component, namely φ Astro 0 , can be fixed at its best-fit value. This ansatz corresponds to a conservative approach where one leaves to the additional DM component the smallest room possible.
Concerning the DM component, we only consider two different decay channels of DM particles (hereafter denoted as χ), which can be seen as benchmarks. In particular, we focus on the hadronic decay channel χ → tt, and on the leptonic one χ → τ + τ − . Since the neutrino spectra produced in channels with hadronic or leptonic final-states have very different shape, considering these two decay channels practically covers all the possible phenomenological scenarios. The neutrino energy spectra provided by DM decays at the production are obtained by using the tables of Ref. [63]. Moreover, we consider the Navarro-Frenk-White (NFW) distribution [64] for the DM halo density profile of the Milky Way and the ΛCDM parameters according to the Planck analysis [65]. In the present analysis, in order to account for the neutrino oscillations we adopt the same approach used in Ref. [32].
In order to quantify how the fit is statistically improved by adding a DM neutrino component of a given mass m DM on top of an astrophysical power-law, we remake the analysis given in Ref. [32] with the latest 6-year HESE data [45]. In particular, we perform a likelihood-ratio statistical analysis on the neutrino energy spectrum, where the Test Statistics (TS) is defined as Following the above discussion, the Test Statistics is eval- to its best-fit value, and by scanning over the only remaining free parameter, namely the normalization of the DM flux. According to Wilks [66] and Chernoff [67] theorems, the TS follows the distribution 1 2 δ (TS) + 1 2 χ 2 (TS). Therefore, we evaluate the preference of the data for a twocomponent flux with respect to the single astrophysical power-law, in number of standard deviations σ. The likelihood function L adopted in Eq. (2) is a binned multi-Poisson likelihood [68], whose expression is equal to where the quantity n i is the observed number of neutrinos in the energy bin i, whereas N i is the expected number (1) with the IceCube effective area [1] and considering an exposure time of 2078 days (for more details see Ref. [32]).
The main results of the present analysis are presented in Figures 2 and 3. The plots display the number of standard deviations σ in the m DM -τ DM plane for the decay channels considered, namely χ → tt and χ → τ + τ − . In particular, the darker the color, the larger the significance in σ of the DM neutrino component. The upper and lower panels of both Figures refer to an astrophysical power-law with spectral index 2.0 and 2.2, respectively. In the plots, the best-fit values (maximum significance) is represented by white stars (the capital letter "H" refers to 6-year HESE analysis) and they are compared to the previous results of Ref. [32] represented here with black dots (the capital letter "M" refers to 2-year MESE analysis). The white solid (dashed) contours enclose the regions in the m DM -τ DM plane where the statistical significance is larger than 3 σ (2 σ). As can be seen from the plots, the maximum value of √ TS depends on the spectral index only, while it is almost independent of the decay channel considered. In particular, the statistical significance at the best-fit is 3.75 σ and 2.60 σ in case of spectral index 2.0 and 2.2, respectively.
Moreover, the present constraints on decaying DM models from IceCube observations are presented by the red lines. In particular, the solid red lines bound from above the regions that are excluded by the 6-year HESE data, while the dashed ones correspond to the same limit deduced by the 2-year MESE data. It is worth noting that the 6-year HESE data bound the possible DM models in a region with m DM ≥ 100 TeV. This feature depends on two effects. On one side the different energy thresholds for HESE data set (20 TeV) and MESE sample (1 TeV) provide different sensitivity of data and hence of TS for light χ. On the other side as can be seen from Figure 1, the second energy bin from left, corresponding to almost 25-40 TeV, shows a defect in the number of events and thus it disfavors any additional second component contributing to this energy. This pushes the possible DM models to higher masses.
Furthermore, the almost horizontal black lines, instead, correspond to the gamma-ray constraints on DM models deduced by the Fermi-LAT measurements of the isotropic diffuse gamma-ray background (IGRB) spectrum [69]. Such limits have been obtained by considering the total electromagnetic energy density ω γ of the IGRB integrated in the energy range from 0.1 GeV to the maximum energy corresponding to m DM /2. We report the constraints for DM contribution ω DM γ equal to 1%, 10% and 100% of ω γ , respectively. In particular, the solid black lines (ω DM γ = ω exp γ ) bound from below the allowed region in the m DM -τ DM plane. However, since it is quite reasonable to assume that the majority of the IGRB spectrum is accounted for by standard astrophysical sources, we consider the limit ω DM γ ≤ 0.1 ω exp γ as a realistic constraint for the DM contribution to the gamma-ray flux. Therefore, only the regions above the short dashed black lines correspond to viable choices of parameters (m DM , τ DM ) for DM models compatible with both neutrino and gamma-ray observations. Such multi-messenger constraints are affected by an uncertainty of about 20%, as discussed in Ref. [32] (see also Ref.s [18,[70][71][72] for different analyses about gamma-rays constraints on DM models).
We note that hadronic channel requires smaller values for the lifetime τ DM and larger DM masses m DM with respect to the channel with leptons in the finalstates in order to account for the IceCube data. This implies that the DM models with quarks as final-states are more in tension with Fermi-LAT data with respect to the models involving leptons, confirming the previous results discussed in Ref. [32]. In particular, in case of a DM decaying into top quarks the best-fit values are m DM 500 TeV and τ DM 2.77 × 10 27 sec, while for the leptonic case we have m DM 400 TeV and τ DM 1.65 × 10 28 sec. In Figure 4 we report the predic-atm. ν atm. μ astro. γ=2.0 DM χ → τ + τ - tions of the two-component neutrino flux for the latter case.
In summary, in this Letter we have focused on the latest IceCube 6-year HESE data that improve the statistics with two more years of data-taking. The IceCube analysis performed on this data sample, with the assumption of a single astrophysical power-law flux, yields a bestfit spectral index γ 6yr = 2.92. This large value for the spectral index (even larger with respect to γ 4yr = 2.58 of IceCube 4-year HESE) is in tension with the theoretical expectation of a hard power-law behavior for the neutrino flux (p-p astrophysical sources) and with the analysis performed on the 6-year up-going muon neutrinos data. Hence, it suggests the presence of a second component at low energy (E ν ≤ 200 TeV) added to a hard power-law flux. Here, following the same approach of Ref. [32], we have scrutinized the interpretation of such an additional component to diffuse neutrino flux in terms of a decaying DM signal.
The 6-year HESE data show an excess between 40 and 200 TeV, with a maximum local significance of 2.6 σ, improving the previous estimate of 2.3 σ corresponding to 2-year MESE [32]. The likelihood-ratio statistical test shows that the statistical improvement of the fit with the two-component flux reaches 3.75 σ (2.60 σ) at the best-fit for the case of γ = 2.0 (γ = 2.2), confirming the previous results on 2-year MESE data [32]. However, differently from 2-year MESE, the best-fit for 6-year HESE seem to prefer larger values for m DM even though this could be partially explained in terms of a lower sensitivity of such data set to low energy events.