Can we measure the neutrino mass hierarchy in the sky?

Cosmological probes are steadily reducing the total neutrino mass window, resulting in constraints on the neutrino-mass degeneracy as the most significant outcome. In this work we explore the discovery potential of cosmological probes to constrain the neutrino hierarchy, and point out some subtleties that could yield spurious claims of detection. This has an important implication for next generation of double beta decay experiments, that will be able to achieve a positive signal in the case of degenerate or inverted hierarchy of Majorana neutrinos. We find that cosmological experiments that nearly cover the whole sky could in principle distinguish the neutrino hierarchy by yielding `substantial' evidence for one scenario over the another, via precise measurements of the shape of the matter power spectrum from large scale structure and weak gravitational lensing.


Introduction
In the past decade, there has been great progress in neutrino physics. It has been shown that atmospheric neutrinos exhibit a large up-down asymmetry in the SuperKamiokande (SK) experiment. This was the first significant evidence for a finite neutrino mass [1] and hence the incompleteness of the Standard Model of particle physics. Accelerator experiments [2,3] have confirmed this evidence and improved the determination of the neutrino mass splitting required to explain the observations. The Sudbury Neutrino Observatory (SNO) experiment has shown that the solar neutrinos change their flavors from the electron type to other active types (muon and tau neutrinos) [4]. Finally, the KamLAND reactor anti-neutrino oscillation experiments reported a significant deficit in reactor anti-neutrino flux over approximately 180 km of propagation [5]. Combining results from the pioneering Homestake experiment [6] and Gallium-based experiments [7], the decades-long solar neutrino problem [8] has been solved by the electron neutrinos produced at Sun's core propagating adiabatically to a heavier mass eigenstate due to the matter effect [9].
As a summary, two hierarchical neutrino mass splittings and two large mixing angles have been measured, while only a bound on a third mixing angle has been established. Furthermore the standard model has three neutrinos and the motivation for considering deviations from the standard model in the form of extra neutrino species has now disappeared [11,12].
New neutrino experiments aim to determine the remaining parameters of the neutrino mass matrix and the nature of the neutrino mass. Meanwhile, relic neutrinos produced in the early universe are hardly detectable by weak interactions but new cosmological probes offer the opportunity to detect relic neutrinos and determine neutrino mass parameters.
It is very relevant that the maximal mixing of the solar mixing angle is excluded at a high significance. The exclusion of the maximal mixing by SNO [4] has an important impact on a deep question in neutrino physics: "are neutrinos their own anti-particle?". If the answer is yes, then neutrinos are Majorana fermions; if not, they are Dirac. If neutrinos and anti-neutrinos are identical, there could have been a process in the Early Universe that affected the balance between particles and anti-particles, leading to the matter anti-matter asymmetry we need to exist [13]. This question can, in principle, be resolved if neutrinoless double beta decay is observed. Because such a phenomenon will violate the lepton number by two units, it cannot be caused if the neutrino is different from the anti-neutrino (see [10] and references therein). Many experimental proposals exist that will increase the sensitivity to such a phenomenon dramatically over the next ten years (e.g., [14] and references therein). The crucial question we want to address is if a negative result from such experiments can lead to a definitive statement about the nature of neutrinos. Within three generations of neutrinos, and given all neutrino oscillation data, there are three possible mass spectra: a) degenerate, with mass splitting smaller than the neutrino masses, and two non-degenerate cases, b) normal hierarchy, with the larger mass splitting between the two more massive neutrinos and c) inverted hierarchy, with the smaller spitting between the two higher mass neutrinos. For the inverted hierarchy, a lower bound can be derived on the effective neutrino mass [10]. The bound for the degenerate spectrum is stronger than for inverted hierarchy. Unfortunately, for the normal hierarchy, one cannot obtain a similar rigorous lower limit.

JCAP05(2010)035
Neutrino oscillation data have measured the neutrino squared mass differences, which are hierarchical. Given the smallness of neutrino masses and the hierarchy in mass splittings, we can characterize the impact of neutrino masses on cosmological observables and in particular on the matter power spectrum by two parameters: the total mass Σ and the ratio of the largest mass splitting to the total mass, ∆. As we will show, one can safely neglect the impact of the solar mass splitting in cosmology. In this approach, two masses characterize the neutrino mass spectrum, the lightest one, m, and the heaviest one, M .
Neutrino oscillation data are unable to resolve whether the mass spectrum consists in two light states with mass m and a heavy one with mass M , named normal hierarchy (NH) or two heavy states with mass M and a light one with mass m, named inverted hierarchy (IH). Near future neutrino oscillation data may resolve the neutrino mass hierarchy if one of the still unknown parameters that relates flavor with mass states is not too small. On the contrary, if that mixing angle is too small, oscillation data may be unable to solve this issue. Analogously, a total neutrino mass determination from cosmology will be able to determine the hierarchy only if the underlying model is normal hierarchy and Σ < 0.1 eV (see e.g., figure 1). If neutrinos exist in either an inverted hierarchy or are denegerate, (and if the neutrinoless double beta decay signal is not seen within the bounds determined by neutrino oscillation data), then the three light neutrino mass eigenstates (only) will be found to be Dirac particles.
In this paper, we investigate whether cosmological data may positively establish the degenerate spectrum from the inverted hierarchy (or vice versa). Our approach is to take cosmic variance limited surveys, rather than specifically planned experiments, so that we can determine if (even in the ideal case) cosmology can make any impact on this question.

JCAP05(2010)035
2 Massive neutrinos and the power spectrum Massive neutrinos affect cosmological observations in a variety of different ways. For example, cosmic microwave background (CMB) data alone constrain the total neutrino mass Σ < 1.3 eV at the 95% confidence level [15]. Neutrinos with mass 1eV become nonrelativistic after the epoch of recombination probed by the CMB, thus massive neutrinos alter matter-radiation equality for a fixed Ω m h 2 . After neutrinos become non-relativistic, their free streaming damps the small-scale power and modifies the shape of the matter power spectrum below the free-streaming length. Combining large-scale structure and CMB data, at present the sum of masses is constrained to be Σ < 0.3 eV [16]. Forthcoming large-scale structure data promise to determine the small-scale (0.1 k 1 h/Mpc) matter power spectrum exquisitely well and to yield errors on Σ well below 0.1 eV (e.g., [17,18,33]). Here, we assume the standard ΛCDM model and explore the changes in the matter power spectra due to the neutrino properties (mass and hierarchy).
The effect of neutrino mass on the CMB is related to the physical density of neutrinos, and therefore the mass difference between eigenstates can be neglected. However individual neutrino masses can have an effect on the large-scale shape of the matter power spectrum. In fact, neutrinos of different masses have different transition redshifts from relativistic to nonrelativistic behavior, and their individual masses and their mass splitting change the details of the radiation-domination to matter-domination regime. As a result the detailed shape of the matter power spectrum on scales k ∼ 0.01 h/Mpc is affected. In principle therefore a precise measurement of the matter power spectrum shape can give information on both the sum of the masses and individual masses (and thus the hierarchy), even if the second effect is much smaller than the first.
We define the relation between the neutrino masses m and M and the parameters Σ and ∆ as NH : (recall that m denotes the lightest neutrino mass and M the heaviest).
In figure 1 we show the current constraints on neutrino mass properties in the m-∆ and Σ-∆ planes. While many different parameterizations have been proposed in the literature to account for neutrino mass splitting in a cosmological context, and performed forecasts, [19][20][21][22] here we advocate using the ∆ parameterization for the following reasons. ∆ changes continuously through normal, degenerate and inverted hierarchies; ∆ is positive for NH and negative for IH. Finally, as we will show, cosmological data are sensitive to ∆ in an easily understood way through the largest mass splitting (i.e., the absolute value of ∆), while the direction of the splitting (the sign of ∆) introduces a sub-dominant correction to the main effect. This parameterisation is strictly only applicable for Σ > 0, but oscillations experiments already set Σ > M 0.05eV.
It is important to note that not the entire parameter space in the Σ-∆ plane (or of any other parameterization of the hierarchy used in the literature) is allowed by particle physics constraints and should be explored: only the regions around the normal and inverted hierarchies allowed by neutrino oscillation experiments should be considered (see figure 1).
To gain a physical intuition on the effect of neutrino properties on cosmological observables, such as the shape of the matter power spectrum, it is useful to adopt the following ana-JCAP05(2010)035 lytical approximation, as described in ref. [21]. The matter power spectrum can be written as: where ∆ 2 R is the primordial amplitude of the fluctuations, n s is the primordial power spectrum spectral slope, T (k) denotes the matter transfer function and D ν (k, z) is the scaledependent linear growth function, which encloses the dependence of P (k) on non-relativistic neutrino species.
Each of the three neutrinos contributes to the neutrino mass fraction f ν,i where i runs from 1 to 3, and has a free-streaming scale k fs,i , Analogously, one can define the corresponding quantities for the combined effect of all species, by using Σ instead of m ν i . Since we will only distinguish between a light and a heavy eigenstate we will have e.g., f ν,m , f ν,Σ , k fs,m , k fs,Σ etc., where in the expression for f ν,m one should use the mass of the eigenstate (which is the mass of the individual neutrino, or twice as much depending on the hierarchy) while in k fs,m one should use the mass of the individual neutrino. The dependence of P (k) on non-relativistic neutrino species is in D ν (k, z), given by where k ≫ k fs,i (z) and p i = (5 − 25 − 24f ν i )/4. The standard linear growth function D(z) fitting formula is taken from [23].
In summary there are three qualitatively different regimes in k-space that are introduced by the neutrino mass splitting where the subscript m refers to the light neutrino eigenstate and Σ to the sum of all masses. This description is, however, incomplete: the transitions between the three regimes is done sharply in k while in reality the change is very smooth. In addition, the individual masses change the details of the matter-radiation transition which (keeping all other parameters fixed) adds an additional effect at scales k > k fs,Σ .
In order to explore what constraints can be placed on ∆ and Σ for a given survey set-up we can use a Fisher matrix approach. The elements of F, the Fisher information matrix [24], are given by computed as σ 2 (θ) = (F −1 ) θθ . We can also calculate expected Bayesian evidence for cosmological parameters using the approach of ref. [25,34]. In the case that we are considering we use the formula from [34] for the expectation value of the evidence, in this case the expected Bayes factor is simply the log of ratio of the Fisher determinants. Following ref. [26] the Fisher matrix for the galaxy power spectrum is with N = [nP (k, µ)/(nP (k, µ) + 1)] 2 and V s is the volume of the survey. The integration over the projected angle along the light of sight 1 µ is analytical and the maximum and minimum wavenumbers allowed depend on the survey characteristics with the constraint that k max must be in the linear regime. The derivatives are computed at the fiducial model chosen. Throughout this paper we assume a fiducial model given by basic parameters of the standard LCDM cosmology [15] and the fiducial values for Σ and ∆ are then further specified in each case. Despite its limitations, the analytic description of the neutrino effect described above is extremely useful when performing an order-of magnitude calculation of an effect. Its corresponding Fisher matrix-approach forecasted errors indicate that while Σ can be constrained tightly, nearly-ideal, full-sky, cosmic-variance dominated surveys will be needed to obtain promising errors on ∆. However, we find that the analytical approximation above overestimates the neutrino effects on the P (k) and therefore under-estimates the forecasted errors, by factors of ∼ few in the regime of interest (Σ < 0.3 eV, ∆ along the NI and IH lines). In what follows we therefore use the publicly available CAMB code [27] to compute the matter power spectrum. Figure 3. ∆χ 2 as a function of the degeneracy parameter ∆ for a fixed total neutrino mass Σ (and fixed cosmology). This is a section along a Σ=constant line of figure 1 of the quantity −2 ln L as it would be seen by a Fisher matrix approach for a IH fiducial model. The vertical lines show the location of the normal and inverted hierarchy. Note the bimodal distribution of the ln L surface, which makes the determination of the hierarchy from measurements of the shape of the power spectrum extremely challenging. The ∆χ 2 normalization matches that achievable from an ideal weak lensing survey as described in the text.

JCAP05(2010)035
In figure 2 we show the dependence of P (k) on the parameter ∆ at z = 0 for fixed Σ and fixed cosmological parameters. The dependence is shown as the fractional change of the matter power spectrum for a unit change of the parameter ∆. This quantity is then fed directly to the Fisher matrix (see eq. (2.11)). In order to compute reliably the above derivatives, CAMB needs to be run at the highest precision settings, with fine k sampling and taking care that interpolations procedures in-built in the code do not introduce a spurious signal.
For the value of the total mass Σ = 0.25 eV, adopted in figure 2, normal (inverted) hierarchy correspond to ∆ ∼ 0.05 (∆ ∼ −0.05), indicating that the effect of the neutrino mass splitting on the P (k) is at the ∼ 0.2% level. The dependence of P (k) on ∆ at k > 0.1 h/Mpc arises because even for a fixed Σ the individual masses affect the tail of the energy distribution of the relativistic species and thus matter-radiation equality. Note that ∂ ln P/∂∆ changes sign with ∆ and there is a location, ∆ = 0, the degenerate case, where P (k) shows no dependence on ∆.
To understand the meaning and implications of this let us consider that the error on ∆ is directly proportional to ∆χ 2 = −2(ln L − ln L fiducial ) where L fiducial denotes the Likelihood for the fiducial model. In addition we can write This quantity is shown in figure 3, where the normalization has been chosen to match the constraints achievable from an ideal full sky weak lensing survey as the one considered in section 3. Figure 3 shows that the likelihood surface is bimodal: for example, for a fixed cosmology and fixed Σ, if the fiducial model is the inverted hierarchy ∆ < 0, there is a corresponding JCAP05(2010)035 value of ∆ > 0 (normal hierarchy), consistent with the neutrino oscillations in the allowed region, with a P (k) virtually indistinguishable from the fiducial model.
The bimodality of the likelihood surface also implies that the Fisher matrix approach to forecasting errors need to be applied with care before it can interpreted in terms of distinguishing the hierarchy. The curvature of the likelihood around the fiducial model gives the formal error on ∆, and this error may be much smaller than the distance between IH and NH ∆ values. But this could be interpreted as a determination of the hierarchy if and only if the likelihood had a unique maximum, which is not the case here. This subtlety has not be noticed in the literature before where Fisher error-estimates for parameters describing neutrino hierarchy were presented. In general, errors have been computed around one or more fiducial models and were sometimes found to be smaller than the distance between normal and inverted. We point out here that this cannot be directly interpreted as being sufficient to distinguish the hierarchy (it is a necessary but not sufficient condition if the likelihood is multi-peaked).
A more detailed inspection of figure 3 also indicates that the ∆χ 2 between the two minima (maxima of ln L) is not exactly zero, but it is very small, and that the location of the second minimum (assuming a fiducial IH) does not coincide with the central value of the oscillations-experiments regions. The evidence ratio can the be used to quantify wether a given survey set up could distinguish the two cases. Since this seems to be a very challenging measurement, we will put ourselves in the most optimistic situation where the neutrinos parameters are described as deviations from the minimal LCDM model and only afterwards consider relaxing this assumption.
The philosophy of the rest of the paper is therefore: "can cosmology in the cosmicvariance limit, and for an ideal experiment, distinguish the neutrino heirarchy?" or in other words, "is there enough information in the sky to measure the neutrino hierarchy? "

JCAP05(2010)035 3 Forecasted constraints from large scale structure
Here we explore what constraints can be placed on ∆ and Σ from ideal, cosmic variancedominated future surveys probing the shape of the matter power spectrum. The two probes of large-scale structure (LSS) we consider are the matter power spectrum itself and weak lensing.
We also compute the Fisher matrix of a CMB experiment like Planck 2 in order to help break degeneracies in the cosmological parameters when determined only by the power spectrum, or weak lensing. Therefore our final Fisher matrix is F = F P (k),W L + F CMB . We compute the combined Fisher matrix for variations in the following cosmological parameters: n s , α s , Ω ν h 2 , ∆, Z, Ω b h 2 , Ω c h 2 , h, A s , where α s denotes the running of the power spectrum spectral slope and Z is related to the optical depth to the last scattering surface via Z = exp(−2τ ); Ω ν is related to Σ via Σ = 94Ω ν h 2 (eV). The reported errors on Σ and ∆ are marginalized over the other cosmological parameters. The marginal errors for ∆ and Σ are shown in figure 4; the left panel is for a direct P (k) measurement approach and the right panel is the weak lensing approach.
Because we are interested in answering the question: "is there enough information in the sky to measure the neutrino hierarchy?" we have chosen survey parameters that are ambitous cosmic variance-limited surveys. For the parameter points shown in the left panel of figure 4 we have assumed a survey that covers the full sky 40,000 square degrees and maps the positions of galaxies up to z = 2 corresponding to about 600Gpc 3 comoving volume and maps the 21cm-HI up to z = 5, corresponding to about 2000 Gpc 3 comoving volume. We also assumed a high number density of galaxies so that we work in the cosmic variancedominated regime (nP ≫ 1). Galaxies are expected to be a biased tracer of the dark matter distribution, here we assume the bias to be scale and redshift-independent and thus not to affect the recovery of the shape of the matter power spectrum.
HI surveys [28] target the hyperfine transition in the hydrogen atom, which in the restframe emits a photon in the radio wavelenghts (21 cm). Therefore they survey the amount of neutral hydrogen in the universe. Because most galaxies and dark matter overdensities contain neutral hydrogen, such surveys provide the most un-biased indirect tracer of the dark matter distribution in the Universe. Further, in this frequency band, the radio spectrum is featureless with the only line being the 21 cm one, its observed frequency yielding a redshift and thus the radial distance of the emitter. Thus, an imaging survey automatically gives a three dimensional map of the HI distribution. The main challenge facing the HI surveys is the contamination by foregrounds [28]. For the characteristics of the survey we have followed the numbers given in [29] which yield to bias of 1 and negligible shot noise.
The survey considered is certainly a challenging one, but our calculations indicate that a cosmic variance-limited galaxy and HI survey can provide enough information to determine the neutrino hierarchy. We find that such a survey could constrain the total sum of neutrino mass with extreme accuracy O(10 −5 ). We also find that if the total neutrino mass is smaller than 0.15 eV, then the IH could be distinguished from the NH through an evidence criteria centered on each peak in ∆.
Weak lensing is the effect where the path of photons propogating from a galaxy are distorted by intervening mass concentrations. The amount of distortion depends on the density and distribution of the mass. For an individual galaxy image the weak lensing effect is to induce a change in ellipticity or 'shear'. By using redshift and shear measurements from every galaxy, information on the growth of structure and the geometry of the Universe

JCAP05(2010)035
can be extracted from 3D cosmic shear observations. Here we will use the 3D cosmic shear approach [30][31][32] where the full 3D shear field is characterised using 3D spherical harmonics and the Fisher matrix methodology of [32]. In line with the cosmic-variance limited approach of this article, we assume a large, cosmic-variance limited weak lensing survey covering 40,000 square degrees, to a median redshift of 3 with 50 galaxies per square arcminute. On the right panel of figure 4, we show the marginalised constraints ∆ and Σ, for this cosmic-variance limited survey. We find that the sum of neutrino mass is constrained to extreme accuracy O(10 −6 ). As the neutrino mass decreases the constraints on the IH and NH become smaller and for massess below ∼ 0.15 eV the evidence ratio for the IH and NH constraints (lower panels of figure 4) would become substantial (in a Jeffrey's scale), allowing for the neutrino heirachy to be distinguished. This is again a very challenging survey, and acts to show highlight how demanding the measurement of neutrino mass-splitting can be; however by using shear measurements from Euclid [33] or LSST [17] we may hope to approach this regime.
The degeneracies between Σ and ∆ are small, and the very small constraint on Σ results in the constraints being effectively un-correlated in the Σ-∆ plane. We note that the constraints on ∆ around the IH and NH peaks are tighter for weak lensing than LSS, this is due to lensing providing constraints on both the geometry and the growth of structure, which provides a smaller raw constraint and a more orthogonal constraint to the CMB resulting in smaller errors. Interestingly, even though the weak lensing constraints on ∆ are smaller than for the power spectrum, the evidence ratio is comparable, because, due to the multidimensional degeneracy directions, a naive correspondence between error-bars and evidence is not applicable (it is to a first approximation the difference between the two error bars that is important).
Note that the evidence calculation explicitly assumes two isolated peaks, and so is only applicable when the fiducial points are seperated by multiple-sigma. As a result of this, the evidence calculations may be slightly optimisic for large masses. However, for Σ < 0.2 eV, the χ 2 difference between the two minima becomes noticeable as well as the shift between the location of one of the two minima and the central ∆ value for the oscillations experiments (which induces an additional χ 2 difference). While this information is not fully accounted for in a Bayesian approach to forecasting the evidence, it may be included at the moment of analyzing the data, using different approaches such as the likelihood ratio, and may slightly improve the significance for the hierarchy determination.
While we have used the oscillation results to center the Fisher and evidence calculations on the NH and IH, combining the oscillation experiments constraints will not improve the evidence; in fact, oscillation experiments give symmetric errors on ∆ (i.e. they do not depend on the sign of ∆).

Conclusions
The shape of the matter power spectrum contains information, in order of decreasing sensitivity, about the sum of neutrino masses, the amplitude of the mass splitting and the hierarchy (i.e., the mass splitting order). We have introduced a novel parameterization of the neutrino mass hierarchy, ∆, that has the advantage of changing continuously between normal, degenerate and inverted hierarchies and whose sign changes between normal and inverted. The absolute value of ∆ describes the maximum mass difference between the eigenstates. We stress that, current constraints from neutrino oscillations have ruled out large part of the parameter space given by the sum of the masses and the ∆ parameter, leaving two narrow

JCAP05(2010)035
Are neutrinos their own anti-particle? ( Determine COSMOLOGY COSMOLOGY unknown >0.25eV Figure 5. Role of cosmology in determining the nature of neutrino mass. Future neutrinoless double beta decay (0νβbeta) experiments and future cosmological surveys will be highly complementary in addressing the question of whether neutrinos are Dirac or Majorana particles. Next generation means near future experiments whose goal is to reach a sensitivity to the neutrinoless double beta decay effective mass of 0.01 eV. We can still find two small windows where this combination of experiments will not be able to give a definite answer, but this region is much reduced by combining 0νβbeta and cosmological observations. regions: for a fixed value of the total mass, the value of ∆ for the normal hierarchy is related to that of the inverted one and ∆ NH ≃ −∆ IH (but, in detail, ∆ NH ≡ |∆ IH |). It is the allowed region that cosmology should explore. We found that the information about ∆ accessible from the power spectrum shape yields a degeneracy: parameters values ∆ and −∆ yield nearly identical power spectra and therefore that the likelihood surface in ∆ is bimodal. This was not noted in the literature before and not taking this into account when using the Fisher matrix-approach to forecast future surveys performance may lead to spurious indications of a surveys ability to determine the hierarchy.
Detecting the signature of the hierarchy in the sky is therefore extremely challenging, and therefore we asked: "can cosmology in the cosmic-variance limit, and for an ideal experiment, distinguish the neutrino heirarchy?" or in other words, "is there enough information in the sky to measure the neutrino hierarchy?" To address these questions we have considered ideal, full-sky, cosmic variance-limited surveys and found that substantial Bayesian evidence (ln B ≥ 1) can be achieved assuming a minimal LCDM base model. Increasing the parameter space by including e.g. effective number of neutrino species, non-inflationary motivated shape of the primordial power spectrum, etc. weakens the constraints and makes a determination of the hierarchy not possible, even though the constraint on the total neutrino mass are less JCAP05(2010)035 affected and remain still at an interesting level. Are such a surveys feasible in the next 5-10 years? There are two candidates for such surveys: a full extragalactic survey in the optical/infrared like Euclid 3 [33] and a full 21cm survey by the SKA. 4 Each of these surveys is scheduled to start operations by 2018. Euclid will make an all sky Hubble-quality map for weak lensing and will directly trace the dark matter using this technique; whilst the cosmic variance limited survey we consider here is ambitous with respect to this survey these result serve as a qualitative measure of this surveys expected performance (costraints should be only a factor ≤ 1.5 larger at worst). Euclid will also target emission line galaxies up to z ∼ 3 (therefore these galaxies will have bias parameter close to 1) however nP , quantifying the the ratio of the signal to shot noise, will be only slightly above 1. The 21cm surveys provide the most un-biased indirect tracer of the dark matter distribution in the Universe and have negligible shot noise.
For the degenerate and inverted mass spectra, the next generation neutrinoless double beta decay experiments can determine if neutrinos are their own anti-particle. For the normal hierarchy, the effective electron-neutrino mass may even vanish. However, if the large-scale structure cosmological data, improved data on the tritium beta decay, or the long-baseline neutrino oscillation experiments establish the degenerate or inverted mass spectrum, the null result from such double-beta decay experiments will lead to a definitive result pointing to the Dirac nature of the neutrino mass. This is summarized in figure 5.
If the small mixing in the neutrino mixing matrix is negligible, cosmology might be the most promising arena to help in this puzzle. Our work shows that depending on the total neutrino mass, there might be substantial evidence by cosmological data to infer the neutrino hierarchy.