Status of neutrino properties and future prospects - Cosmological and astrophysical constraints

Cosmological observations are a powerful probe of neutrino properties, and in particular of their mass. In this review, we first discuss the role of neutrinos in shaping the cosmological evolution at both the background and perturbation level, and describe their effects on cosmological observables such as the cosmic microwave background and the distribution of matter at large scale. We then present the state of the art concerning the constraints on neutrino masses from those observables, and also review the prospects for future experiments. We also briefly discuss the prospects for determining the neutrino hierarchy from cosmology, the complementarity with laboratory experiments, and the constraints on neutrino properties beyond their mass.


Introduction
Flavour oscillation experiments have by now firmly established that neutrinos have a mass. Current experiments measure with great accuracy the three mixing angles, as well as the two mass-squared splittings between the three active neutrinos. In the framework of the standard model (SM) of particle physics neutrinos are massless, and consequently do not mix, since it is not possible to build a mass term for them using the particle content of the SM. Therefore, flavour oscillations represent the only laboratory evidence for physics beyond the SM. Several unknowns in the neutrino sector still remain, confirming these particles as being the most elusive within the SM. In particular, the absolute related to the mechanism of mass generation. Even though they are not the focus of this review, we will briefly touch some of this aspects in the final sections of the review.
Cosmological data have reached a very good level of maturity over the last decades. Measurements of the CMB anisotropies from the Planck satellite have put the tightest constraints ever on cosmological parameters from a single experiment [ ], dramatically improving the constraints from the predecessor satellite WMAP [ ]. From the ground, the Atacama Cosmology Telescope (ACT) polarization-sensitive receiver and the South Pole Telescope (SPT) have been measuring with incredible accuracy CMB anisotropies at the smallest scales in temperature and polarization [ , ]. At degree and sub-degree scales, the BICEP/Keck collaboration [ , ] and the POLARBEAR telescope [ ] are looking at the faint CMB "B-mode" signal, containing information about both the early stages of the Universe (primordial B-modes) and the late time evolution (lensing B-modes). The Cosmology Large Angular Scale Surveyor [ ] is working at mapping the CMB polarization field over % of the sky. The SPIDER balloon [ ] successfully completed its first flight and is in preparation for the second launch likely at the end of . The CMB data are often complemented with information from the large-scale structure of the universe. The SDSS III-BOSS galaxy survey have recently released its last season of data [ ]. Extended catalogues of galaxy clusters have been completed from several surveys (see e.g. [ ] and references therein). In addition, weak lensing surveys (Canada-France-Hawaii Telescope Lensing Survey [ ], Kilo-Degree Survey [ ], Dark Energy Survey [ ]) are mature enough to provide constraints on cosmological parameters that are competitive with those from other observables. They also allow to test the validity of the standard cosmological paradigm by comparing results obtained from high redshift observables to those coming from measurements of the low redshift universe.
The current scenario is just a taste of the constraining power of cosmological observables that will be available with the next generation of experiments, that will be taking measurements in the next decade. ]-will test the universe over a wide range of scales with unprecedented accuracy. The same accuracy will enable the reconstruction of the weak lensing signal from the CMB maps down to the smallest scales and with high sensitivity, providing an additional probe of the distribution and evolution of structures in the universe. On the other hand, the new generation of large scale structure surveys -including the Dark Energy Spectroscopic Instrument [ ], the Large Synoptic Survey Telescope [ ], Euclid [ ] and the Wide Frequency InfraRed Spectroscopic Telescope [ ] -will also probe the late-time universe with the ultimate goal of shedding light on the biggest unknown of our times, namely the nature of dark energy and dark matter.
The aim of this review is to provide the state of the art of the current knowledge of neutrino masses from cosmological probes and give an overview of future prospects. The review is organized as follows: in Sec. , we outline the framework of this review, introducing some useful notation and briefly reviewing the basics of neutrino cosmology. Sec. is devoted to discussing, from a broad perspective, cosmological effects induced by massive neutrinos. In Sec. , we will describe in detail how the effects introduced in Sec. affect cosmological observables, such as the CMB anisotropies, large scale structures and cosmological distances. Sec. and Sec. present a detailed collection of the current and future limits on Σm ν from the measurements of the cosmological observables discussed in Sec. , mostly derived in the context of the ΛCDM cosmological model. Constraints derived in more extended scenarios are summarized in Sec. . Sec. briefly deals with the issue of whether cosmological probes are able to provide information not only on Σm ν , but also on its distribution among the mass eigenstates, i.e. about the neutrino hierarchy. In Sec. , we will briefly go through the complementarity between cosmology and laboratory searches in the quest for constraining neutrino properties. Finally, Sec. offers a summary of the additional information about neutrino properties beyond their mass scale that we can extract from cosmological observables. We derive our conclusions in Sec. . The impatient reader can access the summary of current and future limits from Tables , , and .

Notation and preliminaries . Basic equations
Inferences from cosmological observations are made under the assumption that the Universe is homogeneous and isotropic, and as such it is well described, in the context of general relativity, by a Friedmann-Lemaitre-Robertson-Walker (FLRW) metric. Small deviations from homogeneity and isotropy are modelled as perturbations over of a FLRW background.
In a FLRW Universe, expansion is described by the Friedmann equation for the Hubble parameter H: where G is the gravitational constant, K parameterizes the spatial curvature , a is the cosmic scale factor and the ρ is the total energy density. This is given by the sum of the energy densities of the various components of the cosmic fluid. Considering cold dark matter (c), baryons (b), photons (γ), dark energy (DE) and massive neutrinos (ν), and introducing the redshift 1 + z = a −1 , the Friedmann equation can be recast as: where we have introduced the present value of the critical density required for flat spatial geometry ρ crit,0 ≡ 3H 2 0 /8πG (in general, we use a subscript 0 to denote quantities evaluated today), and the present-day density parameters Ω i = ρ i,0 /ρ crit,0 (since we will be always referring to the density parameters today, we omit the subscript 0 in this case). The scalings with (1 + z) come from the fact that the energy densities of nonrelativistic matter and radiation scale with a −3 and a −4 , respectively. For DE, in writing Eq. ( ) we have left open the possibility for an arbitrary (albeit constant) equation-of-state parameter w. In the case of neutrinos, since the parameter of their equation of state is not constant, we could not write a simple scaling with redshift, although this is possible in limiting regimes (see Sec. . ). We use ρ ν,tot to denote the total neutrino density, i.e., summed over all mass eigenstates. Finally, we have defined a "curvature density parameter" Ω k = −K/H 2 0 . From Eq. ( ) evaluated at z = 0 it is clear that the density parameters, including curvature, satisfy the constrain i Ω i = 1.
Let us also introduce some extra notation and jargon that will be useful in the following. We will use Ω m to refer to the total density of nonrelativistic matter today. Thus, this in general includes dark matter, baryons and those neutrinos species that are heavy enough to be nonrelativistic today. In such a way we have that Ω m + Ω DE = 1 in a flat Universe (or Ω m + Ω DE = 1 − Ω k in general), since the present density of photons and other relativistic species is negligible. Since many times we will have to consider the density of matter that is nonrelativistic at all the redshifts that are probed by cosmological observables, i.e. dark matter and baryons but not neutrinos, we also introduce Ω c+b , with obvious meaning. When we consider dark energy in the form of a cosmological constant (w = −1) we use Ω Λ in place of Ω DE to make this fact clear. Finally, we also use the physical density parameters ω i ≡ Ω i h 2 , with h being the present value of the Hubble parameter in units of km s −1 Mpc −1 . As we shall discuss in more detail in the following, cosmological observables often carry the imprint of particular length scales, related to specific physical effects. For this reason we recall some definitions that will be useful in the following. The causal horizon r h at time t is defined as the distance traveled by a photon from the Big Bang (t = 0) until time t. This is given by: Note that this is actually the comoving causal horizon; in the following, unless otherwise noted, we will always use comoving distances. We also note that the comoving horizon is equal to the conformal time η(t) (defined through dt = adη and η(t = 0) = 0). In a Friedmann Universe (i.e., one composed only by matter and radiation), the physical causal horizon is proportional, by a factor of order unity, to the Hubble length d H (t) ≡ H(t) −1 . For this reason, we shall sometimes indulge in the habit of calling the latter the Hubble horizon, even though this is, technically, a misnomer. A related quantity is the sound horizon r s (t), i.e., the distance traveled in a certain time by an acoustic wave in the baryon-photon plasma. The expression for r s is very similar to the one for the causal horizon, just with the speed of light (equal to in our units) replaced by the speed of sound c s in the plasma: ( ) The speed of sound is given by c s = 1/ 3(1 + R), with R = (p b + ρ b )/(p γ + ρ γ ) being the baryon-to-photon momentum density ratio. When the baryon density is negligible relative to the photons, c s 1/ √ 3 and r s r h / √ 3 = η/ √ 3. Imprints on the cosmological observables of several physical processes usually depend on the value of those scales at some particular time. For example, the spacing of acoustic peaks in the CMB spectrum is reminescent of the sound horizon at the time of hydrogen recombination; the suppression of small-scale matter fluctuations due to neutrino freestreaming is set by the causal horizon at the time neutrinos become nonrelativistic; and so on. Moreover, since today we see those scales through their projection on the sky, what we observe is actually a combination of the scale itself and the distance to the object that we are observing. We find then useful also to recall some notions related to cosmological distances. The comoving distance χ between us and an object at redshift z is , ( ) and this is also equal to η 0 − η(z). The comoving angular diameter distance d A (z) is given by for Ω k = 0 .
( ) The angular size θ of an object is related to its comoving linear size λ through θ = λ/d A (z). This justifies the definition of an object of known linear size as a standard ruler for cosmology. In fact, knowing λ, we can use a measure of θ to get d A and make inferences on the cosmological parameters that determine its value through the integral in Eq. ( ). Another measure of distance is given by the luminosity distance d L (z), that relates the observed flux F to the intrinsic luminosity L of an object at redshift z: Similarly to what happened for the angular diameter distance, this allows to use standard candles -objects of known intrinsic luminosity -as a mean to infer the values of cosmological parameters, after their flux has been measured.

. Neutrino mass parameters
According to the standard theory of neutrino oscillations, the observed neutrino flavours ν α (α = e, µ, τ ) are quantum superpositions of three mass eigenstates ν i (i = 1, 2, 3): where U is the Pontecorvo-Maki-Nakagawa-Sasaka (PMNS) mixing matrix. The PMNS matrix is parameterized by three mixing angles θ 12 , θ 13 , θ 23 , and three CP-violating phases: one Dirac, δ, and two Majorana phases, α 21 and α 31 . The Majorana phases are non-zero only if neutrinos are Majorana particles. They do not affect oscillation phenomena, but enter lepton number-violating processes like 0ν2β decay. The actual form of the PMNS matrix is: In addition to the element of the mixing matrix, the other parameters of the neutrino sectors are the mass eigenvalues m i (i = 1, 2, 3). Oscillation experiments have measured with unprecedented accuracy the three mixing angles and the two mass squared differences relevant for the solar and atmospheric transitions, namely the solar splitting ∆m 2 sol = ∆m 2 21 ≡ m 2 2 − m 2 1 7.6 × 10 −5 eV 2 , and the atmospheric splitting ∆m 2 atm = |∆m 2 31 | ≡ |m 2 3 − m 2 1 | 2.5 × 10 −3 eV 2 (see [ , , ] for a global fit of the neutrino mixing parameters and mass splittings). We know, because of matter effects in the Sun, that, of the two eigenstates involved, the one with the smaller mass has the largest electron fraction. By convention, we identify this with eigenstate " ", so that the solar splitting is positive. On the other hand, we do not know the sign of the atmospheric mass splitting, so this leaves open two possibilities: the normal hierarchy (NH), where ∆m 2 31 > 0 and m 1 < m 2 < m 3 , or the inverted hierarchy, where ∆m 2 31 < 0 and m 3 < m 1 < m 2 .
Oscillation experiments are unfortunately insensitive to the absolute scale of neutrino masses. In this review, we will mainly focus on cosmological observations as a probe of the absolute neutrino mass scale. To a very good approximation, cosmological observables are mainly sensitive to the sum of neutrino masses Σm ν , defined simply as Absolute neutrino masses can also be probed by laboratory experiments. These will be reviewed in more detail in Sec. , where their complementarity with cosmology will be also discussed. For the moment, we just recall the definition of the mass parameters probed by laboratory experiments. The effective (electron) neutrino mass m β can be constrained by kinematic measurements like those exploiting the β decay of nuclei. The effective Majorana mass of the electron neutrino m ββ : can instead be probed by searching for 0ν2β decay.

. The standard cosmological model
Our best description of the Universe is currently provided by the spatially flat ΛCDM model with adiabatic, nearly scale-invariant initial conditions for scalar perturbations.
With the exception of some mild (at the ∼ 2σ level) discrepancies that will be discussed in the part devoted to observational limits, all the available data can be fit in this model, that in its simplest ("base") version is described by just six parameters. In the base ΛCDM model, the Universe is spatially flat (Ω k = 0), and the matter and radiation content is provided by cold dark matter, baryons, photons and neutrinos, while dark energy is in the form of a cosmological constant (w = −1). The energy density of photons is fixed by measurements of the CMB temperature, while neutrinos are assumed to be very light, usually fixing the sum of the masses to Σm ν = 0.06 eV, the minimum value allowed by oscillation experiments. In this way, the energy density of neutrinos is also fixed at all stages of the cosmological evolution (see Sec. . ). From Eq. ( ), and taking into account the flatness constraint, it is clear then that the background evolution in such a model is described by three parameters, for example h, ω c and ω b , with Ω Λ given by 1 − Ω m . The initial scalar fluctuations are adiabatic and have a power-law, nearly scale invariant, spectrum, that is thus parameterized by two parameters, an amplitude A s and a logarithimc slope n s − 1 (with n s = 1 thus corresponding to scale invariance). Finally, the optical depth to reionization τ parameterises the ionization history of the Universe. This simple, yet very successful, model can be extended in several ways. The extension that we will be most interested in, given the topic of this review, is a one-parameter extension in which the sum of neutrino masses is considered as a free parameter. We call this seven-parameter model ΛCDM+Σm ν . This is also in some sense the bestmotivated extension of ΛCDM, as we actually know from oscillation experiments that neutrinos have a mass, and from β-decay experiments that this can be as large as eV. In addition to this minimal extension, we will also discuss how relaxing some of the assumptions of the ΛCDM model affects estimates of the neutrino mass. Among the possibilities that we will consider, there are those of varying the curvature (Ω k ), the equation-of-state parameter of dark energy (w) or the density of radiation in the early Universe (N eff , defined in sec. . ).
There are many relevant extensions to the ΛCDM model that however we will not consider here (or just mention briefly). The most important one concerns the possibility of non-vanishing tensor perturbations, i.e. primordial gravitational waves, that, if detected, would provide a smoking gun for inflation. This scenario is parameterized through an additional parameter, the tensor-to-scalar ratio r. In the following, we will always assume r = 0. In any case, this assumption will not affect the estimates reported here, as the effect of finite neutrino mass and of tensor modes on the cosmological observables are quite distinct. Similarly, we will not consider the possibility of non-adiabatic initial perturbations, nor of more complicated initial spectra for the scalar perturbations, including those with a possible running of the scalar spectral index.

