First direct neutrino-mass measurement with sub-eV sensitivity

We report the results of the second measurement campaign of the Karlsruhe Tritium Neutrino (KATRIN) experiment. KATRIN probes the effective electron anti-neutrino mass, $m_{\nu}$, via a high-precision measurement of the tritium $\beta$-decay spectrum close to its endpoint at $18.6\,\mathrm{keV}$. In the second physics run presented here, the source activity was increased by a factor of 3.8 and the background was reduced by $25\,\%$ with respect to the first campaign. A sensitivity on $m_{\nu}$ of $0.7\,\mathrm{eV/c^2}$ at $90\,\%$ confidence level (CL) was reached. This is the first sub-eV sensitivity from a direct neutrino-mass experiment. The best fit to the spectral data yields $m_{\nu}^2 = (0.26\pm0.34)\,\mathrm{eV^4/c^4}$, resulting in an upper limit of $m_{\nu}<0.9\,\mathrm{eV/c^2}$ ($90\,\%$ CL). By combining this result with the first neutrino mass campaign, we find an upper limit of $m_{\nu}<0.8\,\mathrm{eV/c^2}$ ($90\,\%$ CL).

We report the results of the second measurement campaign of the Karlsruhe Tritium Neutrino (KATRIN) experiment. KATRIN probes the effective electron anti-neutrino mass, mν , via a highprecision measurement of the tritium β-decay spectrum close to its endpoint at 18.6 keV. In the second physics run presented here, the source activity was increased by a factor of 3.8 and the background was reduced by 25 % with respect to the first campaign. A sensitivity on mν of 0.7 eV/c 2 at 90 % confidence level (CL) was reached. This is the first sub-eV sensitivity from a direct neutrinomass experiment. The best fit to the spectral data yields m 2 ν = (0.26 ± 0.34) eV 2 /c 4 , resulting in an upper limit of mν < 0.9 eV/c 2 (90 % CL). By combining this result with the first neutrino mass campaign, we find an upper limit of mν < 0.8 eV/c 2 (90 % CL).

I. INTRODUCTION
The discovery of neutrino flavour oscillations [1,2] proves that neutrinos must have a mass, unlike originally assumed in the Standard Model (SM) of Particle Physics. Neutrino oscillation experiments have shown that the weakly interacting neutrino flavour eigenstates, ν f with f ∈ {e, µ, τ}, are admixtures of the three neutrino mass eigenstates ν i with mass eigenvalues m i . While neutrino-oscillation experiments can probe the differences of squared neutrino mass eigenvalues ∆m 2 ij , 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 m 2 ν = i |U ei | 2 m 2 i where U ei are elements of the Pontecorvo-Maki-Nakagawa-Sakata (PMNS) matrix, which describes the mixing of the neutrino states.
The neutrino masses are at least five orders of magnitude smaller than the mass of any other fermion of the SM, which may point to a different underlying masscreation mechanism [3]. 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 universe [4,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 stringent limits on the sum of neutrino masses 1 of m i < 0.12 eV (95 % CL) [6,7]. However, these limits strongly rely on the underlying cosmological assumptions [8,9]. An independent measurement of the neutrino mass could help to break parameter degeneracies of the cosmological models. A powerful way to probe this neutrino property in the laboratory is via a search for neutrinoless double β-decay, providing limits at the level of m ββ < 0.08 −0.18 eV (90 % CL) [10,11], depending on the nuclear matrix element calculation. In contrast to m ν , the effective mass in double-beta decay is given by m ββ = i U 2 ei m i . The limit is only valid under the assumption that neutrinos are their own anti-particle (Majorana particle) and that light neutrinos mediate the decay.
The most direct way to assess the neutrino mass is via the kinematics of a single β-decay. This method is independent of any cosmological model and of the mass nature of the neutrino, i.e. it may be a lepton of Majorana or Dirac type. The neutrino masses m i lead to a reduction of the maximal observed energy of the decay and a small spectral shape distortion close to the kinematic endpoint of the β-spectrum. In the quasi-degenerate 1 In this paper we use the convention c = 1. mass regime, where m i > 0.2 eV, the mass splittings are negligible with respect to the masses, m i , and the observable can be approximated as m 2 ν = i |U ei | 2 m 2 i . The Karlsruhe Tritium Neutrino (KATRIN) experiment [12,13] exploits the single β-decay of molecular tritium T 2 → 3 HeT + + e − +ν 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) [14,15]. 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 1000 days [12]. In this work, we present the second neutrino-mass result of KATRIN, reaching an unprecedented sub-eV sensitivity and limit in m ν from a direct measurement.

