Search for Sources of Astrophysical Neutrinos Using Seven Years of IceCube Cascade Events

TL;DR

This paper reports on a seven-year search for astrophysical neutrino sources using IceCube cascade events, improving sensitivity and angular resolution, but finding no significant sources, thus setting new constraints in neutrino astronomy.

Contribution

The study extends previous work by adding five years of data, enhancing angular resolution, and testing for Galactic emission, making it the most sensitive cascade analysis to date.

Findings

No significant neutrino source clustering detected.

Constraints on astrophysical neutrino sources are the most stringent to date.

Analysis demonstrates improved sensitivity over previous IceCube cascade studies.

Abstract

Low background searches for astrophysical neutrino sources anywhere in the sky can be performed using cascade events induced by neutrinos of all flavors interacting in IceCube with energies as low as ~1 TeV. Previously, we showed that even with just two years of data, the resulting sensitivity to sources in the southern sky is competitive with IceCube and ANTARES analyses using muon tracks induced by charge current muon neutrino interactions - especially if the neutrino emission follows a soft energy spectrum or originates from an extended angular region. Here, we extend that work by adding five more years of data, significantly improving the cascade angular resolution, and including tests for point-like or diffuse Galactic emission to which this dataset is particularly well-suited. For many of the signal candidates considered, this analysis is the most sensitive of any experiment. No significant clustering was observed, and thus many of the resulting constraints are the most stringent to date. In this paper we will describe the improvements introduced in this analysis and discuss our results in the context of other recent work in neutrino astronomy.

Tables

Table 1. Table 1: Sensitivity and results of the SNR stacking analyses, compared to the previous analysis with tracks (Aartsen et al., 2017b ) . Sensitivity and ULs are given as E 2 ⋅ ( E / 100 ​ TeV ) 0.5 ⋅ d ​ N / d ​ E ⋅ superscript 𝐸 2 superscript 𝐸 100 TeV 0.5 𝑑 𝑁 𝑑 𝐸 E^{2}\cdot(E/100\,\mathrm{TeV})^{0.5}\cdot dN/dE in units 10 − 12 ​ TeV ​ cm − 2 ​ s superscript 10 12 TeV superscript cm 2 s 10^{-12}\,\mathrm{TeV\,cm^{-2}\,s} .

	7yr Cascades					7yr Tracks
Catalog	Sensitivity	p-value	$n_s$	$\gamma$	UL	p-value	$n_s$	$\gamma$	UL
SNR with mol. cloud	9.9	0.12	17.2	3.76	24	0.25	16.5	3.95	2.23
SNR with PWN	6.3	1	0	—	6.3	0.34	9.36	3.95	11.7
SNR alone	7.5	0.082	8.2	2.42	15	0.42	3.82	2.25	2.06

Table 2. Table 2: Sensitivity and results of the diffuse Galactic template analyses, compared to latest previous work: a joint IceCube-ANTARES (Albert et al., 2018 ) for KRA γ subscript KRA 𝛾 \text{KRA}_{\gamma}  models, and seven years of IceCube tracks (Aartsen et al., 2017b ) for Fermi -LAT π 0 superscript 𝜋 0 \pi^{0} decay. Sensitivity, fitted flux, and ULs are given as multiples of the model prediction for KRA γ subscript KRA 𝛾 \text{KRA}_{\gamma}  models, and as E 2 ⋅ ( E / 100 ​ TeV ) 0.5 ⋅ d ​ N / d ​ E ⋅ superscript 𝐸 2 superscript 𝐸 100 TeV 0.5 𝑑 𝑁 𝑑 𝐸 E^{2}\cdot(E/100\,\mathrm{TeV})^{0.5}\cdot dN/dE in units 10 − 11 ​ TeV ​ cm − 2 ​ s − 1 superscript 10 11 TeV superscript cm 2 superscript s 1 10^{-11}\,\mathrm{TeV\,cm^{-2}\,s^{-1}} for Fermi -LAT  π 0 superscript 𝜋 0 \pi^{0} decay.

	7yr Cascades				Previous Work
Template	p-value	Sensitivity	Fitted Flux	UL	p-value	Sensitivity	Fitted Flux	UL
${\text{KRA}}_{\gamma }^5$	0.021	0.58	0.85	1.7	0.29	0.81	0.47	1.19
${\text{KRA}}_{\gamma }^{50}$	0.022	0.35	0.65	0.97	0.26	0.57	0.37	0.90
Fermi-LAT ${\pi }^0$	0.030	2.5	3.3	6.6	0.37	2.97	1.28	3.83

Table 3. Table 3: Summary of the source catalog search. The type, common name, and equatorial coordinates (J2000) are shown for each object. Where non-null ( n ^ s > 0 subscript ^ 𝑛 𝑠 0 \hat{n}_{s}>0 ) results are found, the pre-trials significance p pre subscript 𝑝 pre p_{\text{pre}} and best-fit n ^ s subscript ^ 𝑛 𝑠 \hat{n}_{s} and γ ^ ^ 𝛾 \hat{\gamma} are given. ULs are expressed as E 2 ⋅ d ​ N / d ​ E ⋅ superscript 𝐸 2 𝑑 𝑁 𝑑 𝐸 E^{2}\cdot dN/dE , in units 10 − 12 ​ TeV superscript 10 12 TeV 10^{-12}\,\mathrm{TeV} , at E = 100 ​ TeV 𝐸 100 TeV E=100\,\mathrm{TeV} for unbroken E − 2 superscript 𝐸 2 E^{-2} and E − 3 superscript 𝐸 3 E^{-3} spectra ( Φ 2 subscript Φ 2 \Phi_{2} and Φ 3 subscript Φ 3 \Phi_{3} respectively) as well as at E ≪ 100 ​ TeV much-less-than 𝐸 100 TeV E\ll 100\,\mathrm{TeV} for a spectrum with d ​ N / d ​ E ∝ E − 2 ⋅ exp ⁡ ( E / 100 ​ TeV ) proportional-to 𝑑 𝑁 𝑑 𝐸 ⋅ superscript 𝐸 2 𝐸 100 TeV dN/dE\propto E^{-2}\cdot\exp(E/100\,\mathrm{TeV}) ( Φ 2 ​ C subscript Φ 2 𝐶 \Phi_{2C} ).

