Dual-Baseline Search for Active-to-Sterile Neutrino Oscillations in NOvA American Physical Society 02/01/2025 10588186 10.1103/PhysRevLett.134.081804 Physical Review Letters 0031-9007 134 8 M A AceroB AcharyaP AdamsonN AnfimovA AntoshkinE Arrieta-DiazL AsquithA AurisanoA BackN BalashovP BaldiB A BambahE F BannisterA BarrosA BatK BaysR BernsteinT_J C BezerraV BhatnagarD BhattaraiB BhuyanJ BianA C BoothR BowlesB BrahmaC BrombergN BuchananA ButkevichS CalvezT J CarrollE Catano-MurJ P CesarA ChatlaR ChircoB C ChoudharyA ChristensenM F CicalaT E CoanA CooleybeckC Cortes-ParraD CoveyouL CremonesiG S DaviesP F DerwentP DingZ DjurcicK DobbsM DolceD DoyleD Dueñas TonguinoE C DukesA DyeR EhrlichE EwartP FilipM J FrankH R GallagherF GaoA GiriR A GomesM C GoodmanM GrohR GroupA HabigF HaklJ HartnellR HatcherH HausnerM HeK HellerV HewesA HimmelT HorohoA IvanovaB JargowskyJ JaroszM JudahI KakorinA KalitkinaD M KaplanB Kirezli-OzdemirJ KleykampO KlimovL W KoernerL KolupaevaR KralikA KumarV KusT LackeyK LangJ LesmeisterA ListerJ LiuJ A LockM LokajicekM MacMahonS MagillW A MannM T ManoharanM Manrique PlataM L MarshakM Martinez-CasalesV MatveevB MehtaM D MessierH MeyerT MiaoV MikolaW H MillerS MishraS R MishraA MislivecR MohantaA MorenA MorozovaW MuL MualemM MuetherD MyersD NaplesA NathS NelleriJ K NelsonR NicholE NinerA NormanA NorrickH OhA OlshevskiyT OlsonM OzkaynakA PalJ PaleyL PandaR B PattersonG PawloskiR PettiR K PlunkettL R PraisM RabelhoferA RafiqueV RajM RajaoalisoaB RamsonB RebelP RoyO SamoylovM C SanchezS Sánchez FaleroP ShanahanP SharmaA SheshukovA ShmakovShivamW ShorrockS ShuklaD K SinghaI SinghP SinghV SinghE SmithJ SmolikP SnopokN SolomeyA SousaK SoustruznikM StraitL SuterA SuttonK SuttonS SwainC SweeneyA SztucB Tapia OreguiN TalukdarP TasT ThakoreJ ThomasE TirasM TitusY TorunD TranJ TripathiJ Trokan-TenorioJ UrheimP VahleZ VallariJ D VillamilK J VockerodtM WallbankC WeberM WetsteinD WhittingtonD A WickremasingheT WieberJ WolcottM WrobelS WuW WuW WuY XiaoB YaeggyA YahayaA YankelevichK YoneharaS ZadorozhnyyJ ZalesakR ZwaskaNOvA_Collaboration We report a search for neutrino oscillations to sterile neutrinos under a model with three active and one sterile neutrinos (3+1model). This analysis uses the NOvA detectors exposed to the NuMI beam, running in neutrino mode. The data exposure,13.6×1020protons on target, doubles that previously analyzed by NOvA, and the analysis is the first to useνμcharged-current interactions in conjunction with neutral-current interactions. Neutrino samples in the near and far detectors are fitted simultaneously, enabling the search to be carried out over aΔm412range extending 2 (3)orders of magnitude above (below)1eV2. NOvA finds no evidence for active-to-sterile neutrino oscillations under the3+1model at 90% confidence level. New limits are reported in multiple regions of parameter space, excluding some regions currently allowed by IceCube at 90% confidence level. We additionally set the most stringent limits for anomalousντappearance forΔm412≤3eV2.

Published by the American Physical Society2025]]>