. Short thermal history
Given that cosmological observables carry the imprint of different epochs in the history of the Universe, we find it useful to shortly recall some relevant events taking place during the expansion, and their relation to the cosmological parameters. For our purposes, it is enough to start just when the temperature of the Universe was T ∼ 1 MeV, i.e. around the time of Big Bang Nucleosynthesis (BBN) and neutrino decoupling. At these early times (z ∼ 10 10 ), since matter and radiation densities scale as (1 + z) 3 and (1 + z) 4 , respectively, the Universe is radiation-dominated.
• At T ∼ 1 MeV (z ∼ 10 10 ), the active neutrinos decouple from the rest of the cosmological plasma. Before this time, neutrinos were kept in equilibrium by weak interactions with electrons and positrons, that were in turn coupled electromagnetically to the photon bath. After this time, their mean free path becomes much larger the the Hubble length, so they essentially move along geodesics, i.e., they free-stream. Shortly after neutrino decouple, electrons and positrons in the Universe annihilate, heating the photon-electron-baryon plasma, and, to a much lesser extent, the neutrino themselves (in the Sec. . we shall discuss in more detail the neutrino thermal history). After this time, the Universe can essentially be thought as composed of photons, electrons, protons and neutrons (either free, or, after BBN, bound together into the light nuclei), neutrinos, dark matter and dark energy.
• Soon after, at T ∼ 0.1 MeV, primordial nucleosynthesis starts, and nuclear reactions bind nucleons into light nuclei. After this time, nearly all of the baryons in the Universe are in the form of 1 H and 4 He nuclei, with small traces of 2 H and 7 Li. The yields of light elements strongly depend on the density of baryons, on the density and energy spectrum of electron neutrinos and antineutrinos (as those set the equilibrium of the nuclear reactions through which the nuclei are built) and on the total radiation density (as this sets the expansion rate at the time of nucleosynthesis).
• As said above, at early times (high z) the Universe is radiation-dominated, given that the ratio of radiation to matter scales like (1 + z). However, the radiation density decreases faster than that of matter, and, at some redshift z eq , the matter and radiation contents of the Universe will be equal: ρ m (z eq ) = ρ r (z eq ). This is called the epoch of matter-radiation equality, that marks the beginning of the matter-dominated era in the history of the Universe. From the scaling of the two densities, it is easy to see that 1 + z eq = Ω m /(Ω γ + Ω ν ) in a Universe with massless neutrinos (so that their density always scale as (1 + z) 4 ; see Sec. . for further discussion on this point) . Given the current estimates of cosmological parameters, z eq 3400 [ ].
• At T 0.3 eV, electrons and nuclei combine to form neutral hydrogen and helium, that are transparent to radiation. This recombination epoch thus roughly corresponds to the time of decoupling of radiation from matter. This is the time at which the cosmic microwave background (CMB) radiation is emitted. After decoupling, the CMB photons free-stream until the present time (with some caveats, see below). Most of the features that we observe in the CMB anisotropy pattern are created at this time. Given the current estimates of cosmological parameters, z rec 1090 [ ]. Note that in fact the temperature at recombination is basically fixed by thermodynamics, so once the present CMB temperature is determined through observations, z rec = T (z = z rec )/T (z = 0) depends very weakly on the other cosmological parameters.
• Even if photons decoupled from matter shortly after recombination, the large photon-to-baryon ratio keeps baryons coupled to the photon bath for some time after that. The drag epoch z drag is the time at which baryons stop feeling the photon drag. A good fit to numerical results in a CDM cosmology is given by [ ] Given the current estimates of cosmological parameters, z drag 1060 [ ].
• For a long time after recombination, the Universe stays transparent to radiation. These are the so-called "dark ages". However, in the late history of the Universe, the neutral hydrogen gets ionized again due to UV emission of the first stars, that puts an end to the dark ages. This is called the reionization epoch. After reionization, the CMB photons are scattered again by the free electrons. Given the current estimates of cosmological parameters, z re 8 [ ].
• At some point during the recent history of the Universe, that we denote with z Λ , the energy content of the Universe starts to be dominated by the dark energy component. The end of matter domination, and the beginning of this DE domination is set by ρ DE (z Λ ) = ρ m (z Λ ). For a cosmological constant (w = −1), 1 + z Λ = (Ω Λ /Ω m ) 1/3 . Around this time, the cosmological expansion becomes accelerated.

. Evolution of cosmic neutrinos
In this section, we discuss the thermal history of cosmic neutrinos. As anticipated above, in the early Universe neutrinos are kept in equilibrium with the cosmological plasma by weak interactions. The two competing factors that determine if equilibrium holds are the expansion rate, given by the Hubble parameter H(z), and the interaction rate Γ(z) = n σv , where n is the number density of particles, σ is the interaction cross section, and v is the velocity of particles (brackets indicate a thermal average). In fact, neutrino interactions become too weak to keep them in equilibrium once Γ < H. The left-hand side of this inequality is set by the standard model of particle physics, as the interaction rate at a given temperature only depends on the crosssection for weak interactions, and thus, ultimately, on the value of the Fermi constant (σ w ∼ G 2 F T 2 ). The right-hand side is instead set, through Eq. ( ) by the total radiation density (the only relevant component at such early times): H 2 = (8πG/3)(ρ γ +ρ ν ). This is also fixed at any given temperature, in the framework of the minimal ΛCDM model, so that the temperature of neutrino decoupling, defined through Γ(T ν,dec ) = H(T ν,dec ) does not depend on any free parameter in the theory. A quite straightforward calculation shows that T ν,dec 1 MeV [ ].
While they are in equilibrium, the phase-space distribution f (p) of neutrinos is a Fermi-Dirac distribution : where it has been taken into account that at T 1 MeV, the active neutrinos are certainly ultrarelativistic (i.e., T ν m ν ) and thus E(p) p. The distribution does not depend on the spatial coordinate x, nor on the direction of momentump, due to the homogeneity and isotropy of the Universe. Before decoupling, the neutrino temperature T ν is the common temperature of all the species in the cosmological plasma, that we denote generically with T , so that T ν = T . We recall that the temperature of the plasma evolves according to g 1/3 * s aT = const., where g * s counts the effective number of relativistic degrees of freedom that are relevant for entropy [ ].
Since decoupling happens while neutrinos are ultrarelativistic, it can be shown that, as a consequence of the Liouville theorem, the shape of the distribution function is preserved by the expansion. In other words, the distribution function still has the form Eq. ( ), with an effective temperature T ν (z) (that for the sake of simplicity we will continue to refer to as the neutrino temperature) that scales like a −1 (i.e., aT = const). We stress that this means that, when computing integrals over the distribution function, one still neglects the mass term in the exponential of the Fermi-Dirac function, even at times when neutrinos are actually nonrelativistic.
Shortly after neutrino decouple, electrons and positrons annihilate and transfer their entropy to the rest of the plasma, but not to neutrinos. In other words, while the neutrino temperature scales like a −1 , the photon temperature scales like a −1 g −1/3 * s , and thus decreases slightly more slowly during e + e − annihilation, when g * s is decreasing. In fact, applying entropy conservation one finds that the ratio between the neutrino and photon temperatures after electron-positron annihilation is T ν /T = (4/11) 1/3 . The photon temperature has been precisely determined by measuring the frequency spectrum of the CMB radiation: T 0 = (2.725 ± 0.002) K [ , ], so that the present temperature of relic neutrinos should be T ν,0 1.95 K 1.68 × 10 −4 eV.
The number density n ν of a single neutrino species (including both neutrinos and their antiparticles) is thus given by: is the Riemann zeta function of , and in the last equality we have taken into account that g = 2 for neutrinos. This corresponds to a present-day density of roughly 113 particles/cm 3 . The energy density of a single neutrino species is instead This is the quantity that appears, among other things, in the right-hand side of the Friedmann equation (summed over all mass eigenstates). In the ultrarelativistic (T ν m) and nonrelativistic (T ν m) limits, the energy density takes simple analytic forms: These scalings are consistent with the fact that one expects neutrinos to behave as pressureless matter, ρ ν ∝ (1 + z) 3 , in the nonrelativistic regime, and as radiation, ρ ν ∝ (1 + z) 4 , in the ultrarelativistic regime. Given that the present-day neutrino temperature is fixed by measurements of the CMB temperature and by considerations of entropy conservation, it is clear from the above formulas how the present energy density of neutrinos depends only on one free parameter, namely the sum of neutrino masses Σm ν defined in Eq. ( ). Introducing the total density parameter of massive neutrinos Ω ν ≡ i ρ ν i ,0 /ρ crit,0 , one easily finds from Eq. ( ): where we have already included the effects of non-instantaneous neutrino decoupling, see below. In the instantaneous decoupling approximation, the quantity at denominator would be 94.2 eV.
On the other hand, the neutrino energy density in the early Universe only depends on the temperature, and thus it is completely fixed in the framework of ΛCDM model. Using the fact that for photons ρ γ = (π 2 /15)T 4 , together with the relationship between the photon and neutrino temperatures, one can write for the total density in relativistic species in the early Universe, after e + e − annihilation: where N ν is the number of neutrino families. In the framework of the standard model of particle physics, considering the active neutrinos, one has N ν = 3. However, the above formula slightly underestimates the total density at early times; the main reason is that neutrino are still weakly coupled to the plasma when e + e − annihilation occurs, so that they share a small part of the entropy transfer. Moreover, finite temperature QED radiative corrections and flavor oscillations also play a role. This introduces nonthermal distortions at the subpercent level in the neutrino energy spectrum; the integrated effect is that at early times the combined energy densities of the three neutrino species are not exactly equal to 3ρ ν , with ρ ν given by the upper row of Eq. , but instead are given by (3.046ρ ν ) [ , ]. A recent improved calculation, including the full collision integrals for both the diagonal and off-diagonal elements of the neutrino density matrix, has refined this value to (3.045ρ ν ) [ ]. It is then customary to introduce an effective number of neutrino families N eff and rewrite the energy density at early times as: In this review, we will consider N eff = 3.046 as the "standard" value of this parameter in the ΛCDM model, and not the more precise value found in Ref. [ ], since most of the literature still makes use of the former value. This does not make any difference, however, from the practical point of view, given the sensitivity of present and nextgeneration instruments. It is also customary to consider extensions of the minimal ΛCDM model in which one allows for the presence of additional light species in the early Universe ("dark radiation"). In this kind of extension, the total radiation density of the Universe is still given by the right-hand side of Eq. ( ), where now however N eff has become a free parameter. In other words, Eq. ( ) becomes a definition for N eff , that is, just a way to express the total energy density in radiation. The effect on the expansion history of this additional radiation component can be taken into account by the substitution in the rhs of the Hubble equation ( ), with ∆N eff ≡ N eff − 3.046. Note that this substitution fully captures the effect of the additional species only if this is exactly massless, and not just very light (as in the case of a light massive sterile neutrino, for example -see Sec. ).
It is often useful, to understand some of the effects that we will discuss in the following, to have a feeling for the time at which neutrinos of a given mass become nonrelativistic, or, thinking the other way around, for the mass of a neutrino that becomes nonrelativistic at a given redshift. The average momentum of neutrinos at a temperature T ν is p = 3.15T ν . We take as the moment of transition from the relativistic to the nonrelativistic regime the time when p = m ν . Then, using the fact that T ν (z) = (4/11) 1/3 T 0 (1 + z) = 1.68 × 10 −4 (1 + z) eV, one has 1 + z nr 1900 m ν eV .
( ) This relation can be used to show e.g. that neutrinos with mass m 0.6 eV turn nonrelativistic at recombination. In the following, when discussing the effect of neutrino masses on the CMB anisotropies, we will assume that this is the case. Note however that the actual statistical analyses from which bounds on neutrino masses are derived do not make such an assumption. We also note that, given the current measurements of the neutrino mass differences, only the lightest mass eigenstate can still be relativistic today. Thus at least two out of the three active neutrinos become nonrelativistic at some time between recombination and the present.
We conclude this section with a clarification on the role of neutrinos in determining the redshift of matter-radiation equality. Given the present bounds on neutrino masses, we know that equality likely takes place when neutrino are relativistic. In fact, observations of the CMB anisotropies constrain z eq 3400, so that neutrinos with mass m 1.8 eV, just below the current bound from tritium beta-decay, turn nonrelativistic at equality. Thus, for masses sufficiently below the tritium bound, the total density of matter at those times is proportional to Ω c+b . The radiation density is instead provided by photons and by the relativistic neutrinos (and as such does not depend on the neutrino mass), plus any other light species present in the early Universe. So the redshift of equivalence is given by where the last equality makes it clear that, in the framework of the minimal ΛCDM model, the redshift of equivalence only depends on the quantity ω c + ω b , since N eff is fixed and ω γ is determined through observations (it is basically the CMB energy density).