Type	Source	$\alpha {(}^{\circ })$	$\delta {(}^{\circ })$	$p_{\text{pre}}$	${\hat{n}}_s$	$\hat{\gamma }$	${\Phi }_2$	${\Phi }_3$	${\Phi }_{2C}$
BL Lac	PKS 2005-489	302.37	$-48.82$	0.222	7.0	3.8	5.3	4.1	15
	PKS 0537-441	84.71	$-44.09$	$\cdots$	0.0	$\cdots$	3.6	2.6	10
	PKS 0426-380	67.17	$-37.93$	$\cdots$	0.0	$\cdots$	3.6	2.7	10
	PKS 0548-322	87.67	$-32.27$	0.457	0.5	2.4	4.1	3.1	11
	H 2356-309	359.78	$-30.63$	0.452	0.0	$\cdots$	3.8	2.8	11
	PKS 2155-304	329.72	$-30.22$	0.452	0.0	$\cdots$	3.8	2.8	10
	1ES 1101-232	165.91	$-23.49$	0.030	3.6	2.3	9.2	7.4	25
	1ES 0347-121	57.35	$-11.99$	$\cdots$	0.0	$\cdots$	3.8	3.3	10
	PKS 0235+164	39.66	$16.62$	$\cdots$	0.4	3.3	5.6	3.6	11
	1ES 0229+200	38.20	$20.29$	0.459	0.0	$\cdots$	5.8	3.7	12
	W Comae	185.38	$28.23$	0.475	0.0	$\cdots$	6.0	3.4	11
	Mrk 421	166.11	$38.21$	0.373	0.0	$\cdots$	7.0	3.5	13
	Mrk 501	253.47	$39.76$	0.373	0.0	$\cdots$	7.1	3.4	13
	BL Lac	330.68	$42.28$	0.160	6.5	3.4	9.9	5.0	18
	H 1426+428	217.14	$42.67$	0.311	1.1	2.8	7.9	3.8	14
	3C66A	35.67	$43.04$	0.351	0.0	$\cdots$	7.4	3.5	13
	1ES 2344+514	356.77	$51.70$	0.119	7.5	4.0	13	5.5	23
	1ES 1959+650	300.00	$65.15$	0.137	6.1	4.0	20	5.2	30
	S5 0716+71	110.47	$71.34$	0.480	1.5	3.3	13	2.9	20
Flat Spectrum Radio Quasar	PKS 1454-354	224.36	$-35.65$	0.487	0.6	3.4	3.6	2.8	10
	PKS 1622-297	246.52	$-29.86$	0.315	4.2	4.0	4.8	3.7	13
	PKS 0454-234	74.27	$-23.43$	0.483	0.0	$\cdots$	3.6	2.9	9.9
	QSO 1730-130	263.26	$-13.08$	0.162	1.2	1.7	6.5	5.9	19
	PKS 0727-11	112.58	$-11.70$	0.293	11.1	3.6	5.5	4.8	15
	PKS 1406-076	212.23	$-7.87$	$\cdots$	0.0	$\cdots$	3.8	3.4	10
	QSO 2022-077	306.42	$-7.64$	$\cdots$	0.0	$\cdots$	3.8	3.3	10
	3C279	194.05	$-5.79$	$\cdots$	1.1	2.5	3.9	3.4	10
	3C 273	187.28	$2.05$	0.435	2.3	2.5	4.6	3.9	11
	PKS 1502+106	226.10	$10.49$	$\cdots$	2.7	3.8	5.3	3.7	11
	PKS 0528+134	82.73	$13.53$	$\cdots$	0.0	$\cdots$	5.4	3.7	12
	3C 454.3	343.49	$16.15$	0.288	1.9	2.1	8.0	5.3	17
	4C 38.41	248.81	$38.13$	0.373	0.0	$\cdots$	7.1	3.5	13
Galactic Center	Sgr A*	266.42	$-29.01$	0.357	2.2	3.0	4.5	3.5	12
HMXB/mqso	Cir X-1	230.17	$-57.17$	0.400	0.0	$\cdots$	3.7	2.5	11
	GX 339-4	255.70	$-48.79$	0.016	5.9	2.1	9.2	6.6	26
	LS 5039	276.56	$-14.83$	0.459	4.6	3.6	4.3	3.4	11
	SS433	287.96	$4.98$	0.011	30.9	3.1	14	10.0	33
	HESS J0632+057	98.25	$5.80$	$\cdots$	0.0	$\cdots$	4.7	3.4	11
	Cyg X-1	299.59	$35.20$	0.130	8.6	3.0	11	5.4	20
	Cyg X-3	308.11	$40.96$	0.150	7.7	3.2	11	5.0	19
	LSI 303	40.13	$61.23$	$\cdots$	0.0	$\cdots$	10	2.9	16
Massive Star Cluster	HESS J1614-518	63.58	$-51.82$	$\cdots$	0.0	$\cdots$	3.6	2.5	11
Not Identified	HESS J1507-622	226.72	$-62.34$	0.287	0.0	$\cdots$	4.1	2.8	12
	HESS J1503-582	226.46	$-58.74$	0.353	0.0	$\cdots$	3.9	2.7	11
	HESS J1741-302	265.25	$-30.20$	0.201	5.5	3.0	5.8	4.4	16
	HESS J1837-069	98.69	$-8.76$	0.470	4.3	3.4	3.9	3.5	10
	HESS J1834-087	278.69	$-8.76$	0.102	22.3	3.5	7.5	6.6	20
	MGRO J1908+06	286.98	$6.27$	0.018	28.3	3.0	14	9.6	32
Pulsar Wind Nebula	HESS J1356-645	209.00	$-64.50$	0.286	0.0	$\cdots$	3.8	2.8	12
	PSR B1259-63	197.55	$-63.52$	0.287	0.0	$\cdots$	4.0	2.8	12
	HESS J1303-631	195.74	$-63.20$	0.287	0.0	$\cdots$	4.0	2.8	12
	MSH 15-52	228.53	$-59.16$	0.353	0.0	$\cdots$	3.9	2.7	11
	HESS J1023-575	155.83	$-57.76$	0.096	4.7	4.0	5.7	4.4	17
	HESS J1616-508	243.78	$-51.40$	0.146	1.7	1.7	6.1	4.4	18
	HESS J1632-478	248.04	$-47.82$	0.044	3.8	2.0	8.3	6.0	24
	Vela X	128.75	$-45.60$	$\cdots$	0.1	2.0	3.8	2.7	11
	Geminga	98.48	$17.77$	$\cdots$	0.0	$\cdots$	5.5	3.7	11
	Crab Nebula	83.63	$22.01$	0.461	0.0	$\cdots$	6.0	3.6	12
	MGRO J2019+37	305.22	$36.83$	0.182	6.8	3.0	9.8	4.9	18
