title: Constraints on Decaying Sterile Neutrinos from Solar Antineutrinos
authors: Matheus Hostert Maxim Pospelov
text: 
  
  FTPI-MINN-20-28
  
  
  
  
  
  mhostert@umn.edu
  School of Physics and Astronomy, University of Minnesota, Minneapolis, MN 55455, USA
  William I. Fine Theoretical Physics Institute, School of Physics and Astronomy, University of
  Minnesota, Minneapolis, MN 55455, USA
  Perimeter Institute for Theoretical Physics, Waterloo, ON N2J 2W9, Canada
  
  
  School of Physics and Astronomy, University of Minnesota, Minneapolis, MN 55455, USA
  William I. Fine Theoretical Physics Institute, School of Physics and Astronomy, University of
  Minnesota, Minneapolis, MN 55455, USA
  
  
  
  
  
    Solar neutrino experiments are highly sensitive to sources of ν→ν conversions in the ^8B neutrino flux. In this work we adapt these searches to non-minimal sterile neutrino models recently proposed to explain the LSND, MiniBooNE, and reactor anomalies. The production of such sterile neutrinos in the Sun, followed the decay chain ν_4 →νϕ→ννν with a new scalar ϕ results in upper limits for the neutrino mixing |U_e4|^2 at the per mille level. We conclude that a simultaneous explanations of all anomalies is in tension with KamLAND, Super-Kamiokande, and Borexino constraints on the flux of solar antineutrinos. We then present other minimal models that violate parity or lepton number, and discuss the applicability of our constraints in each case. Future improvements can be expected from existing Borexino data as well as from future searches at Super-Kamiokande with added Gd.
  
  
  Constraints on Decaying Sterile Neutrinos from Solar Antineutrinos
      Maxim Pospelov
      December 30, 2023
  ==================================================================
  
  
  
  
  
  
  
  § INTRODUCTION
   
  Beyond neutrino mixing, the study of solar neutrinos provides important input to the standard solar Model (SSM) <cit.> and has been used to search for several phenomena beyond the SM of particle physics. One example is neutrino decay, originally proposed as an alternative solution to the solar neutrino problem <cit.>. Indeed, after precision measurements of the solar neutrino oscillation parameters by KamLAND <cit.>, strong constraints on the lifetimes of ν_2 and ν_3 have been obtained <cit.> as data is consistent with no additional neutrino disappearance. 
  
  Recently, non-minimal neutrino decay models have received interest in the literature, where exotic decays of a relatively heavy (m_ exotic≫ m_ active) and mostly-sterile neutrino are invoked to explain longstanding experimental anomalies at short baselines (SBL). One category of models concerns “visible" sterile neutrino decay, where new sterile states are produced and decay back to visible active neutrinos. Originally proposed in Ref. <cit.> as an explanation to the Liquid Scintillator Neutrino Detector (LSND) anomaly <cit.>, this scenario has now been revisited <cit.> in light of recent data of short-baseline ν_μ→ν_e and ν_μ→ν_e appearance at MiniBooNE <cit.>, as well as ν_e disappearance at reactors <cit.>. Due to the small mixing angles required in this explanation, the effects of attenuation in solar neutrino fluxes is small. Yet, the total number of heavy neutrinos produced is large, and if these states undergo sufficiently distinctive decays or scattering inside a detector, they can be searched for. In this article, we point out that if antineutrinos are produced in the decay of these heavy neutrinos, then they are strongly constrained by existing searches for neutrino-antineutrino transitions in solar neutrino experiments. 
  
  The flux of antineutrinos from the Sun at the MeV energies is negligible <cit.>, which remains an excellent approximation down to tens of keV in energy <cit.>. Combined with the fact that the detection cross section for ν_e is much larger and easier to measure compared to that of ν_e, this makes solar neutrino experiments sensitive to very small fluxes of antineutrinos <cit.>. The current sensitivity reaches fluxes as small as a few times 10^-5 of the ^8B neutrino flux <cit.>. These searches have been discussed in the context of new physics, such as large ν→ν oscillations. This Lepton number (LN) violating process is rather small in most theories, being suppressed by (m_ν/E)^2, but can be enhanced due to spin-flavor precession <cit.>. The latter arises from the coupling of a large neutrino magnetic moment to the solar magnetic field, which induces ν_e →ν_x conversions, followed by flavor transitions into ν_e due to matter effects. Another possibility to generate such LN violating signatures is neutrino decay. For instance, neutrino mass models where LN is a spontaneously broken global symmetry predict the existence of a pseudo-goldstone boson J, the majoron <cit.>. In these models, solar antineutrinos may be produced from the decay ν_2 →ν_1 J, which is enhanced in dense matter <cit.>. This possibility of production from neutrino decay is, in fact, quite general and can be realized in any LN violating model with neutrinos that decay sufficiently fast, be they ν_2, ν_3, or the new mostly-sterile state ν_4 discussed here (see Ref. <cit.> for an early discussion in the context of a 17 keV sterile neutrino). 
  
  In this work, we explore a new possibility where lepton number can, in fact, be conserved but the decay of a new light boson leads to a large flux of antineutrinos. We derive limits on the electron flavor mixing with ν_4, working only with the gauge-invariant and parity-conserving model of Ref. <cit.>. Focusing solely on Dirac neutrinos, we show that our bounds exclude virtually all of the parameter space preferred that can simultaneously explain LSND and MiniBooNE, as well as the region of interest for reactor anomalies. 
  They also disfavor most but not all parameter space suggested as a solution to the MiniBooNE anomaly. For models with Majorana neutrinos, the constraints become even stronger due to ν_4 →νϕ→ννν decays.
  
  The paper is organized as follows. In sec:model we review the benchmark model for decaying sterile neutrinos, and in sec:solantinus we discuss generic aspects of solar antineutrino searches. The resulting constraints, future prospects, and alternative search methods are then discussed in sec:results. We dedicate sec:alternatives to a survey of minimal alternative models for decaying steriles, and conclude in sec:conclusions.
  
  
  
  
  § DECAYING STERILE NEUTRINO
  
  
  
  
  
  The most significant deviations from the three-neutrino paradigm at SBLs are the LSND excess of ν_e events, with a statistical significance of 3.8 σ when interpreted under a ν_μ oscillation hypothesis, and the MiniBooNE excess of ν_e-like events, with a significance of 4.8 σ when interpreted under a ν_μ→ν_e and ν_μ→ν_e oscillation hypothesis. Reactors at very short-baselines also have some evidence of ν_e disappearance <cit.>, but in that case the neutrino flux predictions are highly uncertain and harder to control <cit.>. Despite the large significance of these anomalies, they remain unsolved. Their standard interpretation under oscillations of a eV-scale sterile neutrino leads to strong tensions between different data. This is driven mainly by the absence of anomalous results in ν_μ disappearance experiments <cit.>, as appearance and disappearance channels are strongly correlated in the oscillation scenario <cit.>. In addition, such new sterile states with eV masses are in strong tension with cosmological observations, which has prompted several studies to resolve this by invoking secret interactions <cit.>.  
  
  Visible sterile neutrino decays are, therefore, a natural “next-to-minimal" explanation to SBL anomalies to consider. The advantages of this scenario are that it does not necessarily lead to strong correlations between appearance and disappearance channels, the mass scale of the new sterile state is not fixed by the oscillation length of the experiments, and that it already contains a secret interaction mechanism, possibly alleviating tension with cosmology. In this work, we focus on the decay of steriles with eV to hundreds of keV masses to a new scalar ϕ, as discussed in Refs. <cit.>. In all such visible decay scenarios, heavy neutrinos decay to mostly-active neutrinos via ν_4 →νϕ , where more neutrinos can be produced from ϕ→νν decay if ϕ is massive as in Ref. <cit.>. Here visible refers to the detectability of the decay products, in contrast to models where neutrinos decay to the wrong-helicity states that do not feel the weak interactions (up to tiny helicity-flipping terms proportional to m^2_ν/E_ν^2). 
  
  Such visible decays can explain the anomalous ν_e-like events at SBL experiments by means of a sub-dominant population of ν_4 states in neutrino beams, which often decays to ν_e-like daughters [Constraints on this scenario have been obtained in Ref. <cit.> using the near detector of NOνA and T2K, as well as MINERνA and PS-191. We note that the constraints have been obtained under simplified assumptions, and that a detailed study with total signal efficiency, as well as appropriate uncertainties is needed in order to derive reliable constraints.]. One typical prediction is that the spectrum of daughter ν_e and ν_e neutrinos is softer than the initial flux of ν_4 parents and associated neutrinos, skewing the effective flavor conversion towards lower energies. While this brings a mild improvement over the oscillation fit to the low energy excess observed at MiniBooNE, it leads to less satisfactory energy spectra at LSND, which is compatible with a signal that grows in energy. In addition, the neutrino flux at LSND comes from both π^+ and μ^+ decay at rest, yielding a large and monochromatic ν_μ flux, and a spectrum of ν_μ and ν_e. Since only the ν_e component is detected via the IBD process, the presence of a neutrino-to-antineutrino transition in the decay chain can convert the large ν_μ flux to signal, since ν_4 states can be produced in pion as well as muon decays. This is a crucial point in the study of Dentler:2019dhz, which found improved compatibility between LSND and MiniBooNE regions of preference when this conversion is significant.
  
  For concreteness, we focus on the gauge-invariant and parity-conserving model of Ref. <cit.>, wherein a SM singlet ν_s is introduced and equipped with sizable couplings to a new scalar singlet ϕ. The sterile neutrino can then couple to light and mostly-active neutrinos in a gauge-invariant fashion by means of mixing between the heaviest neutrino state, ν_4, and the active flavors. The relevant Lagrangian is given by
  
      - ℒ = g_ϕν_sν_sϕ + ∑_α,β m_αβν_αν_β,
  
  where the neutrino mass mechanism is left unspecified and assumed to not play a role in the low-energy phenomenology. In the mass basis, the neutrino mass eigenstates are given by ν_i = ∑_α U_α i^* ν_α, with α∈{e,μ,τ,s}, and U a unitary mixing matrix. Under the assumption of parity conservation in the sterile sector, U is identical to the extended Pontecorvo–Maki–Nakagawa–Sakata (PMNS) matrix, now 4×4. We return to this issue in sec:alternatives. In the decays of ν_4 and ϕ, only the three lightest mass states are produced, and so the it is useful to define the low-energy flavor state ν̂_s = ∑_i=1^3 U_α iν_i. For most processes of interest, however, the non-unitarity corrections introduced by working with ν̂_s instead of the full flavor states ν_s is small and appears only at order |U_α4|^4. Unless stated otherwise, we refer to ν̂_s as simply ν from now on, as the mass eigenstates have decohered on their way from Sun. The new scalar does not couple directly to the SM, and loop-induced couplings will ultimately depend on the UV completion of the model and its neutrino mass mechanism (see, for instance, Refs. <cit.>).
  
  
  
  Due to mixing, heavy neutrinos with masses below the MeV scale would be produced in the Sun via the same processes responsible for ν_e production at a rate |U_e4|^2 times smaller. Once produced, the ν_4 mass eigenstates immediately decay to a light neutrinos and the scalar boson. The scalar then decays to a neutrino-antineutrino pair, giving rise to our signal. Overall, the process of interest is
  
      ν_4 (E_ν_4)  → ν(E_1)  + ϕ(E_ϕ)   
         ^↘ν (E_2)  + ν (E_3).
  
  
  For most cases of interest, E_ν_4≫ m_4, so if ν_4 (ν_4) is produced via weak interactions, it will be left-handed (right-handed) polarized to a very good approximation. We then assume all heavy neutrinos to be polarized with a definite helicity h_4=-1 for neutrinos and h_4=+1 for antineutrinos. Nevertheless, due to the assumption of parity conservation for the ϕ interactions with neutrinos, both helicity flipping (HF) and helicity conserving (HC) decay channels are allowed. Assuming all neutrinos to be ultra-relativistic, we find the squared-amplitudes for polarized ν_4^h_4=-1→ν^hϕ decay,
  
  
      |M_ν_4^-→ν^-ϕ|^2    = ∑_i=1^3 g_ϕ^2  |U_si|^2  m_4^2 E_1/E_ν,
      
      |M_ν_4^-→ν^+ϕ|^2    = ∑_i=1^3 g_ϕ^2  |U_si|^2 m_4^2  [(1-r_ϕ^2)-E_1/E_ν]
  
  where r_ϕ = m_ϕ/m_4. Integrated over phase space, both channels contribute identically to a total decay rate of
  
      Γ^ LAB (ν_4 →νϕ) = ∑_ig_ϕ^2/16πm_4^2/E_ν |U_s4U_s i|^2 (1-r_ϕ^2)^2.
  
  Our decay rate is in agreement with Refs <cit.>. Note that helicity conserving decays prefer larger E_1 values, while helicity flipping decays prefer smaller values of E_1. Therefore, for our present application, helicity-flipping decays are important since the antineutrinos from the subsequent scalar decay tend to be more energetic. Also important is the limit r→1, where the scalar particle has most of the ν_4 energy regardless of the helicity structure of the decay. This is the scenario with the most energetic antineutrinos in the final state, for which a simultaneous explanation of MiniBooNE and LSND is most successful.
  
  The scalar decay length in the lab frame to leading order in the small mixing elements is
  
      Γ^ LAB (ϕ→νν) = ∑_i,j=1^3g_ϕ^2/8πm_ϕ^2/E_ϕ |U_s iU_s j|^2.
  
  As expected, the scalar decays are doubly suppressed by small mixing elements, and so it tends to decay more slowly than ν_4. Nevertheless, the decay of both particles can be considered prompt within astrophysical objects. Finally, note that only due to parity conserving nature of the scalar interaction, both left- and right-handed antineutrinos are produced. In this case, only the right-handed antineutrinos (ν^+) are relevant for detection through weak interactions.
  
  
  
  
  § SOLAR ANTINEUTRINOS
  
  
  The flux of MeV antineutrinos from the Sun in the SSM is negligibly small. The largest antineutrino flux at MeV energies comes from small fractions of long-lived radioactive isotopes in the Sun, namely ^232Th, ^238U, and mainly ^40K. This give rise to an antineutrino flux on Earth of about 200 cm^-2 s^-1 with E_ν≲ 3 MeV <cit.>. This component, however, is still 6 orders of magnitude smaller than the geoneutrino flux at the surface of the Earth at these energies, and can be safely neglected. At larger energies, photo-fission reactions produce an even smaller flux of antineutrinos of about 10^-3 cm^-2 s^-1 <cit.>. It is only down at the much lower energies of tens of keV that antineutrinos start being produced in thermal reactions at a similar rate to neutrinos with fluxes as large as 10^9 cm^-2 s^-1 <cit.>.
  
  
  Existing limits on the flux of solar antineutrinos are usually quoted in terms of an energy-independent probability P_ν_e →ν_e of conversion of ^8B neutrinos into antineutrinos. The most stringent limits were obtained by KamLAND in 2011 <cit.>
  
      P^ KamLAND-2011_ν_e →ν_e(E_ν≥ 8.3  MeV) < 5.3 × 10^-5,
  
  which was recently improved in 2021 <cit.>,
  
      P^ KamLAND-2021_ν_e →ν_e(E_ν≥ 8.3  MeV) < 3.5 × 10^-5,
  
  and by Borexino in 2019 <cit.>
  
      P^ Borexino_ν_e →ν_e (E_ν≥1.8  MeV) < 7.2×10^-5,
  
  all at 90 % C.L. In addition, SuperKamiokande (SK) has derived limits on extraterrestrial ν_e sources during phases I, II and III <cit.>, but the high energy thresholds of E_ν > 17.3 MeV make them irrelevant for the study of ^8B neutrinos. For SK phase IV (SK-IV), improvements to the trigger system were implemented and the detection of neutron capture on Hydrogen was made possible, lowering thresholds to E_ν > 13.3 MeV <cit.>. The constraint on solar antineutrino flux was found to be P^ SK-IV-2013_ν_e →ν_e < 4.6× 10^-4. Recently, further improvements to the neutron tagging algorithm lowered this value to E_ν > 8.3 MeV <cit.>, and using the data 2008 to 2018 the limit was improved to
  
      P^ SK-IV-2020_ν_e →ν_e (E_ν≥ 8.3  MeV) < 3.6×10^-4.
  
  A previous preliminary result was shown in Ref. <cit.> and an even more recent update was presented in Ref. <cit.>. Loading of Gd in the SK water tank is expected to greatly improve the neutron tagging efficiency, and would allow for much more stringent limits. With projections on the signal selection efficiency and background reduction, Ref. <cit.> finds that a limit of P^ SK-IV-Gd_ν_e →ν_e (E_ν≥ 8.3  MeV) ≲ 2.2×10^-5 could be achieved with 0.2% Gd loading <cit.>. Lowering the energy threshold of the trigger could further improve these projections.
  
  In addition to these, SNO has also set limits at the level of P^ SNO_ν_e →ν_e (E_ν∈ [4,14.8]  MeV) < 8.3×10^-3 <cit.> at 90% C.L. All limits quoted above assume a total ^8B flux of 5.88 × 10^6 cm^-2 s^-1, except Borexino which assumes 5.46 × 10^6 cm^-2 s^-1, and KamLAND which assumes 5.94 × 10^6 cm^-2 s^-1. At the lowest energies, a bound can also be obtained by noting that the number of elastic ν - e scattering events in solar neutrino experiments decreases if too many ν_4 states are produced, both due to lower ν_e-e cross sections and suppressed ν_e flux. These effects, however, are insensitive to variations of the total ν_e flux below the tens percent level. The predictions from the sterile neutrino decay model are compared with the ^8B flux in fig:B8flux. The independent bounds quoted by KamLAND, Borexino, and SK-IV are shown in fig:fluxlimits as a function of E_ν.
  
  
  The strength of the limits above is mostly due to the large cross section for Inverse beta decay (IBD) on free protons at MeV electron-antineutrino energies. Beyond dominating over the neutrino-electron elastic scattering cross section by about two orders of magnitude (see fig:B8flux), this channel has a distinct signature that drastically reduces backgrounds. After produced, the positron annihilates and the final state neutron is quickly captured by the free protons. This results in a double-bang signal with a positron kinetic energy T_e ≃ E_ν - 1.8 MeV, and a delayed emission of a ≈ 2.2 MeV gamma. The cross section for this process is well understood at high <cit.> and low <cit.> energies, and relatively simple formulae that are valid for all energy regimes have been derived by Ref. <cit.>. In this work we implement the latter calculation, which is provided as machine-friendly data files by Ref. <cit.>.
  
  
  
  
  
    
  Backgrounds For Borexino, reactor neutrinos represent the largest source of backgrounds, but are effectively constrained by DayaBay measurements. Atmospheric neutrino events with genuine IBD scattering or otherwise inherit large uncertainties from the atmospheric neutrino flux and cross sections, but represent only a small contribution (6.5± 3.2 events). Finally, the ^238U and ^232Th geoneutrino fluxes are energetic enough to contribute to the IBD sample, but are only significant up to 3.2 MeV. Borexino omits the contribution of geoneutrinos from the Earth's mantle in their estimation, which is conservative. This component is the most likely explanation for the ≈ 2σ excess seen in the lowest energy bin <cit.> (see also their latest geoneutrino analysis <cit.>). 
  
  The reactor neutrino flux at KamLAND is dominant below 8.3 MeV, but contributes only about 2.2 events above that value. Due to the smaller overburden at KamLAND and SK, they suffer from larger spallation backgrounds, coming mainly from radioactive decays of ^9Li. After muon tagging and fiducial volume cuts, these are reduced to less than 5 events at both locations. The large number of neutrino-electron scattering events presents a background for SK. For this reason, a cut is applied requiring small shower angles with respect to the direction of the Sun, cosθ_⊙ < 0.9. This does not impact IBD events as the positron angle with respect to the incoming neutrino is significantly larger (⟨cosθ⟩≈  0) than in the predominantly forward process of elastic scattering. The observed event spectra and background predictions by the respective collaborations are shown in fig:exp_spectra.
  
  
  
  
   §.§ IBD Rates from Decays
   
  
  
  
  
  
  
  The largest observable flux from sterile neutrinos would come from ν_4 states produced via weak interactions in the decay of ^8B. The number of IBD events in a given experiment can be computed as
  
  
      N_ events/ E_ν_4  E_1   E_3 = 𝒩 Φ^ν_4/ E_ν_4  P^ dec/ E_1  E_3 σ(E_3)  ⟨ P_ν̂_̂ŝ→ν_e⟩(E_3),
  
  
  where 𝒩 stands for total exposure of the experiment, ⟨ P_ν̂_̂ŝ→ν_e⟩ is the flavour transition probability for Solar antineutrinos averaged over the radius of the Sun (see app:Pse), and
  
  
      Φ^ν_4/ E_ν_4    = |U_e4|^2 Φ^ν_e/ E_ν,
      
       P^ dec/ E_1  E_3    = 1/Γ_ν_4Γ_ϕΓ_ν_4→ν̂_̂ŝϕ/ E_1Γ_ϕ→ν̂_̂ŝν̂_̂ŝ/ E_3.
  
  
  Note that eq:nevents is the analogue of Eq. (9) from Ref. <cit.>, and is simpler since we work with very long baselines and under the assumption that the number of initial ν_μ states is negligible. For the flux of ^8B neutrinos, Φ^ν_e/ E_ν, we implement the high-metallicity fluxes in the SSM AGSS09-B16 <cit.>, where the total ^8B neutrino flux is 5.46 × 10^6 cm^-2 s^-1. For low-metallicity models, our constraints on the new physics coupling are weakened by about 20%.
  
  
  With the predicted number of IBD events at each solar neutrino experiment, we implement our statistical test (described in detail in app:statistics) to place limits on the active-heavy mixing angles. Our χ^2 test statistic models solar neutrino flux and experimental backgrounds uncertainties through bin-uncorrelated nuisance parameters with Gaussian errors. Both flux and background uncertainties are fixed at 10%, except for the SK-IV, for which we inflate those to 20%. For the KamLAND 2021 dataset, we combine the data into 2-MeV-wide bins when performing our fit. We have also performed a total rate fit for each of the energy bins to obtain model-independent limits on the solar antineutrino flux. We find similar to those provided by the collaborations, within 50%. Our limits are always weaker, and therefore more conservative, when compared to the official ones. For SK-IV, no model-independent limit was shown in the final article <cit.>, so we show our own result in fig:fluxlimits.
  
  
  
  
  
  § RESULTS
  
  
  
  
  We plot our 90% C.L. limits in fig:limits as a function of |U_e4|^2 and |U_μ4|^2 for m_ϕ/m_4=0.9. On the same axes, we show the preferred (grey-shaded) regions obtained in Ref. <cit.> to explain LSND and MiniBooNE individually, as well as the combined fit to MiniBooNE, LSND and global data (except reactors and cosmology) as “All w/o cosmo". Weaker constraints from the OPERA <cit.> and KARMEN <cit.> neutrino experiments, as well as beta decay kink searches are shown as dashed grey lines.
  We pick two particular cases with the shortest ν_4 and ϕ lifetimes to compare against our limits, corresponding to m_4 Γ_4= 1 eV^2  and m_4Γ_4=10 eV^2. These lifetimes are achieved for couplings close to the perturbativity limit, namely g_ϕ^2 = (1.5)^2 × 10^-2/(|U_e4|^2+|U_μ4|^2) and g_ϕ^2 =  (12)^2 × 10^-2/(|U_e4|^2+|U_μ4|^2), respectively. It is clear that an explanation of LSND is in large tension with solar antineutrino searches for all three experiments we consider. It should also be noted that the region with large |U_μ4|^2 which is not excluded by our curves is excluded by MiniBooNE itself. As expected, KamLAND leads to the strongest bounds despite its large neutrino energy thresholds 8.3 MeV. Borexino and SK-IV bounds are competitive, with the latter performing better for harder antineutrino spectra.
  
  In fig:limits2, we show our constraints for the case m_ϕ/m_4=0.5. A simultaneous explanation of MiniBooNE and LSND is more challenging now, and only a global fit including MiniBooNE but not LSND is available (“All w/o LSND”). In this case, we constrain the region preferred by MiniBooNE significantly. Lower values of m_ϕ/m_4 are even more challenging from the point of view of explaining the SBL results, as ϕ becomes longer lived and helicity-flipping decays of ν_4, which lead to soft daughter spectra, become more important.
  
  
  Changing the heavy neutrino mass but keeping the ratio R = m_ϕ/m_4 fixed leaves our bounds unaltered for all the mass range of interest. Lowering R, on the other hand, weakens our bounds slightly due to the softer solar antineutrino spectrum, although the weakening saturates once below m_ϕ/m_4 ≈ 0.5. We show constraints on the electron mixing angle in fig:limits_ratio. The strongest constraints are obtained for vanishing muon and tau mixing, and for R→1, as in that case ϕ carries most of the energy of its parent particle ν_4. These limits are independent of the absolute scale of m_ϕ and m_4, provided these masses are below the Q-values of the ^8B decays and above the light neutrino masses.
  
  Accounting for perturbativity bounds on g_ϕ and the baselines of LSND and MiniBooNE, a lower bound on m_4 of ≈ 100 eV can be obtained for neutrino mixing of the order of |U_e4|^2≈10^-3. Below this value, ν_4 is too long-lived and does not lead to interesting signatures at short-baseline experiments. The preferred regions for MiniBooNE and LSND shift to larger mixing angles when either ν_4 or ϕ are longer-lived, while our constraints remain unaffected. The parameters used in Figs. <ref> and <ref> are chosen so as to minimize the decay length of new particles. Bounds from kink searches in beta decay become the strongest above m_4 ≳ 5 keV, and peak searches in meson decay preclude an explanation of LSND with the current model at the MeV scale.
  
  
  We note that a massless scalar has also been discussed as an explanation of the MiniBooNE and LSND anomalies <cit.>, although it is disfavored with respect to the best-fit point of Ref. <cit.> at 99% C.L. Interpreting that study in the present lepton-number and parity conserving model, one would obtain no constraint from solar antineutrino searches. In that case, however, light neutrinos will also decay. We comment more on this possibility and others in <ref>.
  
  
  
  
   §.§ Future opportunities
  
  
  Borexino has collected an additional 𝒪(500) days of data on top of the 2771 days already analyzed in Ref. <cit.>. In addition, the improvements made by the collaboration in the latest geoneutrino analysis <cit.>, such as enlarged fiducial volume and improved background rejection, can be implemented in the solar antineutrino search. Just in terms of exposure, this represents an improvement of 40% with respect to the values we use. 
  
  
  
  
    
  Synergy with DSNB Solar antineutrino searches will become even more stringent with upcoming efforts to detect ν_e events from the DSNB <cit.>. The SK detector is expected to detect this neutrino flux with the addition of Gd to its detector volume <cit.>. The large neutron capture cross sections on Gd and the emission of 8 MeV gammas will help reduce backgrounds and lower the ν_e detection threshold to neutrino energies as low as the IBD threshold, provided E_e > 0.8 MeV <cit.>. The increased detection efficiencies at lower energies, and reduced accidental and mis-reconstructed backgrounds will improve on the limits we set, being limited by intrinsic reactor and atmospheric ν_e backgrounds, but also by an exponentially rising spallation background. Future large liquid-scintillator detectors, such as the Jiangmen Underground Neutrino Observatory (JUNO) <cit.>
  and the Jinping neutrino experiment <cit.>, can also improve on current constraints with an expected threshold of E_ν≳ 8.5 MeV <cit.>.  In the far future, observatories capable of accumulating larger numbers of DSNB events, such as the proposed detector THEIA <cit.>, would play an important role in searching for solar antineutrinos. We note that in the event of a detection of the DSNB, one could also constrain the models considered here by requiring small DSNB absorption by relic neutrinos <cit.>. 
  
  
  
    
  Light neutrino decay Even in the parity conserving model discussed so far, one can avoid solar antineutrinos by resorting to a massless ϕ. In that case, however, the light mostly-active neutrinos will decay. For typical parameters relevant for the SBL anomalies, this decay will happen within 1 AU, both visibly and invisibly. For instance, consider ν_2 →ν_1 ϕ and ν_3 →ν_1 ϕ decays with normal ordering and m_1≈0. For a coupling of g_ϕ=1, we find
  
  
      cτ_2^ LAB   ≈ 0.9  AU(10^-5/|U_s1U_s2|^2) (E_2/10 MeV),
          
      
          cτ_3^ LAB   ≈ 0.03  AU(10^-5/|U_s1U_s3|^2) (E_3/10 MeV).
  
  
  In the convention adopted by the literature, these correspond to τ_2^ 0/m_2 = 4.4×10^-5 s/eV and τ_3^ 0/m_3 = 1.3×10^-6 s/eV. For inverted ordering, both ν_1→ν_3ϕ and ν_2→ν_3ϕ decays are of the faster kind in eq:faster_kind. On top of the cosmological issues with such short lifetimes (for recent discussions, see Refs. <cit.>), the largest g_ϕ values relevant for the allowed regions are already excluded by laboratory experiments, such as SNO <cit.>, which is consistent with no disappearance of solar neutrinos (see Refs. <cit.>). Other datasets have also been discussed to constrain the lifetime of light neutrinos, including measurements of the flavor ratios of cosmic neutrinos <cit.> and of the Glashow resonance <cit.>. It should be noted, however, that light-sterile mixing parameters governing light neutrino decay are related to those of SBL anomalies in a model-dependent fashion. In principle, but not without fine-tuning, the correlation between |U_e4|, |U_μ4|, |U_s4^*U_sj|, and |U_si^*U_sj| for i,j<4 may be relaxed. 
  In models with LN violation or LN charged scalars (see below), provided several constraints are satisfied, solar antineutrinos may become relevant again for massless ϕ as the light-neutrino decays ν_2→ν_1ϕ are open.
  
  
  
  
  § ALTERNATIVE MODELS: VIOLATING PARITY AND LEPTON NUMBER
  
  
  Various other possibilities for visible sterile neutrino decay exist, depending on the Dirac or Majorana nature of neutrinos, as well as on the parity structure of the sterile neutrino sector. While we only focused on the parity conserving model discussed above, we would like to dedicate this section to understanding if other minimal extensions of the SM by a singlet sterile neutrino and a scalar are subject to our constraints. For clarity, we focus on SM extensions with a single new sterile neutrino: ν_s=ν_s^L + ν_s^R in the Dirac case and ν^R in the Majorana case. Our findings for the minimal models are summarized in tab:minimal_models.
  
  
   
  
  
  
    
  Dirac L(ϕ)=0 We start with a generalization of eq:Lscalar by writing
  
      ℒ_int^D-0   = g_ϕν^L_sν^R_sϕ + h.c.
         = g_ϕ (V_L)_si^*(V_R)_sjν_iP_Rν_j ϕ + h.c.,
  
  where index summation is understood. For complex g_ϕ, this is the most generic parametrization of the vertex. We implicitly diagonalized the Dirac mass matrix by means of two unitary matrices V_L and V_R, defined by ν^L_s = (V_L)_siν^L_i, ν^R_s = (V_R)_siν_i^R, where ν_i^L,R are the (chiral) mass eigenstates. Note that the enlarged PMNS matrix is defined as U_ PMNS^' = V_L when the charged lepton mass matrix is diagonal. Abandoning the assumption of parity conservation in the sterile sector that was made previously, V_L=V_R, one can have allow for V_L ≠ V_R by choosing different Yukawa couplings for ν_s^L and ν_s^R. By breaking parity at the level of the Dirac mass matrix, it is possible to independently tune the couplings appearing in the operators ν_i P_R ν_4 and ν_iP_Lν_4. In practice, this allows to tune the rate for visible and invisible decays of ν_4^- neutrinos. The same is true for the decay of the scalars, which can be either visible or invisible, depending on V_L and V_R. In these models, a connection to the SBL anomalies through visible decays always predicts visible solar antineutrinos provided ϕ is heavy enough to decay.
  
  
  
    
  Dirac |L(ϕ)|=2 One can also introduce scalars carrying LN. These type of scalars have been usually discussed in the context of Majoron models, but for our current purposes, we assume no particular connection to neutrino masses. We consider a model with a Dirac field ν_s, and a complex scalar ϕ carrying lepton number L(ϕ)=-2. In all generality, we can write
  
      ℒ_int^D-2   = g_L (ν_s)^c P_L ν_s ϕ + g_R (ν_s)^c P_R ν_s ϕ + h.c.
         = g_L (V_L)_si(V_L)_sjν_iP_Lν_j ϕ
              + g_R (V_R)_si (V_R)_sjν_iP_Rν_j ϕ +h.c.,
  
  where again we implicitly diagonalized the Dirac mass matrix with V_L and V_R. In this case, even for parity conserving mass matrices, one can violate parity in the neutrino-ϕ interactions by tuning the arbitrary g_L and g_R couplings. In this case, 
  ν_4^- states produced in the Sun will always decay to visible antineutrinos provided g_R≠0, independently of the decay products of ϕ^*. For this model, one may attempt to explain SBL anomalies with only the decay products of scalar produced in ν_4^- decays by setting g_R → 0. In that case, no visible solar antineutrinos appear. 
  
  
  
    
  Majorana neutrinos A final possibility is to abandon LN and work with Majorana neutrinos. In this case, a minimal model can be built with only ν^R and a scalar ϕ. LN is violated by the ν^R Majorana mass term, and the most general interaction Lagrangian in this case is 
  
      ℒ_int^ M   = g_R (ν^R)^cν^R ϕ + h.c.
         = g_R (V)_si (V)_sjν_iP_Rν_j ϕ +h.c.,
  
  where now we implicitly diagonalized the Majorana mass matrix by means of a single unitary matrix V=U_ PMNS^'. In this case, all light neutrinos as well as antineutrinos are visible due to the reduced number of degrees of freedom. Both HC and HF decays of ν_4 are controlled by the same parameters, and cannot be disentangled as easily. Solar antineutrinos could appear in this case if all other constraints are satisfied.
  
  
  
    
  Simplified models Finally, we note that Refs. <cit.> work with simplified models, and do not specify the origin of the ℒ⊃ g_e ν_eν_4 ϕ vertex. Although this operator may arise from a Lagrangian as simple as eq:Lscalar, it can be considered more generically as a by-product of non-renormalizable operators such as (LH)^2 ϕ and (LH)ν_sϕ. In these effective models, the active-heavy mixing necessary for ν_4 production in most accelerator experiments, |U_μ4|^2, is independent from g_e, which controls the decay rate of ν_4→ν_eϕ. In this case, the only way to generate ν_e appearance at LSND is via muon decays, μ^+→ e^+ν_eν_4. It also follows that the mixing |U_e4|^2 may be parametrically small, turning off ν_4 production in the Sun via mixing. Four-body decays of the type ^8B→^8Be e^+ ν_4 ϕ are negligible as kaon decays constrain g_e^2 |U_μ 4|^2 < 𝒪(10^-7).
  
  
  If a vector particle is introduced instead, the cosmological history is yet even more involved. We do not study this case here, although our solar antineutrino bounds would also apply to parity-conserving scenarios with small modifications. Note that our constraints are not relevant for fully invisible sterile neutrino decays, as invoked to relax the tension between SBL appearance and disappearance tension in Refs <cit.>.
  
  
  
  
  § CONCLUSIONS
  
  
  Puzzling results from some of the short-baseline neutrino experiments will eventually find an explanation with more data coming from the SBN program at Fermilab <cit.> and the π^+ decay-at-rest experiment at J-PARC, JSNS^2 <cit.>. At the moment it is possible to speculate that some form of new physics in the neutrino sector is responsible for the deviation of LSND and MiniBooNE results from theoretical expectations within the minimal three-generation neutrino model. Among such speculations is the class of models where the excess of antineutrinos at LSND and excess of low-energy electron-like events at MiniBooNE is due to a prompt production and decay of dark sector particles. This new sector is likely to comprise a heavier, mostly sterile neutrino, that can be produced via neutrino mixing in meson decays and nuclear reactions. Such heavier neutrino can generate a cascade decay to an unstable bosonic mediator and light neutrino, giving rise to the admixture of electron antineutrinos in the flux at the end of the decay chain. 
  
  We have shown that up to some model dependence one should expect that regular nuclear processes in the Sun create an antineutrino flux. Such flux is stringently constrained by most of the solar neutrino experiments, at a O( few× 10^-5) level owing to a larger cross sections for the IBD processes, and additional structure to the signal that has been exploited to cut on backgrounds. After application of these constraints, our results disfavor large part of the parameter space of the model in Ref. <cit.>, and disfavor this mechanism as an explanation of the LSND excess, while significantly narrowing possible parameter space for the MiniBooNE excess. The simulations used to produce the results in this paper are publicly available on github[ https://github.com/mhostert/SolarDecayingSteriles  github.com/mhostert/SolarDecayingSteriles.].
  
  In general, our limits disfavor large ν→ν transitions that could improve the combined fit of LSND and MiniBooNE data. Such transitions could in principle be avoided if ϕ is lighter than the lightest neutrino state, in which case, mixing angles and CP phases have to be fine-tuned to avoid light neutrino decays. The alternative models with parity violation or apparent LN violation presented in sec:alternatives may avoid ν→ν transitions even for massive ϕ, but would require a case-by-case study of the SBL physics and additional constraints. 
  
  Our constraints add to the existing list of problems of the decaying sterile neutrino solutions to the SBL puzzle. Chiefly among them is cosmology and astrophysics. As is well known, new and relatively strongly interacting states can be populated by the thermal processes leading to the modifications of observed quantities, such as primordial nucleosynthesis yields and/or total amount of energy density carried by neutrinos at late times. In addition, these models are likely to cause strong modifications to the supernovae neutrino spectrum. One reason for such modification is the possibility of the neutrino number-changing processes, such as νν→ννφ→νννν. Given relatively strong couplings in the models of Refs. <cit.>, the underlying cross sections are far greater than weak interaction cross section, meaning that the neutrinos can share energy and maintain their chemical equilibrium immediately after they leave the star. The main physical effect, the degrading of average energy for the SN neutrinos, can be constrained with the observed signal of SN1987A. This has been explored to constrain neutrino self-interactions inside supernovae by the requirement that neutrinos carry sufficient energy to the outer layers of the collapsing star <cit.>. We point out, however, that a more general statement can be made regarding neutrino energy loss outside the dense environment, which is independent of the explosion mechanism. Details of this effect will be addressed in a future publication.
  
  
  
  
  The authors would like to thank Linyan Wan and Sandra Zavatarelli for correspondence on the SK-IV and Borexino experimental capabilities. We also thank Ivan Esteban and Joachim Kopp for discussions. MP is supported in part by U.S. Department of Energy (Grant No. desc0011842). This research was supported in part by Perimeter Institute for Theoretical Physics. Research at Perimeter Institute is supported by the Government of Canada through the Department of Innovation, Science and Economic Development and by the Province of Ontario through the Ministry of Research, Innovation and Science.
  
  
  
  
  
  
  
  
  § SOLAR FLAVOR TRANSITIONS
  
  
  
  When ϕ decays to light neutrinos, it decays into the state |ν̂_s⟩ = ∑_i=1^3 U_s i^* |ν_i⟩. The average transition probability for ν̂_s to exit the Sun and be detected as a ν_e on the surface of the Earth under the (good) approximation of adiabatic flavor conversion is simply 
  
      ⟨ P_ν̂_S→ν_e⟩_N = ⟨∑_i=1^3 |U_si^m|^2|U_ei|^2⟩/|U_e 4|^2+|U_μ 4|^2+|U_τ 4|^2,
  
  which depends on the mixing matrix elements in matter |U_si^m|^2 at the production point. Here, ⟨…⟩ denotes a weighted average over the production region and the subscript N refers to taking the non-canonical normalization of ν̂_s into account. We have neglected Earth matter effects, which for antineutrinos leads to a reduction (increase) of P_2→ e (P_1→ e) below ∼10%, and assumed the unitarity of the 4×4 mixing matrix U. Neglecting all CP phases, we follow Ref. <cit.> and write U = R_23 R_34R_24R_14 R_13 R_12, where R_ij=R(θ_ij) is the usual rotation matrix in the (i,j) plane. Note that if we assume U_τ 4=-c_34s_23s_24+c_23s_34=0, then to leading order in the small angles, |U_e4|^2∼θ_14^2 and |U_μ4|^2∼ c_23^2θ_24^2, with the rest of the mixing matrix elements of the active 3× 3 sub-matrix retaining their usual definition. 
  
  
  In the limit where Δ m^2_21≪Δ m^2_31≪Δ m^2_41, the flavor evolution of solar neutrinos can be described by an effective two-neutrino model, with in-matter modifications only on θ_12. Up to corrections proportional to the new mixing angles (θ_14, θ_24, and θ_34) as well as to θ_13, the new mixing angle in matter is
  
      tan2θ^M_12≃sin2 θ_12Δ m^2_12/Δ m^2_12cos2 θ_12 - A_ CC,
  
  with A_ CC = ± 2√(2) E_ν G_F N_e(r) for neutrinos (antineutrinos) proportional to the electron density at the production region. Since the sterile component does not feel neutral-current interactions, the neutral-current potential also modifies the relation above. Terms proportional to A_ NC = ∓√(2) E_ν G_F N_n(r), where N_n(r) is the neutron number density in the Sun, are small, however. This is because they are proportional to the new small mixing angles, as well as due to the smaller number of neutrons in the Sun, N_ n/N_e≲ 1/2. In our simulation, we follow Ref. <cit.> and keep all corrections in θ_14, θ_24, θ_34, and θ_13 in the transition probabilities. As an approximation, we assume the neutral-current potential to follow the same radial dependence as the charged-current one, with A_ NC/A_ CC = -1/4.
  
  The full transition probabilites using eq:pse in the energy region of the ^8B flux are shown in Fig. <ref> for both neutrinos and antineutrinos. In the limit |U_μ4|,|U_τ4|→0, the scalar decays produce mostly-ν_e states, and P_ν̂_s→ν_e≃ P_ν_e→ν_e. In this case, the flavor evolution is similar to the standard MSW effect, where neutrino (antineutrinos) undergo resonant (non-resonant) adiabatic flavor conversion, in which production of ν_2^m at the center of the Sun is enhanced (suppressed). For |U_μ4|,|U_τ4|≠0, the situation is more complex, but the high energy behavior can be understood by taking the limit θ_12→π/2 (0) in the mixing elements |U_si|^2 for neutrinos (antineutrinos). Note that ν_3^m ∼ν_3, as it should be since Δ m^2_21≪Δ m^2_31.
  
  
  
  § STATISTICAL METHOD
  
  
  When deriving upper limits on the mixing angles, we minimize the following log-likelihood function  
  
  
      ℒ = 2    ∑_i[ μ_i (θ⃗,β⃗) - D_i + D_i lnD_i/μ_i(θ⃗,β⃗)]
      
         + ∑_i,jβ_i,j^2/σ_i,j^2,
  
  
  where θ⃗ stands for the vector of physics parameters (, |U_α|^2), β⃗ the vector of nuisance parameters with individual entries β_j and associated Gaussian errors σ_j. As an approximation, we assume ℒ to follow a χ^2 distribution when estimating our confidence intervals. 
  
  The most important systematics for our study are the uncertainties on the total ^8B solar neutrino flux and total backgrounds numbers. To be conservative, we assign each energy bin two normalisation systematics, one exclusive to the new physics prediction, modelling uncertainties in the solar flux, and one exclusive to backgrounds. All normalization systematics are assumed to be uncorrelated, which is conservative, and are assigned 10% Gaussian errors.
  
  
  
  
  § POLARIZED DECAY RATES
  
  
  To produce tab:minimal_models, we computed the decay rates in each channel explicitly. We collect all results for ν_4 and ϕ decays assuming massless neutrino final states in each one of the models discussed. The total decay rate for ν_i^h_i→(—)ν_j^h_jϕ can be obtained by summing each polarized matrix element as
  
      Γ_4    = ∫_x^ min^x^ max∑_h_i,h_j|M_h_i h_j|^2/m_4 β 16 π x,
  
  where β=|p_4|/E_4 is the velocity of ν_4 in the laboratory frame, and x^ max_ min =  (1±β)(1-r_ϕ^2)/2. Similarly, for ϕ→(—)ν_i^h_i(—)ν_j^h_j decays,
  
      Γ_ϕ   = ∫_x_ϕ^ min^x_ϕ^ max∑_h_i,h_j|M_h_i h_j|^2/ m_ϕ 16π x_ϕ,
      
  where β_ϕ=|p_ϕ|/E_ϕ is the ϕ velocity in the laboratory frame and (x_ϕ)^ max_ min = (1±β_ϕ)/2. 
  
  
  
  
   §.§ Dirac L(ϕ)=0 case
  
  
  Making use of eq:LD0 and neglecting light neutrino masses, the amplitude squared for ν_i^h_i→ν_j^h_jϕ decays is given by
  
      |M_–|^2    = |g_ϕ (V_L)_si (V_R)_sj|^2 m_i^2 C(x)
          
      
          |M_++|^2    = |g_ϕ (V_R)_si (V_L)_sj|^2 m_i^2 C(x),
          
      
          |M_-+|^2    = |g_ϕ (V_R)_si (V_L)_sj|^2 m_i^2 F(x),
          
      
          |M_+-|^2    = |g_ϕ (V_L)_si (V_R)_sj|^2 m_i^2 F(x),
  
  which are also valid for ν_i^h_i→ν_j^h_jϕ decays. We have defined
  
      C(x)    =2 x - (1-β)(1-r_ϕ^2)/2β,
      
          F(x)    =(1+β)(1-r_ϕ^2)-2x/2β,
  
  which apply for helicity conserving and helicity flipping channels, respectively. Note that C(x) grows while F(x) decreases monotonically with x.
  For scalar decays ϕ→ν^h_i_i ν^h_j_j, we compute |M_h_1 h_2|^2 to find
  
      |M_–|^2    = |g_ϕ (V_L)_si (V_R)_sj|^2 m_ϕ^2,
          
      
          |M_++|^2    = |g_ϕ (V_R)_si (V_L)_sj|^2 m_ϕ^2,
  
  with all other combinations vanishing in the limit of massless final states. 
  
  
  
   §.§ Dirac |L(ϕ)|=2 case
  
  
  Now, switching to eq:LD2, the amplitudes for ν_i^h_1→ν_j^h_2ϕ decays are
  
      |M_– |^2    = 4|g_R  (V_R)_si (V_R)_sj|^2 m_i^2 C(x),
          
      
          |M_++|^2    = 4|g_L  (V_L)_si (V_L)_sj|^2 m_i^2 C(x),
          
      
          |M_-+ |^2    = 4|g_L  (V_L)_si (V_L)_sj|^2 m_i^2 F(x),
          
      
          |M_+-|^2    = 4|g_R  (V_R)_si (V_R)_sj|^2 m_i^2 F(x).
  
  The amplitudes for ν_i^h_1→ν_j^h_2ϕ decays can be obtained with the substitution L↔ R. For scalar decays ϕ→ν^h_1_i ν^h_2_j, we compute |M_h_1 h_2|^2 to find
  
      |M_–|^2    = 4 |g_L (V_L)_si (V_L)_sj|^2 m_ϕ^2,
          
      
          |M_++|^2    = 4 |g_R (V_R)_si (V_R)_sj|^2 m_ϕ^2,
  
  where again the amplitudes for decays into antineutrinos can be obtained by L↔ R.
  
  
  
   §.§ Majorana case
  
  
  In the Majorana case, one can make use of the expressions for the Dirac L(ϕ)=0 case, keeping in mind that V_R=V_L, and that an additional overall factor of 2 should be included for the ν_i total decay rates and an overall factor of 2/(1+δ_ij) for the total ϕ→ν_i ν_j decay rate.
  
  apsrev4-1
  
  
  
  