A search for IceCube events in the direction of ANITA neutrino candidates

Ever since the detection of high-energy neutrinos of cosmic origin by IceCube in 2013 (Aartsen et al., 2013a), experiments and theoreticians alike have continued to probe the nonthermal processes in the Universe to understand their origins. The bulk of these astrophysical neutrinos are believed to be created in hadronic interactions between cosmic rays and ambient matter or radiation fields in the vicinity of cosmic accelerators (Gaisser et al., 1995), and their detections can be used to point back to the acceleration sites. Although the first evidence of a neutrino point source, the blazar TXS 0506+056, was reported in 2018 (Aartsen et al., 2018d, c), the overwhelming majority of the measured neutrino flux remains unexplained.

Additionally, another population of neutrinos could exist at extremely high energies. Cosmogenic neutrinos are believed to be the result of interactions between ultra-high-energy (UHE) cosmic rays with the cosmic microwave background (CMB) (Greisen, 1966; Zatsepin & Kuzmin, 1966). This population is expected to manifest as an isotropic flux at Earth, as cosmic ray primaries can travel outside the vicinity of their accelerators before interacting with the CMB.

The Antarctic Impulsive Transient Antenna (ANITA) project is a balloon experiment, designed with the primary purpose of detecting the UHE cosmogenic neutrino flux (Gorham et al., 2009; Hoover et al., 2010; Gorham et al., 2018a). Although this is the project’s primary scientific goal, the experiment is sensitive to a wide array of impulsive radio signals, and ANITA’s first three flights have resulted in a few interesting detections. In this work, we focus on three events observed by ANITA in its searches, all of which have potential neutrino interpretations. Throughout this work, we refer to and explore these events as “neutrino candidates.” In the third flight, one Askaryan neutrino candidate (AAC) event was simultaneously identified in one analysis searching for Askaryan emission (Askar’yan, 1962) and was found to be subthreshold in another. This Earth-skimming event has a signal shape consistent with impulsive broadband emission characteristic of a neutrino origin, and it also came from a location on the continent consistent with simulated distribution of neutrinos of all flavors (Gorham et al., 2018a). However, the detection of one candidate event is consistent with the background level estimates of $0.7^{+0.5}_{-0.3}$ for these analyses. ANITA also reported two additional events, each consistent with an astrophysical $\nu_{\tau}$ emerging from Earth (Gorham et al., 2016, 2018b). In this scenario, a $\nu_{\tau}$ undergoes a charged-current interaction (CC) with a nucleus in Earth. The $\tau$-lepton produced in this interaction subsequently decays in the atmosphere, producing an extensive air shower (EAS). The polarity of the radio signal makes it possible to identify and reject downward moving cosmic-ray–induced EAS, as the radio signals of these EAS acquire a phase reversal (opposite polarity) from reflection off the Antarctic ice, while an upgoing $\tau$-induced EAS does not acquire this phase reversal. For a complete list of details of these events, see Table 1.

The interpretation of these events as extremely high energy upgoing neutrinos poses many challenges under Standard Model assumptions. First, from the observation angles and reconstructed energies of the ANITA events, neutrinos are extremely unlikely to traverse the long chord lengths (Gorham et al., 2016), even after accounting for the probability increase due to $\nu_{\tau}$ regeneration. Second, if these events are of cosmogenic origin, they would imply fluxes that are in severe tension with limits set by multiple experiments (Aab et al., 2015; Zas, 2018; Aartsen et al., 2016a, 2018b) as well as a self-inconsistency from ANITA data alone. For an isotropic flux of cosmogenic neutrinos, ANITA should have detected many more events at other elevation angles than those of the anomalous ANITA events (AAE) as the detector differential acceptance changes with the observation angle (Romero-Wolf et al., 2019).