Seyfert Galaxy	ESO 139-G12	264.41	$-59.94$	0.247	1.6	2.6	4.6	3.3	13
Star Formation Region	Cyg OB2	308.08	$41.51$	0.144	8.0	3.2	11	5.0	19
Starburst/Radio Galaxy	Cen A	201.36	$-43.02$	$\cdots$	0.0	$\cdots$	3.7	2.7	10
	M87	187.71	$12.39$	0.305	3.2	2.4	7.6	5.2	17
	3C 123.0	69.27	$29.67$	0.302	1.0	2.2	8.0	4.7	16
	Cyg A	299.87	$40.73$	0.050	11.2	3.1	13	6.4	24
	NGC 1275	49.95	$41.51$	0.361	0.0	$\cdots$	7.6	3.5	13
	M82	148.97	$69.68$	0.265	3.4	3.2	19	4.2	28
Supernova Remnant	RCW 86	220.68	$-62.48$	0.287	0.0	$\cdots$	4.1	2.8	12
	RX J0852.0-4622	133.00	$-46.37$	$\cdots$	0.0	$\cdots$	3.7	2.5	11
	$\dagger$$\dagger$footnotemark: RX J1713.7-3946	258.25	$-39.75$	0.005	10.8	2.5	11	8.6	32
	W28	270.43	$-23.34$	0.238	0.8	1.6	5.6	4.7	16
	IC443	94.18	$22.53$	0.461	0.0	$\cdots$	6.1	3.7	12
	Cas A	350.85	$58.81$	0.028	12.4	4.0	24	7.0	38
	TYCHO	6.36	$64.18$	0.069	9.5	3.7	22	6.0	34

Equations

$$\mathcal{L}(n_s,\gamma )=\prod _i\left[\frac{n_s}{N}{\mathcal{S}}_i({\vec{x}}_i|\gamma )+\left(1-\frac{n_s}{N}\right){\mathcal{B}}_i({\vec{x}}_i)\right],$$

Full text

astroparticle physics — neutrinos

††thanks: also at Università di Padova, I-35131 Padova, Italy††thanks: also at National Research Nuclear University, Moscow Engineering Physics Institute (MEPhI), Moscow 115409, Russia††thanks: Earthquake Research Institute, University of Tokyo, Bunkyo, Tokyo 113-0032, Japan

1 Introduction

Neutrino astronomy promises to reveal secrets of distant astrophysical objects that likely can never be observed through other messenger particles. Because neutrinos only interact weakly, they can reach us from enormous distances with no attenuation by intervening matter or background radiation and without deflection by magnetic fields. Because they are only produced by hadronic processes, high energy neutrinos are tracers of high energy cosmic ray production (Halzen & Hooper, 2002). While electromagnetic observations can establish that a source candidate provides sufficient energy density for cosmic ray acceleration, direct cosmic ray observation is hindered by magnetic deflection at lower energies and by attenuation at higher energies. Therefore neutrino astronomy may offer our best chance for identifying the sources of high energy cosmic rays (Ahlers & Halzen, 2018).

Neutrino observation is performed by detecting the Cherenkov radiation emitted by relativistic charged particles produced when neutrinos collide with matter in or near a Cherenkov detector. IceCube, the largest such detector to date, consists of an array of photomultiplier tubes (PMTs) spanning one $\mathrm{km^{3}}$ deep in the Antarctic glacial ice near the geographic South Pole. IceCube is sensitive to all neutrino flavors and interaction types. Charged current (CC) muon neutrino interactions yield long-lived muons that can travel several kilometers through the ice (Chirkin & Rhode, 2004), leading to a track signature in the detector. Neutral current (NC) interactions, and CC interactions of most other flavors, yield hadronic and electromagnetic showers that typically range less than 20 m (Aartsen et al., 2014a), with 90% of the light emitted within 4 m of the shower maximum (Radel & Wiebusch, 2013). The small spatial extent of these showers compared to the PMT spacing and the scattering length of light in the ice (Aartsen et al., 2013b) results in a nearly symmetric cascade signature in the detector.

In 2014, we reported the first observation of a flux of neutrinos above $\sim 60\,\mathrm{TeV}$ inconsistent with the expectation from atmospheric backgrounds at greater than $5\sigma$ significance (Aartsen et al., 2014b). While this measurement was dominated by cascade events, the result was soon confirmed using muon tracks above $\sim 300\,\mathrm{TeV}$ originating in the northern sky (Aartsen et al., 2015b, 2016b).

More recently, IceCube data revealed the first direct evidence for high energy neutrino emission associated with a specific astrophysical source, the gamma-ray blazar TXS 0506+056 (Aartsen et al., 2018a, b). Before and since, no other high energy astrophysical neutrino sources have been identified (e.g. Aartsen et al., 2017a). Most source searches have focused on the muon track channel, which gives excellent sensitivity to upgoing muon tracks induced by CC muon neutrino interactions. As viewed by IceCube, upgoing events correspond to sources in the northern celestial hemisphere.

In much of the southern sky, due to larger background rates, the sensitivity of the muon track channel to sources following an $E^{-2}$ spectrum is weaker by an order of magnitude (Aartsen et al., 2017a) — this factor increases to two orders of magnitude or more if the spectrum is as soft as $E^{-3}$ or if it has a cutoff at $E_{\text{cut}}\lesssim 100\,\mathrm{TeV}$ (see e.g. Aartsen et al. (2017d)).

In an initial analysis of two years of data, we demonstrated that the sensitivity of IceCube in the southern sky can be improved significantly by performing complimentary searches using cascade events arising from neutrino interactions of all flavors (Aartsen et al., 2017d). Here, we extend that work in a number of ways. First, we apply similar, though slightly improved, event selection criteria to seven years of data. Second, we obtain significantly improved cascade angular resolution through the use of a specially-designed Deep Neural Network. Finally, we study additional point-like and diffuse Galactic emission scenarios to which this analysis is expected to be especially sensitive. For many of the signal candidates considered, this analysis is the most sensitive of any experiment to date. In this paper, we will begin by describing the IceCube detector and the cascade event selection and reconstruction. Then we will introduce the astrophysical neutrino source candidates considered and the design and performance characteristics of the statistical methods used. Finally, we will present our results and discuss them in the context of other recent work in neutrino astronomy.

2 IceCube