II. EXPERIMENTAL SETUP
The design requirements to detect the small signature of a neutrino mass in the last few eV of a β-decay spectrum are: a high tritium activity (1 × 10 11 Bq), a low background rate ( 0.1 counts per second (cps)), an eVscale 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 type [13]. The experiment is hosted by the Tritium Laboratory Karlsruhe (TLK) allowing the safe supply of tritium at the 10-gramme scale with continuous tritium reprocessing [16,17]. Fig. 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 favourable properties to achieve high count rates near the endpoint. Moreover, the theoretical calculation of the super-allowed decay of tritium is well understood [25][26][27]. The strong gaseous tritium source (nominal activity 1 × 10 11 Bq) is established by the throughput of 40 g/day of molecular tritium through the 10-m long source tube, which is cooled to 30 K to reduce thermal motion of the tritium molecules.
A system of 24 superconducting magnets [28] guides the electrons out of the source towards the spectrometer section and to the detector. Between the source and spectrometer, differential [29] and cryogenic [30] 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 type [31,32]. The spectrometer acts as a sharp electro-static high-pass filter, transmitting only electrons (of charge q = −e) above an adjustable energy threshold qU , where U is the retarding potential applied at the spectrometer. A reduction of the magnetic field strength by about four orders of FIG. 1. Working principle of the KATRIN experiment as explained in the main text. The view into the tritium source depicts three systematic effects: a) molecular excitations during β-decay, b) scattering of electrons off the gas molecules, and c) spatial distribution of the electric potential in the source Usrc(r, z). The view into the spectrometer illustrates the main background processes. Generally, unstable neutral particles (d) can enter the large spectrometer volume and decay. The resulting electrons are accelerated by qUana towards the detector, making them indistinguishable from signal electrons [14,18]. The unstable neutral particles are either 219 Rn and 220 Rn emanated from the getter pumps [19,20] or highly excited atoms sputtered off the inner spectrometer surfaces [21][22][23]. The sputtering mechanisms originate from both α-decays of 210 Po in the structural material and a Penning trap (e), which gives rise to accelerated positive ions from the volume between the two spectrometers [24]. Note that the white pixels in the detector were excluded from analysis as they do not fulfil the 'golden pixel' quality criteria. The magnetic adiabatic collimation is only sketched for the case without electric field for illustration purposes. magnitudes from the entrance (exit) of the spectrometer B source = 2.5 T (B max = 4.2 T) to its center, the analysing plane (B ana = 6.3 × 10 −4 T), collimates the electron momenta. This configuration creates a narrow filter width of ∆E = 18.6 keV · (B ana /B max ) = 2.8 eV and allows for a large angular acceptance, with a maximum pitch-angle 2 for the β-decay electrons of θ max = arcsin (B source /B max ) = 50.4 • in the source. A 12-m diameter magnetic coil system surrounding the spectrometer fine-shapes the magnetic field and compensates the Earth's magnetic field [33,34].
The transmitted electrons are detected by a 148-pixel silicon PIN-diode focal-plane detector installed at the exit of the spectrometer [35]. 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, to reduce the flux of electrons into the main spectrometer. The up-stream 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 2 The pitch angle refers to the angle between the electron's momentum and the direction of the magnetic field at the position of the electron.
rear wall is biased to a voltage (O(100 mV)) to minimise the difference of its surface potential to that of the beam tube, which minimises the inhomogeneity of the source electric potential. The rear section is equipped with an angular-and energy-selective electron gun [36], 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 cm 2 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 the most precise energy loss function for this process to date [37]. Another key calibration source is gaseous krypton, which can be co-circulated with the tritium gas [38]. Monoenergetic conversion electrons from the decay of the metastable state 83m Kr are used to determine spatial and temporal variations of the electric potential in the tritium source. The variations are caused by a weak coldmagnetised plasma, which arises from the high magnetic field (2.5 T) and a large number of ions and low-energy electrons (≈ 1 × 10 12 m −3 ) in the tritium source. The methods of calibration are described in more detail in Sec. VII.
The beamline is equipped with multiple monitoring devices. A laser Raman system continuously monitors the gas composition, providing a measurement at the 0.1 %-precision level each minute. A silicon drift detector system, installed in the transport section, as well as a betainduced X-ray system at the rear section [39], continuously monitor the tritium activity, obtaining a result at the 0.1 %-precision level each minute. The high voltage of the main spectrometer is continuously measured at the ppm level with a high-precision voltage divider system [40,41] and an additional monitoring spectrometer [42]. The magnetic fields are determined with a highprecision magnetic field sensor system [43].
After the successful commissioning of the complete KATRIN beamline in summer 2017 [44], first tritium operation was demonstrated with a small tritium activity of (5 × 10 8 Bq) in mid-2018 [45]. During the first neutrino mass campaign (KNM1) in 2019 [14], the source was operated in a 'burn-in' configuration at a reduced activity of 2.5 × 10 10 Bq, which is required when structural materials are exposed to high amounts of tritium for the first time. Major technical achievements of the second measurement campaign (KNM2) are the operation of the tritium source at its nominal activity of 9.5 × 10 10 Bq and improved vacuum conditions in the spectrometer [46] that led to a reduction of the background by 25 % to 220 mcps. 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 × 10 6 β-electrons. Fig. 2 a) compares the spectra of both neutrino mass campaigns. A direct comparison of the experimental parameters is given in Tab. I.

