The fate of hints: updated global analysis of three-flavor neutrino oscillations

Our herein described combined analysis of the latest neutrino oscillation data presented at the Neutrino2020 conference shows that previous hints for the neutrino mass ordering have significantly decreased, and normal ordering (NO) is favored only at the $1.6\sigma$ level. Combined with the $\chi^2$ map provided by Super-Kamiokande for their atmospheric neutrino data analysis the hint for NO is at $2.7\sigma$. The CP conserving value $\delta_\text{CP} = 180^\circ$ is within $0.6\sigma$ of the global best fit point. Only if we restrict to inverted mass ordering, CP violation is favored at the $\sim 3\sigma$ level. We discuss the origin of these results - which are driven by the new data from the T2K and NOvA long-baseline experiments -, and the relevance of the LBL-reactor oscillation frequency complementarity. The previous $2.2\sigma$ tension in $\Delta m^2_{21}$ preferred by KamLAND and solar experiments is also reduced to the $1.1\sigma$ level after the inclusion of the latest Super-Kamiokande solar neutrino results. Finally we present updated allowed ranges for the oscillation parameters and for the leptonic Jarlskog determinant from the global analysis.


Introduction
Global fits to neutrino oscillation data in the last several years have shown persistent hints for the normal neutrino mass ordering and values of the CP phase δ CP around maximal CP violation [1][2][3][4][5][6]. In this article we are going to re-assess the status of those hints in light of the new data released at the Neutrino2020 conference, in particular by the T2K [7,8] and NOvA [9,10] long-baseline (LBL) experiments. As we are going to discuss in detail, the hints have mostly disappeared or are significantly decreased: both neutrino mass orderings provide fits of comparable quality to the global data from accelerator and reactor experiments, and the CP conserving value δ CP = 180 • is within the 1σ allowed range.
We discuss in detail the origin of this apparent change of trends and trace back the data samples responsible for the change. We are going to compare the latest status with our pre-Neutrino2020 analysis, NuFIT 4.1, available at the NuFIT website [11]. Most relevant for mass ordering and CP phase are the updates of the neutrino samples for T2K [8], from 1.49 to 1.97 × 10 21 POT, and NOvA [10], from 0.885 to 1.36 × 10 21 POT. The T2K and NOvA anti-neutrino exposures are the same as used for NuFIT 4.1, but both collaborations introduced relevant changes in their analysis and hence we have adapted also our antineutrino fits correspondingly. In addition we have updated the reactor experiments Double-Chooz [12,13] from 818/258 to 1276/587 days of far/near detector data and RENO [14,15] from 2200 to 2908 days of exposure.
Another update concerns the solar neutrino oscillation analysis, to include the latest total energy spectrum and the day-night asymmetry of the SK4 2970-day sample presented at Neutrino2020 [16]. As we will show, thanks to these new data the tension on the determination of ∆m 2 21 from KamLAND versus solar experiments has basically disappeared. The outline of the paper is as follows. In Sec. 2 we discuss the status of the neutrino mass ordering and the leptonic CP phase δ CP , focusing on recent updates from T2K, NOvA, as well as the combination of LBL accelerator and reactor experiments. Despite somewhat different tendencies, we will show quantitatively that results from T2K and NOvA as well as reactors are fully statistically compatible. The status of the tension between solar and KamLAND results is presented in Sec. 3. Section 4 contains a selection of the combined results of this global fit, NuFIT 5.0, which updates our previous analyses [1,2,17,18]. In particular we present the ranges of allowed values for the oscillation parameters and of the leptonic Jarlskog determinant. 1 Parametrization conventions and technical details on our global analysis can be found in Ref. [2]. In particular, in what follows we use the definition ∆m 2 3 with = 1 for ∆m 2 3 > 0: normal ordering (NO), = 2 for ∆m 2 3 < 0: inverted ordering (IO). (1.1) We finish by summarizing our results in Sec. 5. A full list of the data used in this analysis is given in appendix A.
2 Fading hints for CP violation and neutrino mass ordering

