A model for neutrino anomalies and IceCube data

We interpret the neutrino anomalies in neutrino oscillation experiments and the high energy neutrino events at IceCube in terms of neutrino oscillations in an extension of the standard model where three sterile neutrinos are introduced so as to make two light neutrinos to be Pseudo-Dirac particles and a light neutrino to be a Majorana particle. Our model is different from the so-called 3 + n model with n sterile neutrinos suggested to interpret short baseline anomalies in terms of neutrino oscillations. While the Pontecorvo-Maki-Nakagawa-Sakata (PMNS) matrix in 3 + n model is simply extended to n × n unitary matrix, the neutrino mixing matrix in our model is parameterized so as to keep the 3 × 3 PMNS mixing matrix for three active neutrinos unitary. There are also no flavor changing neutral current interactions leading to the conversion of active neutrinos to sterile ones or vice versa. We derive new forms of neutrino oscillation probabilities containing the new interference between the active and sterile neutrinos which are characterized by additional new parameters ∆m2 and θ. Based on the new formulae derived, we show how the short baseline neutrino anomalies can be explained in terms of oscillations, and study the implication of the high energy neutrino events detected at IceCube on the probe of pseudo-Dirac neutrinos. New phenomenological effects attributed to the existence of the sterile neutrinos are discussed.


Introduction
The observation of neutrino oscillations in the atmospheric neutrino 1 Super-Kamiokande [3], the solar neutrino 2 SNO [5], and the reactor neutrino Daya Bay [6] and RENO [7] 3 experiments is one of big discoveries in particle physics since 90's. It implies that neutrinos are massive particles and that the three flavor neutrinos ν e , ν µ , ν τ are mixtures of neutrinos with definite masses ν i (with i = 1, 2, . . .). Although neutrino oscillations among three active neutrino flavors have been confirmed by the analysis based on the experiments mentioned above, there exist several anomalies which are unexpected results coming from short baseline (SBL) experiments such as the reactor antineutrino anomaly [9][10][11], the Gallium solar anomaly [12][13][14][15][16][17][18], and the Liquid Scintillator Neutrino Detector (LSND) anomaly [19,20] (including MiniBooNE anomaly [21]). To resolve those neutrino anomalies in terms of neutrino oscillation, it is required to introduce at least one additional squaredmass difference, ∆m 2 SBL , which is much larger than ∆m 2 Sol and ∆m 2 Atm [9][10][11]. This result suggests indication in favor of the possible existence of eV-mass sterile neutrino.

JHEP12(2019)133
Apart from the anomalous results in accelerator and reactor based neutrino experiments favoring the existence of light sterile neutrino, IceCube experiments [22][23][24] announced the observation of vey high energy neutrino events. The study for the track-to shower ratio of the subset with energy above 60 TeV coming from IceCube has shown that the events are consistent with the hypothesis that cosmic neutrinos have been seen even though their origin and propagation are still elusive [25,26]. In order to examine them, we need to discriminate the flavor composition of cosmic neutrinos, which is possible by looking at the topology of the events. But, the current limited statistics does not allow yet to discriminate the initial flavor. It is also widely perceived that the IceCube may serve as astronomical-scale baseline experiment to uncover the oscillation effects due to very tiny mass splitting of the pseudo-Dirac neutrinos. If the oscillation effects induced by pseudo-Dirac neutrinos with very high energy and long trajectory are prominent, then they may affect the observables such as the neutrino flavor composition detected from the ultra-high energy neutrino experiments. Since the current precision on the observation of very high energy neutrino events does not exclude such an oscillation effect, it would be meaningful to confront any model for pseudo-Dirac neutrinos realized by introducing sterile neutrinos with the high energy astrophysical neutrino events.
In this work, we construct a model having three sterile neutrinos which make two light neutrinos to be pseudo-Dirac particles and a light neutrino to be Majorana particle. Then, we investigate if the sterile neutrino that makes a light neutrino to be Majorana particle can play a crucial role in resolving the neutrino anomalies from SBL experiments, and study the implication of the very high energy astrophysical neutrino data from IceCube on the probe of the pseudo-Dirac neutrinos. Thus, the goal of this paper is to interpret both SBL neutrino anomalies and astronomical neutrino data observed at IceCube in terms of neutrino oscillations [27][28][29][30][31] in the context of the model we construct.
Our model is different from the so-called 3+ n model with n sterile neutrinos suggested to interpret SBL anomalies in terms of neutrino oscillations. While the Pontecorvo-Maki-Nakagawa-Sakata (PMNS) matrix in 3 + n model is simply extended to n × n unitary matrix as in refs. [32,33,94], the neutrino mixing matrix in our model is parameterized so as to keep the 3 × 3 PMNS mixing matrix for three active neutrinos unitary. In this model, there are no flavor changing neutral current interactions leading to the conversion of active neutrinos to sterile ones or vice versa. We will present new forms of neutrino oscillation probabilities modified by introducing new sterile neutrinos. The interference between active flavor and sterile neutrinos due to new additional oscillation parameters ∆m 2 and mixing angle θ triggers new oscillation effects which can be responsible for the explanation of SBL neutrino anomalies and very high energy neutrino events at IceCube. Based on the new formulae for neutrino oscillations, we will discuss how SBL neutrino anomalies can be explained or alleviated and consistently accommodate the recent IceCube high energy neutrino events. Constraints on the oscillation parameters coming from solar and atmospheric neutrino data, cosmological observation for the sum of active neutrino masses, effective neutrino mass in β-decay and neutrinoless-double-beta (0νββ)-decay experiments will be discussed. Our study based on the terrestrial, atmospheric and solar experiments is similar to that in the so-called 3 + 1 model in the light that the two very tiny pseudo-Dirac mass JHEP12(2019)133 splittings are not relevant for those experiments. However, this work includes the study for the implication of IceCube data on the probe of the oscillarion effects induced by those two tiny mass splittings while keeping the results for the explanations of SBL neutrino anomalies, atmopsheric and solar neutrino experiments in terms of neutrino oscillations including an eV scale sterile neutrino.
This work is organized as follows. In section 2, we discuss our model and study neutrino masses and mixings in the new framework. In section 3, we study how the new parameters could be constrained through the cosmological data (the sum of active neutrino masses), and the effective neutrino masses in β-decay and 0νββ-decay experiments. In section 4 we develop the new active neutrino oscillation probabilities modified by incorporating new sterile neutrinos. And we study how SBL neutrino anomalies can be explained. In section 5, we investigate how new effects due to sterile neutrino in the solar and atmospheric oscillations can be constrained by the SBL ν e (ν e ) disappearance and ν µ (ν µ ) disappearance channels, and examine whether the high energy neutrino events from IceCube data can be interpreted in terms of neutrino oscillation. In section 6 we examine astronomical neutrino data observed at IceCube to uncover the oscillation effects of tiny mass splittings. In section 7, we state conclusions by summarizing this work.