III. MEASUREMENT OF THE TRITIUM BETA SPECTRUM
The integral β-spectrum is obtained by repeatedly measuring the count rate R data (qU i ) for a set of 39 non-equidistant high voltage (HV) settings U i , creating a retarding energy qU i for the electrons of charge q. The retarding energy is adjusted in a range of qU i ∈ (E 0 − 300 eV, E 0 + 135 eV), where E 0 = 18 574 eV is the approximate spectral endpoint. Note that for the spectral fit, only 28 of those points in the range of qU i ∈ [E 0 − 40 eV, E 0 + 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 HV set-point (called scan step) varies between 17 s and 576 s and is chosen to optimise the neutrino mass sensitivity. The total time to measure the rate at all 39 retarding energies (a so-called scan) is about 2 h. As can be seen in Fig.  2  with a reproducibility at the sub-ppm-level, these scans are later combined to a single spectrum by adding the counts at a given set point.
The focal-plane detector, shown as magnified inset in Fig. 1, is segmented into 148 individual pixels to account for spatial variations of the electromagnetic fields inside the source and spectrometer, and of the background. Removing malfunctioning pixels and those shadowed by hardware upstream, 117 pixels have been selected. For the presented analysis these golden pixels are grouped into 12 concentric rings, resulting in 336  [14] and KNM2 to this work. The total number of β-electrons is counted in the last 40 eV interval of the integral spectrum, which is used for the spectral fit. The β-electron-to-background ratio is given by the ratio of this number and the background counts in the same energy range, i. e. E0 − 40 eV to E0. data points R data (qU i , r j ) for i ∈ {1, . . . , n qU = 28} and j ∈ {1, . . . , n rings = 12} per scan. Fig. 2 c) displays the spectra for each of the 12 detector rings, where all scans have been combined.

IV. DATA ANALYSIS
In order 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.

A. Spectrum calculation
The prediction of the detection rate R(qU i , r j ) is given by a convolution of the differential β-decay spectrum R β (E) with the experimental response function f (E, qU i , r j ), and a background rate R bg (qU i , r j ): Here N T is the signal normalisation calculated from the number of tritium atoms in the source, the maximum acceptance angle, and the detection efficiency. A s ≈ 1 is an additional normalisation 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 distribution (FSD) [47,48] in the differential β-decay spectrum. The Doppler broadening due to the finite thermal motion the tritium molecules in the source, as well as energy broadenings due to spatial and temporal variations of the spectrometer and source electric potential, are emulated by Gaussian broadenings of the FSD.
The response function f (E, qU i , r j ) includes the energy losses due to scattering and synchrotron radiation in the high B-field regions, as well as the transmission of electrons through the main spectrometer. Compared to previous KATRIN analyses, we now use a modified transmission function to account for the non-isotropy of the β-electrons leaving the source 3 . A detailed description of the spectrum calculation can be found in Ref. [49] and in section VII.

B. Unbiased parameter inference
We infer m 2 ν by fitting the experimental spectrum R data (qU i , r j ) with the prediction R(qU i , r j ) by minimising the standard χ 2 = R data C −1 R function, where C contains the variance and co-variance of the data points. In addition to the neutrino mass squared, m 2 ν , the parameters A s (r j ), R bg (r j ), and the effective endpoint E 0 (r j ) are treated as independent parameters for the 12 detector rings, leading to a total number of free parameters of 1 + 3 × 12 = 37 in the fit. The introduction of ringdependent parameters was chosen to allow for possible unaccounted radial effects. In particular, the effective endpoint E 0 (r j ) could account for a possible radial dependence of electric potential in the source. However, the final fit revealed a negligible (<100 meV) radial variation of the endpoint. Another advantage of ring-dependent parameters is to avoid the averaging of the transmission function over all rings, which would introduce an energy broadening, and hence reduce the resolution.
The following analysis procedure was implemented to minimise the potential for human-induced biases. The full analysis is first performed on a Monte-Carlo copy of each scan, simulated based on the true experimental parameters and with the neutrino mass set to zero. After all systematic inputs (e.g. the magnetic field values and uncertainties) are quantified, the fit is performed on the experimental data set, but with a randomly broadened molecular final-state distribution, which imposes an unknown bias on the observable m 2 ν . Only after three independent analysis teams, using different software and strategies, obtained consistent results would the analysis of the true data with unmodified final-state distribution be performed.

