Main

The discovery of neutrino flavour oscillations1,2 proves that neutrinos must have a mass, unlike originally assumed in the standard model of particle physics. Neutrino oscillation experiments have shown that the weakly interacting neutrino flavour eigenstates νf, where f ∈ {e, μ, τ} for electron, muon and tau-neutrino, are admixtures of the three neutrino-mass eigenstates νi with mass eigenvalues mii ∈ {1, 2, 3}. Although neutrino oscillation experiments can probe the differences of squared neutrino-mass eigenvalues \({{\Delta }}{m}_{ij}^{2}\), the absolute neutrino-mass scale remains one of the most pressing open questions in the fields of nuclear, particle and astroparticle physics today. In this paper, we report a measurement of the effective electron anti-neutrino mass defined as \({{\mbox{}}}{m}_{\nu }^{2}{{\mbox{}}}={\sum }_{i}| {U}_{{{{\rm{e}}}}i}{| }^{2}{m}_{i}^{2}\), where Uei are elements of the Pontecorvo–Maki–Nakagawa–Sakata matrix that describes the mixing of neutrino states.

The neutrino masses are at least five orders of magnitude smaller than the mass of any other fermion of the standard model, which may point to a different underlying mass-creation mechanism3. The determination of the neutrino mass would, thus, shed light on the fundamental open question of the origin of particle masses. Despite the smallness of their masses, neutrinos play a crucial role in the evolution of large-scale structures of our cosmos due to their high abundance in the Universe4,5. A direct measurement of the neutrino mass could, hence, provide a key input to cosmological structure formation models. In this respect, cosmological observations themselves provide a stringent limit on the sum of neutrino masses of ∑mi < 0.12 eV (95% confidence level (CL))6,7 (here we use the convention c = 1 for the speed of light). However, these limits strongly rely on the underlying cosmological assumptions8,9. An independent measurement of neutrino mass could help in breaking the parameter degeneracies of cosmological models. A powerful way to probe this neutrino property in the laboratory is via a search for neutrinoless double-beta (double-β) decay. In contrast to mν, the effective mass in double-β decay is given by \({m}_{\beta \beta }=\left|{\sum }_{i}{U}_{{{{\rm{e}}}}i}^{2}{m}_{i}\right|\). This neutrino-mass interpretation is only valid under the assumption that neutrinos are their own anti-particle (Majorana particle) and that light neutrinos mediate the decay. The current most stringent limits derived from different isotopes are mββ < 79−180 meV (76Ge) (ref. 10), mββ < 75−350 meV (130Te) (ref. 11) and mββ < 61−165 meV (136Xe) (ref. 12), and the spread is related to uncertainties in the model-dependent nuclear matrix element calculation.

The most direct way to assess the neutrino mass is via the kinematics of single-β decays or electron capture processes. This method is independent of any cosmological model and of the mass nature of the neutrino, that is, it may be a lepton of the Majorana or Dirac type. The neutrino masses mi lead to a reduction in the maximum observed energy of the decay and a small spectral-shape distortion close to the kinematic endpoint of the β-decay spectrum. In the quasi-degenerate mass regime, where mi > 0.2 eV, the mass splittings are negligible with respect to masses mi, and the observable value can be approximated as \({{\mbox{}}}{m}_{\nu }^{2}{{\mbox{}}}={\sum }_{i}| {U}_{{{{\rm{e}}}}i}{| }^{2}{m}_{i}^{2}\) (refs.13,14).

The Karlsruhe Tritium Neutrino (KATRIN) experiment15,16 exploits the single-β decay of molecular tritium as

$${{{{\rm{T}}}}}_{2}{\to }^{3}{{{{\rm{HeT}}}}}^{+}+{{{{\rm{e}}}}}^{-}+{\bar{\nu }}_{e}$$
(1)

and currently provides the best neutrino-mass sensitivity in the field of direct neutrino-mass measurements with its first published limit of mν < 1.1 eV (90% CL)17,18. KATRIN is designed to determine the neutrino mass with a sensitivity of close to 0.2 eV (90% CL) in a total measurement time of 1,000 days (ref. 15). Another class of experiments is based on the electron capture of 163Ho, where the decay energy is measured with micro-calorimeters19,20. Note that electron capture experiments based on 163Ho measure the mass of the neutrino ν and β− experiments based on tritium measure that of the anti-neutrino \(\bar{\nu }\). New ideas exist to extend the sensitivity of tritium-based neutrino-mass experiments beyond the KATRIN design sensitivity by new technologies, such as cyclotron radiation emission spectroscopy and the development of atomic tritium sources21,22.

In this work, we present the second neutrino-mass result of KATRIN, reaching an unprecedented sub-electronvolt sensitivity and limit in mν from a direct measurement.

The KATRIN experiment

The design requirements to detect the small signature of a neutrino mass in the last few electron-volts of the β-decay spectrum are a high tritium activity (1 × 1011 Bq), a low background rate (≲0.1 counts per second (cps)), an electron-volts-scale energy resolution and an accurate (0.1% level) theoretical prediction of the integral spectrum.

The KATRIN experiment tackles these challenges by combining a high-activity molecular tritium source with a high-resolution spectrometer of the magnetic adiabatic collimation and electrostatic (MAC-E)-filter type16. The experiment is hosted by the Tritium Laboratory Karlsruhe, allowing the safe supply of tritium at the 10 g scale with continuous tritium reprocessing23,24. Figure 1 shows the 70-m-long KATRIN beamline. The isotope tritium has a short half-life of 12.3 years and a low endpoint of 18.6 keV, both of which are favourable properties to achieve high count rates near the endpoint. Moreover, the theoretical calculation of the super-allowed decay of tritium is well understood25,26,27. The strong gaseous tritium source (nominal activity, 1 × 1011 Bq) is established by the throughput of 40 g per day of molecular tritium through the 10-m-long source tube, which is cooled to 30 K to reduce the thermal motion of tritium molecules.

Fig. 1: Illustration of the 70-m-long KATRIN beamline.
figure 1

The main components are labelled. The transport of β-electrons and magnetic adiabatic collimation of their momenta p are illustrated. a–f, View into the tritium source depicts three systematic effects: molecular excitations during β-decay (a), scattering of electrons off the gas molecules (b) and spatial distribution of the electric potential in the source Usrc(r, z) (c). The view into the spectrometer illustrates the main background processes arising from radon decays inside the volume of the spectrometer47,48,49 (d), highly excited Rydberg atoms sputtered off from the structural material via α-decays of 210Po (refs. 44,45,46) (e) and positive ions created in a Penning trap between the two spectrometers57 (f). Low-energy electrons, created in the volume as a consequence of radon decays or Rydberg-atom ionizations, can be accelerated by qUana towards the focal-plane detector, making them indistinguishable from signal electrons66.

