to the standard deviation of the fitted function.

The distribution of $\chi^{2}$/d.o.f. from the fits in all analysis bins is
used as a diagnostic to ensure that the fitted, linear hypersurfaces provide a
good estimate for the expected number of events for the full range of
simulated detector configurations. We find that the means of these
$\chi^{2}$/d.o.f. distributions are all consistent with 1.0 as expected from
good fits for each of the three categories described above (NC + $\nu_{e}$ CC,
$\nu_{\tau}$ CC and $\nu_{\mu}$ CC). Attempts to use higher order polynomial
fits did not yield a significantly improved $\chi^{2}$/d.o.f., and in fact
often rendered the fits less stable.

To produce the histograms for fitting the hypersurfaces, a choice must be made
for the values of flux, cross section and oscillation parameters. We found
that the hypersurface fits are sensitive to the choice of parameters that have
correlations with the effect they encode. Most notably, this effect is
observed between the mass splitting and DOM optical efficiency as demonstrated
in Fig. 20, which shows the difference between fitted hypersurface gradients
for the DOM efficiency dimension for two values of $\Delta m^{2}_{32}$.

This problem arises because we are only fitting the hypersurfaces in
reconstructed phase space, without accounting for the different true energy
and zenith distributions of MC in each analysis bin, which change with each
detector systematic variation. To mitigate this problem, we fit the
hypersurfaces for 20 different values in mass splitting between $1.5\times
10^{-3}\,\mathrm{eV}^{2}$ and $3.5\times 10^{-3}\,\mathrm{eV}^{2}$, and then
apply a piece-wise linear interpolation to all slopes, intercepts and
covariance matrix elements. The oscillation parameter fit can then dynamically
adapt the hypersurfaces for each value of $\Delta m^{2}_{32}$ that is tested
using these interpolated functions. The effects of other parameter choices
were evaluated as well, but none were found to introduce a significant bias.

In this study, a 5-dimensional hypersurface is used to parameterize the most
relevant detection uncertainties, namely the absolute optical efficiency of
the DOMs, the relative angular acceptance of the modules (two parameters, see
Fig. 4), and a global scaling of absorption and scattering lengths for the
bulk of the medium. As motivated in Section II.3, the DOM efficiency is
constrained by a Gaussian prior to the value of 1.0 $\pm$ 0.1. The ice model
parameters are unconstrained in the fit, and allowed to vary within
conservative ranges determined from calibration data. The hole ice model
parameters are bounded within the ranges $-2.0<p_{0}<1.0$ and
$-0.2<p_{1}<0.2$. The bulk ice model parameters are bounded within
$-0.90<\mathrm{Absorption}<1.10$ and $-0.95<\mathrm{Scattering}<1.15$.

We also tested the impact of a new, more detailed ice model [42], depth-
dependent uncertainties of the ice [74], a simulation including the shadow
cast by the cables onto the DOMs, a modification of the quantum efficiency as
a function of wavelength for DeepCore DOMs and variations in the noise rate.
All of these effects were found to be negligible, so they were not included in
the fit.

Figure 20: Difference between the fit hypersurface gradient for DOM optical
efficiency in each analysis bin, comparing the gradient fit assuming
$\Delta\mathrm{m}^{2}_{32}=3.0\times 10^{-3}~{}\mathrm{eV}^{2}$ minus the
gradient fit for $\Delta\mathrm{m}^{2}_{32}=2.5\times
10^{-3}~{}\mathrm{eV}^{2}$.

###  VI.2 Atmospheric neutrino fluxes

The flux of atmospheric leptons depends on the spectrum and composition of
cosmic rays, the atmospheric conditions at the interaction site and the
hadronic interaction model describing the development of the cosmic-ray
showers. The uncertainties of each of these ingredients produces uncertainties
in the lepton fluxes that need to be accounted for in the analysis. Our
studies focus on how these uncertainties impact the neutrino fluxes, since our
data sample is expected to have a 98% neutrino purity. Our starting point is
the atmospheric $\nu$ flux from Honda et al. [47]. We introduce parameters to
encode the sources of uncertainty in the flux, effectively fitting the
atmospheric neutrino spectrum as part of the study.

####  VI.2.1 Cosmic-ray flux uncertainties

Neutrinos of a given energy can be produced in showers initiated by cosmic-ray
primaries with energies up to 100 times higher [75]. This study is focused on
atmospheric $\nu$ between a few and 100 GeV, so the relevant range of the
cosmic-ray spectrum is between 10 GeV and 10 TeV. In this range, the spectrum
is almost entirely composed of protons and helium [76]. The main uncertainty
in this portion of the spectrum can be encoded by a power-law correction
$E^{\Delta\gamma}$ [73, 77], and since the neutrinos follow closely any
modification to the primary spectrum, the same form can be used for modifying
their flux. The corrected flux as used in this study is then given by

$\Phi_{\nu,\mathrm{mod}}=\Phi_{\nu}\left(\frac{E}{E_{\mathrm{pivot}}}\right)^{\Delta\gamma},$
(6)

where $\Delta\gamma$ is the parameter that can be varied. Here
$E_{\mathrm{pivot}}=24$ GeV was chosen to reduce the correlation between
$\Delta\gamma$ and the overall scale of the flux. The estimate of the
uncertainty on the spectral index is $\pm$0.05 [77], but here we increase it
to $\pm$0.1 to account for a similar response that is expected from
independent effects, such as modifications to charged hadron multiplicity in
the final state of DIS interactions [78], variations to higher-twist
parameters and valence quark corrections in GENIE [27] and changes to the
overall optical absorption of the ice.

####  VI.2.2 Hadronic interaction uncertainties

The atmospheric $\nu$ we detect come from the decay of hadrons produced over a
large region of parameter space with few measurements, interpolated by
phenomenological models [73]. The authors of the atmospheric neutrino flux we
use as a baseline have evaluated their impact, showing it as a function of
energy [79]. For this study, however, we require an uncertainty estimate as a
function of energy, direction and neutrino flavor. We therefore use the
MCEq666https://github.com/afedynitch/MCEq package [80], which computes
inclusive distributions of atmospheric leptons, to evaluate the impact of
variations on each of the ingredients that go into the calculation. In
general, we do this by introducing small variations $dB$ on a parameter $B$ of
some model that gives a lepton flux $\Phi_{l}$, compute the change in flux
$\frac{d\Phi_{l}}{dB}$, and introduce them scaled by a value $b$, as

$\Phi_{l,\mathrm{mod}}=\Phi_{l}+\left(b\cdot\frac{\mathrm{d}\Phi_{l}}{\mathrm{d}B}\right).$
(7)

The modifications introduced to MCEq follow the work of Barr et al. [73],
where the parameter space in primary energy $E_{i}$ and the fraction of energy
taken by the secondary meson $x_{lab}=E_{s}/E_{i}$ is divided into regions
with different uncertainties, chosen by the availability of fixed target
experiments data used in hadronic models. Here we adopt the same scheme and
compute the variations of the fluxes $\frac{d\Phi}{dB}$ after modifying the
expected hadronic yields of a model by a constant factor over each region of
$E_{i},x_{lab}$. These modifications were computed using the Sibyll2.3c [75]
hadronic interaction model and the GSF cosmic-ray flux [76] as a starting
point, since the software that produces the baseline model for this analysis
is not openly available [47]. The choice of cosmic-ray spectra and hadronic
interaction model used to compute the modifications $\frac{d\Phi}{dB}$ were
found to have negligible impact.

The variations of the yields of $K^{\pm}$ and $\pi^{\pm}$ were computed
individually and their impact on the flux was evaluated. As the pion ratio is
well-measured, the uncertainty on $\pi^{-}$ is defined by the uncertainty on
$\pi^{+}$ plus an uncertainty on the pion ratio, following [73]. The
uncertainty on $K^{\pm}$ production is kept uncorrelated.

Figure 21: Impact of the variations of the flux parameters used as a function
of reconstructed energy (left figure) and reconstructed zenith angle (right
figure). The parameters where changed by +1$\sigma$, one at a time. By
construction, the variations are symmetric around the nominal expectation. The
phase space that each parameter affects follows that given in [73].