Masses and mixings
Introducing right-handed singlet neutrinos N R , extra neutrino singlet fermions S, and an SU(2) L singlet scalar field Ψ. we construct the renormalizable Lagrangian given in the charged lepton basis as [27] where N R , S are three generations, L stand for SU(2) L left-handed lepton doublet, Φ = (φ + , φ 0 ) T is the SM Higgs doublet andΦ ≡ iτ 2 Φ * . M R and µ are Majorana masses for the N R and S fields, respectively. The above Lagrangian is invariant under U(1) B−L when µ, M R = 0 by assigning quantum numbers L : 1, N R , S : 1, Ψ : −2, and Φ : 0 under the U(1) L (or U(1) B−L ) symmetry. Then, the parameters µ and M R reflect soft symmetry breaking of U(1) L . Nontrivial vacuum expectation value (VEV) of the scalar field Ψ does not break the electroweak symmetry, but spontaneously breaks the U(1) L (or U(1) B−L ) symmetry. Thus the symmetry breaking scale for Ψ can be different from the electroweak scale. Integrating out the heavy Majorana neutrinos in the Lagrangian eq. (2.1), we obtain effective Lagrangian for neutrino sectors given by, where Y D , Y S , Y DS , M R and µ are all 3 × 3 matrices. After the scalar fields Φ and Ψ get VEVs and taking S to be right-handed, the Lagrangian for neutrinos in the charged lepton basis reads

JHEP12(2019)133
where g is the SU (2) coupling constant, θ W is the Weinberg angle, = (e, µ, τ ), ν L = (ν e , ν µ , ν τ ), and S R = (S 1 , S 2 , . . . S n ). The light neutral fermions S α do not take part in the standard weak interaction and thus are not excluded by LEP results, while the number of active neutrinos coupled with the W ± and Z bosons is N ν = 2.984 ± 0.008 [34]. After electroweak symmetry breaking, eq. (2.3) describes 3 + n Majorana neutrinos. In the case of n = 3 sterile neutrinos, the 6 × 6 Majorana neutrino mass matrix is where the 3 × 3 mass matrices M D , M L , and M S are those for Dirac masses, left-and right-handed Majorana masses, respectively, given by Here we take M R m S m D µ, and neutrinos become pseudo-Dirac particles when M D is dominant over M L and M S in eq. (2.4), which reflects m DS (m D m S )/M R . In order to get physical parameters, we perform basis rotations from interaction eigenstates to mass eigenstates [31,35], where ξ L is the mass eigenstate of neutrino, Here, the 6 ×6 unitary neutrino transformation matrix W ν given by, In the above expression, 0 3 is the 3×3 null matrix, U R is an unknown 3×3 unitary matrix, 2 with ϕ i being arbitrary phases. The matrices V 1 and V 2 are presonsible for the maximal mixing between active neutrinos and sterile neutrinos. The angle θ k is introduced thanks to nondegeneracy between M L and M S in eq. (2.4), and thus it is responsible for the deviation of maximal mixing between active JHEP12(2019)133 neutrino and sterile neutrino, which reflects the breaking of degeneracy of a pair of neutrinos in each generations. It is easily see that maximal mixing between ν k and S c k is recovered in the limit of θ k = 0. We also note that the 3 × 3 unitary matrix U L should be the PMNS neutrino mixing matrix responsible for the mixing among three active neutrinos. Then the neutrino mass matrix M ν is diagonalized in the mass eigenstates (ν 1 , ν 2 , ν 3 , S c 1 , S c 2 , S c 3 ) basis as The expression eq. (2.9) represents that the Majorana mass matrices M L and M S , and Dirac neutrino mass matrix, M D in eq. (2.4) are diagonalized by the mixing matrices U L and To get the real and positive mass squared for the 6 neutrino mass eigenstates, we diagonalize the Hermitian matrix M ν M † ν with the help of eq. (2.9) as follows, where δ stands for δ k ≡ (M L ) k +(M † S ) k originated from the left-and right-Majorana masses (k = 1, 2, 3). From the above equation the mixing parameters θ k and φ k in eq. (2.7) can be obtained (2.11) Now, to accommodate both an eV sterile neutrino for a possible solution to the neutrino anomalies [9][10][11][19][20][21] and the high energy neutrino events in the IceCube detector [22][23][24] to be interpreted as new neutrino oscillations, simultaneously, we assume that m ν 3 m s 3 and m ν 1 ≈ m s 1 , m ν 2 ≈ m s 2 . 4 It is equivalent to take the limit of bothM j |(M L ) j | |(M S ) j | (with j = 1, 2) and |(M S ) 3 | |(M L ) 3 |, leading to 5 δ 1(2) (M L ) 1 (2) and |δ 3 | −2m 3 tan 2θ 3 (2.12) which in turn gives θ 1 ≈ θ 2 ≈ 0 and π/4 < θ 3 < π/2. Since the active neutrinos are massive and mixed, the weak eigenstates ν α (with flavor α = e, µ, τ ) produced in a weak gauge interaction are linear combinations of the mass eigenstates with definite masses. The 4 This possibility could theoretically be realized in a non-renomalizable flavor model considering non-Abelian discrete symmetry plus Abelian symmetry, e.g. refs. [31,36,37] where light active neutrino masses are mainly generated by QCD anomalous U(1) symmetry (via Froggatt-Nielson mechanism [38]) while leptonic mixing matrix is produced (through seesaw formula [39]) by non-Abelian discrete symmetry. 5 The Dirac masses of first and second generations are much larger than the left(right)-handed Majorana masses of those, while the Dirac mass of third generation is larger than the left(right)-handed Majorana masses of that for 3π/8 < θ3 < π/2 and smaller for π/4 < θ3 ≤ 3π/8. For |(ML)3| |(MS)3| is not realized due to a requirement m 2