Full size image

A system of 24 superconducting magnets28 guides the electrons out of the source towards the spectrometer section and to the focal-plane detector. Between the source and main spectrometer, differential29 and cryogenic30 pumping systems reduce the flow of tritium by 14 orders of magnitude.

High-precision electron spectroscopy is achieved with the spectrometer of MAC-E-filter type31,32. The spectrometer acts as a sharp electrostatic high-pass filter, transmitting only electrons (with charge q = −e) above an adjustable energy threshold qU, where U is the retarding potential applied at the spectrometer. A reduction in the magnetic-field strength by about four orders of magnitude from the entrance (exit) of the spectrometer Bsrc = 2.5 T (Bmax = 4.2 T) to its centre, the analysing plane (Bana = 6.3 × 10−4 T), collimates the electron momenta. This configuration creates a narrow filter width of ΔE = 18.6 keV × (Bana/Bmax) = 2.8 eV and allows for a large angular acceptance, with a maximum pitch angle for the β-decay electrons of \({\theta }_{\max }=\arcsin \sqrt{({B}_{{{{\rm{src}}}}}/{B}_{\max })}\) = 50.4° in the source. The pitch angle refers to the angle between the electron’s momentum and the direction of magnetic field at the position of the electron. A 12.6-m-diameter magnetic coil system surrounding the spectrometer finely shapes the magnetic field and compensates the Earth’s magnetic field33,34.

The transmitted electrons are detected by a 148 pixel silicon-PIN-diode focal-plane detector installed at the exit of the spectrometer35. By measuring the count rate of transmitted electrons for a set of qU values, the integral β-decay spectrum is obtained. The main spectrometer is preceded by a smaller pre-spectrometer, which operates on the same principle and transmits only electrons above 10 keV, reducing the flux of electrons into the main spectrometer. The upstream end of the beamline is terminated with a gold-plated rear wall, which absorbs the non-transmitted β-electrons and defines the reference electric potential of the source. The rear wall is biased to a voltage \({{{\mathcal{O}}}}\) (100 mV) to minimize the difference in the surface potential to that of the beam tube, which minimizes the inhomogeneity of the source electric potential.

The rear section is equipped with an angular- and energy-selective electron gun36, which is used to precisely determine the scattering probability of electrons with the source gas, governed by the product of column density (number of molecules per square centimetre along the length of the source) and scattering cross section. Furthermore, we use the electron gun to measure the distribution of energy losses for 18.6 keV electrons scattering off the molecular tritium gas, providing one of the most precise energy-loss functions for this process to date37. Another key calibration source is gaseous krypton, which can be co-circulated with the tritium gas38. Mono-energetic conversion electrons from the decay of the metastable state 83mKr are used to determine spatial and temporal variations in the electric potential of the tritium source. These variations are caused by a weak cold-magnetized plasma, which arises from the high magnetic field (2.5 T) and a large number of ions and low-energy electrons (~1 × 1012 m−3) in the tritium source. The methods of calibration are described in more detail in Methods.

The beamline is equipped with multiple monitoring devices. The throughput of tritium gas within the tritium source tube and tritium circulation loop is measured by a flow meter. A laser Raman system continuously monitors the gas composition, providing a measurement at the ≤0.05%-precision level each minute. A silicon drift detector system installed in the transport section and a β-induced X-ray system at the rear section39 continuously monitor the tritium activity, yielding a result at the 0.03%-precision level each minute. The high voltage of the main spectrometer is continuously measured at the parts-per-million level with a high-precision voltage divider system40,41 and an additional monitoring spectrometer42. The magnetic fields are determined with a high-precision magnetic-field sensor system43.

The current background level of KATRIN of about 220 mcps mainly originates from the spectrometer section. The dominant source of background arises from α-decays of 210Po (refs. 44,45,46) in the structural material of the spectrometer vessel. The recoiling 206Pb creates highly excited Rydberg states at the inner spectrometer surfaces, which can be ionized during propagation in the inner volume by thermal radiation. The second source is 219Rn and 220Rn decays in the spectrometer volume, creating primary electrons that are magnetically trapped and slowly cool down by ionizing the residual gas in the spectrometer47,48,49. The resulting low-energy electrons of both sources are accelerated by retarding energy qUana towards the focal-plane detector, making them indistinguishable from signal electrons using the energy information only.

After the successful commissioning of the complete KATRIN beamline in the summer of 2017 (ref. 50), the first tritium operation was demonstrated with a small tritium activity of 5 × 108 Bq in mid-2018 (ref. 51). During the first KATRIN neutrino-mass (KNM1) campaign in 2019 (ref. 17), the source was operated in a ‘burn-in’ configuration at a reduced activity of 2.5 × 1010 Bq, which is required when structural materials are exposed to high amounts of tritium for the first time. Major technical achievements of the second KATRIN neutrino-mass (KNM2) campaign are the operation of the tritium source at its nominal activity of 9.5 × 1010 Bq and improved vacuum conditions by baking of the spectrometer52 to 200 °C for approximately two weeks that led to a reduction in the background by 25% to 220 mcps. The thermal conditioning of the surfaces reduces the coverage of weakly bound atoms (which contribute to the background) and removes water molecules from the cryogenic copper baffles (which improves the capture efficiency for radon emanating from the main getter pumps located behind the baffles53).

We, thus, increased the β-electron-to-background ratio by a factor of 2.7 with respect to the first campaign. In the last 40 eV of the integral spectrum, we collected a total number of 3.7 × 106 β-electrons. Figure 2a compares the spectra of both neutrino-mass campaigns. A direct comparison of the experimental parameters is given in Table 1.

Fig. 2: Measured rate at each retarding energy for KNM1 (refs. 17,18) and KNM2 campaigns.
figure 2

a, Data points with statistical error (multiplied by a factor of 50) and best-fit model (blue and grey lines) individually shown for each campaign. The count rates are summed over all the detector rings. The graph illustrates a reduced background rate, higher signal strength and overall higher statistics (smaller error bars). The fit description is given in equation (2). b, Normalized residuals for the fit (blue line) to the data from KNM2. The shaded areas indicate statistical and total uncertainties. The contribution of systematic uncertainties is derived with the covariance matrix method, explained in the main text. c, Data points with statistical error (multiplied by a factor of 10) for the 12 detector rings in the KNM2 campaign. The turquoise lines show the simultaneous fit to all the data points. d, Normalized residuals for the fit (turquoise line) to the data of the (exemplary) third detector ring. e, Total measurement time at each retarding energy for the KNM1 and KNM2 campaigns.

