Search for Majorana neutrinos exploiting millikelvin cryogenics with CUORE – Nature

Optimum trigger and analysis threshold

The continuous data stream of CUORE is triggered with the optimum trigger, an algorithm based on the optimum filter⁵¹ characterized by a lower threshold than a more standard derivative trigger³². A lower threshold enables us not only to reconstruct the low-energy part of the spectrum, but also yields a higher efficiency in reconstructing the events in coincidence between different calorimeters, and consequently a more precise understanding of the corresponding background components^52,53.

The optimum trigger transfer function of every event is matched to the ideal signal shape, obtained as the average of good-quality pulses, so that frequency components with low signal-to-noise ratio are suppressed. A trigger is fired if the filtered signal amplitude exceeds a fixed multiple of the noise root mean square (RMS), defined separately for each calorimeter and dataset. We evaluate the energy threshold by injecting fake pulses of varying amplitude, calculated by inverting the calibration function, into the data stream. We reconstruct the stabilized amplitude of the fake pulses, fit the ratio of correctly triggered events to generated events with an error function, and use the 90% quantile as a figure of merit for the optimum trigger threshold. This approach enables monitoring of the threshold during data collection, and is based on the assumption that the signal shape is not energy dependent, that is, that the average pulse obtained from high-energy γ events is also a good template for events of a few keV. The distribution of energy threshold at 90% trigger efficiency is shown in Extended Data Fig. 4.

For this work we set a common analysis threshold of 40 keV, which results in >90% trigger efficiency for the majority (97%) of the calorimeters, while at the same time allowing the removal of multi-Compton events from the region of interest through the anti-coincidence cut.

Efficiencies

The total efficiency is the product of the reconstruction, anti-coincidence, pulse shape discrimination (PSD) and containment efficiencies.

The reconstruction efficiency is the probability that a signal event is triggered, has the energy properly reconstructed, and is not rejected by the basic quality cuts requiring a stable pre-trigger voltage and only a single pulse in the signal window. It is evaluated for each calorimeter using externally flagged heater events⁵⁴, which are a good approximation of signal-like events.

The anti-coincidence efficiency is the probability that a true single-crystal event correctly passes our anti-coincidence cut, instead of being wrongly vetoed owing to an accidental coincidence with an unrelated event. It is extracted as the acceptance of fully absorbed γ events at 1,460 keV from the electron capture decays of ⁴⁰K, which provide a reference sample of single-crystal events.

The PSD efficiency is obtained as the average acceptance of events in the ⁶⁰Co, ⁴⁰K and ²⁰⁸Tl γ peaks that already passed the base and anti-coincidence cuts. In principle, the PSD efficiency could be different for each calorimeter, but given the limited statistics in physics data we evaluate it over all channels and over the full dataset. To account for possible variation between individual calorimeters, we compare the PSD efficiency obtained by directly summing their individual spectra with that extracted from an exposure-weighted sum of the calorimeters’ spectra. We find an average ±0.3% discrepancy between the two values and include it as a global systematic uncertainty in the 0νββ fit. This takes a Gaussian prior instead of the uniform prior used in our previous result³², which had its uncertainty come from a discrepancy between two approaches that has since been resolved.

Finally, the containment efficiency is evaluated through Geant4-based Monte Carlo simulations⁵⁵ and accounts for the energy loss due to geometrical effects as well as bremsstrahlung.

Principal component analysis for PSD

In this analysis we use a new algorithm based on principal component analysis (PCA) for pulse shape discrimination. The method has been developed and documented for CUPID-Mo⁵⁶, and has been adapted for use in CUORE⁵⁷. This technique replaces the algorithm used in previous CUORE results, which was based on six pulse shape variables³⁰. The PCA decomposition of signal-like events pulled from γ calibration peaks yields a leading component similar to an average pulse, which on its own captures >90% of the variance between pulses. We choose to treat the average pulse of each calorimeter in a dataset as if it were the leading PCA component, normalizing it like a PCA eigenvector. We can then project any event from the same channel onto this vector and attempt to reconstruct the 10-s waveform using only this leading component. For any waveform x and leading PCA component w with length n, we define the reconstruction error as:

This reconstruction error metric measures how well an event waveform can be reconstructed using only the average pulse treated as a leading PCA component. Events that deviate from the typical expected shape of a signal waveform are poorly reconstructed and have a high reconstruction error. We normalize the reconstruction errors as a second-order polynomial function of energy on a calorimeter-dataset basis (see Extended Data Fig. 5), and cut on the normalized values by optimizing a figure of merit for signal efficiency over expected background in the Q_ββ region of interest. Using this PCA-based method, we obtain an overall efficiency of (96.4 ± 0.2)% compared to the (94.0 ± 0.2)% from the pulse shape analysis used in our previous results, as well as a 50% reduction in the PSD systematic uncertainty from 0.6% to 0.3%.