Neutrino mixing is a well-established phenomenon, with numerous experiments reporting results that agree with a picture including three neutrino mass states (ν 1 , ν 2 , ν 3 ) that mix to form three neutrino flavors (ν μ , ν e , ν τ ) [1][2][3][4][5][6][7][8][9][10][11]. In this three-flavor (3F) framework, oscillations are governed by two mass-squared splittings, Δm 2 21 and Δm 2 32 , where Δm 2 ji ≡ m 2 j -m 2 i , which correspond to the frequency of oscillation for a given neutrino energy (E ν ) and path length (L), three mixing angles, θ 12 , θ 13 , and θ 23 , which drive the magnitude of oscillation, and a CP violating phase, δ CP , which allows for differences between neutrino and antineutrino oscillations. Over the past two decades, a number of anomalous results have been reported in short baseline accelerator neutrino experiments [12,13], radiochemical experiments [14][15][16], and in the reactor neutrino sector [17]. If these anomalies are interpreted as neutrino oscillations, they would require Δm 2 ∼ 1 eV 2 ≫ Δm 2 32 ; Δm 2 21 , necessitating additional neutrino mass states be added to the model. Measurements of the width of the Z boson at the LEP experiments indicate that any additional neutrinos with m ν < m Z 0 =2 must be sterile [18], meaning that they do not interact via the weak force. The global picture of a mass splitting in the ∼1 eV 2 region is complicated by the presence of a number of null results [19][20][21][22][23][24][25][26][27][28][29][30][31][32][33].

NOvA can probe for active-to-sterile neutrino oscillations by searching for disappearance of neutral-current (NC) interactions, which provides a flavor-agnostic measurement of the active neutrino event rate. We can additionally search for active-to-sterile oscillations among ν μ charged-current (CC) interactions by testing for additional sources of disappearance when compared to 3F oscillations. The NOvA experiment consists of two functionally identical detectors placed 14.6 mrad off-axis of Fermilab's NuMI beam [34]. The NuMI beam is extracted over 10 μs approximately every 1.3 s when 120 GeV protons strike a graphite target resulting in a secondary hadron beam. These hadrons are focused using two magnetic horns, and decay to neutrinos as they travel through a 675 m helium-filled decay pipe. The off-axis placement of the detectors results in a narrow neutrino-energy distribution peaked around 2 GeV, with a width of 0.4 GeV and a subdominant highenergy tail.

The NOvA Near Detector (ND) is positioned at Fermilab in Batavia, Illinois, 1 km downstream of the NuMI target and 100 m underground, while the Far Detector (FD) is located 810 km from the target, at Ash River, Minnesota, on the surface with a 3 m water-equivalent overburden. The NOvA detectors are tracking calorimeters with the ND (FD) constructed of 192 (896) planes of highly reflective PVC cells measuring 3.9 cm × 6.6 cm with a length of 3.9 m (15.5 m) [35]. The planes are alternately placed with horizontal and vertical cells to enable three-dimensional reconstruction, and are filled with a blend of mineral oil based liquid scintillator doped with 5% pseudocumene [36]. The ND has additional planes of instrumented cells separated by steel plates at the rear of the detector ("muon catcher") to range out muons. Light produced by charged particles traversing a cell is collected by a single loop of wavelength-shifting fiber that spans the length of the cell and is read out on both ends of the fiber at one end of the cell by one pixel of a 32 pixel avalanche photodiode. Custom readout electronics are used to shape and digitize the signal, and any signal meeting a minimum pulse height requirement within a 550 μs window around the beam pulse is saved for offline analysis. The cosmic background at the FD is sampled by a 10 Hz minimum bias trigger [37].

The simplest extension to 3F mixing is the 3 þ 1 model [38][39][40][41][42], which introduces seven new parameters: Δm 2 41 , θ 14 , θ 24 , θ 34 , δ 14 , δ 24 , and δ 34 . Under this model, the active neutrino survival probability can be approximated as [27] 1

