DUNE sensitivities to the mixing between sterile and tau neutrinos

Light sterile neutrinos can be probed in a number of ways, including electroweak decays, cosmology and neutrino oscillation experiments. At long-baseline experiments, the neutral-current data is directly sensitive to the presence of light sterile neutrinos: once the active neutrinos have oscillated into a sterile state, a depletion in the neutral-current data sample is expected since they do not interact with the $Z$ boson. This channel offers a direct avenue to probe the mixing between a sterile neutrino and the tau neutrino, which remains largely unconstrained by current data. In this work, we study the potential of the DUNE experiment to constrain the mixing angle which parametrizes this mixing, $\theta_{34}$, through the observation of neutral-current events at the far detector. We find that DUNE will be able to improve significantly over current constraints thanks to its large statistics and excellent discrimination between neutral- and charged-current events.


Introduction
In the past decade, a tremendous experimental effort has been carried out in order to constrain scenarios with additional neutrinos with masses below the electroweak scale. LEP data places severe constraints on the invisible decay of the Z. Hence, if there are additional neutrinos below the electroweak scale, they cannot couple to the Standard Model weak bosons (i.e., they should be sterile). Light sterile neutrinos can lead to observable phenomena in a number of electroweak processes through their impact on the unitarity of the leptonic mixing matrix, including meson decays, muon decay, neutrinoless double beta decay and charged lepton flavor violating transitions (see e.g., refs. [1,2] for recent global fits using these observables). Nevertheless, if their masses are light enough so that they are kinematically accessible in these processes, unitarity is effectively restored at low energies and the bounds from electroweak processes fade away. In this case the best limits are derived from oscillation data [3][4][5][6][7][8], see e.g., refs. [9,10] for a detailed discussion of these constraints.
In recent years, the eV-scale has recently been put on the spot due to a set of experimental anomalies independently reported in LSND [11], MiniBooNE [12,13], reactor [14,15] and Gallium experiments [16]. The current and next generation of oscillation experiments will attempt to refute or confirm these hints. The Icecube experiment has recently put impressive limits on the mixing between sterile neutrinos and muon neutrinos U µ4 [17,18],