C. Systematic effects
The analytical description R(qU i , r j ) of the integral β-spectrum contains various experimental and theoretical parameters, such as the magnetic fields, the column density, and the 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 Sec. VII.
The neutrino mass result presented here is dominated by the statistical uncertainty. The largest systematic uncertainties are related to background properties and the source electric potential. First, radon decays in the main spectrometer lead to a non-Poissonian background rate over-dispersion [20] of about 11 %, effectively increasing the statistical uncertainty. Secondly, background generation mechanisms may show a retarding-potential dependence of the background, parametrize by a slope (m qU bg = (0.00 ± 4.74) mcps/keV.) . Thirdly, a removal of stored electrons from a known Penning trap between the spectrometers [24] after each scan step can lead to a slowly increasing background rate (m tscan bg = (3 ± 3) µcps/s) and thus to a scan-step-duration-dependent background contribution. Finally, spatial and temporal variations of 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 systematic effects on the neutrino mass is listed in Tab. II and described in detail in Sec. VII.

V. RESULTS AND DISCUSSION
The χ 2 minimisation 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 squared neutrino mass we find m 2 ν = 0.26 +0.34 −0.34 eV 2 , see Fig. 3. The independent analysis methods agree within about 5 % percent of the total uncertainty. The total uncertainty on the fit is dominated by the statistical error followed by uncertainties of background parameters and the source electric potential. The full breakdown of uncertainties can be found in Tab. II.
Based on the best-fit result we obtain an upper limit of m ν < 0.9 eV at 90 % CL, using the Lokhov-Tkachov method [50]. The Feldman-Cousins technique [51] 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 a ≈ 0.8σ statistical fluctuation assuming a true neutrino mass of zero. The best fit is defined as the 1D weighted median of the distribution. The combined endpoint is calculated using the correlated mean of the 12 ring-wise endpoints per sample.
We also perform a Bayesian analysis of the data set, with a positive flat prior on m 2 ν . The resulting Bayesian limit at 90 % credibility is m ν < 0.85 eV.
The ring-averaged, fitted effective endpoint is E 0 = (18 573.69 ± 0.03) eV. Taking into account the centerof-mass molecular recoil of T 2 (1.72 eV), as well as the absolute electric source potential φ WGTS (σ(φ WGTS ) = 0.6 V) and the work function of the main spectrometer φ MS (σ(φ MS ) = 0.06 eV), we find a Q-value of (18 575.2 ± 0.6) eV, which is consistent with the previous KATRIN neutrino-mass campaign ((18 575.2 ± 0.5) eV [45]) and the calculated Q-value from the 3 He− 3 H atomic mass difference of (18 575.72 ± 0.07) eV [52]. 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) result [45]. A simultaneous fit of both data sets yields m 2 ν = (0.1 ± 0.3) eV 2 and a corresponding upper limit of m ν < 0.8 eV at 90 % CL, based on Lokhov-Tkachov or the Feldman-Cousins technique. The same result is obtained when multiplying the m 2 ν distributions from Monte-Carlo propagation, or adding the χ 2 profile of the individual fits. As both data sets are statisticsdominated, correlated systematic uncertainties between both campaigns are negligible. Furthermore, we investigated a Bayesian combination, where the KNM1 posterior distribution of m 2 ν is used as prior for the KNM2 analysis, yielding consistent results. More details on the combined analyses can be found in the supplementary material.