JHEP12(2019)133
three neutrino active states emitted by weak interactions are described in terms of the mass eigenstates ξ k = (ν k S c k ) (k = 1, 2, 3) as , (2.14) in which the field redefinitions ν j → e i π 4 ν j , S c j → e −i π 4 S c j (with j = 1, 2) and ν 3 → e i π 4 ν 3 , S c 3 → e i π 4 S c 3 are used. In eq. (2.13), U is the 3 × 3 PMNS mixing matrix U PMNS which is expressed in terms of three mixing angles, θ 12 , θ 13 , θ 23 , and three CP-odd phases (one δ CP for the Dirac neutrino and two ϕ 1,2 for the Majorana neutrino) as [34]  where s ij ≡ sin θ ij , c ij ≡ cos θ ij and P ν is a diagonal phase matrix what is that particles are Majorana ones.
And their mass eigenvalues (real and positive) are given as The neutrino masses for the first and second generations lift slightly the degeneracy of masseigenvalues, and we get almost degenerate pairs of eigenstates with tiny mass differences: the mass-squared differences in each pair ∆m 2 k ≡ m 2 ν k − m 2 s k (with k = 1, 2) are so small that the same mass ordering should apply to both eigenmasses, that is, On the other hand, the mass splitting for third generation is given by leading to −1 < cos 2θ 3 < 0, that is, π/4 < θ 3 < 3π/4 due to the requirement of m 2 s 3 m 2 ν 3 . From eqs. (2.12) and (2.18) a possible range of θ 3 can be derived as

JHEP12(2019)133
in which, especially, for 3π/8 < θ 3 < π/2 the third neutrino pair could be Majorana. As is well-known, because of the observed hierarchy |∆m 2 Sol ≡ m 2 ν 2 − m 2 ν 1 > 0, and the requirement of a Mikheyev-Smirnov-Wolfenstein (MSW) resonance for solar neutrinos [40,41], there are two possible neutrino mass spectra: (i) the normal mass ordering (NO) m 2 where ∆m 2 kj ≡ m 2 ν k − m 2 ν j , for normal mass ordering (inverted mass ordering) respectively [42]. Using above eq. (2.12), we can obtain good approximated forms of m 2 ν 3 and m 2 s 3 in terms of m 3 and θ 3 given as The mass parameter m 3 can be derived from eqs. (2.12) and (2.18) as where ∆S 2 kj ≡ m 2 s k − m 2 s j . Clearly, within the range eq. (2.19) a value of θ 3 going around π/4 can realize the inverted mass hierarchy (IH), m ν 2 > m ν 1 m ν 3 , while a value around π/2 favors the degenerate normal mass ordering (DNO), m ν 3 m ν 2 m ν 1 , and degenerate inverted mass ordering (DIO), m ν 2 m ν 1 m ν 3 . Hence, in this picture, neutrino oscillations can be described by ten parameters: six (two independent ∆m 2 Atm , ∆m 2 Sol , three mixing angles θ 12 , θ 13 , θ 23 , and a Dirac CP phase δ CP ) associated with the standard three-active neutrino oscillations [34] and four (∆m 2 1,2 , ∆S 2 31 , θ 3 ) responsible for the new oscillations involving sterile neutrino. 6 Assuming ∆m 2 1(2) ∆m 2 Sol |∆m 2 Atm | |∆m 2 3 |, we expect that the effects of the pseudo-Dirac neutrinos for the first and second generations can be detected through ABL oscillation experiments [27,31], whereas that for the third generation can be measured through SBL oscillation experiments (or possibly long baseline oscillation experiments). The mass splittings will manifest themselves through very long wavelength oscillations characterized by the ∆m 2 1(2) as well as very short wavelength oscillations characterized by the ∆m 2 3 . The mass splitting ∆m 2 3 could be limited by the active neutrino mass orderings with a requirement of |∆m 2 3 | ∆m 2 Atm :

JHEP12(2019)133
where the hierarchical mass orderings (m ν 2 > m ν 1 m ν 3 and m ν 3 m ν 2 > m ν 1 ) are used. And since the mass splittings ∆m 2 1(2) can modify the large mixing angle solution of the solar neutrino oscillations, they should be limited and detailed fits imply a bound [43] Thus, we simply ignore ∆m 2 1(2) in the study of short baseline neutrino oscillations .
Oscillation experiments are unfortunately insensitive to the absolute scale of neutrino masses. Whereas cosmology is mostly sensitive to the total energy density in neutrinos, directly proportional to the sum of the active neutrino masses m ν = m ν 1 + m ν 2 + m ν 3 . We will mainly focus on cosmological observations as a probe of the absolute neutrino mass scale. Using eq. (2.22), the active neutrino masses in eq. (2.16) can be expressed in terms of the new parameters (∆S 2 31 , θ 3 ) and the two known mass squared differences of oscillation experiments (∆m 2 Atm , ∆m 2 Sol ) as Cosmological and astrophysical measurements provide powerful constraints on the sum of neutrino masses complementary to those from accelerators and reactors. There are several upper limits on the sum of active neutrino masses coming from the CMB data and weak lensing data: [45,46] 0.170 eV ,

eV ,
Weak lensing-only [47] (3.2) where a lower limit could be provided by the neutrino oscillation measurements. In order to extract the new physics effects, first we investigate the influence of ∆S 2 31 and sin 2θ 3 on i m ν i by imposing the experimental results on ∆m 2 Atm , ∆m 2 Sol and constraint given in eq.  [64], whereas the red-dotted line to a future sensitivity of mν e 0.20 [66]. For i=1,2,3 m ν , we take the values from eq. (3.2). and 3.3 eV at 95% CL are given by Planck Collaboration [44] and weak lensing-only [47], respectively, in eq. (3.2). In the plot 7 we consider only eV-mass scale of sterile neutrino since too heavy neutrino is conflict with cosmology ∆N eff ν < 0.2 at 95% CL [61,63]. The existence of sterile neutrino with the eV mass can also be constrained by βdecay experiments [48][49][50] and by 0νββ decay experiments [51][52][53][54][55][56]. The two types of mass ordering, discussed above, should be compatible with the existing constraints on the absolute scale of neutrino masses. The most sensitive experiment to search for the new physics effects in β-decay is to use the tritium decay process 3 H → 3 He + e − +ν e . Nonzero neutrino masses distort the measurable spectrum of the emitted electron. The most stringent upper bounds on theν e mass, mν e , have been obtained from direct searches in the Mainz [64] and Troitsk [65] experiments at 95% CL: The SBL anomalies including MiniBooNE data may indicate the existence of eV-mass sterile neutrino if those are interpreted as new oscillation effects, while present cosmological data coming from CMB + largescale structure and big bang nucleosynthesis (BBN) do not prefer extra fully thermalized sterile neutrinos in the eV-mass range since they violate the hot dark matter limit on the neutrino mass [57][58][59][60]. The amount of thermalisation ∆N eff ν as a function of neutrino parameters (mass splitting, mixing, and initial lepton asymmetry) has been quantitatively derived in ref. [61], implying that the parameter space of (∆m 2 3 , θ3) responsible for the existence of an eV-mass sterile neutrino is allowed by requiring such sterile neutrino does not or partially equilibrium at the BBN epoch when the initial lepton asymmetry is large [61,62].

