Constraining light sterile neutrino mass with the BICEP2/Keck Array 2014 B-mode polarization data
Shouvik Roy Choudhury Sandhya Choubey

The SPLASH and  COSMOS Legacy Survey: the cross power between near-infrared and X-ray background fluctuations
    [
    
=============================================================================================================





§ INTRODUCTION


Sterile neutrinos still remain nothing short of an enigma in neutrino physics. Standard model predicts 3 massless neutrinos, while neutrino oscillation experiments have confirmed that neutrinos have mass, albeit very small. However, presence of anomalies in some short-baseline oscillation experiments <cit.> have been explained with an extra species of neutrino, namely a sterile neutrino, of mass ≃ 1 eV, which amply mixes with the active neutrinos but is uncharged under the standard model gauge group. Again, there are  analyses <cit.> which indicate that all the results cannot be explained comfortably with the sterile neutrino hypothesis. A recent result <cit.> from the MiniBooNE collaboration finds present electron neutrino and anti-neutrino appearance data still consistent with an extra sterile neutrino. 


Apart from terrestrial oscillation experiments, in recent years, cosmology has become a very powerful probe of neutrino physics. In a situation where standard model of particle physics is augmented with only an extra sterile neutrino species, there are two parameters of utmost importance. One is the effective number of relativistic neutrino species, N_eff, whose theoretically predicted value, considering only the standard model of particle physics with 3 massless active neutrinos is N_eff^SM = 3.046 <cit.>, but is supposed to increase when contribution from the sterile neutrino is counted. The other is the effective mass of the sterile neutrino, m^eff_s = Δ N_eff^3/4 m_s^ph, where Δ N_eff = N_eff - 3.046 and m_s^ph is the physical mass of the sterile neutrino. Cosmology can provide strong constraints on these two parameters. N_eff, in general, can have contribution from any relativistic species which is not a photon, and hence it is not restricted only to the neutrino sector. Also, in certain scenarios like very low-reheating scenarios with sterile neutrinos <cit.> or self-interacting sterile neutrinos <cit.>, Δ N_eff can be negative. However, we do not consider such scenarios in this work, and consider only a non-interacting extra species of sterile neutrino. 

Provided we are only considering an extension to standard model with neutrino oscillations in a 3+1 scenario, as long as the sterile neutrino is of similar mass to an active neutrino and amply mixes with the active ones, its cosmological implications are identical to the active neutrino. Sufficient mixing will lead to almost complete thermalization <cit.>. However, even if there is partial thermalization, it will, in general, increase N_eff, leading to a delayed matter-radiation equality and a higher value of the Hubble parameter, H(z_dec), at the CMB decoupling (given other parameters are kept fixed). This has two main consequences <cit.> on the CMB anisotropy power spectrum, first being an increase of the first peak of the spectrum due to early Integrated Sachs Wolfe (ISW) effect, and the second being a horizontal shift of the peaks towards higher multipoles. Along with a horizontal shift, there will also be a vertical shift which will decrease the amplitude of the peaks at high multipoles, a phenomenon related to Silk damping. These effects of an additional relativistic sterile neutrino can be partially compensated if other cosmological parameters are simultaneously varied. For example, if the total matter density ω_m is also increased without altering the baryon density, the redshift of matter-radiation equality can be kept fixed. These degeneracies tend to degrade the constraints on N_eff. However, the CMB power spectra won't be exactly the same even after such adjustments with other parameters, especially because of the neutrino anisotropic stress arising from the quadrupole moment of the cosmic neutrino background temperature anisotropies which alters the gravitational potentials <cit.>. Hence constraints can be put on N_eff from CMB power spectra data.

If a light sterile neutrino has a mass ≃ 1 eV, it only starts to become non-relativistic during CMB, and hence the effect of the mass is not strong on CMB power spectra. Sterile neutrinos with masses much smaller than 1 eV will have negligible effect on CMB power sprectra. However, when CMB power spectra data is used with other cosmological observations like constraining the Hubble parameter from direct measurements via a Gaussian prior or using the Baryon Acoustic Oscillation (BAO) data or both, better bounds on the mass of the sterile neutrino can be obtained <cit.>. A plethora of papers <cit.> are available on the effects of neutrino masses on cosmology. Current bounds on sterile neutrinos from cosmological data imply that fully thermalized sterile neutrinos of mass ≃ 1 eV are disfavoured and can only be accommodated with partial thermalization. See previous analyses on constraining sterile neutrino properties with cosmological data <cit.>.

In this paper, we have, for the first time, used the BK14 data, the data on the B-mode polarization of CMB from BICEP2/Keck collaboration, to constrain the parameters associated with sterile neutrinos in an extended ΛCDM model, which can be simply denoted with ΛCDM + r_0.05 + N_eff + m^eff_s. 
BK14 constrains the tensor-to-scalar ratio to r_0.05< 0.07 at 95% C.L, when combined with Planck 2015 and other datasets <cit.>; while exclusion of the BK14 data leads to a significantly less strong bound of r_0.05< 0.12 <cit.>. BK14 data also contains information on gravitational lensing. Thus we expect this data to affect the constraints on the sterile neutrino parameters. We also provide results with N_eff fixed at 4.046 and 3.5 separately, i.e., assuming full and partial thermalization of the sterile neutrinos respectively, and this model is denoted as ΛCDM + r_0.05 + m^eff_s.

This paper has the following structure: in Section <ref> we provide details about our model parameters and other analysis details and briefly describe the datasets used, in Section <ref> we provide the results of our analysis, and we conclude in Section <ref>.





§ COSMOLOGICAL ANALYSIS





 §.§ Model



Below we list the vector of parameters we have varied in this work in two cosmological models.

For ΛCDM + r_0.05 + N_eff + m^eff_s model: 

    θ≡[ω_c,  ω_b,  Θ_s, τ,  n_s,  ln [ 10^10 A_s], r_0.05, N_eff, m^eff_s].
 

For ΛCDM + r_0.05 + m^eff_s model:

    θ≡[ω_c,  ω_b,  Θ_s, τ,  n_s,  ln [ 10^10 A_s], r_0.05, m^eff_s].
 
with N_eff fixed to the value 4.046, which corresponds to full thermalization of the sterile neutrino with active neutrinos and to the value 3.5, which corresponds to partial thermalization. 

The first six parameters correspond to the ΛCDM model. Here ω_c = Ω_c h^2  and  ω_b = Ω_b h^2 are the physical cold dark matter and baryon densities at present, respectively. Θ_s is the angular sound horizon, i.e., the ratio between sound horizon and the angular diameter distance at decoupling. τ is the reionization optical depth. n_s and A_s are the power-law spectral index and power of the inflationary power spectrum, respectively, at the pivot scale of k_* = 0.05 h Mpc^-1. 

r_0.05 is the tensor-to-scalar ratio, also defined at the pivot scale of k_*=0.05h  Mpc^-1. N_eff, effective number of relativistic species which are not photons, is given by,

    ρ_r = π^2/15[ 1+ 7/8(4/11)^4/3 N_eff] T_γ^4,

where T_γ is the temperature of the photons and ρ_r is the radiation density. In our model, apart from photons, we only have 3 active and one sterile neutrino as relativistic species. The sterile neutrino is assumed not to have any self-interactions, or interactions with other particle species. In our work, we have fixed the active neutrino sector to give a contribution of N_eff^SM =3.046 to N_eff, with two massless and one massive neutrino with mass of 0.06 eV. Thus the contribution to N_eff from the sterile species is simply Δ N_eff = N_eff - 3.046. Here it should be mentioned that in certain scenarios like low reheating  Note that while we use the value 3.046 which is predominant in literature, a recent study <cit.>  had found N_eff^SM =3.045.

When the sterile neutrino is relativistic at early times, assuming the only radiation species are photons and neutrinos, contribution of a light sterile neutrino to N_eff is given by <cit.>,

    Δ N_eff =[7/8π^2/15T_ν^4]^-11/π^2∫ dp  p^3  f_s (p),

where T_ν is active neutrino temperature, p is the neutrino momentum, and f_s(p) is momentum distribution function of the sterile neutrino. At late times its energy density is parametrized as an effective mass <cit.>:

    ω_s ≡Ω_s h^2 = m_s^eff/94.1 eV = h^2 m_s^ph/π^2 ρ_c∫ dp  p^2  f_s (p),

where ρ_c is the critical density, Ω_s h^2 is the sterile neutrino energy density. Since sterile neutrinos don't have electroweak interactions and they have mixing with the active neutrinos, they cannot decouple after the decoupling of active neutrinos. Active neutrinos decouple at a temperature T ∼ 1 MeV, when all of them are relativistic. Hence f_s (p) doesn't depend on the physical mass of the sterile neutrino, m_s^ph. However f_s (p) depends on the production mechanism of the light sterile neutrino. If the production is through a thermal process, one can simply write f_s (p) = (e^p/T_s+1)^-1, the usual Fermi-Dirac distribution function, where T_s is the sterile neutrino temperature. In this case, it can be shown that,

    m^eff_s = Δ N_eff^3/4 m_s^ph;        Δ N_eff = (T_s/T_ν)^4.

Non-thermal production, on the other hand, can lead to various possible scenarios. One of the popular scenarios is the Dodelson-Widrow (DW) mechanism <cit.>, for which f_s (p) = β (e^p/T_ν+1)^-1, where β is a normalization factor. In this case, one gets <cit.>,


    m^eff_s = Δ N_eff  m_s^ph;        Δ N_eff =β.


So, the m^eff_s parametrization can accommodate two different scenarios of sterile neutrino production. Also notice that in the ΛCDM + r_0.05 + m^eff_s model, fixing N_eff = 4.046 leads to m^eff_s being same as m_s^ph.

In our work, we conduct a Bayesian analysis to derive constraints on the sterile neutrino parameters.  For all the parameters listed in Eq. (<ref>), and Eq. (<ref>), we impose flat priors. We also limit the physical mass of the sterile neutrino to m_s^ph≤ 10 eV. The prior ranges are provided on the  Table <ref>. We run chains using the November 2016 version of the Markov Chain Monte Carlo (MCMC) sampler CosmoMC <cit.> which incorporates CAMB <cit.> as the Boltzmann code and the Gelman and Rubin statistics <cit.> to estimate the convergence of chains.




 §.§ Datasets



We use separate combinations of the following datasets:


Cosmic Microwave Background: Planck 2015:

Measurements of the CMB temperature and low-l polarization from Planck 2015 <cit.> are used. We consider the high-l (30 ≤ l ≤ 2508) TT likelihood, and also the low-l (2 ≤ l ≤ 29) TT likelihood. We refer to this combination as TT. We also include the Planck polarization data in the low-l (2 ≤ l ≤ 29) likelihood, and denote this as lowP. We also use the Planck lensing potential measurements via reconstruction through the four-point correlation functions of the Planck CMB data <cit.>. We call this simply as lensing. Residual systematics may be present in the the Planck 2015 high-l polarization data <cit.>, so we refrain from using it.



B Mode Polarization data of CMB: 

Considering the B-mode polarization of CMB, we incorporate the recent dataset publicly available from BICEP2/Keck collaboration which includes all data (multipole range: 20 < l < 330) taken up to and including 2014 <cit.>. This dataset is referred to as BK14. 
 


Baryon Acoustic Oscillations (BAO) Measurements and Related Galaxy Cluster data:

In this analysis, we include measurements of the BAO signal obtained from different galaxy surveys. We make use of the SDSS-III BOSS DR12 <cit.> LOWZ and CMASS galaxy samples at z_eff = 0.38, 0.51 and 0.61, the DR7 Main Galaxy Sample (MGS) at z_eff = 0.15 <cit.>, and the 6dFGS survey at z_eff = 0.106 <cit.>. We call this complete combination as BAO. Here z_eff is the effective redshift of a survey. 


Hubble Parameter Measurements: 

We use a Gaussian prior of 73.24 ± 1.74 km/sec/Mpc  on H_0, which is a recent 2.4% determination of the local value of the Hubble parameter by <cit.> which combines the anchor NGC 4258, Milky Way and LMC Cepheids. We denote this prior as HST.






§ RESULTS


For convenience, we have separated the results in two subsections for the the two different models. The description of models and datasets are given at section <ref> and section <ref>, respectively. We have presented the results, first in the ΛCDM + r_0.05 + N_eff + m^eff_s model, and then in the ΛCDM + r_0.05 + m^eff_s model. All the marginalized limits quoted in the text or tables are at 68% C.L. whereas upper limits are quoted at 95% C.L., unless otherwise specified.



 §.§ Results for ΛCDM + r_0.05 + N_eff + m^eff_s model



	

	

In this section, we present the results for the ΛCDM + r_0.05 + N_eff + m^eff_s model. In Table <ref>  we have provided results without BK14 data, whereas, in Table <ref>, the results are with BK14, to compare. We have presented constraints on the three parameters r_0.05, N_eff, and m^eff_s. with which we have extended the ΛCDM model, and also two derived parameters H_0 and σ_8, which are important in constraining the sterile neutrino mass.   

With only TT+lowP, we see that the bound on the sterile mass is relaxed at m^eff_s< 0.78 eV. The bound gets tightened with BAO data, which partially breaks the degeneracy between m^eff_s and H_0 present in the TT+lowP data, by rejecting lower values of H_0 <cit.> and leads to a bound of m^eff_s< 0.53 eV. This effect can be seen pictorially in Figure <ref> where addition of BAO data leads to a significantly smaller magnitude of anti-correlation between m^eff_s and H_0. The HST prior also breaks the degeneracy partially, as can be seen in Figure <ref>. However, the H_0 values preferred by the HST prior are larger than BAO, which leads to a preference to even smaller masses (m^eff_s< 0.34 eV) to keep the comoving distance to the surface of last scattering fixed <cit.>. Adding HST and BAO together with CMB however does not provide better bound than CMB+HST. Also, the lensing data degrades the bound on m^eff_s. We note that CMB and/or BAO data do not allow full thermalization of sterile neutrinos. However, at 95% C.L., with TT+lowP+HST, we obtained a N_eff = 3.63^+0.44_-0.42. Such high values of N_eff disallow the standard model prediction of N_eff^SM = 3.046 at 95% C.L. but allow N_eff = 4.046, i.e., full thermalization. 
On the other hand, it is also imperative to consider recent constraints on N_eff coming from Big Bang Nucleosynthesis (BBN). Planck 2018 results <cit.> have provided bound of N_eff = 2.95^+0.56_-0.52 (95% C.L.) (which is independent of the details of the CMB spectra at high multipoles) by combining the helium, deuterium, and BAO data with an almost model-independent prior on θ_s derived from Planck data. Another recent study on BBN <cit.> provide a tight bound of N_eff = 2.90±0.22 (68% C.L.), which means at at 95% C.L., there will be only a small overlap in the values of N_eff provided by <cit.> and TT+lowP+HST. It is also to be noted that addition of the HST prior leads to a slightly inferior fit to the data, due to the 3.4σ tension present between Planck and HST regarding the value of H_0 <cit.>. We find that in this ΛCDM+r_0.05+N_eff+ m^eff_s model, compared to TT+lowP, the dataset TT+lowP+HST degrades the χ^2-fit by an amount of Δχ^2 = +3.43. 




Akaike information criterion (AIC):
To understand the improvement/worsening of the quality of fit with addition of sterile neutrino parameters (N_eff and m_s^eff) we need to compare the fit to data given by ΛCDM+r_0.05+N_eff+m_s^eff with that of ΛCDM + r_0.05. Since the number of parameters in the two models are not same, a popular method to compare the fit is the Akaike information criterion (AIC) <cit.>.  For a particular model and data, AIC is defined as, 

    AIC = χ^2_best-fit + 2k

 where k is the number of parameters in the model. The model with lower AIC corresponds to the preferred model.

Thus, comparison with another model (with the same data) can be done with ΔAIC = Δχ^2 + 2Δ k.
Usually models with extra parameters provide better fit to the data since they have a larger parameter space. The 2Δ k term penalises models with extra parameters to prevent any over-fitting. Here 2Δ k = 4.
 
We find that for the TT+lowP+HST data:  

    Δχ^2 = χ^2_best-fit (ΛCDM+r_0.05+N_eff+m_s^eff) - χ^2_best-fit (ΛCDM+r_0.05) = -4.3

i.e., the ΛCDM+r_0.05+N_eff+m_s^eff model provides a better χ^2 fit compared to ΛCDM + r_0.05. But due to the 2 extra parameters, ΔAIC = -0.3. Since this difference is small, it implies that the goodness of fits to the TT+lowP+HST data for the two models are similar.


Since the main aim of this paper is to analyze the role of the BK14 data, Table <ref> lists the bounds on the cosmological parameters, now with BK14 data included in each combination. The inclusion of the BK14 data seems to have almost no effect on the bounds of N_eff and H_0, as can be seen by comparing the results of Table <ref> and Table <ref>. However, bounds on m^eff_s improve slightly across all data combinations. The 1-D marginalized posteriors for m^eff_s and r_0.05 for various datasets are shown in Figure <ref>. While for TT+lowP, we had m^eff_s < 0.78 eV, this bound improves to m^eff_s < 0.68 eV with TT+lowP+BK14. Addition of BAO data further improves this bound to m^eff_s < 0.46 eV. Our most aggressive bound in this paper comes with TT+lowP+BK14+HST: m^eff_s < 0.28 eV.
BK14 data significantly constrains the tensor-to-scalar ratio, r_0.05. TT+lowP provides r_0.05<0.127 whereas TT+lowP+BK14 gives a constraint of r_0.05<0.068. However, we found only a very small correlation between r_0.05 and m^eff_s, and that does not explain the decrease in mass. In fact the correlation coefficient (defined as R_ij≡ C_ij/√(C_ii C_jj), where i and j are the two parameters being considered and C is the covariance matrix of cosmological parameters) between r_0.05 and m^eff_s to be R_m^eff_s,r_0.05 = -0.08 with TT+lowP and R_m^eff_s,r_0.05 = +0.02 with TT+lowP+BK14, i.e., there is no significant correlation before addition of BK14 and also no significant change after. However we also find slightly increased values of σ_8 across all data combinations when BK14 is included. For instance, for TT+lowP, we have σ_8 = 0.802^+0.040_-0.029, which increases to σ_8 = 0.814^+0.036_-0.027 with TT+lowP+BK14. Since σ_8 is the normalization of matter power spectrum on scales of 8h^-1 Mpc, a higher σ_8 prefers lower sterile neutrino mass, as larger neutrino masses create larger suppressions in the matter power spectrum <cit.>. Thus σ_8 and m^eff_s, both are strongly anti-correlated. Indeed, we found R_σ_8,m^eff_s = -0.84 with TT+lowP and R_r_0.05,m^eff_s = -0.81 with TT+lowP+BK14, and hence, even such small changes in σ_8 should also create small changes in m^eff_s, which we find is the case here. This has been depicted in Figure <ref>. Again, notice that the lensing data prefers a lower σ_8 value. As in Table <ref>, TT+lowP+BK14+HST+BAO yields σ_8 = 0.841^+0.025_-0.021, whereas adding the lensing data to this combination yields a lower σ_8 =0.820^+0.019_-0.015. Due to the same anti-correlation between σ_8 and m^eff_s, we see that inclusion of lensing data degrades the m^eff_s bounds.

Overall, we can say that the BK14 data makes the case for fully thermalized eV scale sterile neutrinos slightly worse. The parameter to justify this statement is m_s^eff. As we have shown that addition of the BK14 data does not affect the N_eff bounds, BK14 data does not affect the thermalisation situation, as far as cosmological data is concerned. However, short baseline oscillation experiments predict a fully thermalised sterile neutrino of mass ≃ 1 eV. This requires that both N_eff = 4.046 and m^eff_s≃ 1 eV be allowed by the data. Since adding the BK14 data tightens the bounds on m^eff_s for all of the cosmological dataset combinations, it also takes the m^eff_s value further away from the 1 eV value, while N_eff bounds almost remain unchanged. Effect of BK14 data on sum of active neutrino masses (∑ m_ν) was also studied by us recently in <cit.>, in the ΛCDM + r_0.05 + ∑ m_ν model, where we had also found slightly increased σ_8. This is also indirectly confirmed by the recent Planck 2018 results, where they provide a bound of ∑ m_ν < 0.12 eV with Planck TT,TE,EE+lowE+lensing+BAO data in ΛCDM+∑ m_ν model <cit.>, whereas the bound is ∑ m_ν < 0.11 eV with Planck TT,TE,EE+lowE+lensing+BK14+BAO data in the ΛCDM+r+∑ m_ν model <cit.>. This similar effect was seen to persist even in a 12 parameter extended scenario in a recent study with non-phantom dynamical dark energy <cit.>. In this paper we have shown that such an effect is also present in an extended ΛCDM cosmology with light sterile neutrinos. CMB B-mode polarization has two known sources <cit.>. The first one is the inflationary gravitational waves (IGW), i.e., tensors (expected to produce a bump peaked around l≃ 80, the so called 'recombination bump' in the BB-mode CMB spectra) as tensors induce quadruple anisotropies in the CMB within the last scattering surface. The tensor signature cannot be reproduced by scalar perturbations, and the amplitude of the recombination bump depends on the tensor-to-scalar ratio. The second source is gravitational lensing by large scale structure. It leads to deflection of CMB photons at late times, which converts a small part of the E mode power into B mode. This lensing BB spectra is expected to have a peak around l≃ 1000. The BICEP2/Keck experiment has a multipole range 20<l<330 aiming to constrain the tensor-to-scalar ratio. However since r_0.05 and m^eff_s are only weakly correlated, the slightly stronger constraints on the neutrino masses is possibly coming from gravitational lensing information encoded in the BK14 data, and not from measurement of r_0.05. 







 


H_0 and σ_8 tensions:

It is also worth noting that in ΛCDM model, with TT+lowP, Planck collaboration <cit.> found that H_0 = 67.31± 0.96 km/sec/Mpc, whereas in this ΛCDM + r_0.05 + N_eff + m^eff_s model we find H_0 = 68.35^+1.23_-2.50 km/sec/Mpc. This preference to larger values of H_0 decreases the more than 3σ tension present in the ΛCDM model, between Planck 2015 and HST. One of the main reasons is that marginalizing over N_eff, which allows for N_eff>3.046 and higher N_eff values prefer a higher H_0, to keep the acoustic scale parameter θ_s fixed <cit.>, which is very well constrained by Planck data. Thus H_0 and N_eff are strongly correlated.

The ΛCDM + r_0.05 + N_eff + m^eff_s model also helps in reconciling the σ_8 tension present in the σ_8-Ω_m plane in ΛCDM model between Planck 2015 and weak lensing survey, like CFHTLenS <cit.> and KiDS-450 <cit.>. For instance, the KiDS-450 survey constrains the quantity S_8 ≡σ_8 √(Ω_m/0.3) = 0.745 ± 0.039 which has a 2.3σ tension with Planck TT+lowP, which prefers a much higher value of S_8 = 0.851±0.024 <cit.>. Planck data also prefers higher values of σ_8 compared to CFHTLenS. With TT+lowP in base ΛCDM model, one gets σ_8 = 0.829±0.014 <cit.>. However, in this ΛCDM + r_0.05 + N_eff + m^eff_s model, with TT+lowP, we get σ_8 = 0.802^+0.040_-0.029, which is much lower and thereby the conflict is decreased somewhat. We also get S_8 = 0.824^+0.030_-0.027, which is better agreement with KiDS-450 than ΛCDM. However, the BK14 data prefers slightly higher σ_8 values and thereby increases the tension between Planck and these weak gravitational lensing surveys. This can be visualized in Figure <ref>, where we see that the inclusion of BK14 data drives the 2-D contours upwards to a small extent. In Figure <ref>, we have used the CFHTLenS data with conservative cuts as described in <cit.>.

Another important point is that while ΛCDM + r_0.05 + N_eff + m^eff_s helps in relieving the H_0 and σ_8 tensions present in the ΛCDM model, they are not both relieved together in any region of the allowed parameter space. In the right panel of Figure <ref>, we can see that the regions where σ_8 has lower values, H_0 also has lower values (while we need higher values of H_0 to relieve the H_0 tension), and similarly, where H_0 has higher values, σ_8 also has higher values (while we need lower values of σ_8 to relieve the σ_8 tension). This in turn implies that the two conflicts are not resolved together in this model. And BK14 data worsens the conflicts even more. The HST prior also doesn't help the issue here. As we can see from Tables <ref> and <ref>, the inclusion of this Gaussian prior leads to a preference for much higher N_eff values, and higher σ_8 values as well, increasing the conflict. 




 §.§ Results for ΛCDM + r_0.05 + m^eff_s model


	
	

In this section we verify the stability of the results obtained in the previous section, by going to a smaller parameter space. We stop varying N_eff and fix its value to 4.046 and 3.5. The first one corresponds to complete thermalization of sterile neutrinos, while the later one corresponds to partial thermalization. We have restricted ourselves to CMB data only. For N_eff = 4.046 and N_eff = 3.5, the results are given in Tables <ref> and <ref> respectively.

We see that BK14 does help in obtaining better constraint on the sterile mass also in this reduced parameter space. For N_eff = 4.046, with TT+lowP, we get m^eff_s < 0.66 eV, whereas inclusion of BK14 leads to a tighter bound of m^eff_s < 0.50 eV. Similar case of strengthening of mass bound is seen with N_eff = 3.5, although these bounds are more relaxed compared to the case N_eff = 4.046, as a higher N_eff prefers a higher H_0. Again we see that the BK14 data itself does not affect the H_0 constraints much, but heavily constraints the tensor-to-scalar ratio, and also slightly increases the preferred σ_8 values. The main conclusions made in the previous section on the larger parameter space thus remains unchanged in this smaller parameter space. 

It is imperative to note that for sterile neutrinos produced by a thermal process and obeying Eq. <ref>, for N_eff = 4.046, we have m^ph_s = m^eff_s, whereas for N_eff = 3.5, we have m^ph_s = 1.8 m^eff_s. Hence, for N_eff = 3.5 and with TT+lowP+BK14, we have a corresponding bound of m^ph_s<1.13 eV. This implies that CMB data allows sterile neutrinos with mass ≃ 1 eV, but only with partial thermalization with N_eff≃ 3.5. 
When we compare the quality of fit to the TT+lowP+BK14 data between the ΛCDM+r_0.05+m_s^eff model (N_eff = 3.5 and 4.046) and the ΛCDM+r_0.05 (with N_eff = N_eff^SM), we find that, 

for the N_eff = 4.046 case:

    Δχ^2 = χ^2_best-fit (ΛCDM+r_0.05+m_s^eff) - χ^2_best-fit (ΛCDM+r_0.05) = +7.03

whereas, for the  N_eff = 3.5 case:

    Δχ^2 = χ^2_best-fit (ΛCDM+r_0.05+m_s^eff) - χ^2_best-fit (ΛCDM+r_0.05) = -0.22

 
These correspond to ΔAIC = +9.03 (for N_eff = 4.046) and ΔAIC = +1.78 (for N_eff = 3.5). Thus, the model with partial thermalization of N_eff = 3.5 provides only a slightly worse fit to the data compared to the ΛCDM+r_0.05 model (with N_eff = N_eff^SM), and is preferred by the data much more than the full-thermalization case. This is not surprising as in the previous section we had seen that CMB data alone did not allow complete thermalization. 



§ DISCUSSION



Short Baseline (SBL) Oscillation anomalies have hinted towards a fully thermalized sterile neutrino with mass around 1 eV. In this paper we have studied, for the first time, the light eV scale sterile neutrino situation in cosmology in light of the BICEP2/Keck array 2014 CMB B-mode polarization data. We call this dataset BK14. We first considered an extended-ΛCDM scenario with tensor perturbations and sterile neutrino parameters: ΛCDM + r_0.05 + N_eff + m^eff_s model. Apart from BK14, we have used Planck 2015 temperature and low-l polarization data (TT+lowP), latest BAO data and a Gaussian prior on the Hubble constant (HST) from local measurements. We find that inclusion of the BK14 data has almost no effect on the bounds of N_eff and H_0 but it strengthens the bounds on m^eff_s to a small extent by preferring slightly higher values of σ_8, with which m^eff_s is strongly anti-correlated. The BK14 data also tightly constraints the tensor-to-scalar ratio, r_0.05 but we find negligible correlation between r_0.05 and m^eff_s. This makes us think that the effect on mass is coming from the gravitational lensing information encoded in the B-mode polarization and not from the Inflationary Gravitational Waves. The bound of m^eff_s < 0.46 eV (95% C.L.) is found for the combination of Planck 2015, BAO and BK14 datasets, whereas the bound is m^eff_s < 0.53 eV (95% C.L.) without the BK14 data. Our most aggressive bound of m^eff_s < 0.28 eV (95% C.L.) is obtained with Planck 2015, HST and BK14. The HST prior also leads to high N_eff values which allow full thermalization of the sterile neutrino (at 2σ) but such high values are in conflict with bounds from Big Bang Nucleosynthesis. Also, addition of the HST prior to the TT+lowP data leads to a slightly worse χ^2 fit to the data. On the other hand, it is to be noted that as per the Akaike information criterion (AIC) the ΛCDM+r_0.05+N_eff+m_s^eff model provides equally good fit to the data as the ΛCDM+r_0.05 model, for the TT+lowP+HST data combination. Previous studies have indicated that fully thermalized sterile neutrinos with mass ∼ 1 eV (as predicted by SBL experiments) are disfavoured by cosmological data. Our analysis indicates that it becomes slightly more disfavoured with the inclusion of BK14 data, due to tighter mass bounds. The BK14 data also seems to make the agreement between Planck 2015 and CFHTLenS (weak gravitational lensing data) worse due to the higher σ_8 values. 

We would also like to mention that the Planck 2018 results, released during the preparation of this article, indirectly show tightening of bounds on ∑ m_ν   with BK14. They provide a bound of ∑ m_ν < 0.12 eV with Planck TT,TE,EE+lowE+lensing+BAO data in ΛCDM+∑ m_ν model <cit.>, whereas the bound is ∑ m_ν < 0.11 eV with Planck TT,TE,EE + lowE + lensing + BK14 + BAO data in the ΛCDM+r+∑ m_ν model <cit.>. Thus we expect our main conclusion regarding BK14 helping in improving the bound on sterile neutrino mass will remain unchanged if used with the recent Planck 2018 data instead of Planck 2015 that we have used in this paper. 

While this work was still being completed, a new B-mode polarisation data was released publicly, from the same BICEP2/Keck collaboration. This newly released data includes all the measurements upto and including 2015, and thus we call it BK15 <cit.>. To understand the effect of the new data, we performed an MCMC analysis with TT+lowP+HST+BK15 in the ΛCDM+r_0.05+N_eff+ m^eff_s model (with all other settings remaining unchanged). We found the following bounds: m^eff_s<0.27 eV (95% C.L.), r_0.05 < 0.061, and σ_8 = 0.847^+0.026_-0.021. In the same model, when we had used BK14 instead of BK15, we had found (see table 3), m^eff_s<0.28 eV (95% C.L.), r_0.05 < 0.073, and σ_8 = 0.846^+0.026_-0.020. As we can see, that while the bound on r_0.05 changes, the bounds on m^eff_s and σ_8 almost remain unchanged. We also checked that other parameters of interest, like H_0 and N_eff change negligibly. As before, since r_0.05 and m^eff_s have only a very weak correlation, it doesn't affect the mass bound. On the other hand, since TT+lowP+HST+BK15 almost doesn't change the bound on σ_8, the mass bound almost remains the same. Thus, we find that reanalysis with BK15 instead of BK14 will not change the neutrino mass bounds.

This tension between SBL and cosmological datasets has given rise to a number ideas to reconcile the eV-scale sterile neutrinos with cosmology. These include introduction of new "secret interactions" among sterile neutrinos which modifies the background potential and blocks thermalization <cit.>, modifications to the cosmic expansion rate at the time where sterile neutrinos are produced <cit.>, large lepton asymmetry <cit.>, time varying dark energy component <cit.>, very low reheating temperature <cit.>. The recent results that have come from the MiniBooNE collaboration <cit.> have rekindled interest in the sterile neutrinos. Thus, there seems to be interesting and exciting times ahead in the realm of light sterile neutrinos in cosmology.

The authors sincerely thank the anonymous referee for the thoughtful comments and efforts towards
improving our manuscript. SRC thanks the cluster computing facility at HRI (http://cluster.hri.res.in). SRC also thanks Steen Hannestad for useful discussions. The authors would also like to thank the Department of Atomic Energy (DAE) Neutrino Project of HRI. This project has received funding from the European Union's Horizon 2020 research and innovation programme InvisiblesPlus RISE under the Marie Sklodowska-Curie grant agreement No 690575. This project has received funding from the European Union's Horizon 2020 research and innovation programme Elusives ITN under the Marie Sklodowska-Curie grant agreement No 674896. 



jhep