The entire suite of potential variations is encoded in 17 variables, but not
all of them have a significant impact on the neutrino flux in the energy range
of this analysis. Here, only the six parameters listed under the Flux section
in Table 4 were found to be relevant. The uncertainties on the pion yields
resulting from incident parent particles with $E<30$ GeV (regions A to F in
Fig. [73]) have been added in quadrature and grouped into one parameter, as
their modifications have a degenerate effect on the atmospheric neutrino flux.
The impact of these six parameters on the flux of muon neutrinos is shown in
Fig. 21 as function of reconstructed energy and zenith angle.

We further investigated the impact of hadron-air inelastic cross section
uncertainties that drive energy losses during shower development. Using the
Glauber approach [81, 82] to propagate hadron-proton cross sections
uncertainties [83] to hadron-air interactions as a function of the center of
mass energy, we find that their impact on the atmospheric flux model in our
energy range of interest ($<1$ TeV) is below 1% and therefore negligible.

####  VI.2.3 Atmospheric density

The effect of atmospheric density uncertainty in this analysis was studied by
obtaining a variation of atmospheric density profile, perturbing the Earth’s
atmospheric temperature within a prior range given by the NASA Atmospheric
InfraRed Sounder (AIRS) satellite [84] temperature data. The resulting
atmospheric density profiles are injected into MCEq to calculate new fluxes.
This is performed for a variety of cosmic-ray models and hadronic interaction
models available in MCEq. It was found that even the largest variation from
density perturbation has a negligible effect in this energy range, so this
uncertainty was not included in the fit.

###  VI.3 Neutrino-nucleon cross-section

The DeepCore neutrino event samples can span an energy range from GeV to TeV.
For energies above 20 GeV, the interactions are dominated by deep inelastic
scattering (DIS). Below this threshold, interactions with the nucleons as a
whole become more relevant [85]. Since there is no coherent theoretical frame
that explains both regimes, we keep the same approach as event generators in
our study and divide the uncertainties from neutrino interactions into
regions.

####  VI.3.1 Deep Inelastic Scattering

The interactions of neutrinos with individual quarks can be calculated
precisely. The uncertainties on the predictions of these calculations come
from the data set used to describe the probability of finding a quark within
the nucleons being considered. These parton distribution functions (PDFs) are
mainly obtained from electron scattering data, and the predicted neutrino-
nucleon cross section differs depending on which one is used.

Our interaction generator GENIE uses an outdated PDF for the calculation,
GRV98 [86], because it was the only PDF that incorporated the extrapolations
of $Q^{2}$ below those available from measurements at the time that this study
was performed [87], but that are relevant for the energies in question777See
[88] for a new calculation that extends to low $Q^{2}$.. At high energies it
differs from other more sophisticated and newer computations. To address this
point, we have studied the different predictions from multiple DIS
calculations, comparing them to GENIE.

The comparison included GENIE (GRV98) and three additional cross section
calculations: CSMS [89], BGR [90] and a calculation using the CTEQ6 [91] PDF
set. The corresponding cross sections were obtained from the NeutrinoGenerator
(based on [92]) and GENIE-HEDIS [93] codes. We looked at the total and
differential cross section $d\sigma/dE\,dy$ of neutrinos and antineutrinos.
The cross sections seem to divide in two, with CSMS and CTEQ6 giving very
similar results, and BGR and GENIE in rough agreement. The largest differences
are observed in the CC cross section of antineutrinos.

The differences were traced back to the light sea quark contributions to the
cross section with the $s$-quark content being the most different among
calculations. Both CSMS and CTEQ6 have a significantly larger contribution
from these quarks, which results in a noticeable difference for the
antineutrino cross sections. In CSMS the overall interaction rate of
antineutrinos is higher by about 10%, with those events coming mainly from
interactions with high inelasticity, defined as $y=E_{\mathrm{had}}/E_{\nu}$.
Neutrino cross sections, on the other hand, are dominated by valence quarks at
our energies, and thus are largely unaffected by changes to the sea quarks
contribution.

Figure 22: Inclusive total neutrino-nucleon cross section for neutrinos (left)
and antineutrinos (right) on an isoscalar target (black line) from GENIE,
compared to measurements from CCFR [94], NUTEV [95] and NOMAD [96]. The GENIE
cross section with its DIS fraction fully converted to “CSMS mode” (brown
line) is also shown.

The difference between CSMS and GENIE, the largest one observed, was
parameterized as function of energy and inelasticity in a single parameter
added to the analysis that can reweight every event in our simulation to mimic
the predictions of the CSMS calculation. The DIS calculations discussed,
however, are not valid below 100 GeV because of their limited $Q^{2}$ range.
Without further information on how to extrapolate this effect to lower
energies, we decided to avoid discontinuities in the DIS cross section by
applying the same correction derived at 100 GeV for DIS events at lower
energies and thus guaranteeing consistency in our approach. We tested other
extrapolation methods, such as using the correction at 50 GeV or using the
last few energy points to come up with a linear function to go to lower
energies, and found that the analysis presented here is robust against this
choice.

In our parameterization, a value of DIS CSMS equal to zero corresponds to the
cross section in GENIE, while a value of one approximates the prediction from
CSMS. To evaluate the data we set a prior centered at zero and an uncertainty
of 1.0 to ensure the fit is stable while allowing the parameter to choose
between models. The parameter that controls these changes can be converted to
a fractional value to end up with a cross section in between the models.
Comparisons of the resulting total cross section from GENIE and the CSMS-like
modification are shown on Fig. 22, compared to data. The correction brings the
flux-averaged inelasticity from 0.40 at $E_{\nu}$=40 GeV in GENIE to 0.47 when
corrected.

Additional studies were performed to test the impact of the parameters used in
the Bodek-Yang model to extrapolate PDFs in the low $Q^{2}$ region and on the
impact of the hadron multiplicity resulting from these interactions, as
described in [27]. We also considered the impact of nuclear effects, which
have been demonstrated to have an impact on neutrino-nucleon cross sections
[97], in particular on the inelasticity of neutrino interactions with water.
We found these effects to be either negligible or fully degenerate with
existing parameters and therefore are not included explicitly in the fit.

####  VI.3.2 Non-DIS interactions

The two main non-DIS processes that contribute to events in our sample are
resonance production (RES) and charged current quasielastic scattering (CCQE).
Similar to what was done in previous studies [26, 27], two systematic
parameters are included to account for uncertainties in the nucleon form
factors for each interaction. Both these form factors have a dependency on
$Q^{2}$ of the form

$F(Q^{2})\propto\frac{1}{1-(Q^{2}/M_{A}^{2})^{2}},$ (8)

where $M_{A}$ is the axial mass, an effective parameter that can be measured
experimentally. GENIE allows one to re-weight the cross section to different
values for different interaction types, so we have introduced them as free
parameters in the fit: $M_{A}^{CCQE}$ and $M_{A}^{CCRES}$. Altering these
parameters can change the total cross section predictions for CCQE/CCRES
interactions and hence changing the expected number of CCQE/CCRES events seen
in the sample.

####  VI.3.3 $\nu_{\tau}$ charged current cross section

The charged current cross section of tau neutrinos differs from that of other
flavors because the mass of the tau lepton is comparable to the target nucleon
mass and the neutrino energies being studied. In particular, the terms of the
form $\frac{m_{\mathrm{lepton}}^{2}}{2M_{\mathrm{nucleon}}E_{\nu}}$, which are
heavily suppressed for light lepton ($e$ or $\mu$) cases at all energies,
become relevant for tau neutrinos up to tens of GeV [98]. A summary of various
inclusive cross section calculations is discussed in [99], where a proposed
parameterization that can move between calculations is also given. We tested
the impact of these changes in our study, which has a very small component of
$\nu_{\tau}$, and found them to be well below our threshold, so this source of
uncertainty was not included in the final result.

###  VI.4 Atmospheric muon contamination