JHEP12(2019)133
In figure 1 the dotted lines show contour for the neutrino mass in tritium β decay mν e as a function of ∆S 2 31 and sin 2θ 3 with eq. (2.19), where the black-dotted line corresponds to the upper bound mν e < 2.3 eV [64], whereas the red-dotted line to a future sensitivity of mν e 0.20 [66]. As seen in figure 1, the cosmological bounds given in eq. (3.2) are still tighter than the constraints from tritium β decay. The upcoming KATRIN experiment [66] planned to reach the sensitivity of mν e ∼ 0.20 eV will probe the region of the quasidegenerate mass spectrum of the active neutrinos.
On the other hand, the 0νββ-decay rate [67] effectively measures the absolute value of the ee-component of the effective neutrino mass matrix M ν in eq. (2.4). In the basis where the charged lepton mass matrix is real and diagonal, the 0νββ-decay rate can be expressed as Since the two mass eigenstates of first and second generations in each pseudo-Dirac pair have opposite CP parity, the third generation dominantly contributes to the ββ0ν-decay Using eqs. (2.21) and (2.22) one can easily see that ββ0ν-decay rate becomes almost zero, m ββ ≈ 0. Hence if the ββ0ν-decay rate is measured in near future the model would explicitly be excluded. 8

Short baseline neutrino anomalies
Now, let us study how our model can help to resolve the so-called short baseline neutrino anomalies in terms of neutrino oscillations. To see how the new sterile neutrino states cause such new oscillations at short-baselines with neutrino trajectory less than 1 < km, let us bring out a conversion probability of new oscillations with the help of the neutrino mixing matrix eq. (2.7). The conversion probability 9 between the massive neutrinos that a neutrino eigenstate ν a becomes eigenstate ν b follows from the time evolution of mass eigenstates as Note that the claim of observation of 0νββ-decay of 76 32 Ge [68] is strongly disfavored by the recent results of the GERDA experiment [69]. 9 The transition probability of να → νs i between sterile and active neutrinos due to oscillations of active flavor να (with α = e, µ, τ ) with sterile neutrinos νs i (with i = 1, 2, 3), see eq. (2.7), is given For non-vanishingŨµ3, the new mass squared difference ∆m 2 3 could be constrained by the probability of νµ → νs i , while for tiny or vanishingŨµ3 one could not constrain the ∆m 2 3 . In addition, the mass squared difference ∆m 2 3 could be constrained by the disappearance of muon-type neutrinos and antineutrinos produced in the atmosphere, see eq. (5.2) and figure 5, where the mixing θ3 deviated from the maximal mixing π/4 is involved.

JHEP12(2019)133
where a, b = e, µ, τ, s 1 , s 2 , s 3 , L is the distance between the neutrino detector and the neutrino source, E is the neutrino energy, andM ν ≡ W T ν M ν W ν . We are interested in the flavor conversion between the active neutrinos ν e , ν µ , ν τ satisfying the condition of eq. (2.24) which leads to θ 1(2) ≈ 0. From eq. (4.1) the flavor conversion probability between the three-active neutrinos can explicitly be expressed in terms of the oscillation parameters θ, ∆m 2 , L, E, and mixing components U αi of the 3 × 3 PMNS matrix as where ∆Q 2 kj ≡ m 2 ν k − m 2 s j , and δ k3 = 1 for k = 3 and 0 for k = 3. In the model the mixing parameters θ and ∆m 2 are determined by nature, so experiments should choose L and E to be sensitive to oscillations through a given ∆m 2 . As expected, in the limit of m s i → m ν i and θ i → 0 (i = 1, 2, 3), P να→ν β becomes the standard form of conversion probability for three active neutrinos in vacuum, as shown in ref. [34]. The model has interesting features listed below under the assumption of CPT invariance: • From eq. (4.2), we see that P να→να + P να→ν β + P νω→να ≤ 1 with α = β = ω = e, µ, τ , whereas the probabilities for the standard three-active neutrino oscillations satisfy the relation P να→να + P να→ν β + P νω→να = 1.
• The new oscillation effects attributed to ∆m 2 i , ∆S 2 3k and ∆Q 2 k3 with i = 1, 2, 3 and k = 1, 2 in eq. (4.2) can be maximum in the limit of θ 3 → π/2 (favored by the DNO and DIO), which can be relevant to short baseline neutrino experiments.
• In the limit of θ 3 → π/4 (favored by the IH, see eq. (3.1)), the oscillatory term involving ∆m 2 k (with k = 1, 2) in eq. (4.2) can give new oscillation effects only applicable to ABL oscillations. And thus, it is expected that this case could not provide a solution to the SBL anomalies.
To investigate the effects of new oscillations due to the new sterile neutrinos, we present approximated forms of neutrino oscillation probability based on the formula given by eq. We note that new oscillatory terms vanish in the limit of θ 3 = π/4, whereas they can reach maximum in the limit of θ 3 = π/2. In the limit that Similarly, the probabilities Pν µ→νe and Pν µ→νµ are approximately given by