VI. CONCLUSION AND OUTLOOK
The second neutrino-mass measurement campaign of KATRIN, presented here, is the first direct neutrinomass measurement reaching sub-eV 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 quasidegenerate neutrino-mass models and we have provided model-independent information about the neutrino mass, which allows the testing of non-standard cosmological models [64,65]. Fig. 4 displays the evolution of best-fit m 2 ν results from historical neutrino-mass measurements up to the present day.
Compared to 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 1000 days, the total statistics of KATRIN will be increased by another factor of 50. A reduction of the background rate by a factor of two will be achieved by an optimised electromagnetic field design of the spectrometer section. Moreover, by eliminating the radonand Penning-trap-induced background, the backgroundrelated systematic effects are expected to be significantly reduced. Together with a new high-rate krypton calibration scheme and a new method to improve the determination of the magnetic fields we will minimise the dominant systematic uncertainties to reach the target sensitivity of m 2 ν 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 neutrinos [66]. Together with cosmological probes and searches for neutrinoless double β-decay [67], the upcoming KATRIN data will play a key role in measuring the neutrino mass parameters.
In this section we describe the data analysis chain starting from the data processing to the high-level fit and limit setting. Moreover, we provide details on one of the key calibration campaigns, concerning the source electric potential. For a more extensive description of the KATRIN analysis procedure, the reader is referred to a recent publication [15].
A. Data processing, selection, and combination The first step of the analysis chain is the preparation of the data. Raw data are combined into integral spectral data points, which are then fitted with an analytical spectrum prediction including the response of the experiment.
a. Rate determination The electrons which are transmitted through the main spectrometer are further accelerated by a post-acceleration electrode (PAE) with the potential U PAE = 10 keV before they are detected by the focal-plane detector [35]. The latter provides a high detection efficiency (>95 %) and a moderate energy resolution (2.8 keV FWHM at 28 keV). The total rate per pixel at a given retarding potential qU i is determined by integrating the rate in a wide and asymmetric region of interest (ROI) of 14 keV ≤ E ≤ 32 keV. The asymmetric ROI is chosen to account for energy losses in the deadlayer of the Si-PIN diodes [68], for partial energy deposition due to backscattering of electrons off the detector surface and due to charge sharing between pixels.
The detector efficiency slightly depends on qU i . Three effects are considered: 1) The differential energy spectrum of the transmitted electrons at the FPD is shifted (and slightly scaled) according to qU i . Since the same ROI is used for each qU i , 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 of detection efficiency at E 0 − 1 keV, with respect to the efficiency at E 0 , is determined to be δ ROI = 0.002, with an uncertainty of 0.16 %. The data is corrected according to the efficiency at the given qU i value. 2) As the counting rate at the focal plane detector depends on qU i , so does the probability of pile-up. As the energy of most pile-up events is added, they are not covered by the ROI, thereby effectively changing the detector efficiency with qU i . Based on a random-coincidence model assuming a Poisson-distributed signal, the relative reduction of efficiency at E 0 − 1 keV, with respect to the efficiency at E 0 , is estimated to be δ PU = 0.0002, with an uncertainty of 18 %. The data is corrected accordingly. 3) Finally, a significant fraction of about 20 % of all electrons impinging on the detector surface are backscattered. For low retarding potentials and small energy depositions in the detector, these 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 qU i , 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 E 0 − 1 keV. We neglect this effect in this analysis. Systematic uncertainties related to the efficiency corrections are negligibly small for the presented analysis, which only considers the last 40 eV of the tritium spectrum.
b. Selection of golden pixels and scans The 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, higher noise levels, or decreased rate due to misalignment of the beam-line with respect to the magnetic flux tube.
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 and downtime of the Laser Raman system. The golden pixel and run selection is performed prior to the unblinding of the data.
c. 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 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 n data = n qU × n rings = 336.

B. Spectrum calculation
In order to infer the neutrino mass, the integral spectrum is fitted with the analytical spectrum calculation including the experimental response. As can be seen in Eq. (2), the theoretical spectrum is composed of the differential tritium spectrum R β (E) and the experimental response function f (E, qU i ). The differential spectrum is given by with G F denoting the Fermi constant, cos 2 Θ C the Cabibbo angle, |M nucl | 2 the energy-independent nuclear matrix element, and F (E, Z = 2) the Fermi function.
where E 0 denotes the maximum kinetic energy of the electron, in case of zero neutrino mass, and V f describe the molecular excitation energies, which are populated with the probabilities ζ f . E and m e denote the kinetic energy and mass of the βelectron, respectively.
Beyond the molecular effects, further theoretical corrections arise on the atomic and nuclear level [69]. Relevant for this analysis are only the radiative corrections to the differential spectrum, which are included in the analytical description, but not shown here. 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.
The experimental response function depicts the probability of an electron with a starting energy E to reach the detector. It combines the transmission function T (E − , θ, U ) of the main spectrometer and the 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 P s (θ), which depend on the path length through the source and hence on the pitch angle θ, and the energy-loss function f s ( ) for a given number of scatterings s. The energy-loss function is determined experimentally with the electron gun installed at the rear of the source. A typical response function and corresponding energy loss function is shown in Fig. 5.
The integrated transmission function for an isotropic source of electrons is given by It is governed by the magnetic fields at the starting position B source of the electron, the maximum field B max in the beam line, and the magnetic field in the spectrometer's analysing plane B ana . Synchrotron energy losses 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).