The event selection was designed to eliminate atmospheric $\mu$ with high
efficiency. This is achieved for this subset of the data, with an estimated
contamination from simulation studies of $\sim 2$%. However, due to the
efficient rejection and the computational requirements of simulating this
background, few simulated events survive to the analysis binning. Therefore,
even though the atmospheric muon fluxes should be correlated with those of
neutrinos, changes to how they are modeled are mostly negligible in comparison
to the statistical uncertainty of the simulation. Most of their variation is
captured by a scaling factor, so a single parameter that globally changes
their relative contribution to the sample is used (see Atm. $\mu$ scale in
Table 4).

The few atmospheric $\mu$ events that are still present at the final level of
the selection produce histograms that are sparsely populated, with large bin-
to-bin fluctuations. To overcome this issue, we applied a variable bandwidth
Kernel Density Estimator (KDE) [100] to the binned atmospheric $\mu$ data,
which results in smooth expectation histograms for the nominal data set as
well as for variations produced with modified detector conditions. The error
of this procedure was estimated adapting the theoretical upper bound variance
computation described in [100], adjusting it so that their hypersurface fits
would yield a reduced $\chi^{2}=1$. This conservative estimation of the error
for the KDE procedure is included in the analysis as $(\sigma^{\mu}_{i})^{2}$
in Eq. 4.

###  VI.5 Oscillation parameters

Oscillation probabilities are computed using a three-flavor scheme that
includes the coherent forward scattering experienced by neutrinos as they
cross Earth’s matter. This is implemented in a custom python code [101] that
follows the description in [102]. The computation requires all oscillation
parameters to be defined (three mixing angles, two mass differences squared,
and an imaginary phase), as well as the electron density profile of the Earth.

We approximate the Earth as a collection of 12 radial layers of constant
matter density, following the Preliminary Earth Reference Model (PREM) [21].
The electron-to-nucleon fraction can change depending on the chemical
composition of the layers. Here we use the values proposed in [103], with
0.4656 for the inner and outer core, and 0.4957 for the mantle. Varying these
fractions, as well as including a finer radial description of the Earth,
results in negligible changes, so they are kept fixed in the fit to the data.

We also find that this analysis is insensitive to the values of $\theta_{12}$,
$\Delta m^{2}_{21}$ and $\theta_{13}$, within their current uncertainties.
Therefore, we fix their values to recent global fit results from [19], with
$\theta_{12}=33.82^{\circ}$, $\Delta m^{2}_{21}=0.739\times 10^{-5}$ eV2 and
$\theta_{13}=8.61^{\circ}$. These are the global fit results obtained under
the normal neutrino mass ordering hypothesis ($m_{1}<m_{2}<m_{3}$), which is
assumed throughout this analysis as well. A separate, more sensitive analysis
to determine the neutrino mass ordering using the full DeepCore data set is
currently underway. The analysis is also insensitive to the choice of a CP-
violating phase, so we fix $\delta_{CP}=0$.

##  VII Results

This analysis was performed in a “blind” manner, such that all choices
regarding the event selection criteria, analysis binning and implementation of
systematic uncertainties are made prior to fitting the real detector data in
order to avoid biasing the results. After finalizing all analysis choices and
procedures, we perform a fit to real data by minimizing the modified
$\chi^{2}$ test statistic defined in Eq. 4 over 200 bins. This yields the
best-fit values of

$\begin{split}\sin^{2}\theta_{23}&=0.51\pm 0.05,\;\mathrm{and}\\\ \Delta
m^{2}_{32}&=(2.41\pm 0.07)\times 10^{-3}\,\mathrm{eV}^{2},\end{split}$

assuming normal neutrino mass ordering (NO). The 68% C.L. for each parameter
is derived following the Feldman-Cousins prescription [104]. The best fit
nuisance parameter values are reported in Table 4, while Table 5 shows the
observed and expected number of events at the best-fit point.

Figure 23: Observed $\Delta\chi^{2}_{\mathrm{mod}}$ (solid) compared to the expectation (dashed) and the distribution of 1000 pseudo-data trials (yellow and green bands) produced at the best fit point of the analysis for the atmospheric (top) mixing angle and (bottom) mass splitting. The red lines indicate the 68% C.L. derived using the Feldman-Cousins method [104]. Table 5: Best-fit number of events with 7.5 years of livetime for each neutrino flavor and interaction type, as well as atmospheric $\mu$, along with the observed counts from the data. The rate is also given for comparison to other experiments. Type | Events | Rates [1/$10^{6}$s]
---|---|---
$\nu_{\mu}+\bar{\nu}_{\mu}$ CC | 17656 | 75.03
$\nu_{e}+\bar{\nu}_{e}$ CC | 1820 | 7.74
$\nu_{\tau}+\bar{\nu}_{\tau}$ CC | 603 | 2.56
$\nu_{all}+\bar{\nu}_{all}$ NC | 1222 | 5.19
Atmospheric $\mu$ | 711 | 3.02
Total (best-fit) | 22012 | 93.54
Observed | 21914 | 93.08

To assess the goodness-of-fit, we perform 1000 fits to pseudo-data trials that
are generated by Poisson-fluctuating the expectation for neutrinos and
atmospheric muons in the analysis binning, given the best fit values for all
parameters. Using the resulting distribution of test statistics from these
1000 trials, we find that our observed $\chi^{2}_{\mathrm{mod}}$ has a p-value
of 26.1%, indicating good agreement between simulation and data. Figure 23
shows the expected and observed $\Delta\chi^{2}_{\mathrm{mod}}$ in
$\sin^{2}\theta_{23}$ and in $\Delta m^{2}_{32}$, overlaid with the
distribution of 1000 pseudo-data trials. The observed contours are fully
contained within the $1\sigma$ fluctuations of the trials.

Figure 24: Distributions of reconstructed (top) energy, (middle) cosine zenith
and (bottom) PID score for data compared to the best-fit simulation.
Background includes atmospheric $\mu$ and all neutrino types besides
$\nu_{\mu}+\bar{\nu}_{\mu}$ CC events for which this analysis was optimized.

The distributions of reconstructed neutrino energy, reconstructed zenith angle
and the PID score for data compared with the best-fit MC simulation are shown
in Fig. 24. These 1D projections show the data binned more finely than what
was used in the fit by a factor of 2 for reconstructed energy and zenith, and
a factor of 10 for the PID. Figure 25 shows the reconstructed energy and
cosine zenith projected into L/E space for all events. The best-fit
expectation shows good agreement with the data for all observables.

We find a slight excess of events compared to MC in the very highest energy
bin, which is nevertheless consistent within statistical fluctuations.
However, this analysis extends to higher energies than previous oscillation
analyses using DeepCore data in order to better constrain systematic
uncertainties in this off-signal region. Therefore as a cross-check, we fit
the data again with this last, highest energy bin removed and observe a
negligible shift of $\sim$1% in $\Delta m^{2}_{32}$ and $\sim$0.2∘ in
$\theta_{23}$. We further performed fits to data from each season of data-
taking independently, and find that the resulting best fit atmospheric
oscillation parameters and nuisance parameters are all statistically
compatible across each year.

Figure 25: The L/E distribution for the best-fit expectations overlaid with
the observed data. Background includes atmospheric $\mu$ and all neutrino
types besides $\nu_{\mu}+\bar{\nu}_{\mu}$ CC events. The expectation at the
best fit but without oscillations is shown as a dashed line.

Additional cross-checks are performed with simulation to assess the robustness
of the result against perturbations to the detector model that are not
parameterized by the hypersurface treatment described in Section VI.1. These
include a modification to the wavelength dependence of the HQE PMTs to match
alternative laboratory measurements, an implementation of the cable that
shadows part of the photocathode for particular azimuthal angles [105], and
several different bulk ice models that were developed after this analysis was
finalized, which were found to fit equally-well or better to the LED
calibration data.

The simulation generated for each of these perturbations is processed through
the standard event selection and used to generate pseudo-data that is weighted
using the best-fit values from the original fit to data, and fit with the
standard analysis procedure. No significant bias in the fitted oscillation
parameters is observed for any of these perturbations. We note that the
largest shift observed in the best fit point with a significance of
0.3$\sigma$ was found when using simulation produced with the birefringent ice
properties incorporated in the bulk ice model [106]. While insignificant for
this analysis, this points to the potential need for an improved treatment of
systematic uncertainties related to the glacial ice optical properties for
analyses with higher statistics and more scattered photons used in the
reconstruction of neutrino properties.

