Strong constraints on neutrino nonstandard interactions from TeV-scale ν µ disappearance at IceCube

,


INTRODUCTION
Neutrino oscillations are a phenomenon indicating mechanisms beyond the current Standard Model (SM) of particle physics.Experiments have measured the mixing parameters of neutrino states to excellent precision and confirm that at least two states have non-zero mass [1][2][3][4].Neutrino masses are orders of magnitude lighter than the other SM fermion masses, further suggesting the existence of beyond-Standard-Model (BSM) physics [5,6].
When the SM is treated as an effective field theory, neutrino masses can be introduced through the addition of a dimension-5 operator to the SM Lagrangian, with further BSM physics expected through the addition of dimension-6 operators required for renormalizability [7][8][9][10].One class of these dimension-6 operators introduces neutrino non-standard interactions (NSI), which are comprised of new neutral-current (NC) and charged-current (CC) neutrino interactions with charged fermions [11][12][13][14][15].
This paper presents IceCube's latest constraints on the NC NSI parameter µτ using eight years of muonneutrino-1 induced up-going track data 2 , with the highest range of event energies (500 GeV to ∼10 TeV) employed for an NSI analysis to date.A likelihood analysis is performed on the binned neutrino event counts to search for evidence of NSI via modified coherent forward scattering.The analysis uses the same sample of neutrino events and techniques as used in the recent IceCube search for sterile FIG. 1. Muon neutrino oscillogram.An example of how NSI modify predicted neutrino fluxes.Shown here is the probability ratio of NSI-modified oscillations to the SM prediction for atmospheric neutrinos (chosen value is Re( µτ ) = 0.0031, Im( µτ ) = 0.0, the 90% CL bound on positive Re( µτ )).Effects include flux disappearance at energies of 1 TeV and above for events crossing the largest Earth baselines ( cos(θ) = −1) and flux enhancement at ∼ 100 GeV.Note that the neutrino true energy range corresponds to the stated muon proxy energy range, and that the maximum disappearance for this value of µτ is ∼ 3.4%.