The IceCube detector (Aartsen et al. (2017c)) is composed of 5160 Digital Optical Modules (DOMs) buried at depths of 1450 m to 2450 m in the glacial ice near the geographic South Pole. Each DOM includes a 10” photomultiplier tube (PMT) and custom supporting electronics (Abbasi et al., 2010). The DOMs are mounted on 86 vertical strings holding 60 DOMs each, arranged in an approximately hexagonal grid. Seventy-eight of the strings forming the bulk of the array are spaced 125 m apart horizontally, with uniform vertical DOM spacing of $\sim 17\,\mathrm{m}$. The remaining 8 strings, which are concentrated near the center of the detector with $30-60\,\mathrm{m}$ horizontal spacing, constitute a denser in-fill array called DeepCore (Abbasi et al., 2012). On each of the DeepCore strings, 50 of the DOMs are located in the exceptionally clear ice at depths of 2100 m to 2450 m, with vertical spacing of $7\,\mathrm{m}$. The strings were deployed during the Austral summers of 2004–2011.

Digital readouts are triggered when at least eight DOMs observe a signal above 1/4 of the mean expected voltage from a single photoelectron (PE), each in coincidence with such a signal on a nearest or next-nearest neighboring DOM, within a $6.4\,\mathrm{\mu s}$ time window. When this criterion is met, the data acquisition system (DAQ, Abbasi et al. (2009)) collects the data from all DOMs into an event and initiates a first round of processing. Each waveform is decomposed into series of pulse arrival times and PE counts for use by event reconstruction algorithms (Ahrens et al., 2004; Aartsen et al., 2014a). Simple selection criteria are applied to reject the most unambiguous cosmic ray-induced muon backgrounds, reducing the data rate from $\sim 2.7\,\mathrm{kHz}$ at trigger level to $\sim 40\,\mathrm{Hz}$ at filter level. The filtered dataset is compressed and transmitted via satellite to a data center in the north for further processing.

3 Dataset

After the initial selection applied at the South Pole, the remaining dataset is still dominated by atmospheric muons. In order to search for neutrino sources, neutrino candidates are selected, and their properties are reconstructed based on the light arrival pattern observed in the DOMs. In the following, we discuss a re-optimized method for selecting neutrino-induced cascades and a novel machine learning-based approach to reconstructing their arrival directions and energies.

3.1 Event Selection

The procedure for rejecting the atmospheric muon background depends on the event topology of interest. Neutrino-induced muon tracks with energies $\gtrsim 1\,\mathrm{TeV}$ originating in the northern sky can be selected with high efficiency and low atmospheric muon contamination by identifying events reconstructed at declinations $\delta\gtrsim 5^{\circ}$ with high confidence, as only neutrinos can travel through so much intervening earth and/or ice before producing muons that pass through the detector. Neutrino- and cosmic ray-induced muon tracks originating in the southern sky and entering the detector from above can only be distinguished probabilistically, and only under the assumption that the neutrino spectrum is harder than the atmospheric muon spectrum. Thus the energy threshold increases to $\sim 100\,\mathrm{TeV}$ in the southern sky, resulting in weaker sensitivity especially for a soft neutrino spectrum.

In this work, we instead turn our attention to cascade events produced when the neutrino interaction vertex, and hence first observed light, occurs inside the detector. With this approach we accept all neutrino flavors and most interaction types, approximately independent of declination, while efficiently rejecting downgoing atmospheric muons. An added benefit for astrophysical neutrino searches is that for declinations $\lesssim-30^{\circ}$ the atmospheric neutrino background is naturally suppressed because many are accompanied by incoming atmospheric muons originally produced in the same cosmic ray shower in the upper atmosphere (Schonert et al., 2009).

Most Cherenkov light from a muon traveling through ice is radiated through stochastic processes, resulting in a dense, linear series of cascade-like signatures that may be observed in our detector. The mean distance between these energy deposits decreases with increasing energy. For energies $\gtrsim 60\,\mathrm{TeV}$, incoming muons can be rejected with high confidence using a veto region consisting of just the outermost DOMs, reserving the majority of the instrumented region as a fiducial volume for neutrino detection (Aartsen et al., 2014b). To lower the threshold to $\sim 1\,\mathrm{TeV}$ while holding the incoming muon rejection rate constant, the thickness of the veto region must be increased. Below we summarize this method, which is used as described in Wandkowsky & Weaver (2018) and which further optimizes the approach first introduced in Aartsen et al. (2015a).

We begin with all events passing one or more basic filters at the South Pole. A splitting algorithm is applied to each event, identifying $\sim 75\%$ of unrelated but temporally coincident physical events initially merged in the DAQ output by clustering causally connected sets of pulses. We reject any event in which the first $\geq 3$ pulses appear in the outer layer veto region as described in Aartsen et al. (2013a). An additional veto is applied to reject events in which two or more PE are observed consistent with a downgoing track passing through the interaction vertex or a major energy deposition. Finally, a cut is applied on the interaction vertex location, scaling with observed charge as described in Aartsen et al. (2015a) such that at 100 PE the fiducial volume is reduced to just the DeepCore sub-array, while at $\geq 6000\,\mathrm{PE}$ the fiducial volume consists of all but the outermost layer of DOMs. This final cut enables efficient background rejection down to $\sim 1\,\mathrm{TeV}$ by keeping the probability of observing veto photons approximately independent of energy.

We rely on a traditional maximum likelihood method (Aartsen et al., 2014a) to obtain initial reconstructions used for cascade/track discrimination. The goal of this reconstruction is to unfold the spatial and temporal pattern of energy depositions for each event. Two fits are performed: one which is constrained to find a single dominant cascade-like energy deposition, and one which finds a linear combination of such energy depositions distributed along a possible muon track. Events in which at least $6000\,\mathrm{PE}$ were collected are classified as tracks if the free track fit finds at least two non-negligible depositions more than 500 m apart, or if the free track fit is associated with more charge than the single cascade fit. Events with less total collected light are classified as tracks if at least 1.5 PE are consistent with an outgoing muon track originating at the reconstructed interaction vertex (Wandkowsky & Weaver, 2018). All other events are classified as cascades and are used in the present analysis.

The selection criteria described above were applied to data taken from May 2010 to May 2017 as well as to neutrino and atmospheric muon Monte Carlo (MC) simulations used for performance estimates. The first year of data comes from the nearly-complete 79-string configuration while the remaining six years make use of the complete 86-string detector. In a total of 2428 days of IceCube livetime, 10422 events survive until cascade/track discrimination; of these, 1980 are identified as cascades. Note that while the dominant improvement in this dataset is the increase from two to seven years of data, the neutrino effective area is also enhanced by applying coincident event splitting and veto criteria to data from every initial South Pole filter. This increases the acceptance by 23% (67%) for a signal following an $E^{-2}$ ($E^{-3}$) spectrum.