On the other hand, if the origin of the AAE is considered to be from individual cosmic accelerators, there is no inconsistency with diffuse extremely high energy flux limits. This is especially true for accelerators with short characteristic timescales of emission, as many current limits on neutrino point sources are for integrated emission over various experiments’ live times (Aartsen et al., 2019a) and also as the acceptance of ANITA to a specific location in the sky changes throughout the detector’s flight. If we assume that ANITA detected single events of 1 EeV from a cosmic accelerator with an $E^{-\gamma}$ emission power-law spectrum, then one should expect also a larger flux of neutrinos at TeV–PeV energies, where IceCube will be sensitive. Significant correlation between IceCube and ANITA data would not only provide evidence for a neutrino point source, it would also eliminate nonastrophysical explanations of the AAE, such as background and systematics, or nonastrophysical models, which invoke physics beyond the Standard Model.

The focus of this work is to use IceCube to investigate the hypothesis that the ANITA events were from neutrino point sources, considering several neutrino emission time profiles. In section 2, we discuss the IceCube Neutrino Observatory and the event samples used for these analyses. In section 3, we describe the analysis techniques, and we summarize the results in section 4. In sections 5 and 6, we investigate neutrino propagation through large Earth chord lengths to discuss the implications of our results.

IceCube is a cubic-kilometer neutrino detector with 5160 digital optical modules (DOMs) instrumented on 86 cable strings in the clear glacial ice at the geographic South Pole, at depths between 1450 m and 2450 m (Achterberg et al., 2006; Aartsen et al., 2017b). Neutrinos are detected through the Cherenkov radiation emitted by secondary particles produced by neutrino interactions in the surrounding ice or bedrock. Each DOM consists of a 10-inch photomultiplier tube, onboard readout electronics, and a high-voltage board, all contained in a pressurized spherical glass container (Abbasi et al., 2009, 2010). Parameterization of the scattering and absorption of the glacial ice allows accurate energy and directional reconstruction of neutrino events (Aartsen et al., 2013b).

The improved reconstruction techniques adopted to create the event selection (Carver, 2019; Aartsen et al., 2019b) include updates in the direction reconstruction (Ahrens et al., 2004; Aartsen et al., 2014a) to use information on the deposited event energy in the detector. The median angular resolution benefits from a 10% improvement above 10 TeV (where it is smaller than 0.60^∘) compared to previous selections (Aartsen et al., 2017a).

While in the southern sky the trigger rate is dominated by atmospheric muons from cosmic-ray air showers, all of the ANITA candidates have best-fit directions in the northern sky. Here, Earth attenuates the majority of the atmospheric muon signal, and the background at final selection level in the northern sky is dominated by atmospheric muon neutrinos from cosmic-ray air showers (Haack & Wiebusch, 2018). Poorly reconstructed atmospheric muons from the southern sky as well as neutrino-induced cascades are also non-negligible backgrounds in this region of the sky and are removed using a multivariate boosted decision tree trained to distinguish between neutrino-induced muon tracks, atmospheric muons, and cascades, which is described in (Carver, 2019; Aartsen et al., 2019b).

For the analyses presented here, we focus on the full detector configuration of 86 strings, spanning a time window from 2011 to 2018. Approximately 900,000 events from 2532 days are analyzed.

Many previous IceCube analyses searching for neutrino point sources relied on significant spatial clustering of IceCube data alone or of significant association with known populations of astrophysical objects (Aartsen et al., 2019a, 2016b, 2017a, 2013c; Abbasi et al., 2011; Aartsen et al., 2014b). Here, we adopt the procedure described in (Schumacher, 2019) to search for counterparts to ANITA events. Namely, we perform three separate analyses to test different temporal hypotheses in the neutrino emission. Each of these analyses incorporates the information from the localization of the ANITA events through a joint likelihood. The sky is divided into grid positions, $\mathbf{x}_{s}$, and at each point we maximize the likelihood, $\mathcal{L}$, with respect to the expected number of signal events, $n_{s}$, and other signal parameters contained in the variable $\mathbf{\alpha}$ depending on the different signal hypotheses tested as described in sections 3.2 and 3.3. This likelihood is given by