The hypersurface treatment was also validated by generating simulation using
the best-fit nuisance parameter values for detector systematics. There were no
significant chnages obtained in our results after including this set in the
analysis.

Systematic group |  $\delta(\Delta m^{2}_{32})$ [%] |  $\delta(\sin^{2}\theta_{23})$ [%]
---|---|---
Detector |  -33.6 |  -10.6
Flux |  -5.4 |  -1.4
Cross section |  -6.8 |  -0.3
$A_{eff}$ scale |  -1.0 |  -0.4
Atm. $\mu$ scale |  -1.8 |  -1.1
Table 6: Relative change in $1\sigma$ uncertainty assuming perfect knowledge
of each group of systematic uncertainties. Relative change is calculated from
the width of the 68% C.L. assuming Wilks’ theorem. Figure 26: Contours showing
the 90% C.L. allowed region for $\Delta m^{2}_{32}$ and $\sin^{2}\theta_{23}$
from this study (blue) compared to previous IceCube DeepCore results[26, 27].
All confidence intervals shown are derived assuming Wilks’ theorem for a
consistent comparison.  Figure 27: Contours showing the 90% C.L. allowed
region for the atmospheric neutrino oscillation parameters from this study
(blue) compared to results from MINOS [107], NOvA [108], Super-Kamiokande
[109] and T2K [110]. Daya Bay also measures $\Delta m^{2}_{32}$ in conjunction
to $\theta_{13}$, but the results cannot be displayed in the format above [6].
The DeepCore confidence interval is derived assuming Wilks’ theorem.

Table 6 shows the relative contributions of each group of systematic
uncertainties that were considered in this study to the total error. To
determine the contributions, we assume perfect knowledge of each group of
systematic uncertainties in the fit and calculate the relative change in the
width of the 68% CL for each oscillation parameter. The detector systematic
uncertainties contribute most to the total uncertainty, with flux and cross-
section systematics contributing far less to the error budget. The atmospheric
$\mu$ normalization term has little effect thanks to the very small
contamination in the final sample, while the overall normalization term for
neutrinos has almost no effect on the error of the measurement. All
sensitivity to the oscillation parameters comes from the shape of the
oscillation pattern, with the total flux being relatively unimportant for this
measurement. This breakdown of uncertainties also demonstrates that the
measurement can still be largely improved by increased statistics.

The 90% C.L. allowed region for the atmospheric oscillation parameters is
shown in Fig. 26 compared to previous measurements using IceCube DeepCore. All
contours are derived assuming Wilks’ theorem, which leads to slight over
coverage due to the physical boundary of $\theta_{23}$ (see Fig. 23). Taking
into account the previously published Feldman Cousins-corrected 1$\sigma$
errors [26], we observe an improvement of 44% and 37% in the measurements of
$\Delta m^{2}_{32}$ and $\sin^{2}\theta_{23}$, respectively. These new results
therefore represent the most precise measurement of these parameters using
atmospheric neutrinos to date.

Figure 27 shows the new IceCube DeepCore result in comparison to measurements
performed by other experiments, using both accelerator and atmospheric
neutrinos. MINOS [107], T2K [110] and NOvA [108] measure these parameters
using neutrinos produced in particle accelerator facilities with energies
between hundreds of MeV to and a few GeV. Super-Kamiokande [109] uses
atmospheric neutrinos, but the bulk of their statistics are around the 1 GeV
region. With IceCube DeepCore we use neutrinos with energies higher than any
of these experiments, interacting mainly via deep inelastic scattering, and
are therefore subject to different interaction uncertainties. While our
neutrino source is also the atmosphere and we use the same nominal flux
calculation as Super-Kamiokande, we see a different region of the atmospheric
neutrino spectrum and we include several flux-related nuisance parameters in
the fit to adjust for discrepancies with data. Given the differences on how
these measurements are obtained, the overlap between the results is
noteworthy, but difficult to rigorously quantify using individually reported
uncertainties without resimulations accounting for both correlated and
uncorrelated uncertainties. This could be followed up with future studies
using external data releases from each experiment.

##  VIII Conclusion

We have presented the most precise measurement of oscillations of atmospheric
neutrinos to date, using a newly calibrated and filtered data sample from
IceCube DeepCore. The measurement was made possible thanks to state-of-the-art
calibrations, improved event filtering, and significant improvements to our
methods for evaluating sources of uncertainty.

The recent calibration efforts allow us to successfully describe our data with
high precision, especially the response of sensors to single photons and
better description of the optical properties of the ice. The data used spans a
period of 8 years, over which we find stable noise rates and detector
response.

The new event selection suppresses the main sources of background by a factor
greater than 104, while keeping over 20% of all neutrinos that interact in the
detector volume, resulting in a comparable rate of atmospheric $\mu$ and $\nu$
at Level 5, which is a common starting point for future DeepCore analyses. We
further developed an analysis-specific selection focused on events detected
with minimally scattered light, which reduced the background further to
achieve a neutrino purity of 98%, with a high fraction of $\nu_{\mu}$ CC
events.

We investigated many possible sources of uncertainty, and developed improved
methods for implementing them in the analysis. The detector related
uncertainties, known to have the largest impact on our results, were studied
using significantly more simulation sets for known effects than in previous
studies, and were parameterized in a way that accounts for correlations with
other parameters. The uncertainties due to the atmospheric $\nu$ flux were
significantly expanded so that possible variations to the flux were corrected
for as part of the fit. Neutrino-nucleon cross sections were also evaluated
with in different frameworks, where the free parameters have a more physical
interpretation compared to previous analyses. The impact of atmospheric $\mu$,
which has remained a challenge throughout these studies, was reduced by
limiting their fraction in the final sample, making any changes in their
prediction negligible.

The results presented here are consistent with measurements performed using
human-made accelerator neutrino experiments [107, 110, 108]. At the same time,
these results are obtained using much higher energy neutrinos, and are
therefore insensitive to $\delta_{CP}$, as well as many of the cross-section
uncertainties that are critical to understand for accelerator neutrino
experiments. Our results are of similar precision to accelerator measurements
thanks to the enormous flux of neutrinos that are provided by cosmic ray
interactions, and also the ability to constrain the oscillation valley across
several distinct bins in L/E. Our measurements therefore provide an important,
complementary probe of the parameters $\theta_{23}$ and $\Delta m^{2}_{32}$
within the context of the global neutrino oscillation landscape.

Further improvements to our oscillation measurements are currently underway on
three fronts. First, future analyses of DeepCore data will start from the
common event sample described herein, and employ more complex and resource-
intensive reconstruction strategies in order to retain more signal events,
enhance the neutrino flavor identification, and improve the directional and
energy resolutions for cascade-like events in particular [29]. Second, we
continue to improve the accuracy of our simulation by testing and including
higher order effects, such as direct simulation of the hole ice and the shadow
that the cable casts on the modules, which are becoming more relevant as our
statistical precision increases. Finally, the IceCube Upgrade detector is
planned for installation in 2025/26. This detector will serve as an in-fill to
the existing IceCube DeepCore array, and significantly enhance our ability to
detect and reconstruct GeV-scale neutrino interactions [111]. With these
improvements we will continue to provide high precision measurements of
neutrino oscillations at the highest energies and over the longest baselines.

###### Acknowledgements.