Source data

Full size image
Table 1 Comparison of key numbers for the KNM campaigns
Full size table

Measurement of tritium β-decay spectrum

The integral β-decay spectrum is obtained by repeatedly measuring the count rate Rdata(qUi) for a set of 39 non-equidistant high-voltage settings Ui, creating retarding energy qUi for electrons with charge q. The retarding energy is adjusted in the range of qUi ∈ (E0 – 300 eV, E0 + 135 eV), where E0 = 18,574 eV is the approximate spectral endpoint. Note that for the spectral fit, only 28 out of those points in the range of qUi ∈ [E0 − 40 eV, E0 + 135 eV] are used. Data points further below the endpoint are used to monitor the activity stability, complementing the other monitoring devices mentioned above. The time spent at each high-voltage set point (called the scan step) varies between 17 and 576 s, which corresponds to a total measurement time of 3.4 h for the shortest scan step at qUi = 18,534.2 eV and 64.7 h for the longest step at qUi = 18,565.2 eV; this is chosen to optimize the neutrino-mass sensitivity. The total time to measure the rate at all the 39 retarding energies (a so-called scan) is about 2 h. As shown in Fig. 2e, the measurement-time distribution compensates the steeply falling count rate towards the endpoint and peaks at approximately 10 eV below the endpoint, where a neutrino-mass signal would be at the maximum. Based on predefined experimental conditions and data-quality criteria, 361 golden scans were chosen for the presented analysis. Since the high-voltage values can be set with a reproducibility at the sub-parts-per-million level, these scans are later combined into a single spectrum by adding the counts at a given set point.

The focal-plane detector (Fig. 1, inset) is segmented into 148 individual pixels to account for spatial variations in the electromagnetic fields inside the source and spectrometer, as well as of the background. Removing malfunctioning pixels and those in the outer rings shadowed by parts of the beamline upstream, 117 pixels have been selected. For the presented analysis, these golden pixels are grouped into 12 concentric rings, resulting in 336 data points Rdata(qUi, rj), where i ∈ {1,…,nqU = 28} and j ∈ {1,…,nrings = 12} per scan. Figure 2c displays the spectra for each of the 12 detector rings, where all the scans have been combined.

Data analysis strategy

To infer the neutrino mass, we fit the spectral data with a spectrum prediction, given by an analytical description of the β-decay spectrum and the experimental response function, described in the following sections.

Spectrum calculation

The prediction of detection rate R(qUi, rj) is given by a convolution of the differential β-decay spectrum Rβ(E) with experimental response function f(E, qUi, rj) and background rate Rbg(qUi, rj):

$$R(q{U}_{i},{r}_{j})={A}_{{{{\rm{s}}}}}{N}_{{{{\rm{T}}}}}\int\nolimits_{q{U}_{i}}^{{E}_{0}}{R}_{\beta }(E)f(E,q{U}_{i},{r}_{j})\mathrm{d}E+{R}_{{{{\rm{bg}}}}}(q{U}_{i},{r}_{j}).$$
(2)

Here NT is the signal normalization calculated from the number of tritium atoms in the source, maximum acceptance angle and detection efficiency. Also, As ≈ 1 is an additional normalization factor, which is used as a free parameter in the fit. The differential β-decay spectrum Rβ(E) is given by Fermi’s theory. In the analysis, we include radiative corrections and the molecular final-state distribution54,55 in the differential β-decay spectrum. Doppler broadening due to the finite thermal motion of tritium molecules in the source, as well as energy broadenings due to spatial and temporal variations in the spectrometer and source electric potential, is emulated by Gaussian broadenings of the final-state distribution.

The response function f(E, qUi, rj) includes energy losses due to scattering and synchrotron radiation in the high-magnetic-field regions, as well as the transmission of electrons through the main spectrometer. Compared with previous KATRIN analyses, we now use a modified transmission function to account for the non-isotropy of β-electrons leaving the source. The angular distribution of electrons leaving the source is slightly non-isotropic due to the effective path length of electrons through the source, which depends on the pitch angle. A detailed description of the spectrum calculation can be found in ref. 56 and Methods.

Unbiased parameter inference

We infer \({m}_{\nu }^{2}\) by fitting the experimental spectrum Rdata(qUi, rj) with prediction R(qUi, rj) by minimizing the standard χ2 = RdataC–1R function, where C contains the variance and covariance of the data points. In addition to the neutrino mass squared (\({{\mbox{}}}{m}_{\nu }^{2}{{\mbox{}}}\)), the parameters As(rj), Rbg(rj) and effective endpoint E0(rj) are treated as independent parameters for the 12 detector rings, leading to a total number of (1 + (3 × 12)) = 37 free parameters in the fit. The introduction of ring-dependent parameters was chosen to allow for possible unaccounted radial effects. In particular, the effective endpoint E0(rj) could account for a possible radial dependence of the electric potential in the source. However, the final fit revealed a negligible (<100 meV) radial variation in the endpoint. Another advantage of ring-dependent parameters is to avoid the averaging of transmission function over all the rings, which would introduce energy broadening and hence reduce the resolution.

The following analysis procedure was implemented to minimize the potential for human-induced biases. As the first step, the full analysis is performed on a Monte Carlo copy of each scan, simulated based on the actual experimental parameters (such as magnetic fields and column density), obtained by external calibration measurements. In the next step, the fit is performed on the experimental dataset, but with a randomly broadened molecular final-state distribution, which imposes an unknown bias on the observable \({{\mbox{}}}{m}_{\nu }^{2}{{\mbox{}}}\). In the final step, the analysis of experimental data with unmodified final-state distribution is executed. To prevent human-induced errors, each analysis step was performed by four independent analysis teams, using different software and strategies. Consistent results were obtained at each analysis step before proceeding to the next stage.

Systematic effects

The analytical description R(qUi, rj) of the integral β-decay spectrum contains various experimental and theoretical parameters, such as magnetic fields, column density and probability for given molecular excitations, which are known with a certain accuracy. Different techniques (based on covariance matrices, Monte Carlo propagation, nuisance parameters and Bayesian priors) are applied to propagate these systematic uncertainties in the final result. A detailed description of these methods can be found in Methods.