where $n_{b}$ is the expected number of observed background events and $N$ is the total number of observed events in the time window. The vector $\mathbf{x}_{i}$ contains the event observables such as its reconstructed energy, direction, and reconstruction uncertainty. $P_{A}$ is the spatial probability distribution function (PDF) of ANITA events, which are included in Table 1. $B$ describes the energy and declination PDF of our background, which is parameterized from data and is the same among all analyses. Temporal terms in $B$ are described in sections 3.1 and 3.2. While the signal PDF $S$ describes the signal hypothesis, the parameter $\lambda$ modifies the likelihood formalism in order to take into account low-statistics problems in some of the analyses. In general, the signal PDF, $S$, is defined as

These three terms reflect the spatial, energy, and time PDFs, respectively, of our signal hypothesis. The spatial term, $S^{space}$, expresses the probability for an event with best-fit reconstructed direction $\mathbf{x}_{i}$ to originate from a source at the direction $\mathbf{x}_{s}$, according to a two-dimensional Gaussian function with angular resolution $\sigma_{i}$. The energy PDF, $S^{energy}$, describes the probability of obtaining an event with reconstructed energy $E_{i}$ given a declination $\delta_{i}$ under the hypothesis of an $E^{-\gamma}$ power-law energy spectrum, which helps differentiate signal from the known atmospheric backgrounds in our event selection. The time term, $S^{time}$, describes the time PDF of events observed from the source. While the spatial term is shared between all analyses, the energy and temporal terms are unique to each individual analysis. This joint likelihood procedure is carried out in three complementary search strategies: prompt, rolling, and steady.

No significant correlation is found in any of the analyses above the expectation from background. In order to calculate p-values, results are compared against pseudo-experiments from time-scrambled data (Aartsen et al., 2015b). The most significant observation results from the steady search for AAE-141220, with a p-value of 0.08 before trials correction.

Figure 1 displays the sky maps for the prompt, rolling, and steady analyses from left to right in the top panels for AAE-141220. Bottom panels of Figure 1 show the comparison of the observed TS values for each analysis, at the position of the red lines, to their respective TS distributions from pseudo-experiments using time-scrambled data. Similar plots for AAE-061228 and AAC-150108 are displayed in Figure 6.

In the absence of a significant signal, upper limits (90% confidence level) for the time-integrated $\nu_{\mu}+\bar{\nu}_{\mu}$ flux are set for each ANITA event where possible using the prompt and steady analyses (Figure 2). To calculate upper limits, locations are sampled according to the per-event PDFs reported by ANITA, injecting the same level of flux at each sampled location, and running each iteration through the full analysis procedure, which maximizes the joint likelihood at all locations in the sky. This allows us to place upper limits on point sources whose locations are distributed according to the per-event PDF reported by ANITA. We set these limits for an assumed spectrum given by

where $\Phi_{0}$ is a normalization constant on a point-source flux, which carries units of $\rm{GeV}^{-1}\rm{cm}^{-2}\rm{s}^{-1}$. We constrain the time-integrated muon neutrino flux, $E^{2}F$, where

All of the limits we calculate are provided in Table 2. In the case that an upper limit fluctuates below the sensitivity, we conservatively set the upper limit to the sensitivity value. Prompt limits are placed at the specified time windows for emission centered on the ANITA event times, whereas limits from the steady analysis are for emission over the live time of our data sample. This hard spectrum was chosen conservatively because with the observation of EeV events by ANITA, if the underlying spectrum is softer, then the expected number of observable neutrinos for IceCube would increase. As the time-integrated flux sensitivity for the triggered analysis begins to worsen past $10^{5}$ s, upper limits for $\Delta t>10^{5}$ s are only set using the time-integrated approach.