The authors gratefully acknowledge the support from the following agencies and
institutions: USA – U.S. National Science Foundation-Office of Polar Programs,
U.S. National Science Foundation-Physics Division, U.S. National Science
Foundation-EPSCoR, Wisconsin Alumni Research Foundation, Center for High
Throughput Computing (CHTC) at the University of Wisconsin–Madison, Open
Science Grid (OSG), Advanced Cyberinfrastructure Coordination Ecosystem:
Services & Support (ACCESS), Frontera computing project at the Texas Advanced
Computing Center, U.S. Department of Energy-National Energy Research
Scientific Computing Center, Particle astrophysics research computing center
at the University of Maryland, Institute for Cyber-Enabled Research at
Michigan State University, and Astroparticle physics computational facility at
Marquette University; Belgium – Funds for Scientific Research (FRS-FNRS and
FWO), FWO Odysseus and Big Science programmes, and Belgian Federal Science
Policy Office (Belspo); Germany – Bundesministerium für Bildung und Forschung
(BMBF), Deutsche Forschungsgemeinschaft (DFG), Helmholtz Alliance for
Astroparticle Physics (HAP), Initiative and Networking Fund of the Helmholtz
Association, Deutsches Elektronen Synchrotron (DESY), and High Performance
Computing cluster of the RWTH Aachen; Sweden – Swedish Research Council,
Swedish Polar Research Secretariat, Swedish National Infrastructure for
Computing (SNIC), and Knut and Alice Wallenberg Foundation; European Union –
EGI Advanced Computing for research; Australia – Australian Research Council;
Canada – Natural Sciences and Engineering Research Council of Canada, Calcul
Québec, Compute Ontario, Canada Foundation for Innovation, WestGrid, and
Compute Canada; Denmark – Villum Fonden, Carlsberg Foundation, and European
Commission; New Zealand – Marsden Fund; Japan – Japan Society for Promotion of
Science (JSPS) and Institute for Global Prominent Research (IGPR) of Chiba
University; Korea – National Research Foundation of Korea (NRF); Switzerland –
Swiss National Science Foundation (SNSF); United Kingdom – Department of
Physics, University of Oxford.