Cosmological effects of neutrino masses
The impact of neutrino masses -and in general of neutrino properties -on the cosmological evolution can be divided in two broad categories: background effects, and perturbation effects. The former class refers to modifications in the expansion history, i.e. in changes to the evolution of the FLRW background. The latter class refers instead to modifications in the evolution of perturbations in the gravitational potentials and in the different components of the cosmological fluid. We shall now briefly review both classes; we refer the reader who is interested in a more detailed analysis to the excellent review by Lesgourgues & Pastor [ ].
To start, we shall consider a spatially flat Universe, i.e. Ω k = 0, in which dark energy is in the form of a cosmological constant (w = −1) and there are no extra radiation components (N eff = 3.046). Let us also consider a particular realization of this scenario, that we refer to as our reference model, in which the sum of neutrino masses is very small; for definiteness, we can think that Σm ν is equal to the minimum value allowed by oscillations, Σm ν = 0.06 eV. When needed, we will take the other parameters as fixed to their ΛCDM best-fit values from Planck [ ] . Our aim is to understand what happens when we change the value of Σm ν . Increasing the sum of neutrino masses Σm ν will increase ω ν = Ω ν h 2 according to Eq. ( ). Remember that the sum of the density parameters i Ω i = 1; this constraint can be recast in the form: Since ω γ is constrained by observations, and ω k is zero by assumption, we have four degrees of freedom that we can use to compensate for the change in ω ν , namely: increase h, or decrease any of ω c , ω b or ω Λ . For the moment, for simplicity, we will not distinguish between baryons and cold dark matter, pretending that as nonrelativistic components they have the same effect on cosmological observables. This is of course not the case, but we will come back to this later. Then we are left with three independent degrees of freedom that we can use to compensate for the change in ω ν : h, ω b+c , and ω Λ . We prefer to use Ω Λ in place of ω Λ , so that in the end our parameter basis for this discussion will be {h, ω c+b , Ω Λ }. The first option, increasing the present value of the Hubble constant while keeping Ω Λ and ω b+c has the effect of making the Hubble parameter at any given redshift after neutrinos become nonrelativistic larger with respect to the reference model. This can be understood by looking at Eq. ( ), that we rewrite here in this particular case With respect to the reference model, the first three terms in the RHS are unchanged, while the fourth increases because Ω Λ is fixed but h is larger. The last term does not depend on h (because the factor H 2 0 in front of the square brackets cancel the one in the critical density) but yet increases because ρ ν = Σm ν n ν is larger as long as neutrinos are in the nonrelativistic regime. On the other hand, before neutrino become nonrelativistic, ρ ν is the same in the two models, and the change in the Ω Λ h 2 term is irrilevant, because the DE density is only important at very low redshift. So we can conclude that at z z nr , the two models share the same expansion history, while for z z nr the model with "large" neutrino mass is always expanding faster (larger H), or equivalently, is always younger, at those redshifts. In terms of the length scales and of the distance measures introduced in Sec. . , it is easily seen that the causal and sound horizons at both equality and recombination (as well as at the drag epoch) are unchanged, because the expansion history between z = ∞ and z z nr is unchanged. On the other hand, distances between us and objects at any redshift -for example, the angular diameter distance to recombination -are always smaller than in the reference model, because H is always larger between z z nr and z = 0. H increases with the extra neutrino density, so this effect increases with larger neutrino masses (and moreover, z nr also gets larger for larger masses). Given this, we expect for example the angle subtended by the sound horizon at recombination, θ s = r s (z rec )/d A (z rec ) to become smaller when we increase Σm ν . We conclude this part of the discussion that in this case the redshift of equality z eq does not change, since ω b+c is being kept constant, and neutrinos contribute to the radiation density at early times (see discussion at the end of the previous section).
If we instead choose to pursue the second option, i.e., we keep h and Ω Λ constant while lowering ω c+b , we are again changing the expansion history, but this time on a different range of redshifts. In fact, when neutrinos are nonrelativistic, the RHS of Eq. ( ) is unchanged, because the changes in the present-day densities of neutrinos and nonrelativistic matter perfectly compensate; this continues to hold as long as both densities scale as (1 + z) 3 , i.e., roughly for z < z nr . On the other hand, at z > z nr the neutrino density is the same as in the reference model, while the matter density is smaller, so H(z) is smaller as well. Finally, when the Universe is radiation dominated, the two models share again the same expansion history. Then in this scenario we change the expansion history, decreasing H, for z nr z z eq . The sound horizon at recombination increases, and so does the angular diameter distance, so one cannot immediately guess how their ratio varies. However, a direct numerical calculation shows that, starting from the Planck best-fit model, the net effect is to increase θ s , meaning that the sound horizon will subtend a larger angular scale on the sky when Σm ν increases. For what concerns instead the redshift of matter-radiation equality, it is immediate to see that it decreases proportionally to ω c+b , i.e., equality happens later in the model with larger Σm ν .
Finally, when Ω Λ is decreased, the main effect is to delay the onset of acceleration and make the matter-dominated era last longer. This has some effect on the evolution of perturbations, as we shall see in the following. For what concerns the expansion history, since the model under consideration and the reference model only differ in the neutrino mass and in the DE density, they are identical when neutrinos are relativistic and DE is negligible, i.e., at z > z nr . For z < z nr , instead, starting as usual from Eq. ( ) one finds, with some little algebra, that H(z) is always larger in the model with smaller Ω Λ and larger Σm ν . As in the previous case, both r s and r A at recombination vary in the same direction (decreasing in this case); the net effect is again that θ s becomes larger with Σm ν . Also, since the matter density at early times is not changing in this case, the redshift of equivalence is the same in the two models.
We know comment briefly about ω b . One could choose to modify ω b in place of ω c in order to compensate for the change in ω ν . From the point of view of the background expansion, both choices are equivalent, since the baryon and cold dark matter density only enter through their sum ω b+c in the RHS of Eq. ( ). However, changing the baryon density also produces some peculiar effects, mainly related to the fact that i) it determines the BBN yields, and ii) it affects the evolution of photon perturbations prior to recombination. Thus the density of baryons is quite well constrained by the observed abundances of light elements and by the relative ratio between the heights of odd and even peaks in the CMB, (see Sec. . ) and there is little room for changing it without spoiling the agreement with observations.
Let us now turn to discuss the effects on the evolution of perturbations. Given that we have observational access to the fluctuations in the radiation and matter fields, it is useful to discuss separately these two components. The photon perturbations are sensitive to time variations in the gravitational potentials along the line of sight from us up to the last-scattering surface; this is the so-called integrated Sachs-Wolfe (ISW) effect. The gravitational potentials are constant in a purely matter-dominated Universe, so that the observed ISW gets an early contribution right after recombination, when the radiation component is not yet negligible, and a late contribution, when the dark energy density begins to be important. Coming back to our previous discussion, it is clear to see how delaying the time of equality will increase the amount of early ISW, while anticipating dark energy domination will increase the late ISW, and viceversa. For what concerns matter inhomogeneities, a first effect is again related to the time of matter-radiation equality. Changing z eq affects the growth of perturbations, since most of the growth happens during the matter dominated era. Apart from that, a very peculiar effect is related to the clustering properties of neutrinos. In fact, while neutrinos are relativistic, they tend to free stream out of overdense regions, damping out all perturbations below the horizon scale. The net effect is that neutrino clustering is exponentially suppressed below a certain critical scale, the free-streaming scale, that corresponds to the size of the horizon at the time of the transition from the ultrarelativistic to the nonrelativstic regime. If the transition happens during matter domination, this is given by: ( ) On the contrary, above the free-streaming scale neutrinos cluster as dark matter and baryons do. Thus, increasing the neutrino mass and consequently the neutrino energy density will suppress small-scale matter fluctuations relative to the large scales. It will also make small-scale perturbations in the other components grow slower, since neutrino do not source the gravitational potentials at those scales. It should also be noted that the free-streaming scale depends itself on the neutrino mass -specifically, heavier neutrinos will become nonrelativistic earlier and the free-streaming scale will be correspondingly smaller. Moreover, there is actually a free-streaming scale for each neutrino species, each depending on the individual neutrino mass. In principle one could think to go beyond observing just the small-scale suppression and try to access instead the scales around the nonrelativistic transition(s), in order to get more leverage on the mass and perhaps also on the mass splitting. We shall see however in the following that this is not the case.
The suppression of matter fluctuations due to neutrino free-streaming also affects the path of photons coming from distant sources, since those photons will be deflected by the gravitational potentials along the line of sight, resulting in a gravitational lensing effect. This is relevant for the CMB, as it modifies the anisotropy pattern by mixing photons that come from different directions. Another application of this effect, of particular importance for estimates of neutrino masses, is to use the distortions of the shape of distant galaxies due to lensing, to reconstruct the intervening matter distribution.

Cosmological observables
In this section we review the various cosmological observables, and explain how the effects described in the previous section propagate to the observables.

. CMB anisotropies
The CMB consists of polarized photons that, for the most part, have been free-streaming from the time of recombination to the present. The pattern of anisotropies in both temperature (i.e., intensity) and polarization thus encodes a wealth of information about the early Universe, down to z = z rec 1100. Moreover, given that the propagation of photons from decoupling to the present is also affected by the cosmic environment, the CMB also has some sensitivity to physics at z < z rec . Two relevant examples for the topic under consideration are the CMB sensitivity to the redshift of reionization (because the CMB photons are re-scattered by the new population of free electrons) and to the integrated matter distribution along the line of sight (because clustering at low redshifts modifies the geodesics with respect to an unperturbed FLRW Universe, resulting in a gravitational lensing of the CMB, see next section). However, the CMB sensitivity to these processes is limited due to the fact these are integrated effects.
The information in the CMB anisotropies is encoded in the power spectrum coefficients C T T , i.e., the coefficients of the expansion in Legendre polynomials of the two-point correlation function. In the case of the temperature angular fluctuations ∆T (n)/T : Gaussian fluctuations, all the information contained in the anisotropies can be compressed without loss in the two-point function, or equivalently in its harmonic counterpart, the power spectrum. A similar expression holds for the polarization field and for its cross-correlation with temperature. In detail, the polarization field can be decomposed into two independent components, known as E− (parity-even and curl-free) and B− (parity-odd and divergence-free) modes. Given that, it is clear that we can build a total of six spectra C XY with X, Y = T, E, B; however, if parity is not violated in the early Universe, the T B and EB correlations are bound to vanish. Let us also recall that, in linear perturbation theory, B modes are not sourced by scalar fluctuations. Thus, in the framework of the standard inflationary paradigm, primordial B modes can only be sourced in the presence of tensor modes, i.e., gravitational waves. The shape of the observed power spectra is the result of the processes taking place in the primordial plasma around the time of recombination. In brief, in the early Universe, standing, temporally coherent acoustic waves set in the coupled baryon-photon fluid, as a result of the opposite action of gravity and radiation pressure [ ]. Once the photons decouple after hydrogen recombination, the waves are "frozen" and thus we observe a series of peaks and throughs in the temperature power spectrum, corresponding to oscillation modes that were caught at an extreme of compression or rarefaction (the peaks), or exactly in phase with the background (the throughs). The typical scale of the oscillations is set by the sound horizon at recombination r s (z rec ), i.e. the distance travelled by an acoustic wave from some very early time until recombination, see Eq. . The position of the first peak in the CMB spectrum is set by the value of this quantity and corresponds to a perturbation wavenumber that had exactly the time to fully compress once. The second peak corresponds to the mode with half the wavelength, that had exactly the time to go through one full cycle of compression and rarefaction, and so on. Thus, smaller scales (larger multipoles) than the first peak correspond to scales that could go beyond one full compression, while larger scales (smaller multipoles) did not have the time to do so. In fact, scales much above the sound horizon are effectively frozen to their initial conditions, provided by inflation. This picture is complicated a little bit by the presence of baryons, that shift the zero of the oscillations, introducing an asymmetry between even and odd peaks. Finally, the peak structure is further modulated by an exponential suppression, due to the Silk damping of photon perturbations (further related to the fact that the tight coupling approximation breaks down at very small scales). This description also holds for polarization pertubations, with some differences, like the fact that the polarization perturbations have opposite phase with respect to temperature perturbations.
As noted above, the large-scale temperature fluctuations, that have entered the horizon very late and did not have time to evolve, trace the power spectrum of primordial fluctuations, supposedly generated during inflation. On the contrary, since there are no primordial polarization fluctuations, but those are instead generated at the time of recombination and then again at the time of reionization, the polarization spectra at large scales are expected to vanish, with the exception of the so-called reionization peak.
We can now understand how the CMB power spectra are shaped by the cosmological parameters, in a minimal model with fixed neutrino mass. The overall amplitude and slope of the spectra are determined by A s and n s , since these set the initial conditions for the evolution of perturbations. The height of the first peak strongly depends on the redshift of equivalence z eq (that sets the enhancement in power due to the early ISW), while its position is determined by the angle θ s subtended by the sound horizon at recombination. As we have discussed before, z eq and θ s are in turn set by the values of the background densities and of the Hubble constant. The baryon density further affects the relative heights of odd and even peaks, and also the amount of damping at small scales, through its effect on the Silk scale. The ratio of the densities of matter and dark energy fixes the redshift of dark energy domination and the amount of enhancement of large-scale power due to the late ISW. Finally, the optical depth at reionization τ induces an overall power suppression, proportional to e −2τ , in all spectra, at all but the largest scales. This can be easily understood as the effect of the new scatterings effectively destroying the information about the fluctuation pattern at recombination, at the scales that are inside the horizon at reionization. Reionization also generates the large-scale peak in the polarization spectra, described above. Measuring the power spectra gives a precise determination of all these parameters: simplifying a little bit, the overall amplitude and slope give A s e −2τ and n s (the latter especially if we can measure a large range of scales), the ratio of the peak heights and the amount of small-scale damping fix ω b , while the position and height of the first peak fix θ s and z eq , and thus h and ω b+c . The polarization spectra further help in that they are sensitive to τ directly, allowing to break the A s − τ degeneracy, and that the peaks in polarization are sharper and thus allow, in principle, for a better determination of their position [ ]. It is clear that adding one more degree of freedom to this picture, for example considering curvature, the equation of state parameter of dark energy, or the neutrino mass as a free parameter, will introduce parameter degeneracies and degrade the constraints.
Coming to massive neutrinos, as we have discussed in Sec , there is a combination of the following effects when Σm ν , and consequently ω ν , is increased, depending on how we are changing the other parameters to keep Ω = 1: i) an increase in θ s ; ii) a smaller z eq and thus a longer radiation-dominated era; iii) a delay of the time of dark energy domination. These changes will in turn result in: i) a shift towards the left of the position of the peaks; ii) an increased height of the first peak, that is set by the amount of early ISW; iii) less power at the largest scales, due to the smaller amount of late ISW. A more quantitative assessment of this effects can be obtained using a Boltzmann code, like CAMB [ ] or CLASS [ ], to get a theoretical prediction for the CMB power spectra in presence of massive neutrinos. These are shown in Fig. . In the left panel we plot the unlensed CMB temperature power spectra for a reference model with Σm ν = 0.06 eV (ω ν 6.4 × 10 −4 ) (the other parameters are fixed to their best-fit values from Planck ) and for three models with Σm ν = 1.8 eV (ω ν 1.9 × 10 −2 ), where either h, ω c or Ω Λ are changed to keep Ω = 1. We consider three degenerate neutrinos with m = 0.6 eV each, so that they become nonrelativistic around recombination. We also show the ratio between these spectra and the reference spectrum in the right panel of the same figure.
These imprints are in principle detectable in the CMB, especially the first two, since the position and height of the first peak are very well measured; much less so the redshift of DE domination, due to the large cosmic variance at small 's. However, following the above discussion, it is quite easy to convince oneself that these effects can be pretty much canceled due to parameter degeneracies. In fact, simplifying again a little bit, in standard ΛCDM we use the very precise determinations of the height and position of the first peak to determine θ s and z eq , and from them ω c+b and h. In an extension with massive neutrinos, we still have the same determination of θ s and z eq , but we have to use them to fix three parameters, namely ω c+b , h and ω ν , so that the system is underdetermined. One could argue that the amount of late ISW, as measured by the large-scale power, could be used to break this degeneracy, as it would provide a further constraint on the matter density (given that the DE density is fixed by the flatness condition). Unfortunately, measurements of the large-scale CMB power are plagued by large uncertainties, due to cosmic variance, so they are of little help in solving this degeneracy. Given the experimental uncertainties, then, it is clear that, when trying to fit a theory to the data, there will be a strong degeneracy direction corresponding to models having the same θ s and z eq , and thus with identical predictions for the first peak, and slightly different values of z Λ , with very low statistical weight due to the large , Ω Λ (yellow, always below the green apart from the lowest 's) or ω c (blue). The model in blue has a smaller z eq with respect to the reference; the models in yellow and green have a larger θ s ; in addition, the yellow model also has a smaller z Λ . Right: Ratio between the models with Σm ν = 1.8 eV and the reference model. uncertainties in the corresponding region of the spectrum. In other words, the effects of neutrino masses will be effectively "buried' in the small-plateau, where experimental uncertainties are large. The situation is even worse in extended models, for example if we allow the spatial curvature or the equation of state of dark energy to vary [ ].
In any case, the degeneracy between h and ω c+b is not completely exact, so that the unlensed CMB still has some degree of sensitivity to neutrinos that were relativistic at recombination. For example, the Planck temperature data, in combination with high-resolution observations from ACT and SPT, were able to constrain Σm ν < 1.1 eV after marginalizing over the effects of lensing.