C. Systematic uncertainties
The analytical description R(qU i , r j ) of the integral βspectrum, shown in Eq. (2), contains various signal-and 14 − 32 keV -Detector efficiency qU reduction (ROI) 2 × 10 −3 (1 keV below E0) 0.16 % Detector efficiency qU reduction (pile up) 2 × 10 −4 (1 keV below E0) 18 % 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.
a. Signal-related systematic effects The spectrum prediction includes uncertainties of the nine parameters of the empirical energy loss function (the individual relative uncertainties of the parameters σ eloss,k are between 0.016 % and 3.8 %), a relative uncertainty of the product of column density and scattering cross-section (σ ρd×σ = 0.25 %), relative uncertainties of the theoretical description of the molecular final-state distribution (σ FSD = O(1 %)), and relative uncertainties of the magnetic field in the source (σ Bsource = 1.7 %), in the analysing plane (σ Bana = 1 %), and the maximal field (σ Bmax = 0.1 %). The variations of the β-decay activity during a scan (σ scan = 0.03 %) were negligibly small.
The spatial and temporal variations of the source electric potential were not included in the first neutrino mass campaign. With the increase of the source column density from 1.11 × 10 17 cm −2 to 4.23 × 10 17 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 relevant [12,70,71]. A detailed description of the plasma calibration is given in Sec. VII D.
b. Background-related systematic effects Electrons created in a radon decay in the main spectrometer volume can initially be magnetically trapped. These trapped electrons can create a cluster of 10 -100 secondary electrons [20] by scattering off the residual gas in the main spectrometer, which is operated at a pressure of 10 −11 mbar. These secondary electrons arrive at the detector within a time window of about 1000 s, hence leading to a non-Poissonian rate distribution. The observed background rate can be modelled by a Gaussian distribution, with a width 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 depen- dence of the background on the retarding potential and use dedicated test measurements to constrain the possible slope to m qU bg = (0.00 ± 4.74) mcps/keV. Finally, the pre-and main spectrometers, being both operated at 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 system [72] after each scan step. However, a potentially small increase of 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 scan steps we find a slope of m tscan bg = (3 ± 3) µcps/s, which is included in the βspectrum fit 4 .
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 of the effective endpoint mostly leads to a shift of the spectrum and has a negligible effect on the spectral shape. Accordingly, an unknown radial variation of the electric potential is in good approximation absorbed by the ring-wise endpoint parameters. Consequently, these effect have a minor impact on the neutrino mass analysis. However, any temporal or longitudinal variations of the source potential can lead to spectral distortions, which are parameterised by a Gaussian broadening σ P and the parameter ∆ P which quantifies the longitudinal rear-tofront asymmetry of the electric potential of the source. This asymmetry of the potential results in a shift of the energy spectrum associated with scattered electrons (predominantly originating from the rear of the source tube) compared to the spectrum of the unscattered electrons (predominantly originating from the front of the source tube).
Both parameters are assessed with the help of cocirculating 83m Kr gas. The spectroscopy of its monoenergetic conversion electron lines reveals information about the broadening σ P of the line shape, from which an upper limit of ∆ P is derived.
The calibration was performed at an elevated source temperature of T = 80 K to prevent condensation of the Kr gas (as compared to the T = 30 K set-point used during neutrino-mass measurements). The tritium circulation loop of the tritium source [73] can be operated in two modes, a) in a mode with direct-recycling of the krypton-tritium mixture which is limited to a maximum column density of 40 % (2.08 × 10 17 cm −2 ), but delivers a high krypton rate, and b) in a mode of fractional direct-recycling which can be operated up to 75 % (3.75 × 10 17 cm −2 ) of the nominal column density, but with only about 0.5 % of the maximum krypton activity as compared to mode a).
The plasma parameters are inferred by combining the measurements of the internal conversion lines N 2,3 -32 and L 3 -32. The energy of the conversion electrons, emitted from a particular subshell of the 83m Kr atom, is 30 472.6 eV for the L 3 -32 singlet and 32 137.1 eV (32 137.8 eV) for the N 2,3 -32 doublet [74]. The L 3 -32line, on the one hand, has a high intensity, but its natural line width is not known precisely [75]. The N 2,3 -32-lines, on the other hand, have a low intensity, but their natural line width is negligible compared to a spectral broadening caused by variations of the electric potential. Thus, any broadening of the N 2,3 -32-lines beyond the known spectrometer resolution and the thermal Doppler broadening can be assigned to variations of the electric potential within the source. Based on the (See Fig. 6). Combined with L 3 -32-line measurements both at 40 % and 75 % of the nominal column density, we assess the relative change of the broadening and find σ 2 P (75 %) = (8.0 ± 8.2) 10 −3 eV 2 . Finally, we conservatively extrapolate with exponential scaling this value to 84 % of the nominal column density, leading to σ 2 P (84 %) = (12.4 ± 16.1) 10 −3 eV 2 . 5 In the analysis, we limit the 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 [76]. For many samples of σ 2 P we draw a value for ∆ P from a uniform distribution in the range −σ 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 ultra-high intensity krypton source (with about 10 GBq of the 83 Rb mother isotope), which will make the N 2,3 -32-line measurement possible at nominal column density (mode b). 5 The plasma parameters and uncertainties based on the calibration measurements at 80 K cover the conditions during the neutrino mass measurements at 30 K. The effect by the different temperatures is smaller than the assumed uncertainty. For future campaigns the tritium source is operated at 80 K.