## References

  * Fukuda  _et al._ [1998] Y. Fukuda _et al._ (Super-Kamiokande), Evidence for oscillation of atmospheric neutrinos, Phys. Rev. Lett. 81, 1562 (1998), arXiv:hep-ex/9807003 .
  * Ahmad  _et al._ [2002] Q. R. Ahmad _et al._ (SNO), Direct evidence for neutrino flavor transformation from neutral current interactions in the Sudbury Neutrino Observatory, Phys. Rev. Lett. 89, 011301 (2002), arXiv:nucl-ex/0204008 .
  * Eguchi  _et al._ [2003] K. Eguchi _et al._ (KamLAND), First results from KamLAND: Evidence for reactor anti-neutrino disappearance, Phys. Rev. Lett. 90, 021802 (2003), arXiv:hep-ex/0212021 .
  * Anselmann  _et al._ [1994] P. Anselmann _et al._ (GALLEX), GALLEX results from the first 30 solar neutrino runs, Phys. Lett. B 327, 377 (1994).
  * Abdurashitov  _et al._ [1994] J. N. Abdurashitov _et al._ (SAGE), Results from SAGE, Phys. Lett. B 328, 234 (1994).
  * An  _et al._ [2012] F. P. An _et al._ (Daya Bay), Observation of electron-antineutrino disappearance at Daya Bay, Phys. Rev. Lett. 108, 171803 (2012), arXiv:1203.1669 [hep-ex] .
  * Pontecorvo [1957] B. Pontecorvo, Mesonium and anti-mesonium, Sov. Phys. JETP 6, 429 (1957), [Zh. Eksp. Teor. Fiz.33,549(1957)].
  * Maki  _et al._ [1962] Z. Maki, M. Nakagawa, and S. Sakata, Remarks on the unified model of elementary particles, Prog. Theor. Phys. 28, 870 (1962).
  * Capozzi  _et al._ [2021] F. Capozzi, E. Di Valentino, E. Lisi, A. Marrone, A. Melchiorri, and A. Palazzo, Unfinished fabric of the three neutrino paradigm, Phys. Rev. D 104, 083031 (2021), arXiv:2107.00532 [hep-ph] .
  * de Salas  _et al._ [2021] P. F. de Salas, D. V. Forero, S. Gariazzo, P. Martínez-Miravé, O. Mena, C. A. Ternes, M. Tórtola, and J. W. F. Valle, 2020 global reassessment of the neutrino oscillation picture, JHEP 02, 071, arXiv:2006.11237 [hep-ph] .
  * Esteban  _et al._ [2020] I. Esteban, M. C. Gonzalez-Garcia, M. Maltoni, T. Schwetz, and A. Zhou, The fate of hints: updated global analysis of three-flavor neutrino oscillations, JHEP 09, 178, arXiv:2007.14792 [hep-ph] .
  * Cabibbo [1963] N. Cabibbo, Unitary Symmetry and Leptonic Decays, Phys. Rev. Lett. 10, 531 (1963).
  * Kobayashi and Maskawa [1973] M. Kobayashi and T. Maskawa, CP Violation in the Renormalizable Theory of Weak Interaction, Prog. Theor. Phys. 49, 652 (1973).
  * Feruglio [2015] F. Feruglio, Pieces of the Flavour Puzzle, Eur. Phys. J. C 75, 373 (2015), arXiv:1503.04071 [hep-ph] .
  * Athanassopoulos  _et al._ [1996] C. Athanassopoulos _et al._ (LSND), Evidence for anti-muon-neutrino —$>$ anti-electron-neutrino oscillations from the LSND experiment at LAMPF, Phys. Rev. Lett. 77, 3082 (1996), arXiv:nucl-ex/9605003 .
  * Aguilar-Arevalo  _et al._ [2018] A. A. Aguilar-Arevalo _et al._ (MiniBooNE), Significant Excess of ElectronLike Events in the MiniBooNE Short-Baseline Neutrino Experiment, Phys. Rev. Lett. 121, 221801 (2018), arXiv:1805.12028 [hep-ex] .
  * Kaether  _et al._ [2010] F. Kaether, W. Hampel, G. Heusser, J. Kiko, and T. Kirsten, Reanalysis of the GALLEX solar neutrino flux and source experiments, Phys. Lett. B 685, 47 (2010), arXiv:1001.2731 [hep-ex] .
  * Abdurashitov  _et al._ [2009] J. N. Abdurashitov _et al._ (SAGE), Measurement of the solar neutrino capture rate with gallium metal. III: Results for the 2002–2007 data-taking period, Phys. Rev. C 80, 015807 (2009), arXiv:0901.2200 [nucl-ex] .
  * Esteban  _et al._ [2019] I. Esteban, M. C. Gonzalez-Garcia, A. Hernandez-Cabezudo, M. Maltoni, and T. Schwetz, JHEP 2019 (1), arXiv:1811.05487 [hep-ph] .
  * Esteban  _et al._ [2018] I. Esteban, M. C. Gonzalez-Garcia, A. Hernandez-Cabezudo, M. Maltoni, and T. Schwetz, Nufit 4.0, www.nu-fit.org (2018).
  * Dziewonski and Anderson [1981] A. M. Dziewonski and D. L. Anderson, Preliminary reference Earth model, Phys. Earth Planet. Inter. 25, 297 (1981).
  * Akhmedov [1999] E. K. Akhmedov, Parametric resonance of neutrino oscillations and passage of solar and atmospheric neutrinos through the earth, Nucl. Phys. B 538, 25 (1999), arXiv:hep-ph/9805272 .
  * Wolfenstein [1978] L. Wolfenstein, Neutrino Oscillations in Matter, Phys. Rev. D 17, 2369 (1978).
  * Mikheyev and Smirnov [1985] S. P. Mikheyev and A. Y. Smirnov, Resonance Amplification of Oscillations in Matter and Spectroscopy of Solar Neutrinos, Sov. J. Nucl. Phys. 42, 913 (1985).
  * Aartsen  _et al._ [2015] M. G. Aartsen _et al._ (IceCube), Determining neutrino oscillation parameters from atmospheric muon neutrino disappearance with three years of IceCube DeepCore data, Phys. Rev. D 91, 072004 (2015), arXiv:1410.7227 [hep-ex] .
  * Aartsen  _et al._ [2018] M. G. Aartsen _et al._ (IceCube), Measurement of Atmospheric Neutrino Oscillations at 6–56 GeV with IceCube DeepCore, Phys. Rev. Lett. 120, 071801 (2018), arXiv:1707.07081 [hep-ex] .
  * Aartsen  _et al._ [2019a] M. G. Aartsen _et al._ (IceCube), Measurement of atmospheric tau neutrino appearance with icecube deepcore, Phys. Rev. D 99, 032007 (2019a), arXiv:1901.05366 [hep-ex] .
  * Aartsen  _et al._ [2020a] M. G. Aartsen _et al._ (IceCube), In-situ calibration of the single-photoelectron charge response of the IceCube photomultiplier tubes, JINST 15 (06), P06032, arXiv:2002.00997 [physics.ins-det] .
  * Abbasi  _et al._ [2022a] R. Abbasi _et al._ (IceCube), Low energy event reconstruction in IceCube DeepCore, Eur. Phys. J. C 82, 807 (2022a), arXiv:2203.02303 [hep-ex] .
  * Aartsen  _et al._ [2017a] M. G. Aartsen _et al._ (IceCube), The IceCube Neutrino Observatory: Instrumentation and Online Systems, JINST 12 (03), P03012, arXiv:1612.05093 [astro-ph.IM] .
  * Benson  _et al._ [2014] T. Benson, J. Cherwinka, M. Duvernois, A. Elcheikh, F. Feyzi, L. Greenler, J. Haugen, A. Karle, M. Mulligan, R. Paulos, _et al._ , Icecube enhanced hot water drill functional description, Ann. of Glaciology 55, 105–114 (2014).
  * Abbasi  _et al._ [2012a] R. Abbasi _et al._ (IceCube), The Design and Performance of IceCube DeepCore, Astropart. Phys. 35, 615 (2012a), arXiv:1109.6096 [astro-ph.IM] .
  * Abbasi  _et al._ [2010a] R. Abbasi _et al._ (IceCube), Calibration and Characterization of the IceCube Photomultiplier Tube, Nucl. Instrum. Meth. A 618, 139 (2010a), arXiv:1002.2442 [astro-ph.IM] .
  * Abbasi  _et al._ [2009] R. Abbasi _et al._ (IceCube), The IceCube Data Acquisition System: Signal Capture, Digitization, and Timestamping, Nucl. Instrum. Meth. A 601, 294 (2009), arXiv:0810.4930 [physics.ins-det] .
  * Kelley [2014] J. L. Kelley (IceCube), Event triggering in the IceCube data acquisition system, AIP Conf. Proc. 1630, 154 (2014).
  * Kowalski [2004] M. P. Kowalski,  _Search for neutrino-induced cascades with the AMANDA-II detector_ , Ph.D. thesis, Humboldt-Universität zu Berlin, Mathematisch-Naturwissenschaftliche Fakultät I (2004).
  * Feintzeig [2014] J. Feintzeig,  _Searches for point-like sources of astrophysical neutrinos with the IceCube Neutrino Observatory_ , Ph.D. thesis, Wisconsin U., Madison (2014).
  * Kulacz [2019] N. Kulacz,  _In Situ Measurement of the IceCube DOM Efficiency Factor Using Atmospheric Minimum Ionizing Muons_ , Master’s thesis, University of Alberta (2019).
  * Aartsen  _et al._ [2013] M. G. Aartsen _et al._ (IceCube), Measurement of South Pole ice transparency with the IceCube LED calibration system, Nucl. Instrum. Meth. A 711, 73 (2013), arXiv:1301.5361 [astro-ph.IM] .
  * Chirkin [2013a] D. Chirkin (IceCube), Evidence of optical anisotropy of the South Pole ice, in  _Proceedings, 33rd International Cosmic Ray Conference (ICRC2013): Rio de Janeiro, Brazil, July 2-9, 2013_ (2013) p. 0580.
  * Rongen  _et al._ [2020] M. Rongen, R. C. Bay, and S. Blot, Observation of an optical anisotropy in the deep glacial ice at the geographic south pole using a laser dust logger, The Cryosphere 14, 2537 (2020).
  * Abbasi  _et al._ [2022b] R. Abbasi _et al._ (IceCube), In-situ estimation of ice crystal properties at the South Pole using LED calibration data from the IceCube Neutrino Observatory, The Cryosphere Discussions 2022, 1 (2022b).
  * Chirkin [2013b] D. Chirkin (IceCube), Photon tracking with GPUs in IceCube, Nucl. Instrum. Meth. A 725, 141 (2013b).
  * Rongen, M. [2016] Rongen, M. (IceCube), Measuring the optical properties of IceCube drill holes, EPJ Web of Conf. 116, 06011 (2016).
  * Pearson [1901] K. Pearson, Liii. on lines and planes of closest fit to systems of points in space, The London, Edinburgh, and Dublin Philosophical Magazine and Journal of Science. 2, 559 (1901).
  * Andreopoulos  _et al._ [2015] C. Andreopoulos, C. Barry, S. Dytman, H. Gallagher, T. Golan, R. Hatcher, G. Perdue, and J. Yarba, The GENIE Neutrino Monte Carlo Generator: Physics and User Manual, (2015), arXiv:1510.05494 [hep-ph] .
  * Honda  _et al._ [2015] M. Honda, M. Sajjad Athar, T. Kajita, K. Kasahara, and S. Midorikawa, Atmospheric neutrino flux calculation using the NRLMSISE-00 atmospheric model, Phys. Rev. D 92, 023004 (2015), arXiv:1502.03916 [astro-ph.HE] .
  * Heck  _et al._ [1998] D. Heck, J. Knapp, J. N. Capdevielle, G. Schatz, and T. Thouw,  _CORSIKA: A Monte Carlo code to simulate extensive air showers_, Tech. Rep. FZKA-6019 (1998).
  * Becherini  _et al._ [2006] Y. Becherini, A. Margiotta, M. Sioli, and M. Spurio, A Parameterisation of single and multiple muons in the deep water or ice, Astropart. Phys. 25, 1 (2006), arXiv:hep-ph/0507228 .
  * Gaisser [2012] T. K. Gaisser, Spectrum of cosmic-ray nucleons, kaon production, and the atmospheric muon charge ratio, Astropart. Phys. 35, 801 (2012), arXiv:1111.6675 [astro-ph.HE] .
  * Ahn  _et al._ [2009] E.-J. Ahn, R. Engel, T. K. Gaisser, P. Lipari, and T. Stanev, Cosmic ray interaction event generator SIBYLL 2.1, Phys. Rev. D 80 (2009), arXiv:0906.4113 [hep-ph] .
  * Leonard DeHolton [2022] K. Leonard DeHolton,  _Measuring the Atmospheric Neutrino Oscillation Parameters with IceCube DeepCore_ , Ph.D. thesis, Wisconsin U., Madison (2022).
  * Koehne  _et al._ [2013] J. H. Koehne, K. Frantzen, M. Schmitz, T. Fuchs, W. Rhode, D. Chirkin, and J. Becker Tjus, PROPOSAL: A tool for propagation of charged leptons, Comput. Phys. Commun. 184, 2070 (2013).
  * Agostinelli  _et al._ [2003] S. Agostinelli _et al._ , Geant4—a simulation toolkit, Nucl. Instrum. Meth. A 506, 250 (2003).
  * Rädel and Wiebusch [2012] L. Rädel and C. Wiebusch, Calculation of the cherenkov light yield from low energetic secondary particles accompanying high-energy muons in ice and water with geant4 simulations, Astropart. Phys. 38, 53 (2012).
  * [56] C. Kopper _et al._ , https://github.com/claudiok/clsim.
  * Abbasi  _et al._ [2010b] R. Abbasi _et al._ (IceCube), Calibration and characterization of the IceCube photomultiplier tube, Nucl. Instrum. Meth. A 618, 139 (2010b).
  * Abbasi  _et al._ [2012b] R. Abbasi _et al._ (IceCube), The Design and Performance of IceCube DeepCore, Astropart. Phys. 35, 615 (2012b), arXiv:1109.6096 [astro-ph.IM] .
  * Ke  _et al._ [2017] G. Ke, Q. Meng, T. Finley, T. Wang, W. Chen, W. Ma, Q. Ye, and T.-Y. Liu, Lightgbm: A highly efficient gradient boosting decision tree, Adv. in neural inf. and proc. syst. 30, 3146 (2017).
  * Aartsen  _et al._ [2014] M. G. Aartsen _et al._ (IceCube), Improvement in Fast Particle Track Reconstruction with Robust Statistics, Nucl. Instrum. Meth. A 736, 143 (2014), arXiv:1308.5501 [astro-ph.IM] .
  * Ahrens  _et al._ [2004] J. Ahrens _et al._ , Muon track reconstruction and data selection techniques in AMANDA, Nucl. Instrum. Meth. A 524, 169 (2004).
  * Aartsen  _et al._ [2017b] M. G. Aartsen _et al._ (IceCube), Search for sterile neutrino mixing using three years of IceCube DeepCore data, Phys. Rev. D 95 (2017b).
  * Yáñez [2014] J. P. Yáñez,  _Measurement of neutrino oscillations in atmospheric neutrinos with the IceCube DeepCore detector_ , Ph.D. thesis, Humboldt-Universität zu Berlin, Mathematisch-Naturwissenschaftliche Fakultät I (2014).
  * Aguilar  _et al._ [2011] J. A. Aguilar _et al._ (ANTARES), A fast algorithm for muon track reconstruction and its application to the ANTARES neutrino telescope, Astropart. Phys. 34, 652 (2011), arXiv:1105.4116 [astro-ph.IM] .
  * Tanabashi  _et al._ [2018] M. Tanabashi _et al._ (Particle Data Group), Review of particle physics, Phys. Rev. D 98, 030001 (2018).
  * Terliuk [2018] A. Terliuk,  _Measurement of atmospheric neutrino oscillations and search for sterile neutrino mixing with IceCube DeepCore_ , Ph.D. thesis, Humboldt-Universität zu Berlin, Mathematisch-Naturwissenschaftliche Fakultät I (2018).
  * Euler [2014] S. Euler,  _Observation of oscillations of atmospheric neutrinos with the IceCube Neutrino Observatory_ , Ph.D. thesis, Rheinisch-Westflischen Technischen Hochschule, Aachen, Germany (2014).
  * Whitehorn  _et al._ [2013] N. Whitehorn, J. van Santen, and S. Lafebre, Penalized splines for smooth representation of high-dimensional monte carlo datasets, Comp. Phys. Comm Physics Commun. 184, 2214 (2013).
  * Groom  _et al._ [2001] D. E. Groom, N. V. Mokhov, and S. I. Striganov, Atomic Data and Nuclear Data Tables 78, 183 (2001).
  * Friedman [2001] J. H. Friedman, Greedy function approximation: A gradient boosting machine., Ann. Statist. 29, 1189 (2001).
  * Pedregosa  _et al._ [2011] F. Pedregosa _et al._ , Scikit-learn: Machine learning in Python, JMLR 12, 2825 (2011), arXiv:1201.0490 [cs.LG] .
  * Kolmogorov [1933] A. Kolmogorov, Sulla determinazione empirica di una legge di distribuzione, Inst. Ital. Attuari, Giorn. 4, 83 (1933).
  * Barr  _et al._ [2006] G. D. Barr, T. K. Gaisser, S. Robbins, and T. Stanev, Uncertainties in Atmospheric Neutrino Fluxes, Phys. Rev. D 74, 094009 (2006), arXiv:astro-ph/0611266 .
  * Aartsen  _et al._ [2019b] M. G. Aartsen _et al._ (IceCube), Efficient propagation of systematic uncertainties from calibration to analysis with the SnowStorm method in IceCube, JCAP 10, 048, arXiv:1909.01530 [hep-ex] .
  * Riehn  _et al._ [2020] F. Riehn, R. Engel, A. Fedynitch, T. K. Gaisser, and T. Stanev, Hadronic interaction model Sibyll 2.3d and extensive air showers, Phys. Rev. D 102, 063002 (2020), arXiv:1912.03300 [hep-ph] .
  * Dembinski  _et al._ [2018] H. P. Dembinski, R. Engel, A. Fedynitch, T. Gaisser, F. Riehn, and T. Stanev, Data-driven model of the cosmic-ray flux and mass composition from 10 GeV to $10^{11}$ GeV, PoS ICRC2017, 533 (2018), arXiv:1711.11432 [astro-ph.HE] .
  * Evans  _et al._ [2017] J. Evans, D. G. Gamez, S. D. Porzio, S. Söldner-Rembold, and S. Wren, Uncertainties in Atmospheric Muon-Neutrino Fluxes Arising from Cosmic-Ray Primaries, Phys. Rev. D 95, 023012 (2017), arXiv:1612.03219 [astro-ph.HE] .
  * Katori and Mandalia [2015] T. Katori and S. Mandalia, PYTHIA hadronization process tuning in the GENIE neutrino interaction generator, J. Phys. G 42, 115004 (2015), arXiv:1412.4301 [hep-ex] .
  * Honda  _et al._ [2007] M. Honda, T. Kajita, K. Kasahara, S. Midorikawa, and T. Sanuki, Calculation of atmospheric neutrino flux using the interaction model calibrated with atmospheric muon data, Phys. Rev. D 75, 043006 (2007), arXiv:astro-ph/0611418 .
  * Fedynitch  _et al._ [2019] A. Fedynitch, F. Riehn, R. Engel, T. K. Gaisser, and T. Stanev, Hadronic interaction model sibyll 2.3c and inclusive lepton fluxes, Phys. Rev. D 100, 103018 (2019), arXiv:1806.04140 [hep-ph] .
  * Glauber [1955] R. J. Glauber, Cross-sections in deuterium at high-energies, Phys. Rev. 100, 242 (1955).
  * Glauber and Matthiae [1970] R. J. Glauber and G. Matthiae, High-energy scattering of protons by nuclei, Nucl. Phys. B 21, 135 (1970).
  * Halzen  _et al._ [2012] F. Halzen, K. Igi, M. Ishida, and C. S. Kim, Total Hadronic Cross Sections and $\pi^{\mp}\pi^{+}$ Scattering, Phys. Rev. D 85, 074020 (2012), arXiv:1110.1479 [hep-ph] .
  * [84] J. P. Laboratory, https://airs.jpl.nasa.gov/data/get-data/standard-data/.
  * Formaggio and Zeller [2012] J. A. Formaggio and G. P. Zeller, From eV to EeV: Neutrino Cross Sections Across Energy Scales, Rev. Mod. Phys. 84, 1307 (2012), arXiv:1305.7513 [hep-ex] .
  * Gluck  _et al._ [1998] M. Gluck, E. Reya, and A. Vogt, Dynamical parton distributions revisited, Eur. Phys. J. C 5, 461 (1998), arXiv:hep-ph/9806404 .
  * Bodek and Yang [2003] A. Bodek and U. K. Yang, Higher twist, xi(omega) scaling, and effective LO PDFs for lepton scattering in the few GeV region, J. Phys. G 29, 1899 (2003), arXiv:hep-ex/0210024 .
  * Candido  _et al._ [2023] A. Candido, A. Garcia, G. Magni, T. Rabemananjara, J. Rojo, and R. Stegeman, Neutrino Structure Functions from GeV to EeV Energies, (2023), arXiv:2302.08527 [hep-ph] .
  * Cooper-Sarkar  _et al._ [2011] A. Cooper-Sarkar, P. Mertsch, and S. Sarkar, The high energy neutrino cross-section in the Standard Model and its uncertainty, JHEP 08, 042, arXiv:1106.3723 [hep-ph] .
  * Bertone  _et al._ [2019] V. Bertone, R. Gauld, and J. Rojo, Neutrino Telescopes as QCD Microscopes, JHEP 01, 217, arXiv:1808.02034 [hep-ph] .
  * Pumplin  _et al._ [2002] J. Pumplin, D. Stump, J. Huston, H. Lai, P. M. Nadolsky, and W. Tung, New generation of parton distributions with uncertainties from global QCD analysis, JHEP 07, 012, arXiv:hep-ph/0201195 .
  * Gazizov and Kowalski [2005] A. Gazizov and M. P. Kowalski, ANIS: High energy neutrino generator for neutrino telescopes, Comput. Phys. Commun. 172, 203 (2005), arXiv:astro-ph/0406439 .
  * Garcia and Heijboer [2020] A. Garcia and A. Heijboer (KM3NeT), High-energy neutrino event simulation at NLO in Genie for KM3NeT and other observatories, PoS ICRC2019, 895 (2020), arXiv:1908.10077 [hep-ex] .
  * Seligman [1997] W. G. Seligman,  _A Next-to-Leading Order QCD Analysis of Neutrino - Iron Structure Functions at the Tevatron_ , Ph.D. thesis, Nevis Labs, Columbia U. (1997).
  * Tzanov [2005] M. M. Tzanov,  _Precise measurement of neutrino and anti-neutrino differential cross sections on iron_ , Ph.D. thesis, Pittsburgh U. (2005).
  * Wu  _et al._ [2008] Q. Wu _et al._ (NOMAD), A Precise measurement of the muon neutrino-nucleon inclusive charged current cross-section off an isoscalar target in the energy range 2.5 ¡ E(nu) ¡ 40-GeV by NOMAD, Phys. Lett. B 660, 19 (2008), arXiv:0711.1183 [hep-ex] .
  * Klein  _et al._ [2020] S. R. Klein, S. A. Robertson, and R. Vogt, Nuclear effects in high-energy neutrino interactions, Phys. Rev. C 102, 015808 (2020), arXiv:2001.03677 [hep-ph] .
  * Jeong and Reno [2010] Y. S. Jeong and M. H. Reno, Tau neutrino and antineutrino cross sections, Phys. Rev. D 82, 033010 (2010), arXiv:1007.1966 [hep-ph] .
  * Conrad  _et al._ [2010] J. Conrad, A. de Gouvea, S. Shalgar, and J. Spitz, Atmospheric Tau Neutrinos in a Multi-kiloton Liquid Argon Detector, Phys. Rev. D 82, 093012 (2010), arXiv:1008.2984 [hep-ph] .
  * Li and Racine [2007] Q. Li and J. S. Racine,  _Nonparametric Econometrics: Theory and Practice_ (Princeton University Press, 2007).
  * Aartsen  _et al._ [2020b] M. G. Aartsen _et al._ (IceCube), Computational techniques for the analysis of small signals in high-statistics neutrino oscillation experiments, Nucl. Instrum. Meth. A 977, 164332 (2020b), arXiv:1803.05390 [physics.data-an] .
  * Barger  _et al._ [1980] V. Barger, K. Whisnant, S. Pakvasa, and R. J. N. Phillips, Matter effects on three-neutrino oscillations, Phys. Rev. D 22, 2718 (1980).
  * Rott  _et al._ [2015] C. Rott, A. Taketa, and D. Bose, Spectrometry of the Earth using Neutrino Oscillations, Sci. Rep. 5, 15225 (2015), arXiv:1502.04930 [physics.geo-ph] .
  * Feldman and Cousins [1998] G. J. Feldman and R. D. Cousins, A Unified approach to the classical statistical analysis of small signals, Phys. Rev. D 57, 3873 (1998), arXiv:physics/9711021 .
  * Abbasi  _et al._ [2021] R. Abbasi _et al._ (IceCube), A calibration study of local ice and optical sensor properties in IceCube, PoS ICRC2021, 1023 (2021), arXiv:2107.10435 [astro-ph.HE] .
  * Chirkin and Rongen [2020] D. Chirkin and M. Rongen (IceCube), Light diffusion in birefringent polycrystals and the IceCube ice anisotropy, PoS ICRC2019, 854 (2020), arXiv:1908.07608 [astro-ph.HE] .
  * Adamson  _et al._ [2020] P. Adamson _et al._ (MINOS+), Precision Constraints for Three-Flavor Neutrino Oscillations from the Full MINOS+ and MINOS Dataset, Phys. Rev. Lett. 125, 131802 (2020), arXiv:2006.15208 [hep-ex] .
  * Acero  _et al._ [2022] M. A. Acero _et al._ (NOvA), Improved measurement of neutrino oscillation parameters by the NOvA experiment, Phys. Rev. D 106, 032004 (2022), arXiv:2108.08219 [hep-ex] .
  * Fernandez [2021] P. Fernandez (Super-Kamiokande), Atmospheric neutrino oscillations with Super- Kamiokande and prospects for SuperK-Gd (2021), ” XIX International Workshop on Neutrino Telescopes”.
  * Abe  _et al._ [2023] K. Abe _et al._ (T2K), Measurements of neutrino oscillation parameters from the T2K experiment using $3.6\times 10^{21}$ protons on target, Submitted to Eur. Phys. J (2023), arXiv:2303.03222 [hep-ex] .
  * Ishihara [2021] A. Ishihara (IceCube), The IceCube Upgrade - Design and Science Goals, PoS ICRC2019, 1031 (2021), arXiv:1908.09441 [astro-ph.HE] .