For many astrophysical sources, power-law spectra in photons are common over finite energy ranges. Additionally, diffusive shock acceleration models suggest that the neutrino spectrum, as well as gamma rays from pion decay, should follow a power-law spectrum, justifying the choice of testing power laws for corresponding neutrino spectra. However, for the AAE, interpolating a power law between the energy range at which IceCube is sensitive to the best-fit ANITA event energies could pose a problem. For soft spectra, events detected by ANITA would suggest that many events would be detectable at IceCube. For hard spectra, extrapolating between IceCube and ANITA would imply dramatic bolometric neutrino luminosities for any point source.

However even in the case of non-power-law neutrino emission, the limits we can set on muon neutrinos in the TeV–PeV energy range can constrain generic fluxes of incident tau neutrinos with EeV energies. As has been shown in (Safa et al., 2020), any incident flux with an EeV $\nu_{\tau}$ component that traverses large Earth chord lengths will result in a secondary flux of lower energy neutrinos, to which IceCube would be sensitive. We use the same prescription here to analyze how constraining our limits are on a generic point source flux that includes EeV neutrinos.

For any incident flux of neutrinos from the northern sky, $\Phi(E_{\nu},t)$, the number of expected detected tau-neutrino–induced muon events at IceCube is given by

where the first contribution is from emerging tau-leptons that would decay to muons and then pass an IceCube event selection. The second contribution is from the remaining $\nu_{\tau}$ flux, the majority of which has cascaded down in energy. $N^{p}(E_{\nu})$ is the number of targets effectively seen by an incident neutrino with energy $E_{\nu}$. The effective area of this event selection to muons incident on the detector is displayed in Figure 3. $P_{\tau}^{surv}\left(E_{\nu}\right)$ and $P_{\nu}(E_{\nu})$ represent the survival probability of a $\tau$-lepton and $\nu_{\tau}$, given an incident neutrino energy, respectively, and $\Gamma_{\tau\rightarrow\mu}\big{/}\Gamma_{\rm{total}}$ represents the branching ratio for the tau-decay to muon channel, which is approximately 18%.

Similarly, for ANITA, the number of expected events from upgoing $\tau$-leptons is given by

where $\xi_{acc}$ represents ANITA’s acceptance to $\tau$-lepton air showers, taken from (Romero-Wolf et al., 2019). Values for the acceptance at angles that would require an incident neutrino to traverse a large column depth are set to the acceptance near the horizon. We take the value at an angle corresponding to the maximum acceptance before absorption effects dominate. This removes absorption effects in the reported acceptance, which is accounted for separately with the code used to propagate these fluxes, TauRunner, described in (Safa et al., 2020, 2019). We focus our analysis on the nonobservation of coincident events in IceCube at $\Delta T=10^{3}$ s. A similar procedure can be applied to longer time windows. Qualitatively, it would result in similar limits up to the lifetime of the ANITA flight. For longer emission timescales, limits from IceCube become even more constraining as the implied normalization on the ANITA flux would have to increase to compensate for the fraction of time during which ANITA was not taking data.

To make as conservative a statement as possible, we inject fluxes described by delta functions in energy, $\Phi(E_{\nu},t)=\Phi_{0}\delta(E_{\nu}-E_{0})$, where now the normalization carries units of $\rm{cm}^{-2}\rm{s}^{-1}$. After propagating these mono-energetic fluxes, we record what fraction of the incident flux results in a detectable signal at ANITA. We repeat this procedure for a variety of injected initial neutrino energies so that we can find the energy that yields the maximum probability of a $\tau$-lepton arriving at ANITA with an energy within the quoted reconstructed energy bounds. We find that the optimal flux for ANITA corresponds to an injected $\nu_{\tau}$ flux with $E_{0}=1$ EeV. Normalized cumulative distributions from secondary $\tau$-leptons are shown in Figure 4 for injected neutrinos at angles corresponding to the best-fit reconstructed direction of AAE-141220.

We next inject a flux of EeV tau neutrinos and find the spectral shape of the secondary $\nu_{\tau}$ flux that would be incident on IceCube. As we observed zero coincident events in the time window of $10^{3}$ s around AAE-141220 in the prompt analysis, we calculate the maximum allowed flux normalization (at 90% confidence level) on the primary flux that would evade this nonobservation. The results are displayed in Figure 5.