The neutrino-mass result presented here is dominated by statistical uncertainty. The largest systematic uncertainties are related to background properties and source electric potential. First, radon decays in the main spectrometer lead to a non-Poissonian background-rate overdispersion49 of about 11%, effectively increasing the statistical uncertainty. Second, mechanisms for background generation may show a retarding-potential dependence of the background, parametrized by a slope (\({m}_{{{{\rm{bg}}}}}^{qU}\) = ( 0 ± 4.74) mcps keV–1). Third, a removal of stored electrons from a known Penning trap between the spectrometers57 after each scan step can lead to a slowly increasing background rate during each scan (\({m}_{{{{\rm{bg}}}}}^{{{{{{t}}}}}_{{{{\rm{scan}}}}}}\) = (3 ± 3) µcps s–1) and thus to a scan-step-duration-dependent background contribution. Finally, spatial and temporal variations in the source electric potential modify the spectral shape and therefore lead to a relevant systematic uncertainty for the neutrino-mass measurement. The impact of all the systematic effects on the neutrino mass is listed in Table 2 and described in detail in Methods.

Table 2 Breakdown of uncertainties of the neutrino-mass-squared best fit as sorted by size
Full size table

Best-fit result and upper limit

The χ2 minimization reveals an excellent goodness of fit with a χ2 per degree of freedom of 0.9, corresponding to a p value of 0.8. For the best fit of the neutrino mass squared, we find \({{\mbox{}}}{m}_{\nu }^{2}{{\mbox{}}}=0.2{6}_{-0.34}^{+0.34}\) eV2 based on the Monte Carlo propagation technique (Fig. 3). The independent analysis methods agree within about 5% of the total uncertainty. The total uncertainty on the fit is dominated by the statistical error followed by uncertainties in the background parameters and source electric potential. The full breakdown of uncertainties can be found in Table 2.

Fig. 3: Distribution of fitted \({{\mbox{}}}{m}_{\nu }^{2}{{\mbox{}}}\) and E0 values.
figure 3

The distribution is obtained by the Monte Carlo propagation technique (Methods). The best fit is defined as the one-dimensional weighted median of the distribution. The endpoint is calculated using the weighted mean of the 12 ring-wise effective endpoints, taking into account the correlations between the 12 values for each sample. The 1σ (black line) and 2σ (blue line) contours are indicated.

Source data

Full size image

Based on the best-fit result, we obtain an upper limit of mν < 0.9 eV at 90% CL using the Lokhov–Tkachov method58. The Feldman–Cousins technique59 yields the same limit for the obtained best fit. This result is slightly higher than the sensitivity of 0.7 eV due to the positive fit value, which is consistent with ~0.8σ statistical fluctuation assuming a true neutrino mass of zero. We also perform a Bayesian analysis of the dataset, with a positive flat prior on \({m}_{\nu }^{2}\). The resulting Bayesian limit at 90% credibility is mν < 0.85 eV.

The ring-averaged fitted effective endpoint is E0 = (18,573.69 ± 0.03) eV. The Q value defines the energy release in a nuclear reaction. Taking into account the centre-of-mass molecular recoil of T2 (1.72 eV), as well as the absolute electric source potential ϕsrc (σ(ϕsrc) = 0.6 V) and work function of the main spectrometer ϕMS (σ(ϕMS) = 0.06 eV), we find a Q value of (18,575.20 ± 0.60) eV, which is consistent with the previous KATRIN neutrino-mass campaign, namely, (18,575.20 ± 0.50) eV (ref. 51), and the calculated Q value from 3He−3H atomic-mass difference of (18,575.72 ± 0.07) eV (ref. 60). The good agreement underlines the stability and accuracy of the absolute energy scale of the apparatus.

We combined the neutrino-mass results from this work with the previously published KATRIN (KNM1) results17,18. A simultaneous fit of both datasets yields \({{\mbox{}}}{m}_{\nu }^{2}{{\mbox{}}}\) = (0.1 ± 0.3) eV2 and a corresponding upper limit of mν < 0.8 eV at 90% CL, based on the Lokhov–Tkachov or Feldman–Cousins technique. The same result is obtained when multiplying the \({m}_{\nu }^{2}\) distributions from Monte Carlo propagation or adding the χ2 profile of the individual fits. As both datasets are dominated by statistics, correlated systematic uncertainties between both campaigns are negligible. Furthermore, we investigated a Bayesian combination, where the KNM1 posterior distribution of \({m}_{\nu }^{2}\) is used as the prior for the KNM2 analysis, yielding consistent results. More details on the combined analyses can be found in the Supplementary Information.

Conclusion and outlook

The second neutrino-mass measurement campaign of KATRIN, presented here, reached sub-electronvolt sensitivity (0.7 eV at 90% CL). Combined with the first campaign, we set an improved upper limit of mν < 0.8 eV (90% CL). We, therefore, have narrowed down the allowed range of quasi-degenerate neutrino-mass models and we have provided model-independent information about the neutrino mass, which allows the testing of non-standard cosmological models61,62. Figure 4 displays the evolution of the best-fit \({m}_{\nu }^{2}\) results from historical neutrino-mass measurements up to the present day.

Fig. 4: Comparison of best-fit values and total uncertainties with previous neutrino-mass experiments.
figure 4

The error bars are generated from combined statistical and systematic uncertainties. References: Los Alamos (1991)67, Tokyo (1991)68, Zürich (1992)69, Mainz (1993)70, Beijing (1993)71, Livermore (1995)72, Troitsk (1995)73, Mainz (1999)13, Troitsk (1999)74, Mainz (2005)75, Troitsk (2011)76, KATRIN (2019)17,18 and KATRIN (2021); this work, KATRIN (combined): KATRIN (2019) combined with KATRIN (2021). Note that the published gaseous tritium results from Los Alamos and Livermore were analysed using different molecular final-state distributions compared with current state-of-the-art final-state distributions. These earlier distributions have been shown to contribute to the reported negative \({{\mbox{}}}{m}_{\nu }^{2}{{\mbox{}}}\) central values77.

Source data

Full size image

Compared with its previous measurement campaign, the KATRIN experiment has decreased the statistical and systematic uncertainties by about a factor of three and two, respectively. With the total planned measurement time of 1,000 days, the total statistics of KATRIN will be increased by another factor of 50. A reduction in the background rate by a factor of two will be achieved by an optimized electromagnetic-field design of the spectrometer section45,63. Moreover, by eliminating the radon- and Penning-trap-induced background, the background-related systematic effects are expected to be substancially reduced. Together with a high-rate 83mKr calibration scheme and a method to improve the determination of magnetic fields, we will minimize the dominant systematic uncertainties to reach the target sensitivity of mν in the vicinity of 0.2 eV.