E. Parameter inference
We infer the parameters of interest via a minimisation of the χ 2 function where R data (q U , r) gives the measured count rates at a retarding potential qU i for the detector ring r j , R(q U , r) gives the predictions of these rate, and C is the covariance matrix that includes the statistical uncertainties and can be used to also describe systematic uncertainties. The usage of a χ 2 minimisation is justified, as the numbers of electrons per scan-step qU i and detector ring r j 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 2 ν , and 12 ring-wise values for each of the three spectral parameters: the normalisation factor A s , the background rate R bg , and the effective endpoint E 0 . In addition, the spectrum depends on systematic parameters η (as found in Tab. III), such as the column density, tritium isotopologue concentrations, magnetic fields, etc. 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 KA-TRIN analysis: a. Pull method In the 'pull method' the systematic parameters η i are treated as free parameters in the fit and introduced as nuisance terms in the χ 2 function The nuisance terms allow the parameter to vary around its best estimationη i according to its uncertainty σ ηi as determined from external measurements This method is computationally intensive due to the complexity in calculating the tritium spectrum and the minimisation 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.
b. Covariance matrix method As can be seen in Eq. (VII E), the standard χ 2 estimator includes a covariance matrix C which can describe both the statistical and systematic model uncertainties. The diagonal entries describe the uncertainties which are uncorrelated for each R(qU i , r j ), while the off-diagonal terms describe the correlated uncertainties between the R(qU i , r j ). The covariance matrix C is computed by simulating O(10 4 ) β-spectra, with the 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(qU i , r j ) 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.
c. Monte-Carlo propagation method Generally, in the Monte-Carlo (MC) propagation technique, the uncertainties of the parameters η i are propagated by repeating the fit about 10 5 times, with the systematic parameters varied according to their probability density functions each time. The method then returns the distributions of the fit parameters Θ, which reflect the uncertainty of the systematic parameters.
More precisely, when assessing the systematic effects alone, a MC spectrum without statistical fluctuations is created (based on the best fit parameters), which is then fitted 10 5 times with a model, which systematic parameters of interest are varied each time. Contrarily, when evaluating the statistical uncertainty alone, 10 5 statistically randomized MC spectra are created and fitted with a constant model. Accordingly, for obtaining the total uncertainty, both steps are combined.
To extract information on the parameters η from the data itself, each entry in the resulting histogram of 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 the distribution) and uncertainty (integration of the distribution up to 16 % from both side).
The advantage of this method is that the number of free parameters is kept at a minimum, which facilitates the minimisation procedure. The larger number of fits, however, is time consuming and requires the usage of large computing clusters.
d. 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 prior knowledge on the neutrino mass can naturally be applied via a corresponding informative prior. For example, the prior of m 2 ν can be chosen to be flat and positive, which restricts the posterior distribution to the physically allowed regime and changes the credibility interval accordingly.
Ideally, all systematic effects η i would be included as free parameters constrained with our prior knowledge on the parameter. However, due to the computation-ally 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 other systematic uncertainties are included with a covariance matrix in the likelihood or by model variation (MV). For the MV, a large set of Markov Chains is started with randomized but fixed model variations. The randomisations are drawn from the systematic uncertainty distributions. The resulting set of posterior distributions is averaged.

F. 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-Cousins [51] and the Lokhov-Tkachov [50] belt constructions depicted in Fig. 7. In the Feldman-Cousins technique the acceptance region for m 2 ν is determined by ordering this estimator according to the likelihood ra- ) for a given best-fit neutrino mass squared m 2 ν . This method leads to more stringent upper limits for increasingly negative best-fit values. The method of Lokhov-Tkachov avoids this feature by using the standard Neyman belt-construction for positive values of m 2 ν and defining the experimental sensitivity as the upper limit in the non-physical regime of m 2 ν < 0. The Bayesian 90 % credibility interval is obtained by integrating the posterior distribution of m 2 ν from zero to m 2 ν limit , such that the total probability is 90 % as demonstrated in Fig. 8. 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.

G. Results of individual strategies
The frequentist analyses are performed by three independent teams, which use differing implementations of the spectral calculation and different strategies to propagate systematic uncertainties. A Bayesian analysis is performed, which interfaces to one of the spectrum calculation softwares. 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 data set of KATRIN.
• Strategy 1 is based on a C++ framework [77] using the Minuit minimiser. It employs mainly the pull method for handling of 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 and exclusively uses the covariance matrix approach to propagate systematic uncertainties [78]. 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 using a custom-developed minimiser [79,80]. In this strategy mostly the MC-propagation of uncertainties is 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 Toolkit (BAT) [81]. Here, most of the systematic uncertainties are treated via the model variation (MV) technique.
An overview of the strategies is found in Tab. IV. The resulting best fit and systematic uncertainty breakdown are listed in Tab. V.