. . Secondary anisotropies and the CMB Lensing
As observed above, in addition to the features that are generated at recombination, the so-called primary anisotropies, the CMB spectra also carry the imprint of effects that are generated along the line of sight. We have already given an example of one of these secondary anisotropies when we have mentioned the re-scattering of photons over free electrons at low redshift, that creates the distinctive "reionization bump" in the low-region of the polarization spectra. Another important secondary anisotropy is the gravitational lensing of the CMB (see [ , ]): photon paths are distorted by the presence of matter inhomogeneities along the line of sight. In the context of General Relativity, the deflection angle α for a CMB photon is where χ * is the comoving distance to the last scattering surface, f K (χ) is the angulardiameter distance (Eq. ) thought as a function of the comoving distance, Ψ is the gravitational potential, η 0 − χ is the conformal time at which the photon was along the direction n. If we then define the lensing potential as it is straightforward to see that the deflection angle is the gradient of the lensing potential, α = ∇φ. From the harmonic expansion of the lensing potential, we can build an angular power spectrum as < φ m φ * m >≡ δ δ mm C φφ . The lensing power spectrum C φφ is therefore proportional to the integral along the line of sight of the power spectrum of the gravitational potential P Ψ , which in turn can be expressed in terms of the power spectrum of matter fluctuations P m (see the next section for its definition).
The net effect of lensing on the CMB is that photons coming from different directions are mixed, somehow "blurring" the anisotropy pattern. This effect is mainly sourced by inhomogeneities at z < 5 and has a typical angular scale of 2.5 . In the power spectra, this translates in a several percent level smoothing of the primary peak structure ( 1000), while the lensing effect becomes dominant at 3000. We stress that lensing only We are assuming that the lensing field is isotropic.
alters the spatial distribution of CMB fluctuations, while leaving the total variance unchanged. Lensing, being a non-linear effect, creates some amount of nongaussianity in the anisotropy pattern. Thus, other than through its indirect effect on the temperature and polarization power spectra (i.e. on the two-point correlation functions), lensing can be detected and measured by looking at higher order correlations, in particular at the four-point correlation function. In fact, in such a way it has been possible to directly measure the power spectrum C φφ of the lensing potential φ. Another consequence of the nonlinear nature of lensing is that it is able to source "spurious" B modes by converting some of the power in E polarization, thus effectively creating B polarization also in the absence of a primordial component of this kind. The latter effect represents an additional tool to enable the reconstruction of the lensing potential, especially for future CMB surveys. An alternative reconstruction technique is based on the possibility to cross-correlate the CMB signal with tracers of large-scale structures, such as Cosmic Infrared Background (CIB) maps, therefore leading to an "external" reconstruction [ ] (opposite to the "internal" reconstruction performed with the use of CMB-based only estimators [ , ]). The lensing power spectrum basically carries information about the integrated distribution of matter along the line of sight. Given the peculiar effect of neutrino freestreaming on the evolution of matter fluctuations, CMB lensing offers an important handle for estimates of neutrino masses. Since a larger neutrino mass implies a larger neutrino density and less clustering on small scales, because of neutrino free-streaming, the overall effect of larger neutrino masses is to decrease lensing. In the temperature and polarization power spectra, the result is that the peaks and throughs at high-'s are sharper. Concerning the shape of the lensing power spectrum, for light massive neutrinos the net effect is a rescaling of power at intermediate and small scales (see e.g. [ ]). Thus the lensing power spectrum is a powerful tool for constraining Σm ν and will probably drive even better constraints on Σm ν in future. In fact, it is almost free from systematics coming from poorly understood astrophysical effects, it directly probes the (integral over the line of sight of the) distribution of the total matter fluctuations (as opposed to what galaxy surveys do, as we will see in the next section) at scales that are still in the linear regime. Given a cosmological model, it is quite straightforward, using again CAMB or CLASS, to get a theoretical prediction for the lensing power spectrum, as well as for the lensing BB power spectrum. Note that non-linear corrections (see next section for further details) to the lensing potential are important in this case to get accurate large-scale BB spectrum coefficients [ ]. Additional corrections that take into account modifications to the CMB photon emission angle due to lensing can further modify the large-scale lensing BB spectrum [ ].