From MC simulations, we find that 98% of truly cascade-like events which pass all selection criteria are correctly identified as such. The rate at which CC muon neutrino interactions are successfully classified as track events increases with energy as more light is produced by the outgoing muon. For a conventional atmospheric neutrino spectrum, 30% of the cascade channel consists of misclassified CC muon neutrino interactions; for an astrophysical spectrum following $E^{-2.5}$ or harder, this contribution reduces to 5% or less. This population of misclassified events results in a tolerable background at lower energies as well as a small signal contribution at higher energies.

Because muon track analyses specifically target events with high quality track reconstructions and reject events dominated by individual cascade-like energy depositions, we expect the cascade analysis to be largely statistically independent in spite of the small but nonzero misclassification rate. In fact, the final cascade selection shares just a single $\sim 2\,\mathrm{TeV}$ event in common with the latest muon track selection.

3.2 Event Reconstruction

In past work, we have used a maximum likelihood method to reconstruct neutrino energy and direction of travel from IceCube cascades (Aartsen et al., 2014a). This approach relies on detailed parameterizations of the position- and direction-dependent light absorption and scattering lengths in the ice, neither of which is large compared to the DOM spacing. This results in a complex multi-dimensional likelihood function with many local optima in the right ascension and declination coordinates $(\alpha,\delta)$, such that it is computationally expensive to find the global optimum for any given event and prohibitive to estimate the per-event statistical uncertainties.

In this work, we introduce a novel cascade reconstruction using a deep Neural Network (NN). A NN is a highly flexible function mapping from an input layer to an output layer via a series of hidden layers, where each successive layer consists of a set of values computed based on the values contained in the previous layer. The functional forms of the layer-to-layer connections (the network architecture) must be designed a priori; the numerical parameters of those connections are optimized through a training procedure to yield good results for a given training dataset. NNs are well-suited to problems in high energy physics for which we are typically able to generate high-statistics MC datasets for use in training.

Our NN-based reconstruction draws from recent advances in image recognition and is implemented using Tensorflow (Abadi et al., 2015). The network architecture used here is largely the same as one introduced previously for muon energy reconstruction (Huennefeld, 2018). The method will be described in detail in a separate publication, but here we will outline the main considerations relevant in this analysis.

IceCube data consists of a set of waveforms (represented as a series of pulse arrival times and PE counts) accumulated over time on a number of DOMs distributed throughout the three-dimensional instrumented volume, and thus is in general four-dimensional. Our first step is to compute waveform summary values for use in the input layer. For each DOM, these values consist of the relative time of the first pulse; the time elapsed until 20%, 50%, and 100% of the total charge is collected; the total charge collected; the charge collected within 100 ns and 500 ns of the first pulse; and the charge-weighted mean and standard deviation of relative pulse arrival times.

The detector is divided into three sub-arrays: IceCube, lower DeepCore, and upper DeepCore. Each sub-array is independently well-approximated by a regular spatial grid suitable for processing by several initial convolutional layers, which are able to exploit symmetries in the structure of the input data to facilitate efficient network optimization and usage111Alternative methods are being developed to avoid the reliance on regular detector geometry. (see Huennefeld (2018) for diagrams of the relevant geometry). The output from the convolutional layers is taken as the input for each of two fully-connected neural networks (in which each node in a given layer is connected to every node in the preceding layer). One of these networks is optimized to estimate the physical parameters of interest — the right ascension, declination, and energy $(\alpha,\delta,E)$ — while the other is optimized to estimate the uncertainties on these parameters.

All training was performed using 50% of the relevant signal MC, with the remaining 50% reserved for testing analysis-level performance. Two training passes were performed. The first pass made use of several MC datasets: one with baseline values for key parameters such as DOM quantum efficiency and light absorption and scattering lengths, and several more with modified values within estimated systematic uncertainties. In addition to offering overall increased training statistics, the use of these differing datasets may give the NN some robustness against known systematic uncertainties. The second training pass refined the network to give the smallest errors and, on-average, unbiased reconstructions for the baseline MC. In each pass, a priori per-parameter weighting was applied such that angular resolution is valued over energy resolution by a factor of 5.

The expected performance of the NN angular reconstruction (including systematics; see Section 5.2) is shown as a function of energy in Figure 1. Compared to the reconstructions used in our previous analysis of two years of data (Aartsen et al., 2017d), the NN offers significantly improved angular resolution above 10 TeV (a factor of 2 improvement at 1 PeV). While we do not recover the optimal statistics-limited angular resolution described in Aartsen et al. (2014a), we do obtain performance that improves monotonically with increasing energy up to $\sim 1\,\mathrm{PeV}$. At higher energies, the estimated systematic uncertainty becomes large enough to prevent any further improvement. Note that an additional advantage of the NN angular reconstruction used here is that it naturally provides per-event uncertainty estimates usable in the statistical analysis described in Section 5.1, whereas previous work relied on a parameterization of typical uncertainties derived from signal MC.

The performance of the energy reconstruction is comparable to that used in previous work. The estimated energy is within 60% of the true neutrino energy for 68% of events, averaged over all neutrino flavors and interaction types, and approximately independent of spectrum. This performance estimate, like the sensitivities quoted in Section 5.3, assumes a flavor ratio of 1:1:1 with equal contributions from $\nu$ and $\bar{\nu}$, detected via a mixture of CC and NC interactions.

The energy and declination distributions of cascade events in data are compared with neutrino and atmospheric muon MC in Figure 2. The distributions obtained are similar to those observed in the two year sample (Aartsen et al., 2015a, 2017d).

4 Source Candidates

In this work, we search for neutrino emission from a number of Galactic and extra-Galactic source candidates. Each candidate has been studied previously by IceCube, by ANTARES, a neutrino observatory located deep in the Mediterranean sea (Ageron et al., 2011), or by both, such that direct comparisons can be drawn between the results presented here and past work using IceCube tracks and all interaction flavors in ANTARES. In this section, we outline the neutrino emission scenarios that we have considered.

4.1 Point-like Source Candidates

One way to search for astrophysical neutrino sources with only a minimal set of a priori assumptions about source position is to search the entire sky for the most significant point-like neutrino clustering in excess of the background expectation on a dense grid of pixels that are small compared to the neutrino angular resolution. This approach has most recently been employed by IceCube using tracks (Aartsen et al., 2017a) and cascades (Aartsen et al., 2017d) as well as by ANTARES using tracks and cascades in combination (Albert et al., 2017a), and we include it in the present analysis as well. However, an all-sky scan is subject to a large trial factor and thus is in general less sensitive compared to analyses that use prior information to restrict the set of hypothesis tests.

An alternative approach is to scan only the positions of a modest number of well-motivated source candidates, which substantially reduces the trial factor. In addition, where multiple analyses report results for the same or overlapping catalogs, direct comparisons can be made. Here we scan the same catalog of 74 source candidates that was studied in the previous IceCube cascade paper (Aartsen et al., 2017d).