T2K and NOvA updates
We start by discussing the implications of the latest data from the T2K and NOvA longbaseline accelerator experiments, presented at the Neutrino2020 conference. 2 To obtain a qualitative understanding we follow Refs. [2,20] and expand the oscillation probability relevant for the T2K and NOvA appearance channels in the small parameters sin θ 13 , ∆m 2 21 L/E ν , and A ≡ |2E ν V /∆m 2 3 |, where L is the baseline, E ν the neutrino energy and V the effective matter potential [21]: and we have used |∆m 2 3 | L/4E ν ≈ π/2 for T2K and NOvA. At T2K, the mean neutrino energy gives A ≈ 0.05, whereas for NOvA we find that with the empirical value of A = 0.1 1 Additional figures, ∆χ 2 maps and future updates of this analysis will be made available at the NuFIT website [11].
2 During the preparation of this work Ref. [19] appeared presenting related partial results in qualitative agreement with some of our findings for the LBL analysis. the approximation works best. Assuming that the number of observed appearance events in T2K and NOvA is approximately proportional to the oscillation probability we obtain Taking all the well-determined parameters θ 13 , θ 12 , ∆m 2 21 , |∆m 2 3 | at their global best fit points, we obtain numerically C ≈ 0.28 with negligible dependence on θ 23 . The normalization constants N ν,ν calculated from our re-analysis of T2K and NOvA are given for the various appearance samples in table 1. Hence the expression in eqs. (2.4) and (2.5) serve well to understand the main behaviour under varying the parameters sin 2 θ 23 , δ CP , and the mass ordering.
In table 1 we also show the observed number of events, background subtracted events, as well as the ratio r = (N obs − N bck )/N ν(ν) . In a fit, the values of r have to be accommodated by the expression in the square brackets of eqs. (2.4) and (2.5). In brackets, we give also the r values for the NuFIT 4.1 data set, to illustrate the impact of the latest data.
Similar information is presented graphically in figure 1, showing the predicted number of events for the various appearance event samples as a function of δ CP , changing sin 2 θ 23 as well as the ordering, compared to the observed event number. Here the predictions are calculated using our experiment simulation based on fully numerical oscillation probabilities, while the general behaviour of the curves is well described by eqs. (2.4) and (2.5).
We can clearly observe a number of tendencies. T2K data has r > 1 for neutrinos and r < 1 for anti-neutrinos, implying that the square-bracket in (2.4) [(2.5)] has to be enhanced [suppressed]. If θ 13 is fixed as determined by reactor experiments this can be achieved by choosing NO and δ CP 3π/2 (see Sec. 2.2 for a consistent combination of reactor and LBL data). This has been the driving factor for previous hints for NO and maximal CP violation. We observe from the last row in table 1 that indeed this tendency has become somewhat weaker with the new data, though still clearly present. In this respect an interesting role is played by the CC1π event sample. A value r = 3.6 for NuFIT 4.1 shows a large excess of events in this sample, which has come down to r = 2.4 with the latest data. Figure 1 still shows, that even the most favorable parameter choice cannot  accomodate the observed number of events within 1σ. It seems that part of previous hints can be attributed to a statistical fluctuation in this sub-leading event sample. Let us stress, however, that due to the small CC1π event numbers, statistical uncertainties are large. Indeed, CCQE neutrino and anti-neutrino events consistently point in the same direction and they are both fitted best with NO and maximal CP phase.
Moving now to NOvA, we first observe from figure 1 the larger separation between the NO and IO bands compared to T2K. This is a manifestation of the increased matter effect because of the longer baseline in NOvA. Next, neutrino data have r ≈ 1 which can be accommodated by (NO, δ CP π/2) or (IO, δ CP 3π/2). This behavior is consistent with NOvA anti-neutrinos, however in tension with T2K in the case of NO. We conclude from these considerations that the T2K and NOvA combination can be best fitted by IO and δ CP 3π/2. This is indeed confirmed in figure 2, showing the ∆χ 2 profiles as a function of δ CP . We observe in the lower-right panel that NOvA disfavors (NO, δ CP 3π/2) by about 4 units in χ 2 , whereas in the lower-left panel we see for IO consistent preference of T2K and NOvA for δ CP 3π/2. For the combination this leads to a preferred best fit for IO with ∆χ 2 (NO) ≈ 1.5 (which of course is not significant). We can also see that this effect was less relevant in NuFIT 4.1 ( fig. 2, upper panels) for which we had r = 1.3 -compared to current 1.14 -for NOvA neutrino data. This slightly higher ratio allowed some more enhancement of the square-bracket in eq. (2.4) compared to the present situation, leading to less tension between T2K and NOvA for NO. It also lead to a larger significance of NOvA for NO.
The two-dimensional regions for T2K and NOvA in the (δ CP , sin 2 θ 23 ) plane for fixed θ 13 are shown in figure 3. The better consistency for IO is apparent, while we stress that even for NO the 1σ regions touch each other, indicating that also in this case the two experiments are statistically consistent. We are going to quantify this later in section 2.3.