. Large scale structures . . Clustering
The clustering of matter at large scales is another powerful probe of cosmology. The clustering can be described in terms of the two-point correlation function, or, equivalently, of the power spectrum of matter density fluctuations: is the Fourier transform of the matter density perturbation at redshift z. Note that, contrarily to the CMB, that we are bound to observe at a single redshift (that of recombination), the matter power spectrum can, in principle, be measured at different times in the cosmic history, thus allowing for a tomographic analysis.
As for the CMB, the large-scale (small k's) part of the power spectrum traces the primordial fluctuations generated during inflation, while smaller scales reflects the processing taking place after a given perturbation wavenumber enters the horizon. A relevant distinction in this regard is whether a given mode enters the horizon before or after matter-radiation equality. Since subhorizon perturbations grow faster during matter domination, the matter power spectrum shows a turning point at a characteristic scale, corresponding to the horizon at z eq . Given that perturbations grow less efficiently also during DE domination, increasing z Λ produces a suppression in the power spectrum. Also, increasing h will make the horizon at a given redshift smaller; so the mode k that is entering the horizon at that redshift will be larger.
Varying the sum of neutrino masses has some indirect effects on the shape of matter power spectrum, related to induced changes in background quantities, similarly to what happens for the CMB. As explained in Sec. , increasing Σm ν while keeping the Universe flat has to be compensated by changing (a combination of) ω m , Ω Λ or h. This will in turn result in a shift of the turning point and/or in a change in the global normalization of the spectrum. This can be seen in Fig. , where we show the matter power spectra for the same models considered when discussing the background effects of neutrino masses on the CMB.
As it is for the CMB, these effects can be partly canceled due to parameter degeneracies. Neutrinos, however, have also a peculiar effect on the evolution of matter perturbations. This is due to the fact that neutrinos possess large thermal velocities for a considerable part of the cosmic history, so they can free-stream out of overdense regions, effectively canceling perturbations on small scales. In particular, one can define the free-streaming length at time t as the distance that neutrinos can travel from decoupling until t. The comoving free-streaming length reaches a maximum at the time of the nonrelativistic transition. This corresponds to a critical wavenumber k fs , given in equation ( ) for transitions happening during matter-domination, above which perturbations in the neutrino component are erased.
A first consequence of neutrino free-streaming is that, below the free-streaming scale, there is a smaller amount of matter that can cluster. This results in an overall suppression of the power spectrum at small scales, with respect to the neutrinoless case. Secondly, subhorizon perturbations in the non-relativistic (i.e., cold dark matter and baryons) components grow more slowly. In fact, while in a perfectly matter-dominated Universe, the gravitational potential is constant and the matter perturbation grows linearly with the scale factor, δ m ∝ a, in a mixed matter-radiation Universe the gravitational potential decays slowly inside the horizon. Below the free-streaming scale, neutrinos effectively behave as radiation; then in the limit in which the neutrino fraction f ν = Ω ν /Ω m is small, one has for k k fs δ m (k k fs ) ∝ a 1−(3/5)fν , ( ) while δ m ∝ a for k k fs . These two effects can be qualitatively understood as follows: if one considers a volume with linear size well below the free-streaming scale, this region will resemble a Universe with a smaller Ω m and a larger radiation-to-matter fraction than the "actual" (i.e., averaged over a very large volume) values. This yields a smaller overall normalization of the spectrum, as well as a larger radiation damping; the two effects combine to damp the matter perturbations inside the region. So, looking again at the full power spectrum, the net effect is that, in the presence of free-streaming neutrinos, power at small-scales is suppressed with respect to the case of no neutrinos. A useful approximation is It is useful to stress that since f ν is linear in Σm ν , we have the somehow counterintuitive result that the effects of free-streaming are more evident for heavier, and thus colder, neutrinos. The reason is simply that the asymptotic suppression of the spectrum depends only on the total energy density of neutrinos, as this determines the different amount of non-relativistic matter between small and large scales.
Until now, we have somehow ignored the role of baryons in shaping the matter power spectrum. In fact, on scales that enter the horizon after z drag , the baryons are effectively collisionless and behave exactly like cold dark matter. On the other hand, baryon perturbations at smaller scales, entering the horizon before z drag exhibit acoustic oscillations due to the coupling with photons. This causes the appearance of an oscillatory structure in the matter power spectrum. These wiggles in P m (k), that go under the name of baryon acoustic oscillations (BAO), have a characteristic frequency, related to the value of the sound horizon at z drag . Thus they can serve as a standard ruler and can be used very effectively in order to constrain the expansion history.
In more detail, the acoustic oscillations that set up in the primordial Universe produce a sharp feature in the two-point correlation function of luminous matter at the scale of the sound horizon evaluated at the drag epoch, r s (z d ) ≡ r d ; this sharp feature translates in (damped) oscillations in the Fourier transform of the two-point correlation function, i.e., the power spectrum. Measuring the BAO feature at redshift z allows in principle to separately constrain the combination d A (z)/r d , for measurements in the transverse direction with respect to the line of sight, or r d H(z) for measurements along the line of sight. An isotropic analysis instead measures, approximately, the ratio between the combination called the volume-averaged distance, and the sound horizon r d . Given that the value of the sound horizon is well constrained by CMB observations, measuring the BAO features, possibly at different redshifts, allows to directly constrain the expansion history, as probed by the evolution of the angular diameter distance d A (z) and of the Hubble function H(z), or of their average d V (z). In particular, it is straightforward to see that BAO measurements put tight constraints on the Ω m − H 0 r d plane, along a different degeneracy direction that it is instead probed by CMB [ , ]. Therefore, when estimating neutrino masses, the addition of BAO constraints to CMB data helps breaking the parameter degeneracies discussed in the previous section, yielding in general tighter constraints on this quantity. The linear matter power spectrum for a given cosmological model can be computed using a Boltzmann solver. However, comparison with observations is complicated by the nonlinear evolution of cosmic structures. Note that both CAMB and CLASS are able to handle non-linearities in the evolution of cosmological perturbations with the inclusion of non-linear corrections from the Halofit model [ ] calibrated over numerical simulations. In particular, for cosmological models with massive neutrinos, the preferred prescription is detailed in [ ].
From the observational point of view, P m (k, z) can be probed in different ways. In galaxy surveys, the -D spatial distribution of galaxies is measured, allowing to measure the two-point correlation function and to obtain an estimate of the power spectrum of galaxies P g (k, z). Since in this case one is measuring the distribution of luminous matter only, and not of all matter (including dark matter), this does not necessarily coincide with the quantity for which we have a theoretical prediction, i.e., P m ; in other words, galaxies are a biased tracer of the matter distribution. To take this into account, one relates the two quantities through a bias b(k, z): ( ) The bias is in general both a function of redshift and scale. If it is approximated as a scale-independent factor, then the presence of the bias only amounts to an overall rescaling of the matter power spectrum (at a given redshift). In this case, one marginalizes over the amplitude of the matter spectrum, effectively only using the information contained in its shape. A scale-independent bias is considered to be a safe approximation for the larger scales: as an example, for Luminous Red Galaxies sampled at an efficient redshift of 0.5 (roughly corresponding to the CMASS sample of the SDSS III-BOSS survey), a scale-independent bias is a good approximation up to k 0.
On the other hand, scale-dependent features are expected to appear on smaller scales. In this case, the bias can still be described using a few "nuisance" parameters, that are then marginalized over. In any case the exact functional form of the bias function, the range of scales considered, as well as prior assumptions on the bias parameters, are delicate issues that should be treated carefully. An additional complication arises from the fact that massive neutrinos themselves induce a scale-dependence feature in the bias parameter, due to the scale-dependent growth of structures in cosmologies with massive neutrinos [ , ]. It has to be mentioned that, at any given redshift, there exist a certain scale k NL below which the density contrast approaches the limit δ ∼ 1. In this regime, the evolution of cosmic structures cannot be completely captured by a linear theory of perturbations. The modelling of structures in the non-linear regime relies on numerical N-body simulations that must take into account the astrophysical and hydrodynamical processes at play at those scales. The level of complexity of N-body simulations has been increasing over the years, so that the physical processes included in the simulations and the final results are much closer to the observations than they used to be at the beginning. A recent example is given by the MassiveNuS simulations [ ], based on the Gadgetcode [ ] modified to include the effects of massive neutrinos, and the nuCONCEPT simulations [ ]. . Nevertheless, the uncertainties related to the non-linear evolutions of cosmological structures are still higher than those affecting the linear theory, therefore reducing the constraining power coming from the inclusion of those scales in cosmological analysis. In fact, the conservative choice of not including measurements at k < k NL is usually made when performing cosmological analysis. It is easy to understand that the scale entering the non-linear regime is smaller for higher redshifts.
Additional probes of P m are measurements of Lyman-α (Lyα) forests and -cm fluctuations (see e.g. [ ] and [ ] for reviews). Although they are promising avenues since they can probe the matter distribution at higher redshifts and smaller scales than those usually accessible with typical galaxy samples, they still have to reach the level of maturity required to take full advantage of their constraining power. The observation of high-redshift (z ∼ 2) quasars and in particular the measurement of their flux provides a powerful tool for cosmological studies. Indeed, the absorption of the Lyα emission from quasars by the intervening intergalactic medium -an observational feature known as "Lyα forest" -constitutes a tracer of the total matter density field at higher redshifts and smaller scales than those usually probed by galaxy surveys. Similarly to what is done for galaxy samples, one can compute a correlation function of the measured flux variation, or equivalently its power spectrum P Lyα . The latter is again proportional to the total P m via a bias parameter b Lyα . The Lyα bias factor is in general different from the galaxy bias, as each tracer of the underlying total matter distribution exhibits its own characteristics. The Lyα forest is ideally a powerful cosmological tool, being able to access high redshift. Therefore, at fixed scale k, the physics governing the Lyα spectrum is much closer to the linear regime than that related to the galaxy power spectrum. Furthermore, the redshift window probed by Lyα is complementary to that probed by traditional galaxy surveys, in a sense that at higher redshift the relative impact of dark energy on the cosmic inventory is much less. However, a reliable description of the astrophysics at play in the intergalactic medium is essential for deriving the theoretical model for the Lyα absorption features along the line of sight. This description heavily depends on hydrodynamical simulations that reproduce the behaviour of baryonic gas and on poorly known details of the reionization history. In addition, uncertainties in the theory of non-linear physics of the intergalactic medium at small scales can play a non-negligible role.
Finally, another tracer of the total matter fluctuations is represented by fluctuations in the -cm signal. The -cm line is due to the forbidden transition of neutral hydrogen (HI) between the two hyperfine levels of the ground state (spin flip) of the hydrogen atom. The observational technique resides in the possibility to measure the brightness temperature relative to the CMB temperature. Fluctuations in the -cm brightness are related to fluctuations in HI (or equivalently to the fraction of free electrons x e ), which in turn trace the matter fluctuations. Therefore, one can infer P m observationally by measuring the power spectrum of -cm fluctuations P 21−cm . Apart from the technological challenges associated with the detection of the -cm signal, the main source of systematics come from the difficulties to separate the faint -cm signal from the much brighter foreground contamination, mostly due to synchrotron emission from our own galaxy.

. . Cluster abundances
The variation of the number of galaxy clusters of a certain mass M with redshift dN (z, M )/dz is also a valid source of information about the evolution of the late time Universe (see e.g. [ ] for a review). The expected number of clusters to be observed in a given redshift window is an integral over the redshift bin of the quantity where Ω is the solid angle,χ is the so-called completeness of the survey (a measure of the probability that the survey will detect a cluster of a given mass M at a given redshift z) and dN dM (z, M ) is the mass function giving the number of clusters per unit volume. The latter can be predicted once a cosmological model has been specified. The quantity in Eq.( ) is thus directly sensitive to the matter density Ω m and to the current amplitude of matter overdensities, usually parametrized in terms of σ 8 , the variance of matter fluctuations within a sphere of 8 h −1 Mpc. As a result, this probe can be highly beneficial for putting bounds on Σm ν .
Extended catalogues of galaxy clusters have been published in the last decade by the Atacama Cosmology Telescope (ACT) [ , ], the South Pole Telescope (SPT) [ ] and the Planck [ ] collaborations. CMB experiments are in fact able to perform searches for galaxy clusters by looking for the thermal Sunyaev-Zeldovich (SZ) effect, the characteristic upward shift in frequency of the CMB signal induced by the inverse-Compton scattering of CMB photons off the hot gas in clusters. The redshift of cluster candidates is identified with follow-up observations, whereas their mass is usually inferred with X-ray observations or, more recently, calibrated through weak lensing. Regardless of how it is calibrated, the determination of the cluster mass is the largest source of uncertainty for the cluster count analysis, due to possibly imprecise assumptions about the dynamical state of the cluster and/or survey systematics. A common way to factorise the uncertainties related to the mass calibration is to introduce a mass bias parameter that relates the true cluster mass to the mass inferred with observations.

. . Weak lensing
The weak gravitational lensing effect is the deflection of the light emitted by a source galaxy caused by the foreground large-scale mass distribution (lens). The shape of the source galaxy therefore appears as distorted, i.e. it acquires an apparent ellipticity. The cosmic shear is the weak lensing effect of all the galaxies along the line of sight (see e.g. [ ] for a review). Weak lensing surveys offer the possibility to directly test the distribution of intervening matter at low redshifts, thus providing a powerful tool to investigate the late-time evolution of the Universe. By correlating the apparent shapes of source galaxies at different redshifts, one can compute the shear field γ(n, z) as a function of the angular positionn and redshift z. The shear field is usually decomposed in two components: the curl-free E-modes and the divergence-free B-modes. It can be shown that, in absence of systematics, the B-modes are expected to vanish, whereas the power spectrum of the E-modes is equivalent to the lensing power spectrum C φφ ( ). The integrated lensing potential has been defined in Eq. ( ) for a source located at recombination. The corresponding expression for a source at a generic redshift z can be obtained simply by substituting χ * with the comoving distance of the source.
Thus, the power spectrum of the lensing potential -which is due to intervening matter along the line of sight -is recovered from the measurements of the lensing-induced ellipticity of background galaxies; in a similar way, the lensing power spectrum is recovered from the redistribution of CMB photons due to the forming structures along the line of sight. As we have seen in Sec. . . , the spectrum of the lensing potential is a function of the matter power spectrum integrated along the line of sight. Therefore, it carries information about the distribution and growth of structures, representing a powerful tool for constraining Σm ν . It should be mentioned that the observed shear signal γ obs is a biased tracer of the true shear γ true . This effect, mostly due to noise in the pixels when galaxy ellipticity is measured, is usually taken into account by introducing a multiplicative bias m that relates γ true and γ obs : γ obs = (1 + m)γ true + c, where c is the additional noise bias [ ].
In addition, the shear signal can be cross-correlated with the angular distribution of foreground (lens) galaxies (the so-called galaxy-shear or galaxy-galaxy lensing crosscorrelation). This cross-correlation is a powerful way to overcome the limitations induced in the galaxy-galaxy auto-correlation by the unknown galaxy bias. Indeed, the galaxygalaxy lensing is basically a cross-correlation between the galaxy field and the total matter fluctuation field. Measurements of the galaxy-galaxy lensing cross spectrum can therefore help determine the form of the bias.
Cosmological constraints from weak lensing survey are often summarized in terms of bounds on Ω m and σ 8 . As an additional probe of the large scale structure in the Universe, weak lensing can be profitably used to constrain Σm ν .

. Supernovae Ia and direct measurements of the Hubble constant
Measurements of the distance-redshift relation of Supernovae Ia (SNIa) have provided the compelling evidence of the accelerated Universe [ , ]. SNIa are produced in binary stellar systems in which one of the stars is a white dwarf. Accreting matter from its companion, the white dwarf explodes once it reaches the Chandrasekhar mass limit. Therefore, SNIa are standard candles, because their absolute magnitude can be theoretically inferred from models of stellar evolution. A comparison between the absolute magnitude and the apparent luminosity yields an estimate of their luminosity distance d L (z). The expected value of d L in turn depends on the underlying cosmological model. The constraints coming from SNIa in the Ω m − Ω Λ plane are orthogonal to those obtained from CMB. As a result, the combination of the two probes is extremely efficient in breaking the degeneracy between the two parameters. For this reason, SNIa are very useful for constraining models of dark energy and/or arbitrary curvature. Nonetheless, constraints on Σm ν can benefit from the use of SNIa data, thanks to the improved bounds on Ω m .
As already discussed, the effect of light massive neutrinos on the background evolution of the universe can be also compensated by a change in the value of the Hubble constant H 0 . Therefore, it is clear that any direct measurements of H 0 can be highly beneficial for putting bounds on Σm ν . Direct measurements only rely on local distance indicators (i.e. redshift z ), therefore they are little or not-at-all sensitive to changes in the underlying cosmological model. In contrast, indirect estimates from high-redshift probes, such as primary CMB, can suffer from model dependency.
Direct measurements of H 0 are based on the geometric distance calibration of nearby Cepheids luminosity-period relation and the subsequent calibration of SNIa over Cepheids observed in the same SNIa galaxy hosts (see e.g. [ ] and references therein). The goal is to connect the precise geometric distances measured in the nearby universe (usually referred to as "anchors") with the distant SNIa magnitude-redshift relation in order to extract the estimate of H 0 . The main systematics are of course related to the calibration procedure. Further improvements on the precision of direct measurements of H 0 are expected to come once the precise parallaxes measurements from the Gaia satellite will be available.
Local measurements of H 0 are not directly sensitive to Σm ν . Besides, their results, in combination with cosmological probes, can break the degeneracy between cosmological parameters and improve constraints on Σm ν . The main example is in fact the possibility to break the strong (inverse) degeneracy between H 0 and Σm ν that affects CMB constraints.
Indirect estimates of H 0 can be obtained from CMB and BAO measurements. We have already seen in Sec. that the position and amplitude of the first acoustic peak in the CMB spectrum depends on H 0 in combination with other parameters. In addition, we shall mention that, once the BAO are calibrated with the precise determination of r d from CMB, measurements of d A /r d and Hr d (or d V /r d ) yields bounds on H 0 that are competitive with CMB estimates and direct measurements.
We finally mention an additional independent measurement of H 0 . The detection of gravitational wave (GW) signals emitted by merging compact objects (standard sirens) in combination with the observation of an electromagnetic counterpart has been proposed as a standard siren [ , ]. The GW waveform reconstruction allows for a determination of the luminosity distance to the source. Precise determinations of the source localization can lead to percent accuracy in the luminosity distance estimation. The observation of the electromagnetic counterpart of the GW event is then essential to determine the redshift to the source. The full combination of distance-redshift pair can finally be employed to constrain H 0 . In the absence of the detection of an electromagnetic counterpart, methods to infer the redshift from the GW event have been proposed, see e.g. [ ].

. Summary of the effects of neutrino masses
Before moving to report the current observational constraints, we find it useful to summarize the constraining power of different cosmological observables with respect to the neutrino mass. The discussion is somehow qualitative, also given the high-level complexity of the cosmological models. The purpose is also to underline the importance of combining different cosmological probes.
We start from the CMB. For the present discussion, it is useful to consider separately the information coming from the unlensed CMB (i.e. the primary CMB plus all the secondary effect with the exclusion of lensing) and that coming from the weak lensing of CMB photons. For what concerns the former, the sensitivity of the unlensed CMB to neutrino masses is somehow limited. This is mainly due to a geometrical degeneracy between h and ω ν thanks to which one can simultaneously change the two parameters (decreasing h and increasing ω ν ) to keep θ s constant, thus preserving the position of the first peak, with only limited changes to other part of the spectrum (especially changes in the low-region , where the sensitivity is limited by cosmic variance, induced by variations in Ω Λ ). The height of the first peak is preserved by keeping ω c fixed. Having access to the information contained in the CMB lensing, either through its effect on the temperature and polarization power spectra, or through a direct estimation of the lensing power spectrum, helps because Σm ν also affects the matter distribution and then the amplitude of the lensing potential at small scales. This helps breaking the degeneracy described above.
To illustrate this point, in the upper panel of Fig. we show the parameter correlations derived by an analysis of the Planck observations of the temperature, over a wide range of scale, and large-scale polarization anisotropies. We remember that this dataset contains some information about lensing through the high-part of the temperature power spectrum. The negative degeneracy between Σm ν and H 0 is particularly evident. Given that ω c and ω b are both measured quite well from the CMB, this also translates into a strong degeneracy with Ω m = (ω c + ω b )/h 2 and Ω Λ = 1 − Ω m . Among the other parameters, one can notice mild correlations with A s and τ . These are related to the small-scale effects related to the increased lensing in models with larger Σm ν . The overall amplitude of the spectrum A s e −2τ is very precisely determined by CMB observations. On the other hand, the lensing amplitude depends on A s but not on τ . So, the lensing amplitude can be kept constant by increasing both A s and ω ν . At this point τ has to be increased as well to preserve the scalar amplitude A s e −2τ .
Geometric measurements, like those coming from BAO, SNIa or direct measurements of H 0 , greatly help solving the geometrical degeneracy between H 0 and Σm ν . This is evident by comparing the (H 0 , Σm ν ) square in the lower panel of Fig. , where we show parameter correlations from an analysis of the same dataset as above, with the addition of BAO data, with the corresponding square in the upper panel. Measurements of large scale structures, and especially those that are directly sensitive to the total matter distribution at small scales, are very helpful, in that on the one hand they allow to further constrain Ω m , A s and n s and thus reduce degeneracies with these parameters; on the other hand, they allow to probe the regime in which neutrino free-streaming is important. Finally, it is also clear that a precise measurement of τ from a CMB experiment that is sensitive to the large-scale polarization (meaning that it can access a large fraction of the sky) will be highly beneficial.
We have focused our attention to the ΛCDM+Σm ν model. In extended dark energy models (as well as modified gravity models), for example for arbitrary equations of state of the dark energy fluid, the degeneracy between Σm ν and Ω Λ is amplified. Both massive neutrinos and dark energy-modified gravity affect the late time evolution of the universe, so that the individual effects on cosmological observables (mostly structures) can be reciprocally cancelled.

Current observational constraints on Σm ν
In this section we report current constraints on Σm ν from cosmological and astrophysical observations. These constraints are also summarized in Tab. for the reader's conve- . for the description of these datasets. The darker the color shade, the stronger the degeneracy between the corresponding parameter pair. In both panels, the third row and the third column correspond to the correlation coefficients between Σm ν and the remaining cosmological parameters. From the comparison between the two panels, it is clear that the inclusion of BAO data helps reduce the degeneracy between parameters (see e.g. the correlation between Σm ν and H 0 , Ω Λ ); in a few cases, in fact, the inclusion of BAO reverts the degeneracy (see e.g. the correlation between Σm ν and n s ).
nience. Unless otherwise stated, the results are obtained in the framework of a minimal one-parameter extension of the ΛCDM model with varying neutrino mass, dubbed ΛCDM+Σm ν , in which the three mass eigenstates are degenerate (m i = Σm ν /3). Given the sensitivity of current experiments, the degenerate approximation is appropriate. See Sec. for a more detailed discussion on this point.

. CMB
CMB observations are probably the most mature cosmological measurements. The frequency spectrum is known with great accuracy [ ]. Measurements of the power spectrum of CMB anisotropies in temperature are cosmic-variance limited down to very small scales ( ∼ 1500) and the quality of current CMB data in polarization is already good enough to tighten constraints on cosmological parameters [ , , , , ]. The next generation of CMB experiments will further improve our knowledge of CMB polarization anisotropies [ , , , , ]. The main systematics involved in CMB measurements are due to foreground contamination (atmospheric, galactic, extragalactic), calibration uncertainties and spurious effects induced by an imprecise knowledge of the instrument (see e.g. [ , , , , ] for a sample list of references). The tightest constraints on Σm ν from a single experiment come from the measurements of the Planck satellite [ ]. In the context of a one-parameter extension of the ΛCDM cosmological background, the state of the art after the data release was as follows. The combination of the measurements of the CMB temperature anisotropies up to the multipole 2500 (hereafter, "Planck TT") and the large scale ( < 30) polarization anisotropies (hereafter "lowP") leads to an upper bound of Σm ν < 0.72 eV at % CL. The inclusion of the small scale ( ≥ 30) polarization measurements (which we globally label as "Planck TE,EE") provides a tighter upper bound of Σm ν < 0.49 eV at % CL. This latter bound should be regarded as less conservative, as a small level of residual systematics could still affect the small scale polarization data.
The Planck collaboration also provides the most significant measurements of the CMB lensing potential power spectrum for the multipole range 40 < L < 400 (labeled as "lensing") [ ]. When this dataset is included in the analysis, the constraints on Σm ν become: Σm ν < 0.68 eV for Planck TT+lowP+lensing and Σm ν < 0.59 eV for Planck TT,TE,EE+lowP+lensing [ ]. When combining the lensing reconstruction data from Planck with the measurements of the CMB power spectra, it should be kept in mind that CMB power spectra as measured by Planck prefer a slightly higher lensing amplitude than that estimated with the lensing reconstruction. As a result, the bounds on Σm ν obtained by their combination have less weight for smaller values of Σm ν than the corresponding bounds obtained from CMB power spectra only. Nevertheless, higher values of Σm ν are still disfavoured.
In , new estimates of the reionization optical depth τ have been published by the Planck collaboration [ ], obtained from the analysis of the high frequency CMB maps, in still affected by unexplained systematics effects at large scales. The estimated % credible interval for τ coming from the EE−only low-data is τ = 0.055 ± 0.009. This estimate is lower than the corresponding interval obtained in from the analysis of the low frequency maps (τ = 0.067 ± 0.023), though the two estimates are well in agreement with each other. The lower value of τ has an impact on the constraints on Σm ν , due to the degeneracy between the optical depth and the amplitude of primordial perturbation A S , as they together fix the normalization amplitude A S e −2τ . A lower τ implies a lower A S and thus a lower lensing amplitude, leaving less room for large values of Σm ν (that would further reduce lensing). If the "lowP" dataset is replaced by the new estimate of τ (labeled as "SimLow"), the bounds improve as follows: Σm ν < 0.59 eV for Planck TT+SimLow and Σm ν < 0.34 eV for Planck TT,TE,EE+SimLow [ ]. .

Large-Scale Structure Data
Although the CMB is an extremely powerful dataset, multiple degeneracies between cosmological parameters limit the constraining power on Σm ν from CMB only, as seen in Sec . . Measurements of the large scale structures (LSS) can help solving these degeneracies. LSS surveys map the distribution and clustering properties of matter at later times (or equivalently at lower redshift) than those accessible with CMB data and are directly sensitive to cosmological parameters that CMB data can only constrain indirectly, such as the total matter abundance at late times (see e.g. [ ] for a review).
In this section, we gather constraints on Σm ν from different LSS probes alone and in combination with CMB data.
. . Baryon acoustic oscillations and the full shape of the matter power spectrum from the clustering of galaxies BAO measurements, obtained by mapping the distribution of matter at relatively low redshifts (z < 3) if compared to the redshifts relevant for CMB, constrain the geometry of the expanding universe, providing estimates of the comoving angular diameter distance d A (z) and the Hubble parameter H(z) at different redshifts (or an angle-averaged combination of the two parameters, d V (z) = [zd 2 A (z)/H(z)] 1/3 ). Therefore, BAO constrain cosmological parameters which are relevant for the late-time history of the Universe, helping break the degeneracy between those parameters and Σm ν .
BAO extraction techniques rely on the ability to localise the peak of the two-point correlation function of some tracer of the baryon density, or equivalently the locations of the acoustic peaks in the matter power spectrum, thus neglecting the information coming from the broad band shape of the matter power spectrum itself. In principle, the full shape (FS) of the matter power spectrum is a valuable source of information about clustering properties of the different constituents of the universe and their reciprocal interactions. In particular, full shape measurements of the power spectrum also provide estimates of the growth of structures at low redshifts through the anisotropies induced by the redshift-space distortions (RSD), usually encoded in the parameter f (z)σ 8 (z), where f (z) is the logarithmic growth rate and σ 8 (z) is the normalization amplitude of fluctuations at a given redshift in terms of rms fluctuations in a 8h −1 Mpc sphere.
In , the final galaxy clustering data from the Baryon Oscillation Spectroscopic Survey (BOSS) were released, as part of the Sloan Digital Sky Survey (SDSS) III . Joint consensus constraints on d A (z), H(z) and f (z)σ 8 (z) from BAO and FS measurements at three different effective redshifts (z eff = 0.38, 0.51, 0.61) are employed to derive constraints on Σm ν in combination with Planck TT,TE,EE+lowP [ ]. The % upper bound is Σm ν < 0.16 eV. When relaxing the constraining power coming from CMB weak lensing (through the rescaling of the lensing potential with the lensing amplitude A L ) and the RSD (through the rescaling of the f σ 8 parameter with the amplitude A f σ 8 ), the bound degrades up to Σm ν < 0.25 eV. When using the FS measurements, it has to be noted that the constraining power of this dataset is highly reduced if one considers that ) the majority of the pieces of information encoded in the FS usually comes from the small-scale region of the power spectrum, where the still imprecisely known non-linearities play a non-negligible role; ) the exact shape and scale-dependence of the bias b between the observed galaxy clustering and the underlying total matter distribution is still debated. Therefore, it is useful to disentangle BAO and FS measurements, to gauge the relative importance of the two measurements in constraining Σm ν . For a thorough comparison between the constraining power of the two datasets, we refer the reader to [ ] (see also [ , ] for analyses using older data), where the authors focus on recent BAO measurements and FS measurements. Here, we summarise the conclusion of the paper: "The analysis method commonly adopted [for FS measurements] results in their constraining power still being less powerful than that of the extracted BAO signal".