NEUTRINO OSCILLATIONS IN EARTH WITH NONSTANDARD INTERACTIONS
Neutrino oscillations in matter are influenced by both material density and composition [8-10, 18, 19].For SM CC coherent scattering, the potential in the flavor basis at position x is represented by [20], with H mat → −H mat for antineutrinos, and the SM matter potential , where G F is the Fermi constant and N e (x) is the electron number density [21,22].To include NSI from a mediator of an unknown energy scale, the collection of flavor-violating and flavor-conserving parameters αβ are introduced, with indices α and β corresponding to neutrino flavors e, µ, and τ .These parameters are defined through the contributions of electrons and nucleons: where N e (x) and N n (x) are the particle number densities at matter depth x for electrons and neutrons, respectively.To good approximation, this is constant through Earth, having Y ⊕ n ≈ 1.051 [2].From these generalized NSI parameters, the combined matter+NSI Hamiltonian is where * is the complex conjugate of and the diagonal parameters are real-valued.Past analyses from IceCube have set constraints on each parameter with a maximum reconstructed energy of 100 GeV [23].
In this analysis, the parameter of interest is solely µτ , as the atmospheric neutrino flux is primarily ν µ and ν µ , which for energies ≥ 20 GeV predominantly oscillate to ν τ and ν τ due to ν e decoupling [24,25].As a result, the atmospheric neutrino sample used in this analysis is most suitable for constraining µ − τ flavor-changing NSI.To verify that µτ may be constrained independently, atmospheric fluxes were simulated with each NSI parameter injected at the boundary values presented by Refs.[23,27].Nonµτ parameters, except for τ τ , were found to induce < 0.2% neutrino disappearance at all sample energies and zenith angles, whereas µτ = 0.0031 (analysis 90% CL right bound) induced ∼3.2% neutrino disappearance.While for large µτ the constraints on µτ and τ τ become correlated, strong τ τ IceCube constraints [23] imply the µτ limit generated at τ τ = 0 is accurate over the allowed parameter space.Thus, the results of this work present a standalone constraint on µτ .
The IceCube Neutrino Observatory is a neutrino detector located at the Geographic South Pole, occupying 1 km 3 of ice at depths 1450-2450 m under the Antarctic surface [28].5160 Digital Optical Modules (DOMs) [29], each consisting of a photomultiplier tube encased in a pressurized glass sphere, are distributed in a hexagonal grid along 78 60-DOM strings spaced 125 m laterally with a vertical DOM spacing of 17 m.An 8-string array of high quantum efficiency DOMs called DeepCore [30] is placed near the center of the detector at the depth where the ice is clearest.The DeepCore string lateral spacing ranges from 42 m to 72 m, with a DOM vertical spacing of 7 m.Data from the full array are used for event selection and reconstruction of relevant observables.
Cosmic-ray(CR)-induced air showers produce highenergy muons and neutrinos that comprise the majority of IceCube events.While muons produced in the Southern hemisphere ("down-going", cos(θ true µ ) > 0) often penetrate the detector volume and are a background to muon neutrino signals, muons produced in the Northern hemisphere (cos(θ true µ ) < 0) are absorbed by the Earth, eliminating the muon background to "up-going" muon neutrino signals.A CC ν µ interaction will produce hadronic products and a forward daughter muon with ∼ 50%-80% of the neutrino's energy [31].
As the muon travels it emits Cherenkov photons that are detected by IceCube DOMs, producing a track-like event that can originate either inside the detector or kilometers outside the array [32,33].From analyzing DOM charge and timing data, the zenith angle and energy of the muon are reconstructed, which determines the incident path through Earth and energy of the neutrino.This analysis uses a sample of 305,735 reconstructed muon tracks from neutrino CC interactions detected between May 13th 2011 to May 19th 2019.Events are binned uniformly both in reconstructed muon energy log(E µ reco ) (13 bins, E µ reco ∈ [500 GeV , 9976 GeV]) and cosine of the muon zenith angle (20 bins, cos( NSI signals in IceCube manifest in the form of anomalous neutrino flavour transitions in detected events compared to the SM prediction. When considering a neutrino-only flux (no antineutrinos) and positive values of Re( µτ ), there is appearance of ν µ due to modified ν τ transitions at E true ν 1TeV and −1 ≤ cos(θ) −0.8, whereas for negative Re( µτ ), it is disappearance of ν µ in the same region.This situation is reversed in the antineutrino case as well as in the inverted neutrino mass ordering (IO) case [34].IceCube cannot distinguish between neutrino and antineutrino signals, and thus the exact ν µ /ν µ sample ratio in the analysis sample is un-known3 .For equal rates of neutrinos and antineutrinos, the combined NSI effects result in NSI signals > 50% weaker than what is predicted for a pure-neutrino or antineutrino sample.An example of this is shown in Fig. 1, in which the NSI effect is largely energy-independent disappearance in the up-going direction.The inability of IceCube to discriminate between neutrinos and antineutrinos also requires an independent fit to the inverted neutrino mass ordering (IO) model, which is reported in addition to the normal ordering (NO) results (Fig. 2).

ANALYSIS
This analysis considers a complex-valued µτ parameter with oscillation probabilities calculated for neutrinos crossing the Earth using the nuSQuIDS [36,37] software package.For illustration, we briefly review the origin of the observed parameter degeneracies using an approximate treatment with small deviations present at the lowest energies, though notably these approximations are not used in the analysis but rather the full 3-neutrino mixing model including matter effects.From Ref. [24], the atmospheric neutrino oscillation probability may be approximated for E ν > 100 GeV as where θ 23 and ∆m 2 31 are standard neutrino mixing parameters [38,39], E ν is the neutrino energy, L is the matter baseline, and V d is the constant potential induced by down quarks (fermion contributions to µτ are normalized to the down quark density, with N d ≈ 3N e and N d ≈ N u in Earth [24]).Changing the mass ordering alters the sign of ∆m 2  23 , inverting the result across µτ = 0.For complex µτ there is a degeneracy in the complex plane at all energies so all contours, such as in Fig. 3, are symmetric in the imaginary dimension.Eq. 3 also contains a further degeneracy: CL boundary contours are circular in the µτ complex plane with the center of the circle approaching the origin as E ν → ∞.The final 2D contour including contributions from all energies is also found to closely resemble a circle with a slight offset from Re( µτ ) = Im( µτ ) = 0. 90% CL contours from pseudoexperiments adhered sufficiently to a circular form that accurate results could be obtained by testing hypotheses along the real axis only (201 uniformly distributed points in Re( µτ ) ∈ [−0.01, 0.01] with Im( µτ ) = 0) and extrapolating the circular contour into the complex plane.The results were verified from testing 361 uniformlydistributed hypotheses in the full complex space in addition to the aforementioned set, with Re( µτ ), Im( µτ ) ∈ [−0.01, 0.01].The likelihood threshold for 90% CL contours was evaluated using the Feldman-Cousins prescription [40] and found to be consistent with Wilks' theorem at one degree of freedom, as expected in the presence of the these degeneracies.

SYSTEMATIC UNCERTAINTIES
Systematic uncertainties are incorporated into the analysis through a collection of nuisance parameters that reweight Monte Carlo (MC) event sets through continuous parameterizations.The dominant sources of uncertainty derive from the shape and normalization of the atmospheric and astrophysical neutrino fluxes, optical properties of South Pole glacial ice, the local DOM environment, and neutrino interaction cross-sections.Other sources of systematic uncertainty were investigated and determined to be inconsequential within the overall statistical uncertainty [16,17].
The conventional4 atmospheric ν µ and ν µ flux is modeled through pion and kaon decay in the MCEq cascade equation solver [41,42] with the the SIBYLL2.3chadronic interaction model [43].The spectra of CR primaries relevant to this sample follows an approximate energy dependence of E −2. 65 .CR spectral index uncertainties are implemented via the nuisance parameter ∆γ conv [44][45][46][47].Uncertainties from meson production due to CR-atmosphere and subsequent interactions are accommodated through reweighting fluxes partitioned by incident parent energy and outgoing secondary energy, presented in Ref. [48].The atmospheric density, relevant to cascade formation, is profiled across zenith through temperature data collected by the AIRS satellite [49].The corresponding nuisance parameter, Atm.Density, is introduced through simulated air showers in randomly perturbed density profiles within the provided uncertainty ranges.Kaon energy losses via interaction with atmospheric nuclei are accounted through the total kaonnucleus cross-section uncertainty [50].Uncertainties from charged pion production and interaction are found to be negligible [16,17].Lastly, the total conventional atmo-FIG.2. Real-only result.Top: The −2∆LLH profile from the fit to data.Blue-shaded regions correspond to the the CL regions determined from the −2∆LLH values.Bottom: Comparison of the 90% CL limits from this analysis to IceCube's previous real-only µτ search [25] and the Super-Kamiokande experiment's inaugural constraints [58].
The optical properties of the bulk glacial ice result from depth-dependent impurity concentrations.To minimize the number of relevant parameters and their uncertainties, the absorption and scattering coefficients collected for each 10 m layer are reparameterized into a Fourier series up to a finite cutoff, with modes ordered from the greatest to weakest effects on the propagation of light in the glacial ice.The SnowStorm software implements an efficient method of sampling the Fourier parameter space by perturbing a single central MC set rather than generating multiple MC sets [57].Two energy-dependent basis functions are inferred from correlations between perturbed modes, and the amplitudes of these functions ultimately serve as the nuisance parameters for the bulk ice uncertainties.These nuisance parameters have a bivariate Gaussian prior.
After deployment, the water in the sensor borehole refreezes with optical impurities inhomogeneously distributed relative to the DOM axis, termed "hole ice" [59].
The consequence of hole ice is the effect on the angular sensitivity in photon detection.This contribution to the angular efficiency has been modeled empirically with two additional parameters, p 1 and p 2 .Refs. [16,17] found only one parameter (p 2 ) has a significant contribution to the uncertainty from hole ice such that variations in p 2 cover any effects seen in shifts of the negligible parameter (p 1 ).
The uncertainties associated with the effective sensitivity of DOMs to photons after deployment are characterized by the DOM efficiency nuisance parameter.Factors contributing to the efficiency include those internal to the DOM, such as the photocathode efficiency and wavelength acceptance, and factors external to the DOM, including the aforementioned hole ice and sensor cable shadow [16,17].
The neutrino cross-section determines both the rate of neutrino absorption in Earth [60,61] and of observable interactions [62,63].Uncertainties regarding neutrino interactions at the detector were found by Refs.[64,65] to be negligible while the uncertainties of the neutrino cross-sections on in-Earth absorption are parameterized through linearly scaling cross-sections σ νµ and σ νµ .The corresponding priors are fixed at the largest uncertainties found within the sample energy range [63].
The impact of the systematic uncertainties was determined by calculating the 90% CL sensitivity when selected nuisance parameters were fixed while the others were fit freely.For these tests, the "Asimov" [66] sensitivity was employed, following its validation against the true median sensitivity from 1,000 pseudoexperiments.The most illustrative test fixed categories of parameters organized into three types: hadronic5 , cosmic 6 , and detector 7 .Fixed cosmic nuisance parameters resulted in a ∼ −0.82% relative change in | µτ |, while fixed hadronic parameters have a relative change of ∼ −1.63%.Lastly, the largest uncertainty contribution is from the detector parameters, which have a ∼ −9.80% relative change from the central sensitivity radius.
For a review of the systematic uncertainties treated in this analysis, see Refs.[17] and [16].The prior and posterior widths for the nuisance parameters at the analysis best-fit are listed in the supplementary materials.

RESULTS
The analysis real-valued best-fit Re( µτ ) = −0.0029.The strongest nuisance parameter pull is the cosmic ray spectral index, with a shift of 0.066 (2.2σ), while all other systematic uncertainty best-fit values are within 1σ of their central values.Fig. 2 displays the test statistic profile for the data and the corresponding CL regions in the top panel, followed by a comparison of 90% CL limits derived from other measurements in the bottom panel.The analysis limits and sensitivities are a factor of ∼ 2 improvement beyond the leading constraints from Ref. [25].
In Fig. 3 are the CL regions (68%, 90%, and 95%) in complex µτ space.Fig. 4 compares the analysis result and sensitivity to the next-leading complex µτ limits from Ref. [23], demonstrating an improvement by a factor of ∼ 4. The result is found to be consistent with expected experimental sensitivity.The best-fit LLH is found to be -0.68 standard deviations from the distribution mean, which is consistent with no NSI at a p-value of 25.2% derived from 1000 trial pseudo-experiments.The best-fit µτ was also consistent with the recovered pseudoexperiment best-fit locations when a non-NSI hypothesis assumed.Compared to initial Re( µτ ) constraints placed by Ref. [58] and subsequent measurements8 such as Refs.[25] and [23], this analysis places the best constraints on Re( µτ ) to-date.Further, few analyses constrain complex NSI parameters, such as Ref. [23], and FIG. 4. Global comparison.Comparison of the analysis 90% CL sensitivity and result to the DeepCore 3-year, 5.6-100 GeV result [23].Green and yellow regions represent 90% CL sensitivity envelopes of symmetrically-counted 68% and 95% (respectively) regions calculated from 1,000 pseudoexperiment trials.
this analysis places the strongest constraints on Im( µτ ) to-date (Fig. 4).
To conclude, 305,735 up-going muon-neutrino tracks from 500 GeV to 9976 GeV detected by the IceCube Neutrino Observatory have been analyzed to search for evidence of µτ NSI.The best-fit point value is consistent with the no-NSI hypothesis at a p-value of 25.2%.The 90% CL limits on real-only µτ are −0.0041< µτ < 0.0031, representing the strongest constraints on any NSI parameter in any oscillation channel to date.

FIG. 3 .
FIG. 3. Complex result.Confidence level regions for complex µτ in blue-shaded regions, with the analysis 90% CL sensitivity in green and the red cross marking the data bestfit.

TABLE I .
[16,17]lysis nuisance parameters with corresponding central values and 1σ widths of the prior and posterior distributions.Choices of priors are reported in Refs.[16,17].The ±1σ posterior uncertainties are calculated from distributions obtained from 1000 pseudoexperiments simulated at the model best-fit.Listed parameter constraints are Gaussian distributed except for those with an asterisk, denoting a bivariate Gaussian distribution.