VIII. COMBINED ANALYSIS OF FIRST AND SECOND NEUTRINO MASS CAMPAIGN
The 1000 days of total measurement time of the KA-TRIN experiment will be distributed within 15-20 individual measurement campaigns. The experimental parameters will not be constant over all these campaigns, as the systematic uncertainties are reduced and the experimental performance is improved (e.g. the background rate will be lowered) with time. Nevertheless, all individual campaigns need to be combined eventually to achieve the full statistical power.

A. Combination of individual results
One approach is to simply multiply the m 2 ν distributions obtained with MC-propagation of the individual measurement phases. As the uncertainties of the first and second neutrino mass campaigns are still largely statistically dominated and the main systematic uncertainties are uncorrelated, this method is well justified for the combination presented here. The best fit, obtained from the combined distribution is m 2 ν = (0.11 ± 0.34) eV 2 . Assuming a symmetric Gaussian distribution of m 2 ν , we derive an upper limit of m ν < 0.81 eV ( 90 % C.L). Adding the χ 2 -profiles of the individual fits leads to statistically consistent results.

B. Simultaneous fit of the both neutrino mass campaigns
Another way, to combine both data sets, is by performing a simultaneous fit of both data sets with a common neutrino mass. The first and second neutrino mass campaigns were performed under different experimental conditions, e.g. the column density was increased by a factor of four and the background level was reduced by 25 % in the second campaign. Correspondingly, the systematic parameters and their uncertainties are partly different in the two campaigns. As a consequence, the data of the two campaigns cannot be described with a single effective spectrum prediction. Nevertheless, a simultaneous fit of both data sets can be performed by extending the χ 2 function to where Θ p1 and Θ p2 depict the endpoint, background, and normalisation of the first (p1) and second (p2) measurement period. These parameters can be different in the two measurement phases, as the experimental conditions are changed. η p1 and η p2 correspond to the two 40 20 0  sets of systematic parameters, while η c illustrates the set of common systematic parameters. The neutrino mass squared m 2 ν is a common parameter of both χ 2 functions. For a simultaneous fit the number of free parameters and data points increases by about a factor of two compared to the fit of a single measurement campaign, which makes the minimisation computationally challenging and prone to numerical noise. To simplify the problem, we choose for this analysis to group the detector rings to one effective detector area (uniform fit), reducing the number of free parameters to n free = 33 and the number of data points to n data = 56.
With respect to the stand-alone first neutrino mass analysis as presented in [15], the spectrum calculation has been slightly improved. It now includes the transmission function for non-isotropic electrons, an improved parametrisation of the energy loss function, and a reduced uncertainty on the magnetic fields and column density. Finally, a possible Penning-trap induced timedependent background, as described in the main text, is included. A re-analysis of the first neutrino mass data with these new inputs reveals a best fit of m 2 ν = −1.0 +0.9 ν . The prior distribution for the KNM2 analysis (green) [15] is the posterior distribution from the KNM1 analysis. The resulting posterior distribution for the second campaign is shown in blue.
C.L. These values are in good (10 %) agreement with respect to the published result in [14].
The fit to both data sets is shown in Fig. 9. The χ 2 /ndof = 51.7/49 = 1.1 indicates excellent agreement of the model to the data. The final result of the combined fit reveals m 2 ν = (0.08 ± 0.32) eV 2 and a corresponding upper limit of m ν < 0.76 eV (90 % CL).
In complement to the standard goodness-of-fit, we perform the Parameter-Goodness-of-Fit (PGoF) test [82] to assess the compatibility of KNM1 and KNM2 data sets. This test quantifies the penalty of combining the data sets compared to fitting them independently. We find a PGoF probability of 17 % reflecting a good agreement between the two statistically independent data sets.

C. Bayesian combination of both neutrino mass campaigns
Another way of combining different neutrino mass measurement campaigns is by using the posterior distribution of one campaign as prior information for the other campaign. In this procedure, correlations between the data sets are neglected.
For the analysis of the first neutrino mass measurement campaign, we use the result published in [72], which includes the five major systematic effects. The resulting posterior distribution for m 2 ν , displayed in Fig. 8, is then used as prior distribution of m 2 ν for the analysis of the second neutrino mass measurement campaign, which is otherwise performed with the same procedure as the stand-alone analysis.
The corresponding 90 % credible interval from a positive prior on m 2 ν corresponds to a Bayesian limit of m ν < 0.72 eV, as illustrated in Fig. 10.