. . Weak lensing
The most recent weak lensing datasets have been released by the Kilo-Degree Survey (KiDS [ , ]) and the Dark Energy Survey (DES [ , ]). It is interesting to note that all of the aforementioned datasets provide results in terms of cosmological parameters which are slightly in tension with the corresponding estimates coming from CMB data (which we remind is a high-redshift probe). In particular, the values of Ω m and S 8 = σ 8 (Ω m /0.3) 0.5 inferred from weak lensing data are lower than the best fit obtained with CMB data. The significance of this tension is at ∼ 2σ level for KiDS and more than 1σ level for the -D marginalized constraints on Ω m and S 8 for DES (even though a more careful measure of the consistency between the two datasets in the full parameter space provides "substantial" evidence for consistency, see [ ] for details).
Weak lensing data tend to favour higher values of Σm ν than those constrained by CMB power spectrum data. In fact, lower values of Ω m and S 8 imply a reduced clustering amplitude, an effect that can be obtained by increasing the sum of neutrino masses.
In [ ], the combination of DES shear, galaxy and galaxy-shear spectra with Planck

TT+lowP and other cosmological datasets in agreement with CMB results (i.e. BAO from dFGS [ ], SDSS DR MGS [ ] and BOSS DR [ ], and luminosity distances from the Joint Lightcurve Analysis (JLA) of distant SNIa [ , ]) yields an upper bound at
% CL on the sum of the neutrino masses of Σm ν < 0.29 eV, almost % higher than the corresponding bound obtained dropping DES data (Σm ν < 0.245 eV). Interestingly enough, the DES collaboration shows that a marginal improvement in the agreement between DES and Planck data is obtained when the sum of the neutrino masses is fixed to the minimal mass allowed by oscillation experiments Σm ν = 0.06 eV.
To conclude this section, we also report the upper bound on Σm ν obtained by weak lensing only data from the tomographic weak lensing power spectrum as measured by the KiDS collaboration [ ]. They found Σm ν < 3.3 eV and Σm ν < 4.5 eV at % CL depending on the number of redshift bins retained in the analysis. These bounds are significantly broader than the constraints coming from CMB only data. Nevertheless, they come from independent cosmological measurements and still tighter than the constraints coming from kinematic measurements of β decay.

. . Cluster counts
An additional low-redshift observable is represented by measurements of the number of galaxy clusters as a function of their mass at different redshifts. Cluster number counts provide a tool to infer the present value of the matter density Ω m and the clustering amplitude σ 8 , to be compared with the equivalent quantities probed at higher redshift by the primary CMB anisotropies.
Depending on the prior imposed on the mass bias, cluster counts tend to prefer lower values of Ω m and σ 8 than the corresponding values obtained with primary CMB. The tension between the two datasets can be as high as 3. Recently, [ ] updated constraints on cosmological parameters, including Σm ν , from the SZ clusters in the Planck SZ catalogue, considering cluster count alone and in combination with the angular power spectrum of SZ sources. A comparison with bounds coming from primary CMB anisotropies is also performed. The combination of the two SZ probes (complemented with BAO measurements from Anderson et al to fix the underlying cosmology) confirms the discrepancy in Ω m and σ 8 at the level of 2.1σ and provides an independent upper limit on the sum of the neutrino masses of Σm ν < 1.47 eV at % CL. When combined with primary CMB, the bound reduces to Σm ν < 0.18 eV. This bound is slightly higher than Σm ν < 0.12 eV found by [ ] in absence of SZ data, as we should expect due to the aforementioned tension between SZ and primary CMB estimates of matter density and power.

. . Lyman-α forests
As all the datasets that probe the clustering of matter over cosmological distances, the Lyα power spectrum is sensitive to Σm ν primarly through the power suppression induced by massive neutrinos at small scales. The Lyα spectrum alone can constraint Σm ν at a level of eV (see e.g. [ ]). The constraining power of the Lyα spectrum is evident when it is combined with CMB data. In this case, the Lyα data are used for setting the overall normalization of the spectrum through their sensitivity to Ω m and σ 8 , whereas the CMB fixes the underlying cosmological parameters and helps break degeneracies between Ω m , σ 8 and Σm ν . Recently, Yèche et al [ ] reported constraints on Σm ν from the combination of the one-dimensional (i.e. angle-averaged) Lyα power spectra from the SDSS III-BOSS collaboration and from the VLT/XSHOOTER legacy survey (XQ-). When the power spectra are used alone (complemented with a gaussian prior on H 0 = (67.3 ± 1.0) km s −1 Mpc −1 ), the authors obtain Σm ν < 0.8 eV at % CL. The bounds dramatically improves to Σm ν < 0.14 eV when CMB power spectrum data from Planck TT+lowP are added to the analysis. The tightest bound on Σm ν from Lyα power spectrum comes from [ ], with Σm ν < 0.12 eV from Planck TT+lowP in combination with the Lyα flux power spectrum from BOSS-DR . Interestingly enough, in both analyses, the limit set by Lyα+Planck TT+lowP does not further improve when the Lyα spectra are combined instead with the full set of CMB data from Planck, including small scale CMB polarization (Planck TT,TE,EE+lowP), and with BAO data from dFGS, SDSS MGS, BOSS-DR . The BAO signal can be also extracted from the Lyα spectrum (see [ ] for a pivotal study), providing estimates of the comoving angular diameter distance d A (z) and of the Hubble parameter H(z) at redshift z 2. Recently, the SDSS III-BOSS DR collaboration reported measurements of the BAO signal at z = 2.33 from Lyα forest [ ]. The estimated values of d A and H are in agreement with a ΛCDM model (even though a slight tension with Planck primary CMB is present), although their precision is smaller than the precision obtained with galaxy-derived BAO measurements. Therefore, at present, the impact of Lyα-BAO data on simple extensions of the ΛCDM model is minimal.
We conclude that it is a conservative choice to take the constraints coming from Lyα with some caution (a similar comment applies to constraints coming from aggressive analysis of the broadband shape of the matter power spectrum from galaxy surveys), until this probe will reach the level of maturity comparable with other traditional cosmological probes.

. Local measurements of the Hubble constant and Supernovae Ia
The most recent estimate of the Hubble constant has been reported in [ ]. The authors improved over their previous measurement of H 0 from . % to . % thanks to an increased sample of reliable SNIa in nearby galaxies calibrated over Cepheids.
Their final estimate, based on the combination of three different anchors, is H 0 = (73.24 ± 1.74) km s −1 Mpc −1 , 3.2σ higher than the indirect estimate of H 0 from Planck TT+SimLow (3.4σ higher than Planck TT,TE,EE+SimLow) in the context of a ΛCDM cosmology with Σm ν = 0.06 eV. Previous analyses from the same authors also pointed to a ∼ 2σ tension between direct measurements of H 0 and indirect estimate from primary CMB anisotropies from Planck (although see [ ] for a re-analysis of the same dataset which slightly reduces the discrepancy to within 1σ agreement). A discussion about the possible reasons behind this discrepancy and ways to alleviate it invoking non-standard cosmological scenarios are beyond the scope of this work. We refer the reader to the dedicated works [ , , ] for further reading. Since the Hubble constant and the sum of neutrino masses are anti-correlated, given the tension between the two probes it is clear that the combination of direct measurements of H 0 with CMB data leads to a preference for smaller values of Σm ν with respect to CMB-only constraints. Indeed, several authors have pointed out the tight constraints on Σm ν for such a combination. As an example, [ ] showed that constraints on Σm ν can be as tight as Σm ν < 0.148 eV at % CL when Planck TT+lowP+BAO are complemented with a gaussian prior on H 0 equal to the estimate of the Hubble constant in [ ], to be compared with Σm ν < 0.186 eV from Planck TT+lowP+BAO only. When lowP is replaced by a gaussian prior on τ compatible with the new estimates from Sim-Low, these numbers change to Σm ν < 0.115 eV (Σm ν < 0.151 eV) with (without) the H 0 prior.
For the sake of completeness, we shall also mention that independent estimates of H 0 from BAO measurements conducted by the SDSS III-BOSS DR collaboration [ ] are in agreement with CMB estimates (see also [ ] for a recent discussion). See also Ref.
[ ] for an additional independent estimate of H 0 with a combination of clustering and weak lensing measurements from DES-Y with BAO and BBN data. A discussion about the combination of five independent measurements of H 0 from cosmological probes and local measurements is also reported in [ , ]. Finally, we report that a standard siren measurement of H 0 has been performed after the detection of the neutron star-neutron star merger GW [ , , ]. The Hubble constant has been constrained as H 0 = 70.0 +12.0 −8.0 km s −1 Mpc −1 at % CL. The accuracy of this determination is not comparable with the precise estimates of direct measurements and other cosmological constraints. However, the standard siren approach represents an additional independent estimate of H 0 and appears as a promising avenue as more GW events with electromagnetic counterparts are detected.
Concerning the inclusion of SNIa, the bounds from Planck TT+lowP improve from Σm ν < 0.72 eV to Σm ν < 0.33 eV at % CL when data from the Joint Lightcurve Analysis [ , ] are included . The most relevant systematics that affect SNIa measurements are related to the way in which SNIa light curves are standardized, with issues mostly arising from photometric calibrations and lightcurve fitting procedures. determination of the BAO signal across (d A /r d ) and along (Hr d ) the line-of-sight, respectively, at z = 1.85, and 16% and 9% determination of the same quantities at the highest redshift achievable with Lyα forest z = 3.55 [ ]. Even in the most conservative scenario when DESI BAO only (i.e. without including information from the broadband shape of the matter power spectrum and Lyα forests) are combined with future CMB experiments, the sensitivity on Σm ν greatly improves. It goes down to σ(Σm ν ) = 0.021 eV for CORE TT,TE,EE,PP+DESI BAO, forecasting a ∼ 3σ detection of Σm ν in the minimal mass scenario [ ]. In the case of S +DESI BAO [ ], σ(Σm ν ) is in the range [0.023 − 0.036] eV ( or [0.020 − 0.032] eV) with a prior of τ = 0.06 ± 0.01 (or τ = 0.060 ± 0.006, the expected sensitivity from Planck-HFI [ ]) and f sky = 0.40, depending on the S angular resolution and noise level. For a 1 resolution and a noise level lower than 2.5 µK · arcmin, σ(Σm ν ) could be further improved with a better measurement of τ down to the level of σ(Σm ν ) < 0.015 eV, that would guarantee a > 4σ detection of Σm ν in the minimal mass scenario.
The DESI mission will be complementary to the science goals of the Large Synoptic Survey Telescope (LSST), a Stage-IV ground-based optical telescope. The main science fields in which LSST will mostly operate are [ ]: "Inventory of the Solar System, Mapping the Milky Way, Exploring the Transient Optical Sky, and Probing Dark Energy and Dark Matter". These goals will be achieved by surveying a ∼ 30000 deg 2 area ( / of which in a "deep-wide-fast" survey mode) over years, in six bands (ugrizy), with incredible angular resolution (∼ 0.7 ), producing measurements of roughly billion stars and galaxies. Thanks to its peculiar observational strategy, LSST will provide multiple probes of the late-time evolution of the universe with a single experiment, namely, weak lensing cosmic shear, BAO in the galaxy power spectrum, evolution of the mass function of galaxy clusters, and a compilation of SNIa redshift-distances. The expected sensitivity on Σm ν [ ] is in the range σ(Σm ν ) = [0.030 − 0.070] eV, depending on the fiducial value of Σm ν assumed when performing forecasts (Σm fid ν = [0 − 0.66] eV). Larger fiducial values for the mass yield better sensitivity. These numbers include a marginalization over the uncertainties coming from an extended cosmological scenario, where a number of relativistic species different than 3.046, a non-zero curvature and a dynamical dark energy w 0 − w a component are allowed. They also take into account the combination of the three-dimensional cosmic shear field as measured by a LSST-like survey with Planck-like CMB data and can be improved by a factor of 2 if either BAO or SNIa measurements are also considered, whereas a factor of √ 2 degradation could come from systematic effects. Interestingly enough, the observational strategy of LSST (large and deep survey) could provide the necessary sensitivity to explore the faint effects that the distinct neutrino mass eigenstates have on cosmological probes. This is a highly debated topic and we refer the reader to Sec. for related discussion.
Synergy between these large ground-based observatories and future space missions is expected. We consider here the ESA Euclid satellite and the NASA Wide Field Infrared Survey Telescope (WFIRST) as representative space-borne missions. Euclid will be a wide-field satellite that operates with imaging and spectroscopic instruments for years and covers roughly 15000 deg 2 in the optical and near-infrared bands, observing a billion galaxies and measuring ∼ 100 million galaxy redshifts [ ]. The redshift depth will be up to z ∼ 2 for galaxy clustering and up to z ∼ 3 for cosmic shear. The combination of the galaxy power spectrum measured with Euclid and primary CMB from Planck is expected to give σ(Σm ν ) = 0.04 eV; if instead the weak lensing dataset produced by Euclid is considered in combination with primary CMB, we expect σ(Σm ν ) = 0.05 eV [ ]. Both combinations provide a ∼ 1σ evidence in the minimal mass scenario. Some authors have also pointed out that weak lensing data as measured by Euclid could discriminate between the two neutrino hierarchies if the true value of Σm ν is small enough (i.e. far enough from the degenerate region of the neutrino mass spectrum), see [ ] and references therein.
WFIRST is an infrared telescope with a primary mirror as wide as the Hubble Space Telescope's primary (2.4 m) and will operate for years [ ]. The primary instrument on board, the Wide Field Instrument, will be able to operate both in imaging and spectroscopic mode, observing a billion galaxies. The instrumental characteristics of WFIRST will more than double the surface galaxy density measured by Euclid. With this setup, WFIRST will test the late expansion of the universe with great accuracy employing supernovae, weak lensing, BAO, redshift space distortions (RSD), and clusters as probes. From the BAO and broadband measurements of the matter power spectrum, WFIRST in combination with a Stage-III CMB experiment could provide σ(Σm ν ) < 0.03 eV [ ].
We want to conclude this section by pointing out that the aforementioned missions will be extremely powerful if combined together. Indeed, they are quite complementary [ ]. A significant example concerning the improvement of constraints on massive neutrinos is the combination of all the previously discussed surveys with the lensing reconstruction from CMB. The cross correlation of weak lensing (optical), CMB lensing power spectrum and galaxy clustering (spectroscopic) can highly reduce the systematics affecting each single probe, in particular the multiplicative bias in cosmic shear [ ]. For example, a combination of WFIRST, Euclid, LSST and CMB Stage-III can achieve σ(Σm ν ) < 0.01 eV [ ]. Another example is the calibration of the cluster mass for SZ cluster count analyses. This calibration can be performed through optical surveys such as LSST or through CMB lensing calibration, with comparable results. In Ref. [ ], the authors show that lensing-calibrated SZ cluster counts can provide a detection of the minimal neutrino mass Σm ν at > 3σ level, also in extended cosmological scenarios.