Although IceCube’s sensitivity peaks many orders of magnitude below the reconstructed energies of the ANITA events, the limits set on any potential neutrino source that created AAE-141220 are more constraining by several orders of magnitude than the implied flux by the ANITA observations. If one considers constant emission over the entire live time of the IceCube event selection, then the time-integrated flux limit set by the IceCube nonobservation of AAE-141220 becomes around one order of magnitude less constraining, as is apparent in the steady limits in Figure 2. However, for the implied normalization placed by ANITA observations, this value would increase by approximately two orders of magnitude, due to the limited live time of the ANITA flight. This has the overall effect of increasing the tension between these two normalizations by approximately one more order of magnitude than for the 10³ s followup shown in Figure 5. It is worth noting that the logic for scaling time-integrated limits also applies to AAE-061228, even though we cannot constrain the shorter timescales for this event. However, the emergence angle of this event at ANITA was shallower than that of AAE-141220, which increases the probability of observing such an event at ANITA by approximately one order of magnitude (Fox et al., 2018) for the same assumed initial flux, and thus the limit on assumed long timescale emission would be about one order of magnitude less constraining than the case of AAE-141220.

If the intrinsic spectrum were to contain contributions from energies below 1 EeV, such as the power-law spectra tested in the analyses presented in Section 3, this would introduce a component to which IceCube might be sensitive but which could not produce events at ANITA consistent with the AAE, and thus this additional component would strengthen the constraints displayed in Figure 5. Additionally, if the spectrum consisted of neutrinos of energy greater than 1 EeV, the secondary $\nu_{\tau}$ spectrum would have a similar shape to that shown in Figure 5, as discussed in (Safa et al., 2020), and therefore the limits on the flux normalization would be constant for fluxes of higher energy, while the energy required to produce such a flux would scale with the injected energy. For that reason, these limits are conservative, and severely constrain any incident spectrum which could produce an observable event at ANITA consistent with the AAE.

Recent detections of the AAE are considered anomalous due to the small survival probability of EeV tau neutrinos through long chord lengths. The events are known to be inconsistent with a cosmogenic interpretation but could have been produced by cosmic accelerators, specifically those with short characteristic timescales. We show here that for timescales as small as $10^{3}$ s, assuming AAE-141220 as originating from a neutrino source, limits set using IceCube data are in tension with the point source flux required to detect one event at ANITA by more than four orders of magnitude. These limits are constraining for a variety of flux models, from simple power laws to any generic model that includes a component at or above EeV energies. In addition to the anomalous events, we also find no evidence for a neutrino source in the direction of the neutrino candidate event from a search for Askaryan emission during ANITA-III. As searches for Askaryan emission with ANITA have targeted a diffuse UHE cosmic neutrino flux (Gorham et al., 2018a) and not localized point source fluxes, studies that quantify the acceptance of the ANITA detector (Cremonesi et al., 2019) focus on diffuse acceptances and not effective areas for neutrino fluxes from fixed locations in the sky. For this reason, we do not provide a comparison between the limits we set here and potential implications for point source fluxes based on the observation of the AAC. With knowledge of the effective area of ANITA in the direction of the AAC, and assuming that any astrophysical flux was roughly equal in flavor upon reaching Earth, the same secondary $\nu_{\tau}$ analysis could be performed for the AAC. However, constraints from such a search would be considerably weaker than those for the AAE, as the AAC was Earth-skimming, and thus a greater fraction of any high-energy incident flux would be able to reach the ANITA detector prior to interacting deep within Earth. Therefore, this method of using secondary $\nu_{\tau}$ fluxes from UHE neutrinos in IceCube could be beneficial for future correlation searches with radio detectors and future Cherenkov detectors such as POEMMA (Venters et al., 2019).