Accelerator versus reactor
In the previous section we have discussed the status of the hints on CP violation and neutrino mass ordering in the latest LBL data. In the context of 3ν mixing the relevant oscillation probabilities for the LBL accelerator experiments depend also on θ 13 which is most precisely determined from reactor experiments (and on the θ 12 and ∆m 2 21 parameters which are independently well constrained by solar and KamLAND data). So in our discussion, and also to construct the χ 2 curves and regions shown in figs. 2, 3, and 4 for T2K, NOvA, Minos, and the LBL-combination, those parameters are fixed to their current best fit values. Given the present precision in the determination of θ 13 this yields very similar results to marginalize with respect to θ 13 , taking into account the information from reactor data by adding a Gaussian penalty term to the corresponding χ 2 LBL . Let us stress that such procedure is not the same as making a combined analysis of LBL and reactor data, compare for instance the blue solid versus black/blue dashed curves in fig. 2. This is so because relevant additional information on the mass ordering can be obtained from the comparison of ν µ and ν e disappearance spectral data [22,23]. In brief, the relevant disappearance probabilities are approximately symmetric with respect to the sign of two effective mass-squared differences, usually denoted as ∆m 2 µµ and ∆m 2 ee , respectively. They are two different linear combinations of ∆m 2 31 an ∆m 2 32 . Consequently, the precise determination of the oscillation frequencies in ν µ and ν e disappearance experiments, yields information on the sign of ∆m 2 3 . This effect has been present already in previous data (see, e.g., Ref. [2] for a discussion). We see from the two lower-left panels of figure 4 that the region for |∆m 2 3 | for IO from the LBL combination (blue curve) is somewhat in tension   In the accelerator-reactor combination this leads again to a best fit point for NO, with ∆χ 2 (IO) = 2.7, considerably less than the value 6.2 of NuFIT 4.1. This is explicitly shown, for example, in the LBL-reactor curves in fig. 2. For the NO best fit, a compromise between T2K and NOvA appearance data has to be adopted, avoiding over-shoting the number of neutrino events in NOvA while still being able to accommodate both neutrino and antineutrino data from T2K, see figure 1. This leads to a shift of the allowed region towards δ CP = π and a rather wide allowed range for δ CP for NO, see figures 2 and 3. On the other hand, we see from these figures that for IO, both T2K and NOvA prefer δ CP 270 • . Consequently, if we restrict to this ordering, CP conservation remains disfavored at ∼ 3σ.
The behaviour as a function of sin 2 θ 23 is shown in fig. 3 and the right panels of figure 4. It is mostly driven by the two T2K neutrino samples. As follows from eq. (2.4), their predicted event rate can be enhanced by increasing sin 2 θ 23 . Therefore, in order to compensate for the reduction in IO, a slight preference for the second θ 23 octant emerges for IO. In case of NO, this is less preferrable, since large sin 2 θ 23 would worsen the T2K anti-neutrino fit as well as NOvA neutrino data.

Consistency between T2K, NOvA and reactors
Let us now address the question of whether some data sets are in tension with each other at a worrisome level. A useful method to quantify the consistency of different data sets is the so-called parameter goodness-of-fit (PG) [24]. It makes use of the following test statistic: where i labels different data sets, χ 2 min,i is the χ 2 minimum of each data set individually, and χ 2 min,glob is the χ 2 minimum of the global data, i.e., χ 2 min,glob = min i χ 2 i . Let us denote by n i the number of model parameters on which the data set i depends, and n glob the number of parameters on which the global data depends. Then the test statistic χ 2 PG follows a χ 2 distribution with n degrees of freedom, where [24] n = i n i − n glob . (2.7) We are going to apply this test now to different combination of the three data sets, "T2K", "NOvA", and "React", where "React" is the joint data set of Daya-Bay, RENO and Double-Chooz. 3 The accelerator samples always include appearance and disappearance channels for both neutrinos and anti-neutrinos. In order to study the consistency of the sets under a given hypothesis for the neutrino mass ordering, all minimizations are preformed restricting to a given mass ordering. Furthermore, the solar parameters are kept fixed and hence, we have n T2K = n NOvA = n glob = 4 (namely θ 13 , θ 23 , δ CP , |∆m 2 3 |) and n React = 2 (namely θ 13 , |∆m 2 3 |). The results are shown in table 2. First, we check the pair-wise consistency of two out of the three sets. In all cases we find perfect consistency with p-values well above 10%. The only exception is NOvA vs React for IO which show tension at the 2σ level. A large contribution to this effect comes from the determination of ∆m 2 3 , which agrees better for NO than for IO, see fig. 4 (lower-left panels). The consistency of all three sets (T2K vs NOvA vs React) is excellent for both orderings.
Second, we perform an analysis for fixed sin 2 θ 13 = 0.0224 for all data sets. Since the accelerator experiments provide a comparatively weak constraint on θ 13 we want to remove this freedom from the T2K and NOvA fits and test the consistency under the hypothesis of fixed θ 13 . Under this assumption, all n i as well as n glob quoted above are reduced by 1. The results of this analysis are shown in the lower part of tab. 2. Testing T2K vs NOvA under this assumption, we find better compatibility for IO, consistent with the discussion above and figs. 2 and 3. Let us stress, however, that even for NO the p-value is 9%, indicating consistency at the 1.7σ level. Hence, we find no severe tension between T2K and NOvA. Finally, the joint T2K vs NOvA vs React analysis with fixed θ 13 reveals roughly equal good consistency among the three sets for both orderings, at around 1.5σ. For NO the very slight tension is driven by T2K vs NOvA, whereas for IO the reactor/accelerator complementarity in the determination of ∆m 2 3 provides a few units to χ 2 PG . To conclude this discussion, we find that all involved data sets are perfectly statistically compatible under the hypothesis of three-flavor oscillations.