. -cm surveys
In this section, we will briefly comment about the possibility to use -cm survey data to constrain Σm ν . We refer the reader to the relevant papers for further readings. Measurements of the -cm signal such as those expected from the Square Kilometer Array The combination assumes a gaussian prior on τ = 0.06 ± 0.01 roughly corresponding to the new estimate from [ ]. b The combination assumes σ(τ ) = 0.002 and noise level of 2.5µK · arcmin. c For a fiducial value Σmν = 0 eV and marginalising over w0 − wa dark energy, arbitrary curvature and N eff . Table : Expected sensitivity on Σm ν from different combination of future cosmological data. Unless otherwise stated, the sensitivity σ(Σm ν ) is forecasted assuming a standard cosmological model with Σm ν = 0.06 eV. DESI refers to the simulated DESI-BAO dataset based on expected experimental performances [ ] (see [ , ] for details). FS refers to the use of the (simulated) measurements of the full shape of the matter power spectrum. The last line implies the use of CMB lensing, Euclid and WFIRST to calibrate the multiplicative bias in the shear measurements from LSST [ ].
(SKA) and the Canadian Hydrogen Intensity Mapping Experiment (CHIME) can shed light on the Epoch of Reionization, including a better determination of the reionization optical depth τ . In addition, they map the distribution of neutral hydrogen in the universe, a tracer of the underlying matter distribution. Therefore, constraints on Σm ν can benefit from -cm measurements in two ways: by breaking the degeneracy between Σm ν and τ (see e.g. [ ], where the authors report σ(Σm ν ) = 0.012 eV for a combination of CORE+Euclid lensing and FS+ a prior on τ compatible with expectations from future -cm surveys); by detecting the effect of Σm ν on the evolution of matter perturbations (see e.g. [ , , ]). https://chime-experiment.ca

Constraints on Σm ν in extended cosmological scenarios
The constraints derived so far apply to the simple one-parameter extension of the standard cosmological model, ΛCDM + Σm ν . When derived in the context of more complicated scenarios, such as models that allow arbitrary curvature and/or non-standard dark energy models and/or modified gravity scenarios etc., constraints on Σm ν are expected in general to degrade (although tighter constraints on Σm ν can be also possible in particular extended scenarios) with respect to those obtained in a ΛCDM + Σm ν cosmology. This effect is due to the multiple degeneracies arising between cosmological parameters that describe the cosmological model under scrutiny. In other words, when more degrees of freedom are available -in terms of cosmological parameters that are not fixed by the model -, more variables can be tuned in order to adapt the theoretical model to the data. For example, CMB data measure with incredible accuracy the location (expressed by the angular size of the horizon at recombination θ) and amplitude (basically driven by the exact value of z eq ) of the first acoustic peak. Therefore, we want to preserve this feature in any cosmological model. As explained before, h, Ω m and Σm ν can be varied together in order to do this. Adding other degrees of freedom, like curvature or evolving dark energy, allows for even more freedom, thus making the degeneracy worse. Of course, the addition of different cosmological data, which are usually sensitive to different combinations of the aforementioned parameters, is extremely helpful in tightening the constraints on Σm ν (and, in general, on any other cosmological parameter) in complex scenarios.
In more detail, constraints on the sum of neutrino masses are particularly sensitive to the so-called "geometric degeneracy". This term refers to the possibility of adjusting the parameters in order to keep constant the angle subtended by the sound horizon at last scattering, that controls the position of the first peak of the CMB anisotropy spectrum. The degeneracy is worsened in models with a varying curvature density Ω k or parameter of the equation of state of dark energy w. Constraints on the expansion history, like those provided by BAO or by direct measurements of the Hubble constant, are particularly helpful in breaking the geometric degeneracy. In principle, one could also expect a degeneracy between the effective number of degrees of freedom N eff and Σm ν , but for a different reason: both parameters can be varied in order to keep constant the redshift of matter-radiation equality. However, this can be done only at the expense of changing the CMB damping scale (see Sec. for further details). High-resolution measurements of the CMB anisotropies are therefore a key to partially break the degeneracy. Finally, a nonstandard relation between the matter density distribution and the lensing potential can be modelled by introducing a phenomenological parameter A L , which modulates the amplitude of the lensing signal [ ]. Most of the current constraining power of CMB experiments on Σm ν comes from CMB lensing. Therefore, it is clear that in models with varying A L the limits on neutrino masses are strongly degraded. However, it should also be noted that A L is usually introduced as a proxy for instrumental systematics; if considered as an actual physical parameter, its value is fixed by general relativity to be A L = 1.
To make the discussion more quantitative, we see how this applies to the constraints obtained with present data and future data. In Tab. , we report a comparison of the constraints on Σm ν for some extensions of the ΛCDM model. In the upper part of the table, we report constraints obtained from the PlanckTT+lowP+lensing+BAO dataset combination, described in Sec. . . These are taken from the full grid of results made available by the Planck collaboration and have obtained with the same statistical techniques used for the ΛCDM model. We see that the constraints are degraded by 30% in models with varying N eff , by 50% in models with varying Ω K or w, and by 65% in models with varying A L . This information is also conveyed, for an easier visual comparison, in Fig. , where we show the mass the sum of neutrino masses as a function of the mass m light of the lightest eigenstate. The green and red curves are for normal and inverted hierarchy, respectively. We show % constraints on Σm ν for different models and dataset combinations as horizontal lines. In the lower section of Tab. we instead report a similar comparison, based on the expected sensitivities of future CMB and LSS probes [ ]. The pattern is very similar to that observed for present data, although it should be noted that the increased precision of future experiments will allow to further reduce the degeneracies. In particular, it is found that the constraints on Σm ν are degraded by ∼ 30% in models with varying Ω K or w, and not degraded at all in models with varying N eff (models with varying A L have not been considered in Ref. [ ]). The cases reported in Tab

Cosmology and the neutrino mass hierarchy
Cosmology is mostly sensitive to the total energy density in neutrinos, directly proportional to the sum of the neutrino masses Σm ν ≡ m 1 + m 2 + m 3 . We can express Σm ν in the two hierarchies as a function of the lightest eigenstate m light (either m 1 or m 3 ) and of the squared mass differences ∆m 2 12 and ∆m 2 13 : Σm  ]. BAO refers to simulated data for DESI and Euclid surveys. The fiducial model adopted for the analysis is the following: Σm ν = 0.06 eV, Ω K = 0, w = −1, N eff = 3.046, Y He = 0.24, r = 0.
It has been a long-standing issue whether or not cosmological probes are sensitive to the neutrino mass hierarchy. In principle, we expect physical effects due to the choice of the neutrino hierarchy on cosmological observables. Individual neutrino species that carry a slightly different individual mass exhibit a slightly different free-streaming scale k f s : depending on their individual mass, neutrinos can finish suppressing the matter power at different epochs, leaving three distinct "kinks" in the matter power spectrum. As a consequence, the weak lensing effects on the CMB and on high redshift galaxies can be slightly affected by the choice of the hierarchy. In practice, all of these signatures are at the level of permille effects on the matter and CMB power spectra, well below the current sensitivity [ ].
Given the current sensitivity (roughly Σm ν < 0.2 eV at % CL), it is then a legitimate assumption to approximate the mass spectrum as perfectly degenerate (m i = Σm ν /3) when performing analysis of cosmological data. Very recently, several authors investigated the possibility that such an approximation could fail reproducing the physical behaviour of massive neutrinos when observed with the high sensitivity of future cosmological surveys [ , , , ]. In addition, the issue of whether future survey could unravel the unknown hierarchy has been addressed by several groups [ , , , , , ]. We refer the reader to the relevant papers for a thorough discussion of these issues. Here, we summarise the main results: ) the sensitivity of future experiments will not be enough to clearly separate the effects of different choices of the neutrino hierarchy, for a given value of Σm ν ; therefore the fully-degenerate approximation is still a viable way to model the neutrino mass spectrum in the context of cosmological analysis; ) the possibility to clearly identify the neutrino hierarchy with future cosmological probes is related to the capability of measuring Σm ν < 0.1eV at high statistical significance, in order to exclude the IH scenario. It is clear that the possibility to do this strongly depends on the true value of Σm ν : the closer it is to Σm N H,min ν = 0.06 eV, the larger will be the statistical significance by which we can exclude IH. This is true independently of whether we approach the issue from a frequentist or Bayesian perspective. In the latter case, however, since a detection of the hierarchy would driven by volume effects, this posits the question of what is the correct prior choice for Σm ν . The issue is extensively discussed in [ , , , , , ].