14 cos 2 θ 34 sin 2 2θ 24 sin 2 Δ 41 sin 2 θ 34 sin 2 2θ 23 sin 2 Δ 31 þ 1 2 sin δ 24 sin θ 24 sin 2θ 23 sin Δ 31 ; ð1Þ while the ν μ survival probability can be approximated as Pðν μ → ν μ Þ ≈ 1sin 2 2θ 24 sin 2 Δ 41 þ 2sin 2 2θ 23 sin 2 θ 24 sin 2 Δ 31 sin 2 2θ 23 sin 2 Δ 31 ; ð2Þ where Δ ji ≡ ðΔm 2 ji L=4E ν Þ. Exact oscillation probability calculations are used for the analysis. Equations (1) and (2) show that sterile neutrinos modify the magnitude of oscillations at the atmospheric frequency, Δm 2 31 , and introduce a new sterile frequency driven by Δm 2 41 . At the NOvA FD, the frequency of sterile oscillations at Δm 2 41 > 0.05 eV 2 is too high for NOvA to resolve. These oscillations manifest as a downward normalization shift in the neutrino energy spectrum. At the shorter baseline of the NOvA ND, energy-dependent sterile oscillations arise when Δm 2 41 > 0.5 eV 2 (Fig. 1). FIG. 1. Oscillation probabilities for (a) active neutrino survival probability and (b) muon neutrino survival probability vs E ν and L=E ν with various model parameters. 3F survival probabilities are shown in black. Probabilities under the 3 þ 1 model are shown with Δm 2 41 of 0.05 eV 2 (blue), 0.5 eV 2 (red), 5 eV 2 (green). As the value of Δm 2 41 increases, oscillations happen over shorter baselines, resulting in noticeable oscillations in the ND. The shaded regions approximately correspond to the fraction of neutrinos for each E ν and L=E ν in each detector.

Fitting both detectors simultaneously extends the region of parameter space to which we are sensitive compared to previous analyses [43,44]. In the ν μ CC and NC channels considered, NOvA is sensitive to the atmospheric 3F oscillation parameters, θ 23 and Δm 2 32 , as well as to the sterile-related parameters θ 24 , θ 34 , Δm 2 41 , and δ 24 , but not to θ 14 as this analysis does not consider ν e appearance.

This analysis uses data collected between February 2014 and March 2020, corresponding to 11.0 × 10 20 (13.6 × 10 20 ) protons on target (POT) for the ND (FD). Approximately 0.1 × 10 20 POT of ND NC selected events were removed from the sample for preanalysis validation, meaning this sample corresponds to 10.9 × 10 20 POT.

The neutrino flux at the NOvA detectors is determined using a simulation of particle production and transport through the beamline using GEANT4 9.2p03 [45][46][47]. The simulated flux is then corrected using the Package to Predict the Flux, which modifies the prediction using external hadron-production data [48]. Neutrino interactions are simulated within the NOvA detectors using GENIE 3.0.6 [49,50]. Additional information about the GENIE configuration used can be found in [51]. The outgoing particles are propagated through the detector geometry using GEANT4 10.4p02. Custom routines are used to simulate the capture and transport of scintillation light, as well as the response of the avalanche photodiodes.

The NC interaction candidates are characterized by hadronic activity resulting from the transfer of energy and momentum from the neutrino to the nucleus. The final state lepton is a neutrino and is not detected. All NC candidates are required to have a reconstructed vertex, at least one reconstructed particle, and must cross at least three contiguous planes.

In the ND, the reconstructed vertex is required to be contained within a volume with boundaries 80 cm from the top, bottom, and sides of the detector, 150 cm from the front face, and 260 cm from the rear face excluding the muon catcher. All reconstructed particles are required to be contained within a volume excluding 20 cm from the top, bottom, and sides of the detector, 150 cm from the front face, and 50 cm from the rear face excluding the muon catcher. The more stringent requirements on distance from the front and rear faces of the detector compared with the requirements on the other faces are selected to reject background candidates due to interactions in rock upstream of the detector and CC beam interactions, respectively.

For FD NC candidates, reconstructed particles are required to be fully contained within a fiducial volume with boundaries 100 cm from the top, bottom, and sides of the detector, and 160 cm from the upstream and downstream faces. The cosmic background interaction rate is significantly higher at the FD than the ND due to a shallower overburden. Accordingly, rather than placing explicit requirements on the vertex position, we use this information along with information about the reconstructed shower shapes and energy, number of hits, and the transverse momentum fraction to construct a boosted decision tree focused on rejecting cosmic backgrounds.

Signal events are selected using a convolutional neural network [52,53], which provides probability scores for different neutrino flavors based on energy depositions in the detector. Optimal requirements on the score are determined by using a figure of merit that considers the systematic and statistical uncertainties on the samples. In the FD, the requirements on convolutional neural network score and cosmic rejection score are jointly tuned. Any event passing the 3F ν μ CC or ν e CC selection [51] is additionally removed from the sample of NC interaction candidates.