##  Appendix A Event selection cuts

The full list of cuts used for the common event selection, from Level 2 to
Level 5, is given in Table 7. The description of the algorithms is given in
Sec. III.

Cut name | Keep events if
---|---
Targeting coincident events |
Uncleaned $t$ length | <13000 ns
Cleaned $t$ length | <5000 ns
Targeting noise |
NChannel cleaned | >= 6
Noise engine | Pass
MicroCount | >2
DC Fiducial hits | >2
L4 noise classifier | >0.85
Targeting muons |
NAbove 200 | <10
VertexGuess Z | <-120m
Causal veto Hits | <7
Veto/fiducial hits | <1.5
C2HR6 | >0.37
RTVeto (Nfid <75) | <4
RTVeto (75<Nfid<100) | <5
L4 muon classifier score | >0.9
Vertex position in $\rho$ | <150m
Vertex position in $z$ | -490m<$p_{z}$<-220m
$\cos(\Delta\theta)$ (corridor-fit) |  $\leq$ 0.7 rad
Number of corridor hits |  $\leq$ 2
Table 7: List of variables and cut values used for the common DeepCore events
election. The description of the variables is given in Sec. III as well as in
[27].

##  Appendix B Impact of systematic uncertainties

Figures 28 to 44 show the percent change in the expected number of events in
the analysis histogram when the value of the nuisance parameters included in
the study is changed. The values chosen for this depiction correspond to the
+1$\sigma$ of the posterior distributions of each parameter as obtained from
ensemble fits to pseudo-data and therefore are representative of the scale of
the variations that the study is sensitive to. For example, for Fig. 28 the
DOM optical efficiency was increased by 2%, which increases the expectation of
number of events by up to 6% in the highest energies of the mixed-channel
histogram, while simultaneously reducing the number of events by about 2% in
the track histogram along the region where neutrino oscillations lead to
$\nu_{\mu}$ disappearance.