Resolved tension in the solar sector
The analyses of the solar experiments and of KamLAND give the dominant contribution to the determination of ∆m 2 21 and θ 12 . It has been a result of global analyses for the last decade, that the value of ∆m 2 21 preferred by KamLAND was somewhat higher than the one from solar experiments. The tension appeared due to a combination of two effects: the well-known fact that the 8 B measurements performed by SNO, SK and Borexino showed no evidence of the low energy spectrum turn-up expected in the standard LMA-MSW [21,25] solution for the value of ∆m 2 21 favored by KamLAND, and the observation of a nonvanishing day-night asymmetry in SK, whose size is larger than the one predicted for the ∆m 2 21 value indicated by KamLAND. In our last published analysis [2] we included the energy-zenith spectra or day/night spectra for SK1-3, together with the 2860-day total energy spectrum of SK4 [26]. This last one made the lack of the turn-up effect slightly stronger. As for the day-night variation in SK4, it was included in terms of their quoted day-night asymmetry for SK4 2055-day [27] A D/N,SK4-2055 = [−3.1 ± 1.6(stat.) ± 1.4(syst.)]% . (3.1) Altogether this resulted in slightly over 2σ discrepancy between the best fit ∆m 2 21 value indicated of KamLAND and the solar results. For example the best fit ∆m 2 21 of KamLAND was at ∆χ 2 solar = 4.7 in the analysis with the GS98 fluxes. Here we update the solar analysis to include the latest SK4 2970-day results 4 presented in Neutrino2020 [16] in the form of their total energy spectrum and the updated day-night We show in fig. 5 the present determination of these parameters from the global solar analysis in comparison with that of KamLAND data. The results of the solar neutrino analysis are shown for the two latest versions of the Standard Solar Model, namely the GS98 and the AGSS09 models [29] obtained with two different determinations of the solar abundances [30]. For sake of comparison we also show the corresponding results of the solar analysis with the pre-Neutrino2020 data [2]. As seen in the figure, with the new data the tension between the best fit ∆m 2 21 of KamLAND and that of the solar results has decreased. Quantitatively we now find that the best fit ∆m 2 21 of KamLAND lies at ∆χ 2 solar = 1.3 (1.14σ) in the analysis with the GS98 fluxes. This decrease in the tension is due to both, the smaller day-night asymmetry (which lowers ∆χ 2 solar of the the best fit ∆m 2 21 of KamLAND by 2.4 units) and the slightly more pronounced turn-up in the low energy part of the spectrum which lowers it one extra unit.