The deposited energy of NC interaction candidates is estimated by taking a weighted sum of the hadronic and electromagnetic components of the calorimetric energy in the detector. An additional bias correction is applied as a function of total calorimetric energy. The overall NC energy resolution is 30% [54]. The event selection criteria and neutrino energy estimator used for the ν μ CC samples in the ND and FD are described in [51].

We consider systematic uncertainties on the beam flux, neutrino interactions, and detector modeling [51]. For this analysis, we identified two sources of uncertainty that required custom handling.

Typically for oscillation analyses using NOvA's extrapolation technique [43,44,51], the ND is assumed to have no oscillations, so any differences between simulation and data can be attributed to mismodeling in the simulation. We then tune cross-section models in the ND simulation to the ND data, producing a new central value (CV) and suite of uncertainties. Because sterile neutrinos may induce oscillations in the ND, differences between data and simulation cannot be attributed to cross-section mismodeling, and so we use untuned simulation and uncertainties. Because the meson exchange current component of the simulation is the least well understood, we have developed shape and normalization uncertainties for this component based on the model spread of the València [55], SuSA [56], and GENIE empirical [57] meson exchange current models.

Many NC neutrino candidates selected for this analysis are produced by kaon decays. Because of the lack of available hadron-production data, prior analyses assigned the beam kaon component a 30% normalization uncertainty in addition to Package to Predict the Flux uncertainties. We instead constrain this uncertainty with samples not used in the analyses: a horn-off data sample, which allows us to probe hadron-production uncertainties without the complications of the focusing horns, and a sample of uncontained high-energy muon neutrinos, which gives us access to the focused kaon peak. We fit for the kaon component normalization marginalizing over potential sterile oscillations across the region of parameter space used in this analysis. This technique results in a 10% uncertainty on this component.

PHYSICAL REVIEW LETTERS 134, 081804 (2025)

For each systematic uncertainty we randomly vary model parameters within their uncertainties to produce a new systematically fluctuated "universe," u. A covariance matrix is constructed,

where N iðjÞ represents the number of events in the ith (jth) energy bin and U is the total number of universes. The C i;j for each systematic uncertainty are summed, producing a final systematic covariance matrix. We use two independent analysis techniques for this search. Analysis 1 employs a hybrid test statistic combining a Poisson log-likelihood treatment of statistical uncertainties with a Gaussian multivariate treatment of systematic uncertainties. A covariance matrix encoding only the systematic uncertainties is used to fit for optimal systematic weights, s αi , for each oscillation channel α and analysis bin i. This test statistic is expressed as

where O i are the observed data, S i ¼ P α s αi p αi represents the systematic weights applied to the prediction, p αi , and C αiβj is a covariance matrix encoding the systematic uncertainty in each oscillation channel (α, β) and analysis bin (i, j) and their correlations.

Analysis 2 employs a traditional Gaussian multivariate formalism,

where P iðjÞ (O iðjÞ ) is the number of predicted (observed) events in bin iðjÞ. This analysis adds statistical uncertainties to the diagonal of the covariance matrix using the combined Neyman-Pearson formalism, V ij ¼ f3=½ð1=O i Þþ ð2=P i Þgδ ij , which yields a smaller bias in best fit model parameters than either the Neyman or Pearson construction [58].

For both analyses the 3F atmospheric oscillation parameters θ 23 and Δm 2 32 are varied in the fit, with a loose Gaussian constraint applied to Δm 2 32 to pin the fit to a 3 þ 1 flavor paradigm. This constraint is centered at j2.51j AE 0.15 × 10 -3 eV 2 in both mass orderings, and is derived from a 2020 global fit to data [59] including atmospheric neutrino oscillations. The width of the constraint is conservatively set to double the 3σ range. The sterile parameters Δm 2 41 , θ 24 , θ 34 , and δ 24 are also freely varied when fitting. The other sterile parameters are held fixed at 0 in the fit, due to constraints from solar and reactor experiments [60] and unitarity [61].

We select 2 826 066 (103 109) ν μ CC (NC) candidates from the ND data compared with the 3F prediction of 2 450 000 AE 530 000 (115 000 AE 30 000). Additionally, we select 209 [62] (469) ν μ CC (NC) candidates from the FD data, compared with a 3F prediction of 200 AE 45 (471 AE 116) using prefit parameter values Δm 2 32 ¼ 2.51 × 10 -3 eV 2 and θ 23 ¼ 49.6° [59].