Figure 28: DOM optical efficiency. Figure 29: Hole ice p0. Figure 30: Hole ice
p1. Figure 31: Ice absorption. Figure 32: Ice scattering. Figure 33: Atm. flux
$\Delta\gamma$. Figure 34: Atm. flux A-F parameters. Figure 35: Atm. flux G
parameter. Figure 36: Atm. flux H parameter. Figure 37: Atm. flux W parameter
($K^{+}$). Figure 38: Atm. flux W parameter ($K^{-}$). Figure 39: Atm. flux Y
parameter. Figure 40: Axial mass CCQE. Figure 41: Axial mass RES. Figure 42:
Deep inelastic scattering correction to CSMS. Figure 43: NC normalization.
Figure 44: Atmospheric $\mu$ scale.

##  Appendix C Data long-term stability

We verified the stability of the data by studying the compatibility of several
variables across all years in the study using a Kolmogorov-Smirnov test. The
observables used to bin the data are shown in Fig. 45. Two selection variables
and one control variable are also shown in Fig. 46. The L4 muon classifier
score is of special interest, since it was trained using a subset of data from
2014. As shown in Fig. 46, all of the years have a similar behavior.

Figure 45: Agreement between years for the reconstructed quantities used in
the fit. Top panels show the distribution for every year, compared with the
average. Bottom panels display the Kolmogorov-Smirnov p-values calculated
between each season of data.

Figure 46: Agreement between years for two of the variables used in the event
selection (L4 muon classifier score and depth of the event vertex) and one
low-level control variable that shows the typical charge per hit observed. For
a description of the content of the panels, see Fig. 45.