New interpretation of reactor neutrino results
In order to probe the effects of the new sterile neutrino from the experimental results obtained at the reactor neutrino experiments, Daya Bay, RENO and Double Chooze (see ref. [70]), we interpret them in terms of neutrino oscillations including the new sterile neutrino by using the new reactor antineutrino probability in eq. (4.9). Adopting ∆m 2 31 by the values determined from atmospheric neutrino oscillation, we can obtain new values of θ 13 along with θ 3 by equating eq. (4.9) with the value ofν e survival probability for the three-active neutrino oscillation given by In the right panel of figure 2, we fix Eν e = 3.5 MeV and ∆S 2 31 = 0.6 eV 2 and the cyan and red plots correspond to the cases of (θ 3 , θ 13 ) = (1.28, 0.156) and (1.50, 0.167) , respectively. Note that the value of θ 13 corresponds to the central value evaluated at the given value of θ 3 with the new oscillation formula eq. (4.9), and the red-triangle, red-square (black-square), blue-circle, black-star, and black-circle error bars represent the values of the parameter R defined by the ratio of reactor antineutrino flux to the theoretical prediction obtained from the RENO [70], Double Chooz [70,71], Daya Bay [72,73], Palo Verde [74], and Chooz [75], where the error bars represent the experimental uncertainties. The values of R from RENO, Double Chooz, Daya Bay, Palo Verde, and Chooz have been obtained by subtracting the effects of θ 13 -driven oscillations. And the yellow band stands for the world average of R updated after including Daya Bay result R = 0.945±0.007(exp.) [72,73], compared with the past global average R past = 0.927 [94] indicated as black-dashed horizontal line, where the uncertainty of common reactor model is ±0.023 with respect to the Huber-Mueller model.

Reactor neutrino flux anomaly
The reactor antineutrino anomaly [9][10][11] is the experimental result presenting a deficit of the rate ofν e in several SBL reactor neutrino experiments with L ∼ (10 − 100) m and Eν e ∼ MeV. In reactor neutrino experiments, electron antineutrinos are detected through the inverse neutron decay processν e + p → n + e + in liquid-scintillator detectors. To  [72,73] and the black-dashed horizontal line for the previous world average [94].
interpret the deficit of observed reactor neutrino fluxes relative to the prediction (Huber-Mueller model [10,76]) in terms of neutrino oscillation including the new sterile neutrinos, it is relevant to use the probability given by eq. (4.4).
We note that the reactor neutrino flux anomaly is not clearly explained at L 500 m as in the 3+1 model [33,94]: eq.

LSND anomaly and MiniBooNE data
The LSND experiment [19,20] reported observation of a statistically suggestive excess of ν e events in a beam ofν µ produced by µ + decay at rest, µ + → e + + ν e +ν µ with 3.8σ significance. 12 The MiniBooNE experiments also observed ν e andν e appearance in ν µ and ν µ beams, respectively having the same L/E as in LSND [21]. In the left panel of figure 3, the black error bars stand for the LSND [19,20] data and the blue (red) error bars stand for the MiniBooNE data excess in (anti-)neutrino mode, at the baseline 541 m from the beryllium target, in particular, in the interval of energies 200 < E ν < 475 MeV, which corresponds to L/E range beyond that probed in the LSND experiment [21]. It has been 11 The recent NEOS and DANSS results [77,78] have a tension with the Gallium and reactor anomalies in the 3+1 model [94], while their results show a preference on our model predictions. 12 The similar KARMEN collaboration [79] did not measure any excess ofνe events over the background at a mean distance L 17.7 m with Eν µ = 12 ∼ 52.8 MeV, which did not fully exclude the LSND result. shown that the observed excesses in MiniBooNE experiment are in agreement with the LSND result, and provide a good fit to a large ∆m 2 solution in a two-neutrino oscillation framework, even though the two experiments have completely different neutrino energies, neutrino fluxes, reconstruction, backgrounds, and systematic uncertainties.
The LSND and MiniBooNE excesses could be explained byν µ →ν e oscillation whose probability is given by eq. (4.5). As in the case of 3 + 1 model [33,94], to explain the excess ofν e events, we need ∆S In the numerical estimation, we vary the parameter space (θ 13 , θ 23 ) in the ranges of (0.157−0.175, 0.730−0.885) and (0.148 − 0.164, 0.730 − 0.885) which correspond to i=e,µ,τ m ν i = 2.673 and 0.705 eV, respectively, for two bench mark points. It is worthwhile to note that the parameter space of (∆S 2 31 , θ 3 ), as presented in figure 1, is constrained through i=e,µ,τ m ν i by the cosmological data and the effective neutrino mass in tritium β-decay. However, statistical uncertainties have to be reduced by gaining more data in order to confirm our model in the following two respects: (i) the two data points in the left panel of figure 3 from MiniBooNE at L/E 1.5 (m/MeV) seem to favor m ν 0.705 eV, which is disfavored by Planck Collaboration (TT+lowP) at 95% CL [45,46] while still not excluded by the weak lensing only data [47], (ii) it seems that the ν e (blue-bars) andν e (red-bars) modes in MiniBooNE data, in principle, could be discriminated by considering the CP violating term in eq. (4.2). In the model setup, however, it seems not possible to discriminate between them within the short baseline due to the CP violating terms proportional to −4 sin(∆m 2 21 L/2E).

JHEP12(2019)133
In addition, searching for theν e appearance in LSND and MiniBooNE also implies ν µ disappearance whose oscillation probability is given by eq. (4.6). In the right panel of figure 3, we present how the survival probability Pν µ→νµ behaves along with L/E(m/MeV) for the same cases in the left panel of figure 3. The large deviation from Pν µ→νµ = 1 at large L/E as shown in the right panel of figure 3 is the characteristic feature of this model we consider, which makes our model different from 3+1 model. Accelerator based experiments forν µ (or ν µ ) disappearance, such as MINOS and MINOS+ neutrino experiments [80], may be sensitive to oscillations involving sterile neutrinos [81] for the regions of 10 −2 L/E(km/GeV) 0.75 at near detector (L = 1.04 km) and 8 L/E(km/GeV) 7×10 2 at far detector (L = 735 km). The former region can cover L/E values of order 1m/MeV in the right panel of figure 3, while the latter can cover range of L/E ∼ O(10 ∼ 100) [km/GeV], comparing with the results of the atmospheric muon neutrino events observed 13 in Super-Kamiokande [84] and the left plot in figure 5. To see such an oscillation dip due to new sterile neutrino as shown in the right panel of figure 3, we need the neutrino baseline ∼ 6 km (14 km) for detector having ∼ 3 GeV (7 GeV) peak energy with ∆S 2 31 = 0.6 eV 2 . In addition, such ν µ disappearance oscillation effect could be observed at L/E ∼ 1.2 [m/MeV] with E νµ ∼ 0.5 GeV and L = 600 m in the Short-Baseline Neutrino (SBN) Program experiment [85] for the same model parameters in figure 3. And, at the far detector (L = 735 km) a possibility on searching for signatures of sterile neutrinos in the ν µ disappearance by using eq. (4.6) will be considered in section 5.2.