We consider one source in more detail: the supermassive black hole at the center of the Galaxy, Sagitarius A*. Based on hints from gamma-ray observations (e.g. Herold & Malyshev, 2019), there may be emission up to some unknown high energy cutoff from a spatially extended region centered approximately on this object. Therefore we evaluate constraints on the flux from this region as a function of possible spatial extension and for several possible spectral cutoffs.

The gamma-ray blazar TXS 0506+056 does not appear in the a priori catalog described above. In light of this, and in anticipation of future identifications of unexpectedly promising source candidates based on neutrino observations, we treat this object as a monitored source to be studied separately from the catalog scan described above.

For source classes for which we can predict approximate relative signal strengths, it may be possible to increase the signal-to-background ratio using a source-stacking method (e.g. Abbasi et al., 2011). Because the present analysis offers good sensitivity in the southern sky, roughly independent of possible spatial extension up to a few degrees, we include stacking analyses for three Galactic supernova remnant (SNR) catalogs derived from SNR Cat (Ferrand & Safi-Harb, 2012) and previously studied using IceCube tracks (Aartsen et al., 2017b). These SNRs are categorized based on their environment: those with associated molecular clouds, those with associated pulsar wind nebulae (PWN), and those with neither. The angular extension of these objects reach up to $1.63^{\circ}$, and each catalog comprises a preponderance of objects in the southern sky.

4.2 Diffuse Galactic Emission

Cosmic ray interactions with interstellar gas in the Milky Way are expected to produce neutral and charged pions, where neutral pions would decay to observable gamma rays and charged pions would yield potentially observable neutrinos. The hadronic gamma-ray emission up to 100 GeV has been identified by Fermi-LAT using a multi-component fit (Ackermann et al., 2012). A corresponding neutrino flux prediction can be obtained by extrapolating this measurement to energies above $1\,\mathrm{TeV}$ in the context of Galactic cosmic ray production and propagation models.

The original model fits by Fermi somewhat under-predict the measured gamma ray flux in the Galactic plane, and especially near the Galactic center, above $10\,\mathrm{GeV}$. The $\text{KRA}_{\gamma}$ models obtain better agreement with gamma ray data in this regime by introducing galactocentric cosmic ray diffusion parameter variability and an advective wind (Gaggero et al., 2015, 2017). Model-dependent neutrino flux predictions are provided assuming cosmic ray injection spectra with exponential cutoffs at 5 PeV or 50 PeV per nucleon; we refer to these as $\text{KRA}_{\gamma}^{5}$ and $\text{KRA}_{\gamma}^{50}$, respectively.

The latest constraints on diffuse Galactic neutrino emission depend on the $\text{KRA}_{\gamma}$ models and were obtained in a joint IceCube and ANTARES analysis (Albert et al., 2018) which made use of complimentary features of the IceCube track analysis (Aartsen et al., 2017b) and ANTARES track and cascade combined analysis (Albert et al., 2017b). In this work we search for emission following $\text{KRA}_{\gamma}^{5}$ as the primary diffuse Galactic emission result; we also test for emission following $\text{KRA}_{\gamma}^{50}$. Finally, we test for emission following the spatial profile of the Fermi-LAT $\pi^{0}$-decay measurement, assuming an $E^{-2.5}$ neutrino energy spectrum.

4.3 Fermi Bubbles

The Fermi bubbles consist of a pair of gamma ray emission regions that extend to $\sim 55^{\circ}$ above and below the Galactic center (Su et al., 2010). Most of the Fermi bubble region yields a relatively hard gamma ray spectrum up to $\sim 100\,\mathrm{GeV}$, with some evidence for spectral softening above that energy (Ackermann et al., 2014). The gamma-ray emission has been speculated to be of hadronic origin (Crocker & Aharonian, 2011), powered by cosmic ray acceleration in the vicinity of the Galactic Center; however, the true origin of the Fermi bubbles has not yet been experimentally identified.

We derive constraints on emission from the Fermi bubbles following spectra of the form $dN/dE\propto E^{-2.18}\cdot\exp(-E/E_{\text{cut}})$, for $E_{\text{cut}}\in\left\{50\,\mathrm{TeV},100\,\mathrm{TeV},500\,\mathrm{TeV}\right\}$ — the same spectra tested in recent work by ANTARES (Hallmann & Eberl, 2018). If there is neutrino emission from the Fermi bubbles with a significantly softer spectrum or lower cutoff energy, this analysis would not be sensitive to it.

5 Analysis Methods and Performance

The source searches described in the previous section use established methods from recent IceCube work. In this section, we review the statistical methods and describe the systematic uncertainty treatment applied here. Then we discuss the sensitivity of this analysis to the source candidates under consideration.

5.1 Statistical Methods

In this work we consider two broad categories of source candidates: point-like and extended template, where the latter include diffuse Galactic emission and emission spanning the Fermi Bubbles. Both analysis types are based on the standard likelihood (Braun et al., 2008) given by a product over all events $i$ in the dataset:

[TABLE]

where $N$ is the total number of events; $n_{s}$ is the expected number of signal events; $\gamma$ is the signal spectral index; $\vec{x}_{i}$ represents the event right ascension, declination, angular uncertainty, and energy $\left\{\alpha_{i},\delta_{i},\sigma_{i},E_{i}\right\}$; $\mathcal{S}_{i}(\vec{x}_{i}|\gamma)$ is the probability density function (PDF) assuming event $i$ is part of the signal population; and $\mathcal{B}_{i}(\vec{x}_{i})$ is the PDF assuming event $i$ is part of the atmospheric or unrelated astrophysical background populations. For all source types, $n_{s}$ is free to vary between 0 and $N$. For point-like sources, the signal spectral index $\gamma$ is free to vary between 1 and 4, while for extended templates $\gamma$ is fixed to a source-dependent constant value ($\gamma=2.5$ for diffuse Galactic emission and $\gamma=2.18$ for emission from the Fermi bubbles).

The details of our signal and background likelihoods, $\mathcal{S}_{i}$ and $\mathcal{B}_{i}$, follow established methods applied previously to IceCube tracks for individual (Aartsen et al., 2017a) and stacked (e.g. Abbasi et al., 2011) point-like sources as well as for extended templates (Aartsen et al., 2017b). We do not require a specialized treatment, in contrast to our previous cascade analysis (Aartsen et al., 2017d), thanks to increased statistics in the experimental dataset as well as new per-event angular uncertainty estimates given by the NN reconstruction.