JHEP07(2018)079
while in the electron sector strong bounds on U e4 have been set by the Daya Bay experiment [19]. In the near future, experiments such as SOX [20] or STEREO [21] (among others) will constrain further the mixing with electron neutrinos, while the short-baseline neutrino program at Fermilab will tighten the bounds on the mixing with muon neutrinos [22]. A joint analysis of Bugey-3, Daya Bay and MINOS data has also been performed to constrain the cross-product |U e4 | 2 |U µ4 | 2 [23]. Conversely, placing equally competitive limits on the mixing with tau neutrinos is a much more difficult task, due to the technical challenges associated to the production and detection of a ν τ beam.
Indirect constraints on the mixing with ν τ can be derived from the observation of matter effects in atmospheric neutrino oscillations. For example, the IceCube experiment have set the limit |U τ 4 | 2 < 0.15 (at 90% CL) for an active-sterile mass splitting equal to 1 eV 2 [18] while the Super-Kamiokande experiment have set the bound |U τ 4 | 2 < 0.18 (at 90% CL) for an active-sterile mass splitting above 0.1 eV 2 [4]. 1 A non-zero θ 24 and θ 34 active-sterile mixing produces striking signatures in the zenith and energy distribution of cascade events in IceCube DeepCore, and after some years of data taking it is possible to probe the θ 34 parameter space [25]. On the other hand, a more direct test for the mixing between sterile neutrinos and tau neutrinos can be performed using long-baseline experiments. At long-baseline experiments most of the initial ν µ flux has oscillated into tau neutrinos by the time it reaches the far detector, thanks to ν µ → ν τ oscillations driven by the atmospheric mass-squared splitting. The OPERA experiment has constrained the impact of sterile neutrinos on this oscillation channel, using charged-current ν τ events at the far detector, setting the bound 4|U µ4 | 2 |U τ 4 | 2 < 0.116 (at 90% CL) for an active-sterile masssquared splitting above 0.1 eV 2 [26]. However, their results are severely limited by statistics, since the ν τ charged-current cross section is still low at multi-GeV neutrino energies.
Alternatively, the mixing between sterile neutrinos and tau neutrinos can be tested at long-baseline experiments searching for a depletion in the neutral-current event rates at the far detector. In fact, both the MINOS and the NOvA experiments have provided competitive constraints using this approach [27,28]. Future long-baseline experiments, with larger detectors, more powerful beams and a better control of systematic uncertainties, may be able to push these limits even further. In this work, we focus on the potential of the DUNE experiment [29]. Previous studies of sterile neutrino oscillations using the DUNE far detector data can be found, e.g., in refs. [9,[30][31][32][33][34][35][36]. 2 However, to the best of our knowledge the neutral-current data sample has not been considered in any of these works. The liquid Argon detector technology has excellent particle identification capabilities and therefore a very good discrimination power between charged-and neutral-current events. In addition, the statistics collected at DUNE will exceed considerably (by a rough order of magnitude) the number of events collected at MINOS or NOvA. Thus, DUNE offers an excellent benchmark to conduct a search for sterile neutrino mixing using neutral-current data. 1 An important constraint on the tau-sterile mixing angle θ34 has been obtained by combining IceCube DeepCore data [18] and short baseline data in ref. [24]. However, the constraint is given for a specific value of the sterile mass squared difference larger than 1 eV 2 .
2 For a sensitivity study using the DUNE near detector to probe sterile neutrino oscillations at ∆m 2 41 ∼ 1eV 2 we refer the reader to ref. [37].

JHEP07(2018)079
Although current hints of sterile-active neutrino mixing with ν e and ν µ occurs for a ∆m 2 of 0.1 eV 2 , in this paper we consider a broader range of ∆m 2 's similar to what Daya Bay has performed for the ν e disappearance search for sterile neutrinos, see ref. [19]. If a sterile neutrino only mixes with ν τ , then searches using ν e and ν µ disappearance as well as ν e appearance in a ν µ beam will not constrain such sterile-tau mixing.
The manuscript is organized as follows. In section 2 we derive the oscillation probabilities in the ν µ → ν s andν µ →ν s oscillation channels at the far detector of long-baseline experiments, and discuss the different limits of interest depending on the active-sterile mass-squared splitting. Section 3 summarizes the main features of the DUNE experiment and the details relevant to our numerical simulations. Our results are presented in section 4, and in section 5 we summarize and draw our conclusions. Some useful expressions for the elements of the mixing matrix using our parametrization can be found in appendix A.

Oscillation probabilities in the 3 + 1 framework
In this section we derive approximate expressions for the oscillation probabilities, which will be useful in understanding the results of our numerical simulations later on. The mixing matrix U that changes from the flavor to the mass basis in the 3 + 1 neutrino framework is a 4 × 4 unitary matrix: where α ≡ e, µ, τ, s and i ≡ 1, 2, 3, 4. In this work we are interested in the effect of oscillations into sterile states on the event rates measured at the DUNE far detector.
Assuming that no oscillations have taken place at the near detector, this can be done searching for a depletion in the number of neutral-current (NC) events at the far detector with respect to the prediction obtained using near detector data. For a perfect beam of muon neutrinos with flux φ νµ (i.e., assuming no beam contamination from other neutrino flavors), the number of NC events at the far detector can be expressed as: and is therefore sensitive to oscillations in the ν µ → ν s channel. Here, σ N C ν is the neutralcurrent cross section for the active neutrinos, which is independent of the neutrino flavor. In the absence of a sterile neutrino, the NC event rates should be the same at the far and near detectors up to a known normalization factor coming from the different distance, detector mass, efficiency, and the different geometric acceptance of the beam at the two sites. In fact, the combined fit between near and far detector data should provide a very efficient cancellation of systematic errors associated to the flux and cross section in this channel [28].
In addition to the standard solar and atmospheric mass-squared differences, in the 3+1 framework the oscillation probabilities depend on three new splittings ∆m 2 4k ≡ m 2 4 − m 2 k , with k = 1, 2, 3. Given the values of the neutrino energy and distance corresponding to the far detector at DUNE, for illustration purposes we can effectively neglect the solar JHEP07(2018)079 mass splitting and focus on the effects of the oscillation due to the atmospheric and the sterile mass-squared splittings. 3 Under the approximation ∆ 21 ∆ 31 , ∆ 41 , the oscillation probability in the ν µ → ν s channel is given (in vacuum) by: where we have defined ∆ ij ≡ ∆m 2 ij L/4E. The probability in eq. (2.2) is completely general, but does not allow to see the number of independent parameters which enter the oscillation probability. A 4 × 4 unitary matrix U can be parametrized in terms of six mixing angles and three Dirac CP-violating phases. 4 In the following, we choose to parametrize it as the product of the following consecutive rotations: Here, O ij denotes a real rotation with an angle θ ij affecting the i and j sub-block of the mixing matrix, while V ij denotes a similar rotation but this time including a complex phase. For example: where s ij ≡ sin θ ij and c ij ≡ cos θ ij . In this notation, θ i4 are the new mixing angles with the fourth state, and δ 14 , δ 24 are the two new CP-violating phases. In this parametrization, the complex phase associated with the V 13 rotation corresponds to the standard CP-violating phase in three-families, δ 13 ≡ δ CP , and the 3 × 3 sub-block of the matrix shows only small deviations from a unitary matrix, which at leading order are proportional to s 2 j4 and therefore within current bounds [8].
For simplicity, from now on we consider θ 14 = 0, which is a valid approximation given the strong constraints set by reactor experiments in the range of ∆m 2 41 considered in this work [23]. In this case there is no sensitivity to the δ 14 phase, which disappears from the mixing matrix, and the relevant elements of the mixing matrix read insensitive to these and therefore they will be ignored here. 2. The oscillation maximum due to the active-sterile mass-squared splitting matches the distance to the far detector (i.e., ∆ 41 ≈ ∆ 31 ):  Note that if U * µ4 U s4 + U * µ3 U s3 ≈ 0 there is a significant cancellation in the probability. This will be discussed in more detail later in this section.
3. The oscillations due to the active-sterile mass-splitting are already averaged-out at the far detector 5 (i.e., ∆ 41 ∆ 31 ): (2.9) 5 A similar expression in this limit, but assuming a real mixing matrix, can be found in ref. [38]. As mentioned above, a destructive interference between the standard and non-standard contributions to the oscillation amplitude is possible for certain values of the active-sterile mixing parameters and, in particular, for certain values of the CP phase δ 24 . This is shown in figure 1 for different values of ∆m 2 41 around the atmospheric scale, when the oscillation probability simplifies to eq. (2.8). The solid lines in all panels have been obtained for ∆m 2 41 = ∆m 2 31 : notice that a cancellation of the oscillation amplitude takes place in this case for δ 24 = 0, as shown in the left panel in figure 1. In this case, the contribution from the interference (last term in eq. (2.8)) is negative and cancels almost exactly the two other contributions to the oscillation probability. In fact, it is straightforward to show that, in the limit c 13 = c 24 = c 34 = 1, the amplitude of the oscillation is proportional to c 2 23 |s 24 c 23 − s 34 s 23 e iδ 24 | 2 , which vanishes exactly if δ 24 = 0 and s 24 c 23 = s 34 s 23 . This cancellation is only partial (or negligible) for other values of the CP phase, as expected, and this can be seen from the middle and right panels in the figure. For other values of the active-sterile mass splitting the oscillation pattern is more complex, as shown by the dotted blue and dashed yellow lines in figure 1. In the most general case, the dependence of the probability with the energy becomes non-trivial due to the interference of different terms oscillating at different frequencies. Moreover, as we will see in section 4 the cancellation in the probability can also be severe in the limit ∆m 2 41 ∆m 2 31 . Given the strong limits that have been set on the θ 24 angle by the oscillation experiments looking for oscillations involving a sterile neutrino in the eV scale, it is worth to address explicitly the case when θ 24 → 0. Under this assumption, the probability simplifies considerably with respect to the expression in eq. (2.6):

JHEP07(2018)079
In contrast with eq. (2.6), in this case there is no sensitivity to δ 24 and, most importantly, there is no dependence with the sterile mass-squared splitting. The oscillations in this case are solely driven by the atmospheric mass-squared splitting, and the size of the effect

JHEP07(2018)079
is directly proportional to s 2 34 . Moreover, the dependence with the standard oscillation parameters goes as c 4 13 sin 2 2θ 23 ∼ O(1). Finally, it is worth to mention that matter effects will modify the oscillation probability in eq. (2.6). We have checked that the size of these modifications is relatively small and, therefore, the vacuum probabilities are precise enough to understand the behavior of the numerical simulations in the following sections. However, in our numerical analysis, matter effects have been properly included using a constant matter density of 2.96 g · cm −3 .

Simulation
In contrast to usual analyses searching for signals of sterile neutrino oscillations at short distances, in this work we want to take advantage of the capabilities of the DUNE far detector, located at a distance of L = 1300 km from the source. In particular, we focus on the potential of NC measurements to discriminate between the 3-flavor and 4-flavor scenarios. To this end, we rely on the excellent capabilities of the DUNE far detector to discriminate between charged-current (CC) and NC events. All the simulations in the current work have been performed using a modified version of the GLoBES [39,40] library which includes a new implementation of systematic errors as described in ref. [41]. The neutrino oscillation probabilities in a 3+1 scenario have been implemented using the new physics engine available from ref. [42].
In our simulation of the signal, we have computed separately the contributions to the total number of events coming from ν e , ν µ and ν τ NC interactions at the detector. For simplicity, we have assumed a 90% flat efficiency as a function of the reconstructed visible energy. The experimental observable for a NC event is a hadronic shower with a certain visible energy (energy deposited in the detector in the form of a track and scintillation light). The correspondence between a given incident neutrino energy and the amount of visible energy deposited in the detector has to be obtained from the simulation of neutrino interactions and detector reconstruction of the particles produced in the final state. To this end, we use the migration matrices provided by the authors of ref. [43], which were obtained using the LArSoft simulation software [44] accounting for the far detector geometry, neutrino-argon interactions and propagation of the final state particles in the detector active volume. The authors of ref. [43] used bins in visible energy of 50 MeV for the reconstructed energy of the hadron shower, as opposed to the DUNE CDR studies where wider bins of 125 MeV were considered [29]. In the present work we have considered two sets of matrices: the original set provided by the authors of ref. [43], with 50 MeV bins, and a (more conservative) rebinned version of these matrices where the bin size was increased to 250 MeV. We performed our simulations for the two options (with 50 MeV bins and 250 MeV bins) and found similar results for the two sets of matrices. Therefore, in the following we will adopt the more conservative 250 MeV bin size as our default configuration.
The main backgrounds for this search would be ν e , ν µ and ν τ CC events that might be mis-identified as NC events. We have assumed that the background rejection efficiency for CC events is at the level of 90%. However, this is probably a conservative estimate: for instance, muons leave long tracks in liquid Argon (LAr) that are difficult to misidentify

JHEP07(2018)079
Signal Background  Table 1. Expected total number of events with a (reconstructed) visible energy between 0.5 and 8 GeV at the DUNE far detector. The number of events is shown for the signal and background contributions separately. This corresponds to 7 yrs of data taking (equally split between neutrino and antineutrino running modes) with a 40 kton detector and 1.07 MW beam power, yielding a total of 300 kt·MW·yr. In all cases, signal and background rejection efficiencies have already been accounted for. In the case of N ντ CC , the number of events already includes the branching ratio for hadronic τ decays. Usual oscillations (in the three-family scenario) have been considered in the computation of the backgrounds, setting θ 23 = 42 • and the rest of the oscillation parameters in agreement with their current best-fit values.
as NC events, except when they have very low energies or are not completely contained in the detector. On the other hand, the active neutrino flavors would be affected by standard oscillations. Consequently, the number of ν µ CC events would be largely suppressed since most of the initial muon neutrinos have oscillated to tau neutrinos by the time they reach the detector. Given the energetic neutrino flux at DUNE, some of the oscillated ν τ flux will interact at the detector via CC, producing τ leptons. In most of the cases (≈ 65%), the τ decays hadronically producing a shower: these events constitute an irreducible background and consequently no rejection efficiency has been assumed in this case. We have assumed a Gaussian energy resolution function for the ν µ and ν e background contributions, following the values derived in ref. [43] from LArSoft simulations, while the hadronic showers produced from hadronic tau decays have been smeared using the same migration matrices as for the NC signal.
The expected total number of signal and background events is summarized in table 1, where the different background contributions are shown separately for clarity. As can be seen from this table, the largest background contribution comes from ν µ CC events misidentified as NC, due to the large flux available at the far detector, while the contributions coming from ν e and ν τ CC events are much smaller and approximately of equal size. In all cases, both signal and backgrounds receive contributions from right-and wrong-sign neutrino events due to the intrinsic contamination of the beam. The number of events has been computed for visible energies between 0.5 GeV and 8 GeV, which is the region used in our analysis, using the beam configuration with 80 GeV protons as in ref. [45]. Additional experimental details for the DUNE setup considered in this work can be found in refs. [29,45].
The expected NC event distributions are shown in figure 2, as a function of the (reconstructed) visible energy, for the three-family scenario (white histogram) and for the case when there is a sizable mixing angle with the sterile neutrino (blue/light gray histogram). As expected, a depletion in the number of events can be observed in the 3 + 1 case with JHEP07(2018)079 respect to the three-family scenario. Moreover, the events pile up at low energies due to the energy carried away by the outgoing neutrino in the final state. One can also see that the energy distribution of the background (shown by the green/dark gray histogram) is dictated by the standard oscillations suffered by the active neutrinos as they propagate to the far detector, which is well-known. In this case, all particles in the final state would be observed, and there is practically no pile-up at low energies. Due to this, the sensitivity to oscillations in the ν µ → ν s channel is enhanced when some energy information is included in the fit, as we will see in the next section. This is exploited in our numerical analysis implementing a binned χ 2 in the visible (deposited) energy in the detector, with a χ 2 function defined in eq. (B.1).
The effect of systematic uncertainties is accounted for through the addition of pullterms to the χ 2 , as specified in the appendix B. In addition to an overall normalization uncertainty for the signal and background (which is bin-to-bin correlated), a shape uncertainty for the signal (bin-to-bin uncorrelated) has been included to account for possible systematic uncertainties related to the shape of the event distributions. Moreover, all nuisance parameters are taken to be uncorrelated between the neutrino and antineutrino channels as well as between the different contributions to the signal and/or background events. Unless otherwise stated, the final χ 2 is obtained after marginalization over the nuisance parameters and the relevant standard oscillation parameters (sin 2 2θ 23 , sin 2 2θ 13 , ∆m 2 31 ) within current experimental uncertainties [46][47][48]. Specifically, we consider the following Gaussian priors: σ(sin 2 2θ 13 ) = 0.005, σ(sin 2 2θ 23 ) = 0.05 and σ(∆m 2 31 )/∆m 2 31 = 0.04. Unless otherwise specified, we have assumed a conservative 10% Gaussian prior for all nuisance parameters, included as pull-terms in the χ 2 . In practice, however, the cancellation of JHEP07(2018)079 systematic errors in the NC channels is expected to be extremely efficient, since the near detector can be used to measure the same convolution of the flux and cross section as in the far detector. This contrasts with oscillation measurements in appearance mode (ν α → ν β ) using CC data, where the initial and final neutrino flux spectrum (and flavor) differ due to the impact of standard oscillations, making the cancellation of systematic uncertainties extremely challenging. 6 In spite of these difficulties, the DUNE collaboration expects to reach a precision at the percent level in the ν µ → ν e andν µ →ν e appearance channels. In view of this, we expect the 5%-10% values considered in this work for the NC sample to be conservative.
Before concluding this section, let us comment on the relevance of the near detector data and its possible impact on the fit. In this work, we have not simulated the near detector explicitly: its design is still undecided and its expected performance is therefore unclear yet. A detailed simulation of the near-far detector data combination is beyond the scope of this work and can ultimately be performed only by the experimental collaboration. In this work, instead, we have assumed that the oscillations due to the new state have not developed yet at the near detector. For neutrino energies in the region around 2-3 GeV, and for a near detector located at a distance of L ∼ O(500) m, this is a valid approximation as long as ∆m 2 41 < 1 eV 2 . Under this assumption, the near detector measurements will provide a clean determination of the convolution of the NC cross section and the muon neutrino flux, which can then be extrapolated to the far detector with a small uncertainty. At this point, it should be mentioned that our assumed prior uncertainties for the systematic errors in the fit would correspond to the values used for the analysis of the far detector event rates. Thus, they correspond to estimates on the size of the final systematic errors that have to be propagated to the far detector, once the near detector data has already been accounted for. Finally, it should also be stressed that in the case that θ 14 = θ 24 = 0 there would be no effect on the near detector data regardless of the new mass-squared splitting. The reason is that, as it was shown in eq. (2.10), the dependence with ∆m 2 41 drops from the oscillation probabilities: this guarantees no effect at the near detector, while at the far detector data the oscillation would be driven by the atmospheric scale. Thus, in this case the effect in the oscillation would be observable for large enough θ 34 .

Results
In this section we show our numerical results for the expected sensitivities to the new mixing parameters in the different scenarios discussed in section 2. By the time DUNE starts taking data the constraints on the sterile mixing angles θ 14 and θ 24 might be very tight. Nevertheless DUNE is also sensitive to the θ 34 sterile mixing angle, which is currently the less constrained among the three sterile-active mixing angles. Therefore, we initially consider the simpler case where two of the new mixing angles fixed to zero, θ 14 = θ 24 = 0 and study the sensitivity of the DUNE experiment to θ 34 . Next we proceed to turn on the mixing angle θ 24   Under the assumption θ 24 = θ 14 = 0, the expression for the vacuum sterile neutrino appearance probability is given by eq. (2.10) and does not depend on any of the new oscillation frequencies induced by the sterile, nor any of the CP-violating phases. An interesting question to ask in this case is if DUNE will be able to improve over current constraints on θ 34 , assuming that the experiment will measure event distributions in agreement with the expectation in the three-family scenario. In this case, the "observed" event distributions are simulated setting all θ i4 = 0, and are then fitted using increasing values of θ 34 . The sensitivity to θ 34 is shown in the left panel of figure 3. As seen in the figure, our results show a considerable dependence on the size and implementation of systematic errors. Assuming a (conservative) 10% systematic error on both normalization (σ norm ) and shape (σ shape ), we find that DUNE will be sensitive down to values of sin 2 θ 34 ∼ 0.12, at 90% C.L. (1 d.o.f.). For comparison we also show the limit on this mixing angle obtained from atmospheric neutrino data collected by the Super-Kamiokande (SK) collaboration [4], for ∆m 2 41 > 0.1 eV 2 . If prior uncertainties could be reduced to the 5% level for both normalization and shape errors, we find that DUNE would be able to improve over the SK constraint by more than a factor of two. It should be stressed that the DUNE constraint JHEP07(2018)079 would be valid for any value of ∆m 2 41 , as long as θ 24 , θ 14 0. In the next subsections we will study in detail the phenomenology in case θ 24 = 0.
The lines labeled as "Rate only" in the left panel of figure 3 do not include a binned χ 2 and only consider the total event rates in the computation of the χ 2 . The change in sensitivity can be appreciated from the comparison between the dashed pink and dotdashed red lines, for 10% systematic errors (or between the dot-dot-dashed green and solid blue lines, for 5% systematic errors). As can be seen, the inclusion of energy information leads to a noticeable improvement in the results. Therefore, in the rest of this section we will only consider a binned χ 2 , using equally-sized bins in visible energy, as described in section 3.
Finally, we comment on the analysis shown in the right panel of figure 3. Different to the analysis shown in the left panel, in the right panel we performed a discovery reach analysis for θ 34 taking all systematical errors at the 5% level. In this case, we assume that the experiment will measure event distributions in agreement with the expectation in the four-family scenario. In order to quantify the impact of also having a nonzero θ 14 and θ 24 , for simplicity, the four-flavor parameters where fixed to their 'true' values (except for θ 34 ) with the values in the plot labels. In this case (for θ 14 = 0 or θ 24 = 0) there is a dependence with the sterile mass squared difference, and we have fixed its value to ∆m 2 41 = 0.5 eV 2 which is one of the values considered in our results in figure 5. For ∆m 2 41 = 0.5 eV 2 we are therefore in the regime where the sterile oscillation is averaged-out at the far detector, in relation to eq. (2.9). In fact, is in this regime where constraints in the θ 24 − θ 34 plane are reported by different the collaborations (as will be addressed in section 4.3). It is then worth to mention that θ 14 is tightly constrained by reactor experiments for the ∆m 2 41 considered [19] and therefore its impact (even for θ 24 = 0) is marginal, as shown in the right panel of figure 3. Thus, since by the time DUNE will be running smaller values of θ 14 and θ 24 are expected, the case when θ 14 ≈ θ 24 ≈ 0 is of particular relevance. In this last case DUNE, with the considered configuration, will produce a 'indication' of a nonzero θ 34 (i.e. if it happens to be as large as ∼ 18 • ) with a significance of ∼ 2 σ for the assumed systematical errors. The scenario where θ 24 = 0 leads to a more interesting phenomenology, since in this case the oscillation probability also depends on the active-sterile mass-squared splitting. In this case, assuming as our true hypothesis a 3+1 with nonzero θ 34 and θ 24 , it is relevant to ask if the experiment would be able to reject the three-family hypothesis. This is shown in figure 4, as a function of the possible true values of ∆m 2 41 and sin 2 θ 24 . The true value of θ 34 is set to be nonzero, while θ 14 = 0 is assumed for simplicity. In all panels, the expected events distributions are computed using the indicated values as true input values. The obtained "observed" event distributions are then compared to the expected result in the three-family scenario, i.e., in absence of a sterile neutrino. The contours indicate the sets of true values (θ 24 , ∆m 2 41 ) for which the three-family hypothesis would be successfully rejected at 90% C.L.. The different panels in figure 4 show the dependence of our results , this leads to a decreased sensitivity in this region of the parameter space for δ 24 = 0 with respect to the results obtained for δ 24 = π. The interference has the opposite effect in the region ∆m 2 41 ∆m 2 31 : for negative values of cos δ 24 the second term in eq. (2.7) is negative and suppresses the probability, leading to worse results for δ 24 = π. In fact, it can be easily shown that, in the limit θ 23 = π/4, c 2 13 = 1 and at the first oscillation maximum (sin 2 ∆ 31 = 1) the oscillation probability in eq. (2.7) approximates to P µs ≈ c 2 24 (s 2 34 + 2s 24 s 34 c 34 cos δ + s 2 24 c 2 34 ) ,

JHEP07(2018)079
where the effect of the interference term can be easily appreciated. Conversely, in the limit where the new frequency is averaged-out (∆m 2 41 ∆m 2 31 ) the results show a very mild dependence with the value of δ 24 . This can be easily explained from the expression in eq. (2.9), which shows two terms that depend on the value of δ 24 : the first one is directly proportional to (c 2 13 s 2 23 −1/2) O(δθ 23 −s 2 13 /2), where δθ 23 ≡ θ 23 −π/4, and is therefore very suppressed; while the second term is proportional to sin 2∆ 31 and it is completely off-peak at the first oscillation maximum. In fact, in the same limit (θ 23 = π/4, JHEP07(2018)079 c 2 13 = 1) and at the first oscillation maximum it is easy to show that the term proportional to cos δ 24 in the oscillation probability in eq. (2.9) is additionally suppressed with cos 2θ 23 , which is small for θ 23 near maximal mixing. The central panel in figure 4 shows the dependence of the results with the true value of θ 34 . In this case, all priors for the systematic uncertainties are set at the 10% and we have fixed δ 24 = 0. As shown in the figure, in the region where ∆m 2 41 ∆m 2 31 there is a strong dependence of the results with the true value of θ 34 , while the contours do not show large variations for larger mass splittings. This behavior can again be easily traced back to the approximate oscillation probabilities in section 2.
Finally, the right panel in figure 4 shows the dependence of the results with the assumed priors for the systematic uncertainties. In this panel, the true values of δ 24 and θ 34 have been set as indicated in the top label. The solid line uses our default implementation for the systematic uncertainties, where all priors are set to 10% for both the shape and normalization and for both signal and background. The dot-dashed line, on the other hand, shows the room for improvement if all prior uncertainties can be reduced down to 5%. As can be seen from the figure, the improvement is dramatic and leads to a successful rejection of the three-family hypothesis in practically all the parameter space, with the sole exception of the region around ∆m 2 41 ∆m 2 31 (which is very difficult to reject, since this is the region where significant cancellations can take place for δ 24 = 0).

Expected allowed regions in the θ 24 − θ 34 parameter space
If the observed event distributions show an agreement with the three-family expectation, one would proceed to derive a limit on the mixing angles θ 24 and θ 34 . However, as we saw in section 2 the oscillation probabilities show a large dependence with the new CP-violating phase δ 24 , and strong cancellations between the different contributions may occur. The effect of the cancellations is much more severe in the limit ∆m 2 41 → 0 than for larger values of the active-sterile mass splitting and, therefore, we expect very different results as a function of this parameter. Figure 5 shows the expected allowed regions in the θ 24 and θ 34 plane if the observed event distributions are found to be in agreement with the three-family hypothesis. In this case, the "observed" event distributions are simulated assuming the three-family hypothesis, and fitted in a 3+1 scenario. The value of the χ 2 function, for a given pair of test values θ 24 − θ 34 , is obtained after minimization over the new CP-violating phase δ 24 and over all nuisance parameters. As for the mass splitting ∆m 2 41 , it has been kept fixed during the fit to the test value indicated in each panel to show the difference in the results. For simplicity, we have also kept all the standard parameters fixed during the minimization procedure; however, minimization over the standard parameters is not expected to affect significantly the results shown here.
As shown in figure 5, the resulting allowed regions are very different if the results are tested using ∆m 2 41 ∆m 2 31 or a ∆m 2 41 in the averaged-out regime. In the former case, a strong cancellation in the oscillation probability can always be achieved setting the value of δ 24 ∼ π, as outlined in section 2 and section 4. , for a simulation assuming θ i4 = 0 as true input values. The lines labeled as "10% sys" ("5% sys") have been obtained assuming 10% (5%) prior uncertainties for the signal (both shape and normalization) and 10% for the background (normalization only). For comparison, the right panel shows the latest results from the NOvA experiment from a NC search [28] (gray region), and from atmospheric data from the Super-Kamiokande experiment [4] and IceCube DeepCore data [18] (darker gray regions labeled with red lines), also at the 90% C.L.. large enough to allow for an efficient cancellation in the probability. Thus, in this regime DUNE could disfavor just the upper left and lower right corner of the parameter space. Conversely, in the limit ∆m 2 41 ∆m 2 31 the impact of the new CP-violating phase δ 24 is much milder and does not allow for a cancellation in the oscillation probability. A closed contour is therefore obtained in this case.
NOvA has observed 95 neutral current events at the far detector while 83.5±9.7(stat.)± 9.4(syst.) events where expected in the three-flavor case. Since no evidence for an sterile neutrino oscillation was found, they placed the following constraints (assuming cos 2 θ 14 = 1) for the active-sterile mixings: θ 24 < 20.8 • and θ 34 < 31.2 • at 90% of C.L for a ∆m 2 41 compatible with no oscillation at the near detector (0.05 eV 2 ≤ ∆m 2 41 ≤ 0.5 eV 2 ). This results correspond to an exposure-equivalent of 6.05 × 10 20 POT and a total systematical errors ∼ 12% ref. [28]. Experiments observing atmospheric neutrinos like the Super-Kamiokande experiment have also constrained the tau-sterile mixing angle. SK, after an analysis of 4, 438 live-days of data, found no evidence for sterile neutrinos constraining |U µ4 | 2 < 0.041 and |U τ 4 | 2 < 0.18 for ∆m 2 41 > 0.1 eV 2 at 90% of C.L [4]. Similarly, IceCube, by the use of three years of atmospheric neutrino data from the DeepCore detector, which was consistent with three-flavor neutrinos, placed a bound on the active-sterile mixing: |U µ4 | 2 < 0.11 and |U τ 4 | 2 < 0.15 for ∆m 2 41 = 1 eV 2 at 90% of C.L [18]. For comparison, in the right panel of figure 5 we show the currently allowed NOvA regions from ref. [28] as well as from atmospheric data from the Super-Kamiokande experiment [4] and IceCube DeepCore data [18]. 7 7 It is worth to notice that all constraints shown in the right panel of figure 5 are valid for an sterile mass squared difference ∆m 2 41 > 0.1 eV 2 and therefore they do not apply to the case shown in left panel.

JHEP07(2018)079
As shown in the figure, DUNE is expected to improve over a factor of two with respect to the current allowed region set by NOvA with a good control of systematics below 5%. With 5% systematics DUNE will also improve over the current IceCube constraint for θ 24 < 9 • . Finally, we want to comment on the impact of having a nonzero electron-sterile mixing in the results shown in figure 5. Even in the case where θ 14 current constraint were completely relaxed in the analysis (i.e. 'free'), our result for θ 34 does not change at all. Only the θ 24 bound (for θ 34 1) is affected. However, θ 14 is tightly constrained by Daya Bay experiment [19] for the ∆m 2 41 range where the current limits on θ 34 , shown in the right panel of figure 5, apply. On the other hand, in the left panel of the same figure (with θ 14 unconstrained) DUNE is not able to exclude the small window 24 • θ 24 30 • with 10% systematical errors. This, also reinforces our conclusion about the loose of constraining power for ∆m 2 41 ∆m 2 31 due to cancellations in the probability.

Summary and conclusions
The experimental anomalies independently reported in LSND, MiniBooNE, reactor and Gallium experiments have put the possible existence of an eV-scale sterile neutrino under intense scrutiny. In the near future a new generation of short-baseline experiments will come online to refute or confirm these hints, and will place strong constraints on the mixing of a light sterile neutrino with electron and muon neutrinos. Achieving similar bounds on the mixing with tau neutrinos is a much more difficult task, given the technical challenges associated to the production and detection of ν τ . At long-baseline experiments, however, oscillations in the ν µ → ν τ channel guarantee that most of the beam will have oscillated into ν τ by the time it reaches the far detector, thanks to the atmospheric masssquared splitting. By searching for a depletion in the number of neutral-current (NC) events measured at the far detector, experiments like NOvA or MINOS have been able to probe the mixing between ν τ and a fourth neutrino.
In this work, we have studied the potential of the future DUNE experiment to conduct a search for sterile neutrinos using the NC data expected at the far detector, taking advantage of the excellent capabilities of liquid Argon to discriminate between charged-current and NC events. For simplicity, we have focused on a 3 + 1 scenario, where only one extra sterile neutrino is introduced. In this case, the mixing matrix has to be extended including three additional mixing angles (θ 14 , θ 24 and θ 34 ) and two CP-violating phases δ 14 and δ 24 (our parametrization is given by eq. (2.3)). The oscillation probabilities will generally depend on an additional oscillation frequency dictated by the mass-squared splitting between the active and sterile states, ∆m 2 41 . First, we have derived the oscillation probabilities in different regimes paying particular attention to the dependence with the new CP-violating phases. Unlike in other studies where the mass of the sterile was required to be at (or around) the eV scale, here we have allowed it to vary between 10 −5 eV 2 and 10 −1 eV 2 ; thus, in eqs. (2.7)-(2.9) we provide approximate expressions for the oscillation probabilities in three different regimes, depending on the mass of the sterile state: (i) ∆m 2 41 → 0; (ii) ∆m 2 41 ∆m 2 31 ; and (iii) ∆m 2 41 ∆m 2 31 . We have then proceeded to simulate the expected sensitivity of the DUNE experiment using the expected NC events collected at the far detector. We have studied the variation JHEP07(2018)079 of our results with the implementation and size of the systematic errors. The details of our numerical simulations and the χ 2 implementation can be found in section 3.
First, working under the assumption θ 24 = θ 14 = 0, we have determined the sensitivity of the DUNE experiment to the third mixing angle θ 34 . In this case, the oscillation probability is independent of the new CP-violating phases; furthermore, oscillations are solely driven by ∆m 2 31 , see eq. (2.10). We find that DUNE will be able to improve over current constraints on this parameter set by the SK experiment, and will be sensitive to values of sin 2 θ 34 ∼ 0.12 (at 90% CL) for our default implementation of systematic uncertainties. If systematic errors could be reduced down to 5%, the experimental sensitivity would reach sin 2 θ 34 ∼ 0.07 (at 90% CL).
Next we proceeded to study the case where θ 24 = 0. In this case, the oscillation probabilities depend on the active-sterile mass-squared splitting. The phenomenology becomes more complicated and, in particular, strong cancellations in the probability can take place for certain values of δ 24 and ∆m 2 41 . First, we considered the 3+1 scenario as the true hypothesis, and determined for which values of the mixing parameters DUNE would be able to reject the three-family scenario. Our results are summarized in figure 4, where we show the dependence of the sensitivity with the CP phase δ 24 , the mixing angle θ 34 and the size of the systematic errors. We found that the sensitivity of the experiment to the presence of a sterile neutrino, measured as its ability to reject the three-family scenario, depends heavily on the value of the CP phase. For example, for ∆m 2 41 = 10 −4 eV 2 , sin 2 θ 34 = 0.1 and δ 24 = 0, DUNE would be able to reject the three-family scenario for sin 2 θ 24 4 × 10 −3 eV 2 ; conversely, for δ 24 = π (and assuming the same value for θ 34 and ∆m 2 41 ), θ 24 could be almost two orders of magnitude larger and the three-family scenario would not be rejected by the data. The behavior of our results can be easily understood in terms of the oscillation probabilities, as explained in detail in section 4.2.
Finally, we considered the opposite situation, and assumed that the experiment will find a result that is in agreement with the three-family expectation. In this case, we determined the allowed confidence regions that would turn from the analysis of the simulated data. Our results are shown in figure 5. The simulated data were tested using two very different values of the active-sterile mass-squared splitting. In the averaged-out regime (∆m 2 41 ∆m 2 31 ), a closed contour is obtained; we find that DUNE would be able to improve over NOvA constraints in this place by a factor of two or more, depending on the size of the systematic errors assumed. Conversely, in the case of ∆m 2 41 ∆m 2 31 the experimental results would allow values of θ 24 and θ 34 to be as large as 30 • . The reason is, again, the possibility of having a strong cancellation in the oscillation probability, which could lead to a nonobservable effect in the event distributions even in presence of very large mixing angles.
The DUNE experiment has unprecedented discrimination between neutral-current and charged-current events for a long-baseline experiment: this will allow for a measurement or constraint of the ν τ fraction of a possible sterile neutrino(s). Given the difficulties associated to the production and detection of ν τ 's, measurement or limiting this fraction by other means is very challenging. In this paper, we show that the DUNE experiment can provide an excellent constrain or discover a sterile neutrino that primarily mixes with only the ν τ . Where σ a and σ b are total normalization systematical errors in signal and background, respectively. For simplicity we have assumed σ a = σ b = σ norm . c i are the bin-to-bin uncorrelated systematics with error σ shape . The last four terms in eq. (B.1) are penalties to the χ 2 function due to the systematics included in eq. (B.2), and also due to the standard oscillation parameters λ that are marginalized over assuming they have been measured asλ j ± σ j .
Open Access. This article is distributed under the terms of the Creative Commons Attribution License (CC-BY 4.0), which permits any use, distribution and reproduction in any medium, provided the original author(s) and source are credited.