Earth matter effects
Let us explore the Earth matter effect [40,41,[86][87][88][89] by examining the propagation of atmospheric neutrinos produced in cosmic-ray air showers from the Earth's atmosphere to the inside of the Earth. When muon neutrinos pass through the Earth matter, the MSW effect [40,41] should be taken into account. Path-length ranges from 10 km to 1.27 × 10 4 km depending on arrival zenith angle. The matter density encountered by neutrinos propagating is on average ρ ⊕ ∼ 3 g/cm 3 in the Earth's crust and outer mantle, ∼ 5 g/cm 3 in the inner mantle, and between 10 and ∼ 13 g/cm 3 in the core [90]. Muon neutrino oscillations modified due to matter effects can produce distinctive signatures of sterile neutrinos in the large set of high energy atmospheric neutrino data (both in the TeV energy window from IceCube [91] and at lower energy from DeepCore [82]). For E ν > 100 GeV, three-active neutrino oscillation length larger than the diameter of the Earth and can be neglected. In order to find appropriate physics parameters (θ, ∆m 2 , L, E) for atmospheric neutrino oscillations, we consider the effective Hamiltonian in-matter H m in flavor basis, which has the form of 6 × 6 matrix

JHEP12(2019)133
where k = 1, 2, 3 and α = e, µ, τ . Here the parameters A α = 2E m V α is a measure of the importance of matter effect with the matter-induced effective potential; V e , V µ , V τ , and V s = 0 are the potentials experienced by the electron, muon, tau, and sterile neutrinos respectively, and E m is the neutrino energy in matter. For anti-neutrinos V α → −V α . ν e 's have charged-current (CC) interactions with electrons and neutral-current (NC) interactions with electrons and nucleons, V e = √ 2G F (N e − N n /2), while ν µ 's and ν τ 's have only NC interactions, V µ = V τ = √ 2G F (−N n /2), and any ν s 's have no interactions, V s = 0, where Fermi's constant, G F , and the average electron and neutron densities along the neutrino path, N e and N n , respectively.
The mass matrix of the massive neutrinos in matter can be diagonalized through a new unitary mixing matrix W m , (4.13) The diagonalization of H m by the unitary matrix W m in matter gives a condition (4.14) and its effective mass-squared eigenvalues in matter which are positivẽ For θ 1(2) ∼ 0 corresponding to ∆m 2 1(2) < 10 −12 eV 2 in eq. (2.24) we assume no Earth matter effects occurs, leading to θ m k → θ k . Then from eq. (4.14) the energy of upwardgoing atmospheric muon neutrinos in matter can be derived as and sin 2(θ m 3 − θ 3 ) → sin 2(θ 3 − θ m 3 ) forν µ , indicating that the energy E m varies with the matter mixing angle θ m 3 . Note that as θ m 3 → θ 3 + nπ/2 (with n = 0, 1, 2, 3, . . .), the matter effect gets faded. In figure 4 we plot the effective mass-squaredm 2 s 3 (black-dotted lines), m 2 ν 3 (blue-solid lines) and E m (red-bands) in terms of θ m 3 . Thanks to the positiveness of m 2 s 3 ,m 2 ν 3 and E m , they are predicted for limited regions of θ m 3 as can be seen in figure 4. The lines overlapped with the red bands represent he effective masses affected by Earth matter effect, in the limit of E m → 0,m 2 ν 3 approached to m 2 ν 3 , implying that matter effect becomes negligible for the muon neutrinos with energies 1 GeV passing through the interior of Earth. Hence, it is expected that the ν µ ( orν µ ) disappearance probability for high energy ( 1 GeV) muon neutrinos passing through the Earth interior to search for JHEP12(2019)133 Recently, the high-energy IceCube detector has measured the atmospheric muon neutrino spectrum at energy E νµ = 320 GeV ∼ 20 TeV in hope of finding the oscillation signatures of light sterile neutrinos [91], but no evidence for anomalous ν µ orν µ disappearance is observed. The results coincide with the model prediction on the matter effects at the same energy window, contrary to the models in refs. [92][93][94][95][96][97][98][99][100][101][102][103][104][105][106] which has a tension with ν µ (orν µ ) disappearance experiments. And the IceCube sub-detector DeepCore at energy window of atmospheric muon neutrinos 10 ∼ 100 GeV [82] can also have potential to search for the signature of light sterile neutrinos as discussed in the next section.

New effects in solar and atmospheric oscillations
Now, let us examine the oscillation effects due to the sterile neutrinos on the long baseline experiments such as KamLAND [34], T2K [107], MINOS and MINOS+ Collaboration [80], solar neutrino and atmospheric neutrino experiments [82,84].

New effects in ν e disappearance from KamLAND, T2K and Solar neutrino oscillation
For the long baseline such as KamLAND experiment [108], the survival probability ofν e events in the model we consider is approximately given by

JHEP12(2019)133
where we assume CP invariance and the term sin 2 ∆m 2 3 4E L is averaged out for long baseline (e.g. L 180 km). Applying eq. (5.1) to KamLAND data [109,110], we can extract the value of θ 3 with precise measurements of ∆m 2 21 and θ 12 . However, we see that eq. (5.1) is almost the same as the expression for the standard oscillation probability for three-active neutrinos because the new term concerned with θ 3 is negligible due to the tiny value of sin 4 θ 13 5 × 10 −4 . Thus, new effects due to the sterile neutrinos on KamLAND, T2K experiments and solar neutrino oscillation are negligible.
Since the T2K experiment at the ND280 near detectors covers L/E values of order 1 m/MeV, an information such as θ 3 in eq. (5.1) on sterile neutrinos can also be extracted. It has performed a search for ν e disappearance in a neutrino beam whose ν e component is peaked at an energy of 500 MeV [107]. From eq. (5.1) in the limit L 4πE νe /∆m 2 21 the ν e disappearance probability driven by the sterile neutrino can easily be obtained for the baseline of 280 m, and it shows that its effect on sterile neutrino is negligible due to the tiny value of sin 4 θ 13 .

New effects in ν µ disappearance from SuserK and IceCube
In addition to the TeV muon neutrinos discussed in previous section, the IceCube Collaboration can also have potential to search for the signature of light sterile neutrinos by observing atmospheric neutrinos in the tens-of-GeV range through its sub-detector DeepCore (for reference, see ref. [82]). The atmosphere of the Earth is constantly being bombarded by cosmic rays, primarily made up of protons and helium nuclei produced by atmospheric objects [34]. These cosmic rays have been observed over a wide range of energy, from 1 GeV to 10 11 GeV. We recall atmospheric neutrinos with energy of ∼ a few GeV which are mostly produced by primary cosmic rays with energy of ∼ 100 GeV [34]. According to the three-active neutrino oscillation probability, large deficits of muon neutrino have been observed in upward-going events and at L/E ∼ 230 km/GeV [34] as shown by blue sinusoidal curves in figure 5 which is well consistent with the Super-Kamiokande's atmospheric neutrino data [84]. For the baselines L 4πE/∆m 2 Sol , 4πE/∆m 2 1(2) , which is relevant to the atmospheric neutrino, the survival and conversion probability of muon neutrino are approximately given by  figure 5 shows ν µ survival probability as a function of L/E(km/GeV), where the blue sinusoidal curve and the cyan rapid oscillations stand for the probability predicted from the standard form of probabilty for three active neutrinos and that from the new oscillation probability with parameters θ 3 = 1.28 and ∆S 2 31 = 0.6 eV 2 , and the data points are atmospheric neutrino events observed in Super-Kamiokande [84]. As shown in the plot the atmospheric neutrino data from Super-Kamiokande [84] are well consistent with the new ν µ ↔ ν τ oscillation (cyan curve), showing the first oscillation dip appeared at ∼ 500 km/GeV, and interestingly enough, the error bars with red-circle showing large deficits compared with the results in the three-active neutrino framework are also well consistent with the new oscillation curve, which may be due to the existence of new sterile neutrino.
Since the high energy muon neutrinos (10 ∼ 100 GeV) passing through the interior of the Earth are little affected by matter in this model, it is expected that the IceCube sub-detector DeepCore results [82] are similar to that provided by Super-Kamiokande atmospheric neutrino data [84] representing no evidence for light sterile neutrinos. In addition to the IceCube DeepCore results [82], the track-like muon events selected in the MINOS and MINOS+ data [70] shows the first oscillation dip at ∼ 500 km/GeV mainly due to the muon neutrino oscillation into not sterile neutrinos but tau neutrinos.

JHEP12(2019)133
However, there exists a possibility to probe the existence of sterile neutrinos by comparing the results of P νµ→νµ with those of 1−P νµ→ντ . In fact, both results are almost equivalent for the standard oscillation for the three-active neutrinos. But, they are different for the new oscillation affected by the sterile neutrinos. In the right panel of figure 5, we plot 1−P νµ→ντ , where the blue sinusoidal curve and the cyan oscillation curve correspond to the standard three-active neutrino oscillations and new oscillation with θ 3 = 1.28, ∆S 2 31 = 0.6 eV 2 , respectively. In the ν µ ↔ ν τ oscillation, the oscillation parameters are fixed to the best-fit values 14 for NO given in eq. (2.20). Such discrepancy between P νµ→νµ and 1−P νµ→ντ would be proved through both the ν µ → ν τ appearance experiments and the ν µ disappearance experiments in the future. Recently, the OPERA collaboration has confirmed that muon neutrinos primarily oscillate into tau neutrinos [111]. In the case ofm 2 ν 3 → m 2 ν 3 , as shown in figure 4, an MSW resonance is expected to occur in multi-GeV, similar to the standard oscillation for the three-active neutrino [34], for the atmospheric and accelerator ν e,µ neutrinos traveling in earth matter with the propagation eigenstates of active neutrinos [40,41].

Implications of IceCube data
The astrophysical neutrinos with very high energy fly galactic and extra galactic distances far beyond the earth-sun distance can give us an opportunity to probe pseudo-Dirac neutrinos with very tiny mass splittings as mentioned before. Taking into account astronomicalscale baseline satisfying 4πE/∆m 2 Sol,Atm L ∼ 4πE/∆m 2 k with k = 1, 2 to uncover the oscillation effects of very tiny mass splitting ∆m 2 k , the probability of neutrino flavor conversion from eq. (4.2) reads where the oscillatory terms involving the atmospheric and solar mass-squared differences and the large mass-squared differences ∆m 2 3 , ∆S 2 3k and ∆Q 2 3k with k = 1, 2 are averaged out over such long distances. As shown in [112], the matter effects inside the Gamma Ray Burst (GRB) sources as well as the earth are not significant, which makes us to consider vacuum oscillation only for astrophysical neutrinos. Neutrino telescope such as IceCube [113] observes neutrinos from extragalactic sources located far away from the earth and with neutrino energy 10 5 GeV E 10 7 GeV. Given neutrino trajectory L and energy E, the oscillation effects become prominent when ∆m 2 k ∼ E/4πL, where L ≡ L(z) is a distance-measure with redshift z, which is different from comoving or luminosity distance, given by in which Ω Λ , Ω m , and H 0 stand for the dark energy density of the Universe, the matter density of the Universe, and the present Hubble expansion rate, respectively. The asymptotic value of L(z) is about 2.1 Gpc achieved by large value of z, which means that the smallest ∆m 2 k that can be probed with astrophysical neutrinos with E is 10 −17 eV 2 (E/PeV) [28,114]. If this is the case, in order to observe the oscillation effects the oscillation lengths should not be much larger than the flight length before arriving at neutrino telescopes in earth for given tiny mass splittings, that is, On the other hand, oscillation effects induced by tiny mass splittings ∆m 2 k for the pseudo-Dirac neutrinos can affect the track-to-shower ratio for the number of shower N S and track events N T measured from IceCube experiment, which is given by [25][26][27], with a spectral index ω, the detector effective areas A α (E) and initial flavor composition φ 0 β . Here p T is the probability that an observed event produced by a muon neutrino is a track event, which is mildly dependent on energy and approximately equals to 0.8 [22][23][24]. The prediction for the ratio N T /N S depends on the initial flavor composition  Figure 6. Plots of the track-to-shower ratio N T /N S as a function of L (log 10 [path length/10 kpc]) for NO (red curve) and IO (blue curve) with ∆m 2 1 = 10 −15 eV 2 , ∆m 2 2 = 10 −16 eV 2 and θ 3 = 1.28. Each panel corresponds to the specific initial flavor composition (φ 0 e : φ 0 µ : φ 0 τ ) at the source. For three neutrino mixing angles and Dirac-type CP phase we take the global fit results in eq. (2.20). Red and blue curved lines correspond to normal and inverted neutrino mass orderings, respectively, for ω = 2.2, whereas light red and light blue regions represent the corresponding results for ω = 1.8−2.6. Gray shaded regions are forbidden by N T /N S = 0.18 +0. 13 −0.05 [25,26].
We confront the predictions of N T /N S with experimental results by taking ∆m 2 1 = 10 −15 eV 2 , ∆m 2 2 = 10 −16 eV 2 and θ 3 = 1.28 as a benchmark point as well as the best-fit values in eq. (2.20) for the neutrino mixing angles and CP phase. As can be seen from eqs. (6.5), (6.6), the tiny mass splittings ∆m 2 k(=1,2) can be searched for, looking at high energy cosmic neutrinos by measuring the track-to-shower ratio N T /N S as the function of L (log 10 [path length/10 kpc]).
In the numerical analysis, we use the spectral index given by ω = 2.2 ± 0.4 [22][23][24] and the best-fit values for NO (IO) in eq. (2.20). Figure 6 shows the plots of the track-to-shower ratio N T /N S as a function of L (log 10 [path length/10 kpc]) for the neutrino energy 60 TeV E ν 3 PeV studied in refs. [25,26]. According to four specific assumptions at each panel for the flavor compositions at the source (φ 0 e : φ 0 µ : φ 0 τ ), for ω = 2.2, ∆m 2 1 = 10 −15 eV 2 and ∆m 2 2 = 10 −16 eV 2 the normal (inverted) mass ordering is presented as the red (blue) curved line, whereas light red and light blue regions represent JHEP12(2019)133 the corresponding results for ω = 1.8 − 2.6. Gray shaded regions are forbidden by the measurement N T /N S = 0.18 +0. 13 −0.05 [25,26]. In figure 6, we see that the oscillation effect occurs at distance 1.5 Mpc and it is averagged out at distance beyond 1.5 Mpc. The predictions of N T /N S for the given inputs and the specific initial flavor compositions 1/3 : 2/3 : 0 and 1/2 : 1/2 : 0 are consistent with the measurement, whereas those for the other two initial flavor compositions are disfavored. In the plots, we draw the horizontal dashed lines corresponding to the cases without oscillation effects (i.e. the cases for ∆m 2 1,2 = 0.). The gap between the predictions with and without oscillations is due to the oscillatory term in eq. (6.1). Therefore, substantial reduction of uncertainty in the masurement of N T /N S would test not only the model itself but also the oscillation effects induced by pseudo-Dirac nature of neutrinos.

Conclusion
We have proposed a convincing model containing sterile neutrinos to interpret both SBL neutrino anomalies and high energy neutrino data observed at IceCube in terms of neutrino oscillations. Different from the so-called 3 + 1 model where the PMNS matrix is simply extended to 4 × 4 unitary matrix as in refs. [32,33,94], the new 4 × 4 neutrino mixing matrix in our model is parameterized in a way to keep the 3 × 3 PMNS mixing matrix for three active neutrinos unitary. A characteristic feature of this model is that there are no flavor changing neutral current interactions leading to the conversion of active neutrinos to sterile ones or vice versa. We have presented new forms of neutrino oscillation probabilities modified by introducing new sterile neutrinos. In this scenario, there are new mass squared differences (∆m 2 1,2,3 , ∆S 2 31 ) and new mixing angle θ 3 in addition to the standard oscillation parameters associated with only three active flavor neutrinos. While ∆m 2 1,2 are responsible for astronomical baseline high energy neutrino oscillations, ∆m 2 3 (or ∆S 2 31 ) and θ 3 are usable to interpret SBL neutrino oscillations.
Our model can explain SBL neutrino anomalies in terms of neutrino oscillations at the same level of 3 + 1 model. However, there still exist small tensions (1) in the reactor and Gallium data, flux normalization and 5 MeV bump observed from ν e → ν e disappearance, (2) between the MiniBooNE data in the region L/E 1.5 [m/MeV] and the model predictions at ∆S 2 31 = 0.6 eV 2 for ν µ → ν e appearance for m ν favored by Planck Collaboration (TT+lowP) at 95% CL [45,46], and (3) between the Super-Kamiokandes atmospheric neutrino data and the IceCube DeepCore results including the MINOS and MINOS+ data released in Neutrino 2018 .
We have shown that resolution of the reactor antineutrino flux anomaly is required to reduce statistical uncertainties and/or to understand its underlying physics at baselines L 500 m. In the present model the values of the parameters θ 3 and ∆S 2 31 favored by the Planck data of m ν 0.705 eV are not conflict with the NEOS and DANSS results [77,78] which are in a tension with the Gallium and reactor anomalies as in 3+1 model [94].
We have shown that the LSNDν e appearance data favors probability driven by while the excess of two data at L/E 1.5 [m/MeV] seems to be disfavored by a cosmological bound m ν 0.705 eV. Since the LSND and MiniBooNE data can be interpreted as ν µ →ν e oscillation, experimental search forν µ (or ν µ ) disappearance would test our model in which it could be observed at L/E ∼ 1.2 [m/MeV] with E νµ ∼ 0.5 GeV and L = 600 m for ∆S 2 31 = 0.6 eV 2 with θ 3 = 1.28 in the SBN Program experiment [85]. In addition, we have studied the earth matter effect, and found that it is negligible for muon neutrinos having energies 1 GeV when they pass through the interior of the Earth. Finally, we have found that the existence of light sterile neutrino we consider does not affect solar neutrino oscillation and thus no constraint on new parameters came out from it. It has been shown that the Super-Kamiokande's atmospheric neutrino data are consistent with the new ν µ ↔ ν τ oscillation affected by sterile neutrino, showing the first oscillation dip appeared at ∼ 500km/GeV which is a characteristic feature of this model. The most recent data of DeepCore [82] and MINOS and MINOS+ [70] experiments do not show any signature of light sterile neutrinos. In addition, we have shown that the probabilities of Pν µ→νµ and 1 − P νµ→ντ versus L/E [km/GeV] have a clear discrepancy with different oscillation signatures of light sterile neutrinos, unlike the expected from the three neutrino standard form. Such discrepancy could be probed through both the ν τ appearance and the ν µ disappearance experiments in the future. We have discussed the implications of the very high energy neutrino events detected at IceCube on the probe of the oscillation effcts induced by two pseudo-Dirac mass splittings.