As in previous work, we define the test statistic as the log likelihood ratio $\mathcal{T}=-2\ln\{\mathcal{L}(n_{s}=0)/\mathcal{L}(\hat{n}_{s},\hat{\gamma})\}$, where $\hat{n}_{s}$ and $\hat{\gamma}$ are the values which maximize $\mathcal{L}$, subject to the constraints specified above. This test statistic is used to compute significances, sensitivities, discovery potentials, and upper limits (ULs). For the all-sky (source candidate catalog) scan, we compute a post-trials significance based on the most significant pixel (source candidate) tested, in order to guarantee the reported false positive rates. Sensitivities (90% CL), upper limits (90% CL), and discovery potentials ($5\sigma$) are defined as in our previous analysis (Aartsen et al., 2017d) and are computed using the Neyman construction (Neyman, 1937).

5.2 Systematic Uncertainties

The dominant systematic uncertainties in this analysis include the optical properties of the ice, the quantum efficiency of the DOMs, and the neutrino interaction cross section. These uncertainties affect the angular resolution and the signal acceptance. As in our previous cascade analysis (Aartsen et al., 2017d), we treat these effects as approximately separable. However, we have improved our approach to each consideration; we describe our latest method in the following.

The NN reconstruction is trained to yield optimal performance on baseline MC; the angular resolution for real data events will be somewhat worse. To estimate how much worse, we perform dedicated simulations of events similar to those observed, but using depth-dependent ice model variations intended to cover the uncertainties in the model. By comparing the median resolution from these modified simulations with that from the baseline MC, we obtain a function of energy that quantifies how much worse the resolution may be than expected from the baseline. This factor ranges from 10% at 1 TeV to $\sim 50\%$ at 2 PeV, and is taken as a correction to the angular separation between the reconstructed and true direction for each event in the baseline MC. This factor is similarly applied to the angular uncertainty estimates $\sigma_{i}$ for both MC and data events. In this way, we directly account for systematic uncertainties impacting angular resolution in the quantiles shown in Figure 1 as well as in all p-values and sensitivity flux calculations in the analysis.

The above treatment accounts for the analysis-level impact of systematic uncertainties for each observed event. To address the uncertainties in the detection efficiency, and thus in sensitivity, discovery potential, and upper limit fluxes, we compute the energy-integrated signal acceptance, as a function of declination and for each considered spectrum, based on additional MC datasets produced with varied modeling assumptions (the same modified datasets used in NN training; see Section 3.2). We find that for plausible ice model and detector variations, the signal acceptance variation ranges from $\sim 10\%$ for an unbroken $E^{-2}$ spectrum to $\sim 17\%$ for $E^{-2}$ with an exponential cutoff at 100 TeV, roughly independent of declination. As was done in the previous analysis, we estimate an uncorrelated $4\%$ impact from uncertainties in the neutrino interaction cross section. These values are added in quadrature on a per-spectrum basis to obtain a final estimate of uncertainties via signal acceptance effects. In the remainder of this paper, all sensitivity, discovery potential, and upper limit fluxes include this factor.

5.3 Sensitivity

All sensitivities discussed in the remainder of this paper are per-neutrino flavor (assuming a flavor ratio of 1:1:1 at the detector), but summed over $\nu$ and $\bar{\nu}$. The point source sensitivity flux as a function of source declination is shown in Figure 3 for several spectral scenarios: unbroken power laws following hard ($\gamma=2$) and soft ($\gamma=3$) spectra, and spectral cutoff scenarios $dN/dE\propto E^{-2}\cdot\exp(-E/E_{\text{cut}})$ with $E_{\text{cut}}\in(100\,\mathrm{TeV},\,1\,\mathrm{PeV},\,\gg 1\,\mathrm{PeV})$. Where published values are available for previous IceCube work with tracks (Aartsen et al., 2017a) or cascades (Aartsen et al., 2017d), or for the most recent ANTARES track and cascade combined analysis (Albert et al., 2017a), these are shown for comparison. We find that the present analysis improves upon the previous IceCube work with cascades at all declinations and across the tested spectra, with the largest improvements reaching a factor larger than 4 in the southern sky. Furthermore, we now obtain the best sensitivity of any analysis for hard sources in the southern-most $\sim 30\%$ of the sky ($\sin(\delta)<-0.4$). This search also achieves sensitivity comparable to that of ANTARES for spectra with cutoffs as low as $E_{\text{cut}}=100\,\mathrm{TeV}$, but with much weaker declination dependence.

The sensitivities of the SNR stacking analyses are listed in Table 1. In this work we obtain a sensitivity below previously set ULs (Aartsen et al., 2017b) only for the SNR-with-PWN catalog, which consists of eight southern SNRs and one northern SNR. It is nevertheless interesting to revisit all three catalogs here because, while they all include southern source candidates, in previous work the results necessarily were dominated by northern candidates due to the strongly declination-dependent signal acceptance of the IceCube track selection.

The sensitivities of the diffuse Galactic template analyses are listed in Table. 2. This analysis obtains $\sim 30\%\ (40\%)$ better sensitivity to $\text{KRA}_{\gamma}^{5}$ ($\text{KRA}_{\gamma}^{50}$) than the recent joint IceCube+ANTARES analysis (Albert et al., 2018). Compared to the IceCube analysis using seven years of tracks (Aartsen et al., 2017b), this analysis obtains $\sim 15\%$ better sensitivity to emission following the spatial profile of the Fermi-LAT $\pi^{0}$-decay measurement. These improvements are possible because the expected emission follows a soft ($\gamma\sim 2.5$) spectrum and is concentrated near the Galactic center at $\delta\sim-30^{\circ}$, where IceCube track analyses are subject to a large background of atmospheric muons but the present cascade analysis efficiently rejects this background as well as some of the atmospheric neutrino background; the improvement is larger for the $\text{KRA}_{\gamma}$ models than for the Fermi-LAT $\pi^{0}$ model because the former are specifically tuned to increase the concentration of the expected flux near the Galactic center.

The sensitivity flux for the Fermi Bubble analyses is $\sim 30\%$ below the upper limits shown in Figure 8, approximately independent of spectral cutoff. This analysis obtains sensitivity that is at least one order of magnitude better than the recent ANTARES search (Hallmann & Eberl, 2018), with the improvement increasing with spectral cutoff energy, $E_{\text{cut}}$. Because we assume an even more extended template than ANTARES, covering a total solid angle of about 1.18 sr compared to $\sim 0.66\,\mathrm{sr}$, this factor is even larger if considered in terms of flux per solid angle. Once again, this improvement is due to efficient rejection of atmospheric backgrounds for the cascade dataset used in this work.

6 Results