Complementarity with laboratory searches
Cosmological observables are ideal probes of the neutrino absolute mass scale, though they are not the only probes available. In fact, laboratory avenues such as kinematic measurements in β-decay experiments (see e.g. [ ]) and neutrino-less double-β decay (0ν2β) searches (see e.g. [ , ]) provide complementary pieces of information to those carried by cosmology.
Kinematic measurements are carried on with β-decay experiments mostly involving tritium 3 H. The shape of the decay spectrum close to the end point is sensitive to the (electron) neutrino mass and can be parametrized in terms of constraints on the electron neutrino effective mass defined in Eq. ( ). The current best limits on m β come from the Troitzk and Mainz experiments, with m β < 2.05 eV [ ] and m β < 2.3 eV [ ] at % CL. The new generation 3 H β-decay experiment KATRIN (Karlsruhe Tritium Neutrino ) is expected either to reach a sensitivity of m β < 0.2 eV at % CL, an order of magnitude improvement with respect to current sensitivities, or to detect the neutrino mass if it is higher than m β = 0.35 eV. Note that a detection of non-zero neutrino mass in KATRIN would imply Σm ν 1 eV, and would then be in tension with the cosmological constraints obtained in the framework of the ΛCDM model. This could point to the necessity of revising the standard cosmological model, although it should be noted that none of the simple one-parameter extensions reported in Tab. could accommodate for such a value.
Future improvements in kinematic measurements involve technological challenges, since KATRIN reaches the experimental limitations imposed to an experiment with spectrometers. Future prospects are represented by the possibility of calorimetric measurements of 136 Ho (HOLMES experiment [ ]) and measurements of the 3 H decay spectrum via relativistic shift in the cyclotron frequency of the electrons emitted in the decay (Project experiment [ ]). Although the bounds coming from β-decay experiments are very loose compared to bounds from cosmology, nevertheless they are appealing for the reason that they represent model-independent constraints on the neutrino mass scale, only relying on kinematic measurements. 0ν2β decay is a rare process that is allowed only if neutrinos are Majorana particles. A detection of 0ν2β events thus would solve the issue related to the nature of neutrinos, whether they are Dirac or Majorana particles. Searches for 0ν2β directly probe the number of 0ν2β events, which is related to the half life of the isotope involved in the decay T 1/2 . The latter can be translated in limits on the Majorana mass m ββ (defined in Eq. ) once a nuclear model has been specified. In practice, a bound on T 1/2 is reflected in a range of bounds on m ββ , due to the large uncertainties associated with the exact modelling of the nuclear matrix elements. Additional complications are due to model dependencies: when translating bounds on T 1/2 to bounds on m ββ , a mechanism responsible for the 0ν2β decay has to be specified. This is usually the exchange of light Majorana neutrinos, though alternative mechanisms could be responsible for the lepton number violation that not necessarily allow a direct connection between T 1/2 and m ββ . Finally, it can be shown that in the case of NH, disruptive interference between mixing parameters could prevent a detection of 0ν2β events, regardless of the neutrino nature and the lepton-number violation mechanism.
We report here some of the more recent limits on m ββ from 0ν2β searches. Constraints are reported as a range of % CL upper limits, due to the uncertainty on the nuclear matrix elements. We also specify the isotope used in each experiment. The current bounds are m ββ < 0.120 − 0.270 eV from Gerda Phase-II (  The next generation 0ν2β experiments, such as LEGEND, SuperNEMO, CUPID, SNO+, KamLAND -Zen, nEXO, NEXT, PANDAX-III, aims to cover the entire region of IH, reaching a 3σ discovery sensitivity for m ββ of 20 meV or better, roughly an order of magnitude improvement with respect to the current limits (see Ref. [ ] for a more detailed discussion and for a full list of references).
As outlined above, laboratory searches and cosmology are sensitive to different combinations of neutrino mixing parameters and individual masses. Therefore, it makes sense to compare their performances in terms of constraints on the neutrino mass scale. It is also beneficial to combine these different probes of the mass scale, in order to overcome the limitations of each single probe and increase the overall sensitivity to the neutrino masses [ , , ]. This is possible because, once the elements of the mixing matrix are known, specifying one of three mass parameters among (m β , m ββ , Σm ν ), together with the solar and atmospheric mass splittings, uniquely determines the other two. Oscillation experiments measure precisely the values of the mixing angles and of the squared mass differences, with an ambiguity on the sign of ∆m 2 31 , so that these parameters can be simply fixed to their best-fit values, given the larger uncertainties on the absolute mass parameters. The value of the Dirac phase, on the other hand, is known with lesser precision, and the Majorana phases, relevant for the interpretation of 0ν2β searches, are not probed at all by oscillation experiments. However this ignorance can be folded into the analysis using standard statistical techniques. Finally, the relation between the mass parameters also depends on the mass hierarchy. This can be taken into account either by performing different analysis for NH and IH, or by marginalizing over the hierarchy itself (see e.g. Ref. [ ]).
Combining the different probes of the absolute mass scale, with the support of oscillation results, leads to some interesting considerations. First of all, basically all of the information on the absolute mass scale comes from cosmology and 0ν2β searches. This confirms the naive expectation that can be made by comparing the sensitivity of the different probes. However, we recall again that the robust limits on m β from kinematic experiments represent an invaluable test for the consistency of the more model dependent constraints coming from cosmology and 0ν2β decay experiments. At the moment, cosmology still provides most of the information on the neutrino masses, although the sensitivity of 0ν2β experiments is rapidly approaching that of cosmological observations. A summary of the current limits is reported in Fig. of Ref. [ ]. To better illustrate the complementarity of cosmology and 0ν2β searches, we show in Fig. how they constrain, together with oscillation experiments, the allowed space in the (m ββ , m light ) plane. In more detail, we show the region in that plane that is singled out by oscillation experiments, for normal and inverted hierarchy. The width of the allowed regions traces the uncertainties on the CP-violating phases. We show current upper % bounds on m ββ from 0ν2β searches as horizontal lines, and current % bounds on m light from cosmology as vertical lines. These are translated from the bounds on Σm ν using information from oscillation experiments and assuming normal hierarchy. Assuming inverted hierarchy would however make a barely noticeable difference on the scale of the plot. It can be seen that in general cosmological observations are more constraining than 0ν2β searches.
In the future, however, one can expect that the constraining power of these two probes will be roughly equivalent. This can be seen in Fig. where, similarly to Fig , we show the allowed space in the (m ββ , m light ) plane for future cosmological and 0ν2β probes. As shown in Ref. [ ], the constraining power of 0ν2β searches for Σm ν would also depend crucially on the possibility of reducing the uncertainty on the nuclear matrix elements for the 0ν2β isotopes. In fact, provided that neutrinos are Majorana particles and that the leading mechanism responsible for the decay is a mass mechanism, the combination of cosmological probes and 0ν2β measurements could not only lead to a detection of the mass scale, but could also solve the hierarchy dilemma and provide useful information about (at least one of) the Majorana phases [ , , ].

Constraints on N eff
Until now, we have focused on the capability of cosmological observations to constrain neutrino masses. However, as noted in the introduction, cosmology is also a powerful probe of other neutrino properties. The main example is without any doubt the effective number of neutrino families (also called effective number of relativistic degrees of freedom) N eff , defined in equation Eq. ( ). As it is clear from its definition, N eff is simply a measure of the total cosmological density during the radiation-dominated era. More precisely, it represents the density in relativistic species, other than photons, normalized to the energy density of a massless neutrino that decouples well before electron-positron annihilation (that, we remember, is not actually the case). As explained in Sec. . , the standard framework, in which photons and active neutrinos are the only relativistic degrees of freedom present, and neutrino interactions follow the SM of particle physics, predicts N eff = 3.046 after electron-positron annihilation [ , , ]. Given its meaning, it is clear that a deviation from the expected value of N eff can hint to a broad class of effects -in fact, all those effects that change the density of light species in the early Universe. Those effects are not necessarily related to neutrino physics, as the definition of N eff in terms of the number of relativistic degrees of freedom suggests. For example, the existence of a Goldstone boson that decouples well before the QCD phase transition would appear as an increased number of degrees of freedom, with ∆N eff ≡ N eff − 3.046 = 0.027 [ ]. Speaking however about changes in N eff that are somehow related to neutrino physics, the most notable example is probably the existence of one (or more) additional, sterile light eigenstate, produced through some mechanism in the early Universe. In such a situation, one would have N eff > 3.046, as well as an additional contribution to Σm ν . Note that a light sterile neutrino would not necessarily contribute with ∆N eff = 1, as it does not share the same temperature as the active neutrinos.
In this section we will focus on cosmological constraints on sterile neutrinos. However, for completeness, we mention a few other examples of scenarios in which ∆N eff can possibly be different from zero. One is the presence of primordial lepton asym-  Figure : Majorana mass m ββ of the electron neutrino as a function of the mass m light of the lightest neutrino eigenstate, for normal (green) or inverted (red) hierarchy. The filled regions correspond to the uncertainty related to the CP-violating phases. The horizontal dashed lines show % current upper limits from 0ν2β searches. In particular, we show the tightest and loosest limits among those reported in the text, namely the most stringent from KamLAND-Zen (labeled "KamLAND-Zen, optimistic NME"), and the less stringent from CUORE (labeled "CUORE, pessimistic NME"). NME refers to uncertainty related to the nuclear matrix elements. We also show vertical dashed lines corresponding to % upper limits on Σm ν from cosmological observations, translated to upper limits on m light using the information from oscillation experiments. In particular we show different model and dataset combinations, from right to left: PlanckTT+lowP in the ΛCDM + Σm ν model, PlanckTT+lowP+BAO in the ΛCDM + Σm ν + Ω K model, PlanckTT+lowP+BAO in the ΛCDM + Σm ν model. The vertical lines shown in the plot assume normal hierarchy, but the difference with the case of inverted hierarchy is very small on the scale of the plot. ], in which the latest reheating episode of the Universe happens just before BBN, at temperatures of the order of a few MeV, so that neutrinos do not have time to thermalizeăcompletely. In this case, one has ∆N eff ≤ 0. Finally, non-standard interactions between neutrino and electrons can modify the time of neutrino decoupling [ ], so that the entropy transfer from e + e − annihilation is larger than in the standard picture and N eff is larger. We note that the effects related to these new scenarios are often more complicated that just a change in N eff : for example, both in the case of lepton asymmetries and low reheating, the neutrino distribution function is changed in a non trivial way, affecting also the other moments of the distribution (like the number density, the average velocity, etc). Finally, to mention a possibility that is not related to changes in N eff , cosmology can also probe the free-streaming nature of neutrinos, for example by looking for the effects of non-standard interactions among neutrinos [ , , , ], or between neutrinos and dark matter [ , , ]. Let us briefly recall how N eff is constrained by cosmological observations [ ]. Increasing N eff will make the Universe expand faster (larger H) during the radiationdominated era, and thus be younger at any given redshift. Then the comoving sound horizon at recombination will be smaller, going like 1/H, while the angular diameter distance to recombination stays constant, because H is unchanged after equality, so that θ s is smaller. Also, for fixed matter content, this will make the radiation-dominated era last longer. Recalling our discussion in Sec. . , the effect on the CMB spectrum is that the first peak is enhanced due to the larger early ISW, and all the peaks are moved to the right. However, as we have already learned, these effects can be canceled by acting on other parameters. There is however a more subtle and peculiar of effect of N eff , that is related to the scale of Silk damping. The damping scale roughly scales as 1/ √ H, i.e., as √ t, as expected for a random walk process. Then the ratio between the angle subtended by the sound horizon and that subtended by the damping length scales like H −1 /H −1/2 = H −1/2 . Since θ s is fixed by the position of the first peak, this means that increasing N eff , the damping length is projected on larger angular scales, or, equivalently, that damping at a given scale is larger. In conclusion the net effect is to lower the damping tail of the CMB spectrum. This effect is difficult to mimic with other parameters, at least in the standard framework. The damping length also depends on the density of baryons, so in principle one could think of changing this to compensate for the effect of N eff ; however, the baryon density is very well determined by the ratio of the heights of the first and second peak, so that it is in practice fixed. One possibility, in extended models, is to vary the fraction of primordial helium. Since the mean free path of photons depends on the number of free electrons, and helium recombines slightly before hydrogen, changing the helium to hydrogen ratio alters the Silk scale. However, this requires the assumption of nonstandard BBN, since, in the framework of standard BBN, the helium fraction is fixed by ω b and N eff themselves, so it is not a free parameter.

�����
We first review constraints on N eff in a simple one-parameter extension of ΛCDM, in which N eff is left free to vary, and the mass of active neutrinos is kept fixed to the minimum value allowed by oscillations. This case can be considered as the most agnostic, in some sense, in which one does not make any hypothesis on the new physics that is changing N eff (and thus on any other effects this new physics might produce). Moreover, one can think of these as limits for a very light (massless) sterile neutrino. Finally, constraining N eff is a robustness check for the standard ΛCDM mode. In fact, measuring N eff = 3.046 within the experimental uncertainty can be see as a great success of the standard cosmological model. It can be regarded as an indirect detection of the CνB, or, at least, of some component who has the same density, within errors, as we would expect for the three active neutrinos . From PlanckTT+lowP, one gets N eff = 3.13 ± 0.32; adding BAO gives N eff = 3.15 ± 0.23 [ ]. Both measurements, with a precision of ∼ 10%, are in excellent agreement with the standard prediction. Moreover, according to these results, ∆N eff = 1 is excluded at least at the 3σ level.
Using also information about the full shape of the matter power spectrum, the BOSS collaboration finds N eff = 3.03 ± 0.18 [ ]. We note that adding information from direct measurements of the Hubble constant results in larger values of N eff (N eff = 3.41 ± 0.22 from Planck TT,TE,EE+lowP+lensing+BAO+JLA+H 0 , see [ ]); this is due to the tension with the value of H 0 that is inferred from the CMB, that is alleviated in models with larger N eff . The next generation of cosmological experiments will improve these constraints by roughly an order of magnitude, getting close to the theoretical threshold of ∆N eff = 0.027 discussed at the beginning of this section, corresponding to a Goldstone boson decoupling before the QCD phase transition. Moreover, it will be possible to confirm the effects of non-instantaneous decoupling, since future sensitivities will allow to distinguish, at the 1-σ level, between N eff = 3 and N eff = 3.046. The combination of CORE TT,TE,EE,PP will put an upper bound at % CL of ∆N eff < 0.040 on the presence of extra massless (m 0.01 eV) species [ ] in addition to the three active neutrino families. The CORE collaboration puts limits also on the scenario in which the three active neutrinos have a fixed temperature, but their energy density is rescaled as (N eff /3.046) 3/4 . This scenario can account for an enhanced neutrino density (if N eff > 3.046) and reduced neutrino density (if N eff < 3.046 as for example in the case of lowreheating scenarios). In this case, CORE TT,TE,EE,PP yields N eff = 3.045 ± 0.041. Forecasts from S show that, in order to get closer to the threshold of ∆N eff = 0.027, a sensitivity of 1µK · arcmin and f sky > 50% are needed for a 1 beam size [ ]. Efficient delensing will help improve the limits on N eff : delensed spectra will have sharper acoustic peaks, allowing to constrain N eff not only through the impact on the Silk scale, but also through the phase shift in the acoustic peaks [ ]. Finally, having access to a larger sky fraction -and therefore to a larger number of modes observed -will be beneficial for constraints on N eff [ ]. We conclude this summary about future limits by noticing that the inclusion of LSS data, such as BAO measurements from DESI and Euclid, provides only little improvements with respect to CMB only constraints (e.g., The fact that, when probed, there is no hint for deviations from the free-streaming behaviour should strengthen our belief that we are really observing the CνB. This constraint has been obtained in the context of a ΛCDM + Σmν cosmology, with Σm fid ν = 0.06 eV. The constrain applies to the scenario of extra light relics in addition to the three massive neutrino families, i.e. N eff ≥ 3.046. b The constrain applies to the scenario of three massive neutrinos with energy density rescaled by N eff , i.e. N eff can be either lower or greater than . . c The combination includes delensed CMB spectra and a gaussian prior on the optical depth τ = 0.06 ± 0.01. . For a summary of current and future limits on N eff , we refer to Tab. . Let us now come to the case of a massive sterile neutrino. A sterile neutrino would contribute both to N eff and to ω ν . Its effect on the cosmological observables will thus be related to changes in these two quantities, as explained through this review. In fact, in principle, we should specify the full form of the distribution function of the sterile neutrino, and its effects could not be fully parameterized through N eff and to ω ν . Fortunately, one has that, when the distribution function is proportional to a Fermi-Dirac distribution, all the effects on the perturbation evolution of a light fermion can be mapped into two parameters [ ]: its energy density in the relativistic limit (and thus its contribution to N eff ) and its energy density in the nonrelativistic limit (and thus its density parameter, let us denote it with ω s to distinguish it from the active neutrinos).

Bounds
This covers several physically interesting cases, namely those of a sterile neutrino that either (i) has a thermal distribution with arbitrary temperature T s , or (ii) is distributed proportionally to the active neutrinos, but with a suppression factor χ s (this corresponds to the Dodelson-Widrow (DW) prediction for the non-resonant production scenario [ ]; see also Ref. [ ]). Defining an effective mass m eff s by mimicking Eq. , i.e.: Planck data are consistent with no sterile neutrinos: the % allowed region in parameter space is N eff < 3.7, m eff s < 0.52 eV from PlanckTT + lowP + lensing + BAO. However, it should be noted that they do not exclude a sterile neutrino, provided it contribution to the total energy density is small enough. We recall that a light sterile neutrino has been proposed as an explanation of the anomalies observed in short-baseline (SBL) experiments (see e.g. Ref. [ ] and references therein). However, a sterile neutrino with the mass (m s 1 eV) and coupling required to explain reactor anomalies would rapidly thermalize in the early Universe (see e.g. Refs. [ , ]) and lead to ∆N eff = 1, strongly at variance with cosmological constraints (excluded at more than % confidence considering the above combination of Planck and BAO data). We conclude this section by quoting the forecasts for future cosmological probes. In the context of a ΛCDM + Σm ν model with Σm fid ν = 0.06 eV and m fid s = 0 eV, the combination of CORE TT,TE,EE,PP with BAO measurements from DESI and Euclid will provide ∆N eff < 0.054 and m s < 0.035 eV [ ].

Summary
The absolute scale of neutrino masses is one of the main open questions in physics to date. Measuring the neutrino mass could shed light on the mechanism of mass generation, possibly related to new physics at a high energy scale. From the experimental point of view, neutrino masses can be probed in the laboratory, with β-and double β-decay experiments, and with cosmological observations. In fact, cosmology is at the moment the most sensitive probe of neutrino masses. Upper limits from cosmology on the sum of neutrino masses are possibly based on combinations of different observables. Results from the CMB alone can be regarded as very robust: these are of the order of Σm ν < 0.7 eV ( % CL). The addition of geometrical measurements, like those provided by BAO -also very robust -bring down this limit to Σm ν < 0.2 eV ( % CL). More aggressive analyses can get the bound very close to the minimum value allowed by oscillation experiments in the case of inverted hierarchy, but are based on observations where control of systematics is more difficult and thus should be taken with caution. It should also be borne in mind that cosmological inferences of neutrino masses are somehow model dependent. In extended cosmological models, especially those involving non-vanishing spatial curvature or dark energy, the constraints on Σm ν are degraded, even though they still remain very competitive with those obtained from laboratory experiments. Combination of future CMB and LSS experiments could reach, if systematics are kept under control, a sensitivity of meV in the first half of the next decade, allowing a 4σ detection of neutrino masses if the hierarchy is normal and the lightest eigenstate is massless. In that case, it will also be possible to exclude the inverted hierarchy scenario with a high statistical significance.
Present data are also compatible with the standard description of the neutrino sector, based on the standard model of particle physics. CMB measurements constrain the number of relativistic species at recombination to be N eff = 3.13 ± 0.32 at % CL. The inclusion of LSS data further tightens the constraints to N eff = 3.03 ± 0.18 at % CL. These results exclude the presence of an additional thermalized species at more than 3σ level. Cosmological data are also consistent with no sterile neutrinos. Thus no new physics in the neutrino sector is presently required to interpret cosmological data.
The standard picture will be tested more thoroughly by future experiments, that will allow to probe to an unprecedented level the physics of neutrino decoupling. An example would be the possibility to constrain non-standard neutrino-electron interactions. Future cosmological probes will also possibly reach the sensitivity necessary to detect, at the 1-σ level, the increase in the number of degrees of freedom due to a Goldstone boson that decouples well before the QCD phase transition.