Global fit results
Finally we present a selection of the results of our global analysis NuFIT 5.0 using data available up to July 2020 (see appendix A for the complete list of the used data including references). We show two versions of the analysis which differ in the inclusion of the results of the Super-Kamiokande atmospheric neutrino data (SK-atm). As discussed in Ref. [2] there is not enough information available for us to make an independent analysis comparable in detail to that performed by the collaboration, hence we have been making use of their tabulated χ 2 map which we can combine with our global analysis for the rest of experiments. This table was made available for their analysis of SK1-4 corresponding to 328 kton-years data [31]. The collaboration has presented new oscillation results obtained from the analysis of updated SK4 samples, both by itself [32] and in combination with the SK1-3 phases [16]. They seem to indicate that their hint for ordering discrimination has also decreased. Unfortunately the corresponding χ 2 maps of these analyses have not been made public. Hence in what follows we refer as "with SK-atm" to the analysis including the tabulated SK1-4 328 kiloton years data χ 2 map, i.e., the same as in NuFIT 4.0 and 4.1.
Here we graphically present the results of our global analysis in the form of onedimensional ∆χ 2 curves ( fig. 6) and two-dimensional projections of confidence regions ( fig. 7). The corresponding best fit values as well as 1σ and 3σ confidence intervals for the oscillation parameters are listed in table 3. 5 Defining the 3σ relative precision of the parameter by 2(x up − x low )/(x up + x low ), where x up (x low ) is the upper (lower) bound on a parameter x at the 3σ level, we obtain the following 3σ relative precision (marginalizing over ordering): where the numbers between brackets show the impact of including SK-atm in the precision of the determination of such parameter. The ∆χ 2 profile of δ CP is not gaussian and hence its precision estimation above is only indicative.
In table 3 we give the best fit values and confidence intervals for both mass orderings, relative to the local best fit points in each ordering. The global confidence intervals (marginalizing also over the ordering) are identical to the ones for normal ordering, which have also been used in eq. (4.1). The only exception to this statement is ∆m 2 3 in the analysis without SK-atm: in this case a disconnected interval would appear above 2σ corresponding to negative values of ∆m 2 3 (i.e., inverted ordering). Projecting over the combinations appearing on the elements of the leptonic mixing matrix we derive the following 3σ ranges (see Ref. [33]   Note that there are strong correlations between these allowed ranges due to the unitary constraint. The present status of leptonic CP violation is further illustrated in fig. 8 where we show the determination of the the Jarlskog invariant defined as: CP sin δ CP = cos θ 12 sin θ 12 cos θ 23 sin θ 23 cos 2 θ 13 sin θ 13 sin δ CP .

(4.3)
It provides a convention-independent measure of leptonic CP violation in neutrino propagation in vacuum [34] -analogous to the factor introduced in Ref. [35] for the description of CP violating effects in the quark sector, presently determined to be J quarks CP = (3.18 ± 0.15) × 10 −5 [36]. From the figure we read that the determination of the mixing angles implies a maximal possible value of the Jarlskog invariant of J max CP = 0.0332 ± 0.0008 (±0.0019) (4.4) at 1σ (3σ) for both orderings. Furthermore we see that with the inclusion of the new results, the best fit value J best CP = −0.0089 is only favored over CP conservation J CP = 0 with ∆χ 2 = 0.38, irrespective of SK-atm.

Summary
Let us summarize the main findings resulting from the Neutrino2020 updates in neutrino oscillations.
• The best fit in the global analysis remains for the normal mass ordering, however, with reduced significance. In the global analysis without SK-atm, inverted ordering is disfavored only with a ∆χ 2 = 2.7 (1.6σ) to be compared with ∆χ 2 = 6.2 (2.5σ) in NuFIT 4.1. This change is driven by the new LBL results from T2K and NOvA which indeed by themselves favor IO (with θ 13 as determined by the reactor data and θ 12 and ∆m 2 21 by the solar and KamLAND results). The best fit for NO in the combined global analysis is driven by the better compatibility between the ∆m 2 3 determined in ν µ disappearance at accelerators with that from ν e disappearance at reactors (see left panel in fig. 4).
• Despite slightly different tendencies in some parameter regions, T2K, NOvA and reactor experiments are statistically in very good agreement with each other. We have performed tests of various experiment and analysis combinations, which all show consistency at a CL below 2σ (section 2.3).
• If atmospheric data from Super-Kamiokande is included, inverted ordering is disfavored with a ∆χ 2 = 7.3 (2.7σ) compared to ∆χ 2 = 10.4 (3.2σ) in NuFIT 4.1. Hence, a modest indication for NO remains. Let us note that in the recent Super-Kamiokande update presented at Neutrino2020 [16] (with increased statistic and improved mass ordering sensitivity) the ∆χ 2 for IO is reduced by about 1 unit compared to the analysis we are using in our global fit. Therefore we expect that once the χ 2 map for the new SK analysis becomes available, the combined hint in favor of NO may further decrease.
• The best fit for the complex phase is at δ CP = 195 • . Compared to previous results (e.g., NuFIT 4.1 [11]), the allowed range is pushed towards the CP conserving value of 180 • , which is now allowed at 0.6σ with or without SK-atm. If we restrict to IO, the best fit of δ CP remains close to maximal CP violation, with CP conservation being disfavored at around 3σ.
• New solar neutrino data from Super-Kamiokande lead to an upward shift of the allowed region for ∆m 2 21 , which significantly decreased the tension between solar and KamLAND data. They are now compatible at 1.1σ, compared to about 2.2σ for the pre-Neutrino2020 situation.
Overall we have witnessed decreasing significance of various "hints" present in previous data. This is consistent with the fate of fluctuations which is that of fading away as time goes by.