The result of the unbiased all-sky scan is shown in Figure 5. The most significant source candidate was found at $(\alpha,\delta)=(271.23^{\circ},7.78^{\circ})$ with a pre-trial p-value of $1.8\times 10^{-3}$ ($2.9\sigma$), corresponding to a post-trial p-value of 0.69.

The results of the source candidate catalog scan are tabulated Table 3. The most significant source was RX J1713.7-3946, a well-known SNR that is also included in the SNR-alone catalog. For this source candidate we found a pre-trial p-value of $5.0\times 10^{-3}$ ($2.6\sigma$), corresponding to a post-trial p-value of 0.28. Flux upper limits for each source are plotted, along with the sensitivity and $5\sigma$ discovery potential of this analysis, in Figure 4 as a function of source declination for each of the benchmark point source spectra discussed in the previous section. For the one monitored source, TXS 0506+056, we find $n_{s}=0$. Note that the measured flux for TXS 0506+056 is just $E^{2}\cdot dN/dE\sim 10^{-12}\,\mathrm{TeV\,cm^{-2}\,s^{-1}}$, or about $5\times$ lower than the cascade sensitivity at $\delta=5.69^{\circ}$, and thus the null result we find here is consistent with previous results (Aartsen et al., 2018b).

We set constraints on extended emission in the vicinity of the supermassive black hole at the center of the Galaxy, Sagitarius A∗, in Figure 6. For this object we find a small but non-zero best fit ($p_{\text{pre}}=0.357$). We then compute ULs, assuming a spectrum of the form $dN/dE\propto E^{-2}\cdot\exp(E/E_{\text{cut}})$ for various choices of $E_{\text{cut}}$, as a function of possible Gaussian source extension, $\sigma_{\text{Sgr~{}A}^{*}}\in[0,5^{\circ}]$. In these calculations, we include the source extension only in the signal simulation but not in the likelihood test. The relative independence of this result with respect to assumed source extension underscores the importance of atmospheric background rejection at the event selection level, relative to per-event angular reconstruction, in the overall performance of this analysis.

The results of the SNR stacking analyses are shown in Table 1. We find $n_{s}=0$ for SNR with PWN and mild excesses for the other two catalogs, the most significant of which is an excess with $p=0.082$ for SNR alone. The SNR-with-PWN category is the only one for which this analysis finds a sensitivity flux below the previous UL from the track analysis (Aartsen et al., 2017b); the UL found here constitutes a reduction of $\sim 50\%$.

The results of the diffuse Galactic extended template analyses are shown in Table 2. The primary hypothesis test, for emission following the $\text{KRA}_{\gamma}^{5}$ model, was also the most significant with a p-value of 0.021 ($2.0\sigma$) and a best-fit flux222Note that fitted fluxes, unlike ULs, are central values and are thus not subject to the penalty factors described in 5.2 of $0.85\times\text{KRA}_{\gamma}^{5}$. The best-fit fluxes for each template are consistent with ULs set by previous work (Aartsen et al., 2017b; Albert et al., 2018).

Prior to this analysis, the most significant ($1.5\sigma$) indication for diffuse Galactic emission came from an IceCube analysis using a spatially-binned method and only events originating in the northern sky in order to constrain the spectrum of possible emission following the Fermi-LAT $\pi^{0}$ template (Aartsen et al., 2017b). As an a posteriori test, we extend the template analysis described in Section 5.1 to include the spectral index $\gamma$ as a free parameter. A 2D scan of the resulting likelihood for the Fermi-LAT $\pi^{0}$ model is shown in Figure 7, with contours from the spatially-binned track analysis shown for comparison. In both analyses, the best fit is obtained for a harder spectrum close to $\gamma=2$, with both normalization and spectral index consistent within less than $1\sigma$. These independent results would remain statistically insignificant even under a combined analysis. Nevertheless, they are consistent with each other and with a possible astrophysical signal, potentially imperfectly tracing the spatial dependence prescribed by the $\text{KRA}_{\gamma}$ and Fermi-LAT $\pi^{0}$ models, at a level only starting to approach the reach of existing detectors and methods.

For emission from the Fermi bubbles, we obtain $n_{s}=5.2$, with a p-value of $0.30~{}(0.51\sigma)$. Flux upper limits based on these tests are shown in Figure 8. In the absence of significant emission, we set the most stringent limits to date on possible high energy neutrino emission from this intriguing structure.

7 Conclusion and Outlook

In this work, we apply a novel NN reconstruction to seven years of IceCube cascade data in order to search for high energy neutrino emission from a number of astrophysical source candidates. By improving the angular resolution and time-integrated signal acceptance with respect to our previous analysis using two years of data (Aartsen et al., 2017d), we obtain significant gains in sensitivity, with the best sensitivity of any experiment to date for sources concentrated in the southern sky. Nevertheless, we did not find significant evidence for emission from any of the sources considered.

While we have considered several neutrino source candidates, the ensemble of tests is far from exhaustive. We have begun to revisit multi-wavelength EM data in an effort to identify new catalogs of sources of interest for individual and stacking analyses. Furthermore, as in our previous paper (Aartsen et al., 2017d), we have still used IceCube cascades primarily in just time-integrated analyses. In future work we intend to explore time-dependent source candidates, including e.g. high-variability blazars as well as transients such as gravitational wave candidates reported by Advanced LIGO (Abbott et al., 2016). The NN reconstruction is especially promising for rapid follow-up of transient source candidates because once the NN is trained, compute time for the reconstruction is negligible.

In future work, we plan to revisit the event selection criteria. The selection used in this paper already achieves very good rejection of atmospheric backgrounds using explicit cuts on low-level parameters in the data. However, it is possible to improve the signal acceptance by including machine learning methods not only in the cascade reconstruction but in the event selection as well (e.g. Niederhausen & Xu, 2018).

Finally, we have deliberately attempted to maintain statistical independence between this analysis and others performed using IceCube tracks. We have separately developed multiple throughgoing (e.g. Aartsen et al., 2017a, 2016b) and starting (Aartsen et al., 2016a, 2019) track selections, each with differing energy- and declination-dependent background rates and signal acceptances. Combined analyses using tracks and cascades may offer the best sensitivity achievable using the existing IceCube detector alone. Joint IceCube–ANTARES analyses so far have not included IceCube cascades (Adrian-Martinez et al. (2016), updated results in preparation). All-flavor, multi-detector analysis will likely give the best possible sensitivity in a future analysis.

Acknowledgements: The IceCube collaboration acknowledges the significant contributions to this manuscript from Michael Richman. 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, Particle astrophysics research computing center at the University of Maryland, Institute for Cyber-Enabled Research at Michigan State University, and 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).

\startlongtable