Statistical analysis

The high-level statistical 0νββ decay analysis consists of an unbinned Bayesian fit on the combined data developed with the BAT software package⁵⁸. The model parameters are the 0νββ decay rate (Γ_0ν), a linearly sloped background, and the ⁶⁰Co sum peak amplitude. Γ_0ν and the ⁶⁰Co rate are common to all datasets, with the ⁶⁰Co rate scaled by a preset dataset-dependent factor to account for its expected decay over time. The base background rate, expressed in terms of counts keV⁻¹ kg⁻¹ yr⁻¹, is dataset-dependent, whereas the linear slope to the background is shared among all datasets, because any structure to the shape of the background should not vary between datasets. Γ_0ν, the ⁶⁰Co rate, and the background rate parameters have uniform priors that are constrained to non-negative values, whereas the linear slope to the background has a uniform prior that allows both positive and negative values.

In addition to these statistical parameters, we consider the systematic effects induced by the uncertainty on the energy bias and energy resolution^59,60, the value of Q_ββ, the natural isotopic abundance of ¹³⁰Te, and the reconstruction, anti-coincidence, PSD and containment efficiencies. We evaluate their separate effects on the 0νββ rate by adding nuisance parameters to the fit one at a time and studying both the effect on the posterior global mode ({hat{{Gamma }}}_{0nu }) and the marginalized 90% CI limit on Γ_0ν.

A list of the systematics and priors is reported in Extended Data Table 1. The efficiencies and the isotopic abundance are multiplicative terms on our expected signal, so the effect of each is reported as a relative effect on Γ_0ν. By contrast, the uncertainties on Q_ββ, the energy bias, and the resolution scaling have a non-trivial relation to the signal rate; therefore, we report the absolute effect of each on Γ_0ν. The dominant effect is due to the uncertainty on the energy bias and resolution scaling in physics data. We account for possible correlations between the nuisance parameters by including all of them in the fit simultaneously.

We chose a uniform prior on our physical observable of interest Γ_0ν, which means we treat any number of signal events as equally likely. Other possible uninformative choices could be considered appropriate, as well. Because the result of any Bayesian analysis depends to some extent on the choice of the priors, we checked how other reasonable priors affect our result⁵⁷. We considered: a uniform prior on (sqrt{{{Gamma }}_{0nu }}), equivalent to a uniform prior on m_ββ and also equivalent to using the Jeffreys prior; a scale-invariant uniform prior on logΓ_0ν; and a uniform prior on 1/Γ_0ν, equivalent to a uniform prior on ({T}_{1/2}^{0nu }).

These priors are all undefined at Γ_0ν = 0, so we impose a lower cut-off of Γ_0ν > 10⁻²⁷ yr⁻¹, which with the given exposure corresponds to approximately one signal event. The case with a uniform prior on (sqrt{{{Gamma }}_{0nu }}) gives the smallest effect, and strengthens the limit by 25%, whereas the flat prior on 1/Γ_0ν provides the largest effect, increasing the limit on ({T}_{1/2}^{0nu }) by a factor of 4. In fact, all these priors weigh the small values of Γ_0ν more. Therefore, our choice of a flat prior on Γ_0ν provides the most conservative result.

We compute the 0νββ exclusion sensitivity by generating a set of 10⁴ toy experiments with the background-model, that is, including only the ⁶⁰Co and linear background components. The toys are split into 15 datasets with exposure and background rates obtained from the background-only fits to our actual data. We fit each toy with the signal-plus-background model, and extract the distribution of 90% CI limits, shown in Extended Data Fig. 4.

We perform the frequentist analysis using the Rolke method⁶¹, obtaining a lower limit on the process half-life of ({T}_{1/2}^{0nu } > 2.6times {10}^{25},{rm{yr}}) (90% CI). The profile likelihood function ℒ for Γ₀ν is retrieved from the full Markov chain produced by the Bayesian analysis tool. The non-uniform priors on the systematic effects in the Bayesian fit are thus incorporated into the frequentist result as well. We extract a 90% confidence interval on Γ_0ν by treating −2logℒ as an approximate χ² distribution with one degree of freedom. The lower limit on ({T}_{1/2}^{0nu }) comes from the corresponding upper edge of the confidence interval on Γ_0ν. Applying the same method to the toy experiments, we find a median exclusion sensitivity of ({T}_{1/2}^{0nu } > 2.9times {10}^{25},{rm{yr}}).