The best fits from the two analyses agree well with the data and with each other (Fig. 2). The technique used by Analysis 1 allows us to present the best fit with systematic uncertainty pulls applied (dashed), resulting in improved agreement between data and simulation. This agreement indicates that any discrepancy between data and simulation can be accounted for by our systematic uncertainties.

We present 90% limits for both the Δm 2 41sin 2 θ 24 and Δm For each, a Δχ 2 surface is constructed over a grid of the two parameters that define the space, with a fit performed for the remaining parameters at each point. To correct the critical χ 2 values we use a hybrid Feldman-Cousins technique. Covariance matrix fitting techniques do not result in fitted pull terms for each systematic uncertainty. Accordingly, for systematic uncertainties we use a Highland-Cousins technique [63], where for each new universe a value of each systematic parameter is drawn from its a priori distribution. For oscillation parameters, we use the Profiled Feldman-Cousins technique [64].

NOvA's Δm 2 41sin 2 θ 24 limits [Fig. 3(a)] are competitive at high Δm 2 41 , and exclude new regions of interest in the IceCube 90% allowed region [65], around 6 eV 2 < Δm 2 41 < 11 eV 2 . Analysis 1 excludes slightly more values of sin 2 θ 24 at low Δm 2 41 and Analysis 2 excludes slightly more at high Δm 2 41 . Sensitivity to sin 2 θ 24 primarily comes from muon neutrino disappearance, which is independent of θ 34 [Eq. ( 2)]. For high values of Δm 2 41 sensitivity is driven by the ND data, meaning differences in the limits come from different handling of the systematic uncertainties. Sensitivity at low Δm 2 41 arises primarily from FD data, meaning that the weaker limit in Analysis 2 comes from undercoverage of the combined Neyman-Pearson statistical technique.

Our Δm 2 41sin 2 θ 34 contours [Fig. 3(b)] represent world-leading limits for Δm 2 41 < 0.1 eV 2 . Sensitivity to sin 2 θ 34 comes from our NC samples. For oscillations at the sterile frequency, oscillation probability ∝ cos 2 θ 34 sin 2 θ 24 , resulting in reduced sensitivity to sterile oscillations in the ND [Eq. ( 1)]. For this space our sensitivity comes primarily from oscillations at the atmospheric frequency and therefore FD data, which are statistically limited and without strong dependence on Δm 2 41 . In this space, Analysis 1 excludes slightly more parameter space than Analysis 2 across the full range considered. The differences in the contours in this space are attributed to the different statistical treatments of the two analyses. Finally, in Fig. 3(c), we present our results in terms of the effective mixing parameter, sin 2 2θ μτ ¼ 4jU μ4 j 2 jU τ4 j 2 ¼ sin 2 θ 24 sin 2 θ 34 , which can be thought of as describing anomalous sterile-driven ν τ appearance. Because the analyses are consistent and the Feldman-Cousins procedure is resource intensive, we choose to present this contour using only Analysis 1. NOvA's ND is at a higher L=E ν than other experiments with limits in this space, meaning that we are able to probe to lower values of Δm 2 41 resulting in the NOvA 90% limit being world-leading across large areas below Δm 2 41 ¼ 3 eV 2 . Notably, this limit excludes a new region of phase space around Δm 2 41 ¼ 1 eV 2 , the preferred region of Δm 2 41 for current anomalies. In conclusion, an improved search for sterile neutrino oscillations under the 3 þ 1 oscillation paradigm has been performed using NOvA data. We use two covariance matrix-based techniques that allow us to probe a wider range of Δm 2 41 values than previous NOvA analyses [43,44]. Differences between the limits for the two analyses can be taken as an uncertainty due to analysis choices such as statistical treatment, systematic treatment, binning of the Δχ 2 surface, and fitting technique. We find that the NOvA data are consistent with 3F oscillations at 90% confidence, and our limits agree with sensitivity studies performed using 3F oscillations. Our limits are the first presented in some regions of phase space, while excluding new regions of parameter space currently allowed by IceCube at 90% confidence level. This Letter additionally sets the most stringent limits for anomalous ν τ appearance for Δm 2 41 ≲ 3 eV 2 , including the strongest limits around Δm 2 41 ¼ 1 eV 2 .