High-precision β-spectroscopy with the KATRIN experiment has proven to be a powerful means to probe the neutrino mass with unprecedented sensitivity and to explore physics beyond the standard model, such as sterile neutrinos64. Together with cosmological probes and searches for neutrinoless double-β decay65, the upcoming KATRIN data will play a key role in measuring the neutrino-mass parameters.

Methods

Here we describe the data analysis chain starting from data processing to high-level fit and limit setting. Moreover, we provide details on one of the key calibration campaigns, concerning the source electric potential. A more extensive description of the KATRIN analysis procedure is provided elsewhere18.

Data processing, selection and combination

The first step of the analysis chain is data preparation. Raw data are combined into integral spectral data points, which are then fitted with an analytical spectrum prediction including the response of the experiment.

Rate determination

Electrons transmitted through the main spectrometer are further accelerated by a post-acceleration electrode with potential UPAE = 10 keV before they are detected by the focal-plane detector35. The latter provides a high detection efficiency (>95%) and moderate energy resolution (2.8 keV full-width at half-maximum at 28 keV). The total rate per pixel at a given retarding potential qUi is determined by integrating the rate in a wide and asymmetric region of interest (ROI) of 14 ≤ E ≤ 32 keV. The asymmetric ROI is chosen to account for energy losses in the dead layer of the Si PIN diodes78, for partial energy deposition due to the backscattering of electrons off the detector surface as well as charge sharing between pixels. The detector efficiency slightly depends on qUi. Three effects are considered. (1) The differential energy spectrum of the transmitted electrons at the focal-plane detector is shifted (and slightly scaled) according to qUi. Since the same ROI is used for each qUi, some electrons are no longer covered by the fixed ROI when lowering the potential, which effectively changes the detection efficiency. Based on reference measurements, the relative reduction in detection efficiency at E0 − 300 eV, with respect to the efficiency at E0, is determined to be δROI = 0.00048, with an uncertainty of 0.16%. The data are corrected according to the efficiency at the given qUi value. (2) The count rate at the focal-plane detector as well as the probability of pile-up depends on qUi. The energy of most pile-up events is added such that they are not covered by the ROI, thereby effectively changing the detector efficiency with qUi. Based on a random-coincidence model assuming a Poisson-distributed signal, the relative reduction in efficiency at E0 − 300 eV is estimated to be δPU = 0.00017, with an uncertainty of 18%. The data are corrected accordingly. (3) Finally, a fraction of about 20% of all the electrons impinging on the detector surface are backscattered. The majority of these electrons are back-reflected to the same detector pixel (within the integration time of the energy filter) by the pinch magnetic field or retarding potential. However, for low retarding potentials and small energy depositions in the focal-plane detector, the backscattered electrons have a chance of remaining undetected by overcoming the retarding potential a second time in the direction towards the tritium source. The lower the qUi, the higher is the probability for lost electrons, effectively changing the detection efficiency. Based on Monte Carlo simulations, we estimate the efficiency reduction to be δBS < 0.001 at E0 − 300 eV. We neglect this effect in this analysis. Systematic uncertainties related to efficiency corrections are negligibly small for the presented analysis, which only considers the last 40 eV of the tritium spectrum.

Selection of golden pixels and scans

The focal-plane detector is segmented into 148 pixels of equal area. For the analysis, 117 pixels have been chosen. The rejected pixels show undesirable characteristics such as broadened energy resolution or higher noise levels (a total of 6 pixels), partial shadowing by the silicon-drift-detector-based beam monitor (a single pixel) and decreased rate due to misalignment of the beamline with respect to the magnetic-flux tube (a total of 24 pixels) that stems from tiny deviations in the orientation of the 24 superconducting magnets from the design28.

Out of the 397 scans recorded, 361 passed a strict quality assessment and were included in the neutrino-mass analysis. The other runs were rejected for reasons of failed set points of spectrometer electrodes (9 runs), downtime of the laser Raman system (23 runs) and other operational conditions (4 runs). The golden pixel and run selection is performed before the unblinding of the data.

Data combination

The 117 golden pixels are grouped into 12 detector rings and the detector counts recorded within each ring are summed to obtain 12 independent spectra. All the golden scans are combined by adding the counts recorded at the same retarding energies. This method leads to a total number of data points of ndata = nqU × nrings = 336.

Spectrum calculation

To infer the neutrino mass, the integral spectrum is fitted with the analytical spectrum calculation including the experimental response. As shown in equation (2), the theoretical spectrum is composed of the differential tritium spectrum Rβ(E) and experimental response function f(E, qUi). The differential spectrum is given by

$$\begin{array}{ll}{R}_{{{{\rm{\beta }}}}}(E)&=\frac{{G}_{{{{\rm{F}}}}}^{2}{\cos }^{2}{{{\Theta }}}_{{{{\rm{C}}}}}}{2{\uppi }^{3}}{\left|{M}_{{{{\rm{nucl}}}}}\right|}^{2}F(E,Z^{\prime} =2)\\ &\cdot (E+{m}_{e})\sqrt{{(E+{m}_{e})}^{2}-{m}_{e}^{2}}\\ &\cdot \mathop{\sum}\limits_{f}{\zeta }_{f}\ {\varepsilon }_{f}(E)\sqrt{{{\varepsilon }_{f}(E)}^{2}-{m}_{\nu }^{2}}{{\Theta }}({\varepsilon }_{f}(E)-{m}_{\nu })\end{array},$$
(3)

where GF denotes the Fermi constant, cos2ΘC is the Cabibbo angle, ∣Mnucl∣2 is the energy-independent nuclear matrix element and F(E, Z′ = 2) is the Fermi function. Moreover, εf(E) = E0 − Vf − E, where E0 denotes the maximum kinetic energy of the electron in the case of zero neutrino mass and Vf describes the molecular excitation energies populated with probabilities ζf. Also, E and me denote the kinetic energy and mass of the β-electron, respectively.

Beyond the molecular effects, further theoretical corrections arise on the atomic and nuclear level79. Only radiative corrections to the differential spectrum are relevant for this analysis, which are included in the analytical description, but not shown here56. Furthermore, we emulate the effect of Doppler broadening, as well as spatial and temporal source and spectrometer electric potential variations, by broadening the final-state distribution with a Gaussian distribution (discussed in the main text).

The experimental response function

$$\begin{array}{ll}f(E-qU)&=\int\nolimits_{\epsilon = 0}^{E-qU}\,\int\nolimits_{\theta = 0}^{{\theta }_{\max }}{{{\mathcal{T}}}}(E-\epsilon ,\theta ,U)\sin \theta \\ &\cdot \mathop{\sum}\limits_{s}{P}_{s}(\theta )\ {f}_{s}(\epsilon )\ \mathrm{d}\theta \ \mathrm{d}\epsilon \end{array}$$
(4)