These new limits, in conjunction with the inconsistency of isotropic flux interpretations, leave no room for an astrophysical interpretation of the AAE in the context of the Standard Model for time windows as short as $10^{3}$ s. However, it has been shown that these events can be explained using physics beyond the Standard Model, as many models suggest that the AAE lend support for axionic dark matter, sterile neutrinos, supersymmetry, or heavy dark matter (Cherry & Shoemaker, 2019; Anchordoqui et al., 2018; Huang, 2018; Dudas et al., 2018; Connolly et al., 2018; Fox et al., 2018; Collins et al., 2019; Esteban et al., 2019; Heurtier et al., 2019b, a; Abdullah et al., 2019; Anchordoqui & Antoniadis, 2019; Borah et al., 2019; Chipman et al., 2019; Cline et al., 2019; Esmaili & Farzan, 2019; Hooper et al., 2019; Chauhan & Mohanty, 2019). Many of these models, excluding the axionic dark matter explanation (Esteban et al., 2019) or those heavy dark matter scenarios that are tuned to prevent signatures in IceCube (Hooper et al., 2019), can be constrained by this nonobservation at IceCube. Dedicated tests to quantify these constraints are beyond the scope of this work and may be the focus of a future study. In addition to explanations that point to new physics, it has recently been suggested that the AAE could be explained by downgoing CR-induced EAS that reflected off of subsurface features in the Antarctic ice (Shoemaker et al., 2019). Another possible explanation could be coherent transition radiation from the geomagnetically induced air shower current, which could mimic an upgoing air shower (de Vries & Prohira, 2019; Motloch et al., 2017). Explaining these anomalous events with systematic effects or confirming the need for new physics requires a deeper understanding of ANITA’s detection volume. Efforts such as the HiCal radio frequency pulser, which has flown alongside ANITA in the last two flights (Prohira et al., 2018), are already underway to try to characterize the various properties of the Antarctic ice surface.

The IceCube Collaboration acknowledges the significant contributions to this manuscript from Anastasia Barbano, Alex Pizzuto, and Ibrahim Safa. The authors gratefully acknowledge the support from the following agencies and institutions: USA – U.S. National Science Foundation-Office of Polar Programs, U.S. National Science Foundation-Physics Division, Wisconsin Alumni Research Foundation, Center for High Throughput Computing (CHTC) at the University of Wisconsin–Madison, Open Science Grid (OSG), Extreme Science and Engineering Discovery Environment (XSEDE), U.S. Department of Energy-National Energy Research Scientific Computing Center, the particle astrophysics research computing center at the University of Maryland, Institute for Cyber-Enabled Research at Michigan State University, and the astroparticle physics computational facility at Marquette University; Belgium – Funds for Scientific Research (FRS-FNRS and FWO), FWO Odysseus and Big Science programmes, and Belgian Federal Science Policy Office (Belspo); Germany – Bundesministerium für Bildung und Forschung (BMBF), Deutsche Forschungsgemeinschaft (DFG), Helmholtz Alliance for Astroparticle Physics (HAP), Initiative and Networking Fund of the Helmholtz Association, Deutsches Elektronen Synchrotron (DESY), and High Performance Computing cluster of the RWTH Aachen; Sweden – Swedish Research Council, Swedish Polar Research Secretariat, Swedish National Infrastructure for Computing (SNIC), and Knut and Alice Wallenberg Foundation; Australia – Australian Research Council; Canada – Natural Sciences and Engineering Research Council of Canada, Calcul Québec, Compute Ontario, Canada Foundation for Innovation, WestGrid, and Compute Canada; Denmark – Villum Fonden, Danish National Research Foundation (DNRF), Carlsberg Foundation; New Zealand – Marsden Fund; Japan – Japan Society for Promotion of Science (JSPS) and Institute for Global Prominent Research (IGPR) of Chiba University; Korea – National Research Foundation of Korea (NRF); Switzerland – Swiss National Science Foundation (SNSF); United Kingdom – Department of Physics, University of Oxford.