depicts the probability of an electron with starting energy E to reach the focal-plane detector. It combines the transmission function \({{{\mathcal{T}}}}(E-\epsilon ,\theta ,U)\) of the main spectrometer and electron’s energy losses ϵ due to inelastic scattering with the tritium molecules in the source. The scattering energy losses are described by the product of the s-fold scattering probabilities Ps(θ), which depend on the path length through the source and hence on the pitch angle θ (angle between the energy momentum and magnetic field at the position of the electron), and the energy-loss function fs(ϵ) for a given number of scatterings s. The energy-loss function is experimentally determined with the electron gun installed at the rear of the source. A typical response function and the corresponding energy-loss function is shown in Extended Data Fig. 1.

The integrated transmission function for an isotropic source of electrons is given by

$$T(E,U)=\int\nolimits_{\theta = 0}^{{\theta }_{{{\mbox{max}}}}}\ {{{\mathcal{T}}}}(E,\theta ,U)\cdot \sin \theta \ \mathrm{d}\theta$$
(5)
$$=\left\{\begin{array}{ll}0&,\,E-qU < 0\\ 1-\sqrt{1-\frac{E-qU}{E}\frac{{B}_{{{\mbox{src}}}}}{{B}_{{{{\rm{ana}}}}}}\frac{2}{\gamma \,+\,1}}\ &,\,0\le E-qU\le {{\Delta }}E\\ 1-\sqrt{1-\frac{{B}_{{{\mbox{src}}}}}{{B}_{{{\mbox{max}}}}}}&,\,E-qU > {{\Delta }}E\end{array}\right.\ .$$
(6)

It is governed by the magnetic fields at the starting position Bsrc of the electron, the maximum field Bmax in the beamline and the magnetic field in the spectrometer’s analysing plane Bana. The acceptance angle of the MAC-E-type filter is given by \({\theta }_{\max }=\arcsin \sqrt{({B}_{{{{\rm{src}}}}}/{B}_{\max })}\) = 50.4°. Synchrotron energy losses56 of β-electrons in the high magnetic field in the source and transport systems are included as an analytical correction to the transmission function (not shown here).

Systematic uncertainties

The analytical description R(qUi, rj) of the integral β-decay spectrum (equation (2)) contains various signal- and background-related parameters, which are known with a certain accuracy. In the following, we describe these parameters, their uncertainties and their treatment in the neutrino-mass analysis.

Signal-related systematic effects

Spectrum prediction includes uncertainties in the nine parameters of the empirical energy-loss function (the individual relative uncertainties of parameters σeloss,k are between 0.016% and 3.800%); a relative uncertainty of the product of column density and scattering cross section (σρd×σ = 0.2500%); relative uncertainties of the theoretical description of the molecular final-state distribution (σFSD = \({{{\mathcal{O}}}}\) (1.000%)); and relative uncertainties of the magnetic field in the source (\({\sigma }_{{B}_{{{{\rm{src}}}}}}\) = 1.700%), in the analysing plane (\({\sigma }_{{B}_{{{{\rm{ana}}}}}}\) = 1.000%) and the maximal field (\({\sigma }_{{B}_{\max }}\) = 0.100%). The variation in β-decay activity during a scan was determined from the product of column density (as determined from the tritium throughput) and tritium purity or alternatively from the beam monitor (silicon drift detector system in the transport section). The obtained standard deviation of the relative activity is negligibly small (σscan = 0.035%).

The spatial and temporal variations in the source electric potential were not included in the first neutrino-mass campaign. With the increase in source column density from 1.11 to 4.23 × 1017 cm–2 in this campaign, the creation rates and densities of the charge carriers, as well as the fraction of scattered electrons, increased accordingly, which makes plasma effects more relevant15,80,81. A detailed description of the plasma calibration is given in the main text.

Background-related systematic effects

Electrons created during radon decay in the main spectrometer volume can initially be magnetically trapped. Each of these trapped electrons can create a cluster of 10–100 secondary electrons by scattering off the residual gas in the main spectrometer49,66, which is operated at a pressure of 10−11 mbar. These secondary electrons arrive at the focal-plane detector within a time window of about 1,000 s, hence leading to a non-Poissonian rate distribution49. The observed background rate can be modelled by a Gaussian distribution, with a width that is 11.2% wider than expected from a purely Poisson distribution. This overdispersion is treated as an increased statistical uncertainty in the analysis.

As the transmission conditions for the background electrons slightly depend on the retarding-potential setting, a small retarding-potential dependence of the background can occur. In the analysis, we allow for a linear dependence of the background on the retarding potential and use dedicated test measurements to constrain the possible slope to \({m}_{{{{\rm{bg}}}}}^{qU}\) = (0 ± 4.7) mcps keV–1.

Finally, the pre- and main spectrometers, both being operated at a high voltage, create a Penning trap between them. Stored electrons and the subsequently produced positive ions, which can escape the trap into the main spectrometer, are a source of background, as illustrated in Fig. 1. To mitigate this background, the trap is emptied with an electron-catcher system57 after each scan step. However, a potentially small increase in the background rate within a scan step cannot be excluded, which could lead to a background dependence on the duration of the scan step. By fitting a linear increase to the rate evolution within all the scan steps, we find a slope of \({m}_{{{{\rm{bg}}}}}^{{{{{{t}}}}}_{{{{\rm{scan}}}}}}\) = (3 ± 3) µcps s–1, which is included in the β-decay spectrum fit. This effect was first observed in later physics runs, in which the scan-step duration was significantly increased. It was, therefore, only included after the data were already unblinded to a subgroup of the collaboration. As the evaluation of the size of this uncertainty was provided by an independent task group of the collaboration, the reported result remains free of bias.

Source-potential calibration

The absolute electric potential of the source does not affect the spectral shape of the measured spectrum. An unknown absolute source potential is largely absorbed by the effective endpoint, which is a free parameter in the fit. A change in the effective endpoint mostly leads to a shift in the spectrum and has a negligible effect on the spectral shape. Accordingly, an unknown radial variation in the electric potential is in good approximation absorbed by the ring-wise endpoint parameters. Consequently, these effects have a minor impact on the neutrino-mass analysis. However, any temporal variation, like high-frequency plasma instabilities, slow drifts of the absolute potential or longitudinal variations of the source potential can lead to spectral distortions, which are parameterized by Gaussian broadening σP and parameter ΔP that quantifies the longitudinal rear-to-front asymmetry of the electric potential of the source. This asymmetry of the potential results in a shift in the energy spectrum associated with the scattered electrons (predominantly originating from the rear of the source tube) compared with the spectrum of the unscattered electrons (predominantly originating from the front of the source tube). It enters into the model via equation (4), where the energy-loss function fs=1(ϵ) of singly scattered (s = 1) electrons is shifted by ΔP relative to the energy-loss function fs=0(ϵ) of unscattered (s = 0) electrons.

Both parameters are assessed with the help of co-circulating 83mKr gas, assuming that the possible plasma instabilities or longitudinal plasma profile are not affected by its presence in a minute concentration. The spectroscopy of its mono-energetic conversion electron lines reveals information about the broadening σP of the lineshape, from which an upper limit of ΔP is derived82.

The calibration was performed at an elevated source temperature of T = 80 K to prevent the condensation of Kr gas (compared with the T = 30 K set point used during the neutrino-mass measurements). The tritium circulation loop of the tritium source83 can be operated in two modes: (1) in a mode with the direct recycling of the krypton–tritium mixture, which is limited to a maximum column density of 40% (2.08 × 1017 cm−2), but delivers a high 83mKr rate; (2) in a mode of fractional direct recycling, which can be operated up to 75% (3.75 × 1017 cm−2) of the nominal column density, but with only about 0.5% of the maximum 83mKr activity compared with the mode described in (1).

The plasma parameters are inferred by combining the measurements of internal conversion lines N2,3-32 and L3-32. The energy of the conversion electrons, emitted from a particular subshell of the 83mKr atom, is 30,472.6 eV for the L3-32 singlet and 32,137.1 eV (32,137.8 eV) for the N2-32 (N3-32) doublet84. The L3-32 line, on one hand, has a high intensity, but its natural linewidth is not precisely known85. The N2,3-32 lines, on the other hand, have a low intensity, but their natural linewidth is negligible compared with spectral broadening caused by variations in the electric potential. Thus, any broadening of the N2,3-32 lines beyond the known spectrometer resolution and thermal Doppler broadening can be assigned to variations in the electric potential within the source. Based on the N2,3-32-line doublet measurement at 40% of the nominal column density, we find \({\sigma }_{{{{\rm{P}}}}}^{2}\)(40%) = (1.4 ± 0.3) × 10–3 eV2 (Extended Data Fig. 2). Combined with L3-32 line measurements at both 40% and 75% of the nominal column density, we assess the relative change in broadening and find \({\sigma }_{{{{\rm{P}}}}}^{2}\)(75%) = (8.0 ± 8.2) × 10–3 eV2. Finally, we conservatively extrapolate by exponentially scaling this value to 84% of the nominal column density, leading to \({\sigma }_{{{{\rm{P}}}}}^{2}\)(84%) = (12.4 ± 16.1) × 10–3 eV2. The plasma parameters and uncertainties deduced from the calibration measurements at 80 K cover the conditions during the neutrino-mass measurements at 30 K. The effect of different temperatures is smaller than the assumed uncertainty. For future physics runs, the tritium source will be operated at 80 K.

The contribution to broadening from slow drifts and shifts as determined from the run-wise analysis of the tritium spectra during KNM2 only yields a value of 10−3 eV2. This value is significantly smaller and can be determined with much higher precision than broadening \({\sigma }_{{{{\rm{P}}}}}^{2}\) obtained from the 83mKr measurement.

In the analysis, we limit broadening σP to positive values. We construct the probability density function of the asymmetry parameter ΔP based on phenomenological considerations on the relation of both plasma parameters given by ∣ΔP∣ ≤ σP/1.3 (ref. 82). For many samples of \({\sigma }_{{{{\rm{P}}}}}^{2}\), we draw a value for ΔP from a uniform distribution in the range of −σP/1.3 ≤ ΔP ≤ σP/1.3. The resulting distribution can be approximated by a Gaussian distribution centred around 0 meV with a width of 61 meV.

We expect to reduce the systematic uncertainties in future campaigns by operating the tritium source at the same temperature during the neutrino-mass and krypton measurements and by using an ultrahigh intensity 83mKr source (with about 10 GBq of the 83Rb parent radionuclide), which will make the N2,3-32 line measurement possible at a nominal column density (mode b).

Parameter inference

We infer the parameters of interest via a minimization of the χ2 function

$$\begin{array}{lll}{\chi }^{2}&=&{\left({{\bf{R}}}_{{{{\rm{data}}}}}(q{\bf{U}},{\bf{r}})-{\bf{R}}(q{\bf{U}},{\bf{r}}| {{\bf{\varTheta }}},{\bf{\upeta} })\right)}^{T}\\ &&\cdot {C}^{-1}\left({\bf{R}}_{{{{\rm{data}}}}}(q{\bf{U}},{\bf{r}})-{\bf{R}}(q{\bf{U}},{\bf{r}}| {{\bf{\varTheta }}},{\bf{\upeta} })\right),\end{array}$$
(7)

where Rdata(qU, r) gives the measured count rates at a retarding potential qUi for the detector ring rj, R(qU, r) gives the predictions of these rates and C is the covariance matrix that includes the statistical uncertainties and can also be used to describe systematic uncertainties. The usage of a χ2 minimization is justified, as the numbers of electrons per scan step qUi and detector ring rj are sufficiently large (>700) to be described by a Gaussian distribution instead of a Poisson distribution.

The fit has (1 + (3 × 12)) = 37 free parameters Θ, including a single parameter for the neutrino mass squared \({m}_{\nu }^{2}\), and 12 ring-wise values for each of the three spectral parameters, namely, normalization factor As, background rate Rbg and effective endpoint E0. In addition, the spectrum depends on systematic parameters η (Extended Data Table 1), such as the column density, tritium isotopologue concentrations, magnetic fields and so on. These parameters are known with a certain accuracy, and their uncertainty needs to be propagated to the final neutrino-mass result. Four different methods are used for the KATRIN analysis, as discussed below.

Pull method

In the ‘pull method’, systematic parameters ηi are treated as free parameters in the fit and introduced as nuisance terms in the χ2 function:

$${\chi }_{{{{\rm{tot}}}}}^{2}({{\bf{\varTheta }}},{\bf{\upeta} })={\chi }^{2}({{\bf{\varTheta }}},{\bf{\upeta} })+\mathop{\sum}\limits_{i}{\left(\frac{\hat{{\eta }_{i}}-{\eta }_{i}}{{\sigma }_{{\eta }_{i}}}\right)}^{2}\,.$$
(8)

The nuisance terms allow the parameter to vary around its best estimation \(\hat{{\eta }_{i}}\) according to its uncertainty \({\sigma }_{{\eta }_{i}}\) as determined from external measurements.

This method is computationally intensive due to the complexity in calculating the tritium spectrum and minimization with respect to multiple free parameters. For example, it is not practical to treat the uncertainties of the molecular final-state distribution, which is given as a discrete list of excitation energies and corresponding probabilities, with this method. The advantage of this method is that we make the maximum use of the data. If the spectral data contain information about the systematic parameters η, it is automatically taken into account.

Covariance matrix method

As shown in equation (7), the standard χ2 estimator includes covariance matrix C, which can describe both statistical and systematic model uncertainties. The diagonal entries describe uncertainties that are uncorrelated for each R(qUi, rj), whereas the off-diagonal terms describe the correlated uncertainties between R(qUi, rj). Covariance matrix C is computed by simulating \({{{\mathcal{O}}}}\)(104) β-decay spectra, with systematic parameters ηi varied according to their probability density functions in each spectrum. From the resulting set of spectra, the variance and covariance of the spectral points R(qUi, rj) are determined.

As the covariance matrices for individual or combined systematic effects are computed before fitting, this method is efficient with respect to computational costs. The dimension of the covariance matrix is given by the number of data points. Therefore, the efforts for matrix calculation and inversion can be diminished by reducing the number of data points.

Monte Carlo propagation method

Generally, in the Monte Carlo propagation technique, the uncertainties of parameters ηi are propagated by repeating the fit about 105 times, with the systematic parameters varied according to their probability density functions each time. The method then returns the distributions of the fit parameter Θ, which reflects the uncertainty of the systematic parameters.

More precisely, when assessing the systematic effects alone, a simulated spectrum without statistical fluctuations is created (based on the best-fit parameters), which is then fitted 105 times with a model whose systematic parameters of interest are varied each time. Contrarily, when evaluating the statistical uncertainty alone, 105 statistically randomized Monte Carlo spectra are created and fitted with a constant model. Accordingly, for obtaining the total uncertainty, both steps are combined.

To extract information on parameters η from the data, each entry in the resulting histogram of the best-fit parameters is weighted with the likelihood of the corresponding fit. The final distributions are then used to estimate the best-fit values (mode of distribution) and uncertainty (integration of distribution up to 16% from both sides).

The advantage of this method is that the number of free parameters is kept at a minimum, which facilitates the minimization procedure. The larger number of fits, however, is time consuming and requires the usage of large computing clusters.

Bayesian inference

In Bayesian inference, one computes the posterior probability for the parameters of interest Θ from a prior probability and a likelihood function according to Bayes’ theorem. For the KATRIN analysis, the Bayesian approach has the advantage that the prior knowledge of neutrino mass can naturally be applied via a corresponding informative prior. For example, the prior of \({m}_{\nu }^{2}\) can be chosen to be flat and positive, which restricts the posterior distribution to the physically allowed regime and accordingly changes the credible interval.

Ideally, all systematic effects ηi would be included as free parameters constrained with our prior knowledge on the parameter. However, due to the computationally expensive spectrum calculation and the fact that the Bayesian inference requires a large number of samples in the Markov chain for sampling the posterior distribution, only the qU-dependent background is currently treated in this way. All the other systematic uncertainties are included with a covariance matrix in the likelihood or by model variation. For the model variation, a large set of Markov chains is started with randomized but fixed model variations. The randomizations are drawn from the systematic uncertainty distributions. The resulting set of posterior distributions is averaged.

Limit setting

We present two frequentist methods and one Bayesian method for the construction of an upper limit of the neutrino mass. For the former, we use the classical Feldman–Cousins59 and Lokhov–Tkachov58 belt constructions (Extended Data Fig. 3). In the Feldman–Cousins technique, the acceptance region for \({\widehat{m}}_{\nu }^{2}\) is determined by ordering this estimator according to the likelihood ratio \(\frac{{{{\mathcal{L}}}}\left({\widehat{m}}_{\nu }^{2}\ | \ {m}_{\nu }^{2}\right)}{{{{\mathcal{L}}}}\left({\widehat{m}}_{\nu }^{2}\ | \ \max (0,\ {\widehat{m}}_{\nu }^{2})\right)}\) for a given best-fit neutrino mass squared \({m}_{\nu }^{2}\). This method leads to more stringent upper limits for increasingly negative best-fit values. The Lokhov–Tkachov method avoids this feature by using the standard Neyman belt construction for positive values of \({\widehat{m}}_{\nu }^{2}\) and defining the experimental sensitivity as the upper limit in the non-physical regime of \({\widehat{m}}_{\nu }^{2}\) < 0.

The Bayesian 90% credible interval is obtained by integrating the posterior distribution of \({m}_{\nu }^{2}\) from zero to \({{{\mbox{}}}{m}_{\nu }^{2}{{\mbox{}}}}^{{{{\rm{limit}}}}}\), such that the total probability is 90% (Extended Data Fig. 4). Note that the interpretation of the limit obtained in this way is different from the frequentist confidence limits, and hence, the numerical values may not coincide.

Results of individual strategies

The frequentist analyses are performed by three independent teams, which use differing implementations of spectral calculation and different strategies to propagate systematic uncertainties. A Bayesian analysis is performed, which interfaces to one of the spectrum calculation software. The analysis by independent teams is a powerful means to cross-check individual analyses. The following four data analysis strategies were applied to the second dataset of KATRIN.

  • Strategy 1 is based on a C++ framework (C++14 / Gnu compiler)86 using the Minuit minimizer. It mainly employs the pull method for handling systematics. We note that for the presented analysis, the input value for the ‘scan-step-duration-dependent background’ was corrected after the official unblinding of the data.

  • Strategy 2 is implemented in a MATLAB framework (R2020a) and exclusively uses the covariance matrix approach to propagate systematic uncertainties87. We note that for the presented analysis, the ‘scan-step-duration-dependent background’ systematic was implemented only after the official unblinding of the data.

  • Strategy 3 is based on a C++ framework (C++17 / Gnu compiler) using a custom-developed minimizer88,89. In this strategy, the Monte Carlo propagation of uncertainties is mostly applied.

  • Strategy 4 performs a Bayesian interpretation of the data. For this approach, the spectrum calculation software of ‘Strategy 3’ is interfaced with the Bayesian analysis toolkit90. Here most of the systematic uncertainties are treated via the model variation technique.

An overview of the strategies is found in Extended Data Table 2. The resulting best fit and systematic uncertainty breakdown are listed in Extended Data Table 3.