∎

\thankstext

emailE-Mail: analysis@icecube.wisc.edu \thankstextaEarthquake Research Institute, University of Tokyo, Bunkyo, Tokyo 113-0032, Japan ¹¹institutetext: III. Physikalisches Institut, RWTH Aachen University, D-52056 Aachen, Germany ²²institutetext: Department of Physics, University of Adelaide, Adelaide, 5005, Australia ³³institutetext: Dept. of Physics and Astronomy, University of Alaska Anchorage, 3211 Providence Dr., Anchorage, AK 99508, USA ⁴⁴institutetext: Dept. of Physics, University of Texas at Arlington, 502 Yates St., Science Hall Rm 108, Box 19059, Arlington, TX 76019, USA ⁵⁵institutetext: CTSPS, Clark-Atlanta University, Atlanta, GA 30314, USA ⁶⁶institutetext: School of Physics and Center for Relativistic Astrophysics, Georgia Institute of Technology, Atlanta, GA 30332, USA ⁷⁷institutetext: Dept. of Physics, Southern University, Baton Rouge, LA 70813, USA ⁸⁸institutetext: Dept. of Physics, University of California, Berkeley, CA 94720, USA ⁹⁹institutetext: Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA ¹⁰¹⁰institutetext: Institut für Physik, Humboldt-Universität zu Berlin, D-12489 Berlin, Germany ¹¹¹¹institutetext: Fakultät für Physik & Astronomie, Ruhr-Universität Bochum, D-44780 Bochum, Germany ¹²¹²institutetext: Université Libre de Bruxelles, Science Faculty CP230, B-1050 Brussels, Belgium ¹³¹³institutetext: Vrije Universiteit Brussel (VUB), Dienst ELEM, B-1050 Brussels, Belgium ¹⁴¹⁴institutetext: Dept. of Physics, Massachusetts Institute of Technology, Cambridge, MA 02139, USA ¹⁵¹⁵institutetext: Dept. of Physics and Institute for Global Prominent Research, Chiba University, Chiba 263-8522, Japan ¹⁶¹⁶institutetext: Dept. of Physics and Astronomy, University of Canterbury, Private Bag 4800, Christchurch, New Zealand ¹⁷¹⁷institutetext: Dept. of Physics, University of Maryland, College Park, MD 20742, USA ¹⁸¹⁸institutetext: Dept. of Physics and Center for Cosmology and Astro-Particle Physics, Ohio State University, Columbus, OH 43210, USA ¹⁹¹⁹institutetext: Dept. of Astronomy, Ohio State University, Columbus, OH 43210, USA ²⁰²⁰institutetext: Niels Bohr Institute, University of Copenhagen, DK-2100 Copenhagen, Denmark ²¹²¹institutetext: Dept. of Physics, TU Dortmund University, D-44221 Dortmund, Germany ²²²²institutetext: Dept. of Physics and Astronomy, Michigan State University, East Lansing, MI 48824, USA ²³²³institutetext: Dept. of Physics, University of Alberta, Edmonton, Alberta, Canada T6G 2E1 ²⁴²⁴institutetext: Erlangen Centre for Astroparticle Physics, Friedrich-Alexander-Universität Erlangen-Nürnberg, D-91058 Erlangen, Germany ²⁵²⁵institutetext: Département de physique nucléaire et corpusculaire, Université de Genève, CH-1211 Genève, Switzerland ²⁶²⁶institutetext: Dept. of Physics and Astronomy, University of Gent, B-9000 Gent, Belgium ²⁷²⁷institutetext: Dept. of Physics and Astronomy, University of California, Irvine, CA 92697, USA ²⁸²⁸institutetext: Dept. of Physics and Astronomy, University of Kansas, Lawrence, KS 66045, USA ²⁹²⁹institutetext: SNOLAB, 1039 Regional Road 24, Creighton Mine 9, Lively, ON, Canada P3Y 1N2 ³⁰³⁰institutetext: Department of Physics and Astronomy, UCLA, Los Angeles, CA 90095, USA ³¹³¹institutetext: Dept. of Astronomy, University of Wisconsin, Madison, WI 53706, USA ³²³²institutetext: Dept. of Physics and Wisconsin IceCube Particle Astrophysics Center, University of Wisconsin, Madison, WI 53706, USA ³³³³institutetext: Institute of Physics, University of Mainz, Staudinger Weg 7, D-55099 Mainz, Germany ³⁴³⁴institutetext: Department of Physics, Marquette University, Milwaukee, WI, 53201, USA ³⁵³⁵institutetext: Physik-department, Technische Universität München, D-85748 Garching, Germany ³⁶³⁶institutetext: Institut für Kernphysik, Westfälische Wilhelms-Universität Münster, D-48149 Münster, Germany ³⁷³⁷institutetext: Bartol Research Institute and Dept. of Physics and Astronomy, University of Delaware, Newark, DE 19716, USA ³⁸³⁸institutetext: Dept. of Physics, Yale University, New Haven, CT 06520, USA ³⁹³⁹institutetext: Dept. of Physics, University of Oxford, 1 Keble Road, Oxford OX1 3NP, UK ⁴⁰⁴⁰institutetext: Dept. of Physics, Drexel University, 3141 Chestnut Street, Philadelphia, PA 19104, USA ⁴¹⁴¹institutetext: Physics Department, South Dakota School of Mines and Technology, Rapid City, SD 57701, USA ⁴²⁴²institutetext: Dept. of Physics, University of Wisconsin, River Falls, WI 54022, USA ⁴³⁴³institutetext: Dept. of Physics and Astronomy, University of Rochester, Rochester, NY 14627, USA ⁴⁴⁴⁴institutetext: Oskar Klein Centre and Dept. of Physics, Stockholm University, SE-10691 Stockholm, Sweden ⁴⁵⁴⁵institutetext: Dept. of Physics and Astronomy, Stony Brook University, Stony Brook, NY 11794-3800, USA ⁴⁶⁴⁶institutetext: Dept. of Physics, Sungkyunkwan University, Suwon 440-746, Korea ⁴⁷⁴⁷institutetext: Dept. of Physics and Astronomy, University of Alabama, Tuscaloosa, AL 35487, USA ⁴⁸⁴⁸institutetext: Dept. of Astronomy and Astrophysics, Pennsylvania State University, University Park, PA 16802, USA ⁴⁹⁴⁹institutetext: Dept. of Physics, Pennsylvania State University, University Park, PA 16802, USA ⁵⁰⁵⁰institutetext: Dept. of Physics and Astronomy, Uppsala University, Box 516, S-75120 Uppsala, Sweden ⁵¹⁵¹institutetext: Dept. of Physics, University of Wuppertal, D-42119 Wuppertal, Germany ⁵²⁵²institutetext: DESY, D-15738 Zeuthen, Germany

Search for steady point-like sources in the astrophysical muon neutrino flux with 8 years of IceCube data

IceCube Collaboration\thanksrefemail: M. G. Aartsen\thanksrefChristchurch M. Ackermann\thanksrefZeuthen J. Adams\thanksrefChristchurch J. A. Aguilar\thanksrefBrusselsLibre M. Ahlers\thanksrefCopenhagen M. Ahrens\thanksrefStockholmOKC D. Altmann\thanksrefErlangen K. Andeen\thanksrefMarquette T. Anderson\thanksrefPennPhys I. Ansseau\thanksrefBrusselsLibre G. Anton\thanksrefErlangen C. Argüelles\thanksrefMIT J. Auffenberg\thanksrefAachen S. Axani\thanksrefMIT P. Backes\thanksrefAachen H. Bagherpour\thanksrefChristchurch X. Bai\thanksrefSouthDakota A. Barbano\thanksrefGeneva J. P. Barron\thanksrefEdmonton S. W. Barwick\thanksrefIrvine V. Baum\thanksrefMainz R. Bay\thanksrefBerkeley J. J. Beatty\thanksrefOhio,OhioAstro J. Becker Tjus\thanksrefBochum K.-H. Becker\thanksrefWuppertal S. BenZvi\thanksrefRochester D. Berley\thanksrefMaryland E. Bernardini\thanksrefZeuthen D. Z. Besson\thanksrefKansas G. Binder\thanksrefLBNL,Berkeley D. Bindig\thanksrefWuppertal E. Blaufuss\thanksrefMaryland S. Blot\thanksrefZeuthen C. Bohm\thanksrefStockholmOKC M. Börner\thanksrefDortmund F. Bos\thanksrefBochum S. Böser\thanksrefMainz O. Botner\thanksrefUppsala E. Bourbeau\thanksrefCopenhagen J. Bourbeau\thanksrefMadisonPAC F. Bradascio\thanksrefZeuthen J. Braun\thanksrefMadisonPAC H.-P. Bretz\thanksrefZeuthen S. Bron\thanksrefGeneva J. Brostean-Kaiser\thanksrefZeuthen A. Burgman\thanksrefUppsala R. S. Busse\thanksrefMadisonPAC T. Carver\thanksrefGeneva C. Chen\thanksrefGeorgia E. Cheung\thanksrefMaryland D. Chirkin\thanksrefMadisonPAC K. Clark\thanksrefSNOLAB L. Classen\thanksrefMunster G. H. Collin\thanksrefMIT J. M. Conrad\thanksrefMIT P. Coppin\thanksrefBrusselsVrije P. Correa\thanksrefBrusselsVrije D. F. Cowen\thanksrefPennPhys,PennAstro R. Cross\thanksrefRochester P. Dave\thanksrefGeorgia M. Day\thanksrefMadisonPAC J. P. A. M. de André\thanksrefMichigan C. De Clercq\thanksrefBrusselsVrije J. J. DeLaunay\thanksrefPennPhys H. Dembinski\thanksrefBartol K. Deoskar\thanksrefStockholmOKC S. De Ridder\thanksrefGent P. Desiati\thanksrefMadisonPAC K. D. de Vries\thanksrefBrusselsVrije G. de Wasseige\thanksrefBrusselsVrije M. de With\thanksrefBerlin T. DeYoung\thanksrefMichigan J. C. Díaz-Vélez\thanksrefMadisonPAC H. Dujmovic\thanksrefSKKU M. Dunkman\thanksrefPennPhys E. Dvorak\thanksrefSouthDakota B. Eberhardt\thanksrefMainz T. Ehrhardt\thanksrefMainz B. Eichmann\thanksrefBochum P. Eller\thanksrefPennPhys P. A. Evenson\thanksrefBartol S. Fahey\thanksrefMadisonPAC A. R. Fazely\thanksrefSouthern J. Felde\thanksrefMaryland K. Filimonov\thanksrefBerkeley C. Finley\thanksrefStockholmOKC A. Franckowiak\thanksrefZeuthen E. Friedman\thanksrefMaryland A. Fritz\thanksrefMainz T. K. Gaisser\thanksrefBartol J. Gallagher\thanksrefMadisonAstro E. Ganster\thanksrefAachen S. Garrappa\thanksrefZeuthen L. Gerhardt\thanksrefLBNL K. Ghorbani\thanksrefMadisonPAC W. Giang\thanksrefEdmonton T. Glauch\thanksrefMunich T. Glüsenkamp\thanksrefErlangen A. Goldschmidt\thanksrefLBNL J. G. Gonzalez\thanksrefBartol D. Grant\thanksrefEdmonton Z. Griffith\thanksrefMadisonPAC C. Haack\thanksrefAachen A. Hallgren\thanksrefUppsala L. Halve\thanksrefAachen F. Halzen\thanksrefMadisonPAC K. Hanson\thanksrefMadisonPAC D. Hebecker\thanksrefBerlin D. Heereman\thanksrefBrusselsLibre K. Helbing\thanksrefWuppertal R. Hellauer\thanksrefMaryland S. Hickford\thanksrefWuppertal J. Hignight\thanksrefMichigan G. C. Hill\thanksrefAdelaide K. D. Hoffman\thanksrefMaryland R. Hoffmann\thanksrefWuppertal T. Hoinka\thanksrefDortmund B. Hokanson-Fasig\thanksrefMadisonPAC K. Hoshina\thanksrefMadisonPAC,a F. Huang\thanksrefPennPhys M. Huber\thanksrefMunich K. Hultqvist\thanksrefStockholmOKC M. Hünnefeld\thanksrefDortmund R. Hussain\thanksrefMadisonPAC S. In\thanksrefSKKU N. Iovine\thanksrefBrusselsLibre A. Ishihara\thanksrefChiba E. Jacobi\thanksrefZeuthen G. S. Japaridze\thanksrefAtlanta M. Jeong\thanksrefSKKU K. Jero\thanksrefMadisonPAC B. J. P. Jones\thanksrefArlington P. Kalaczynski\thanksrefAachen W. Kang\thanksrefSKKU A. Kappes\thanksrefMunster D. Kappesser\thanksrefMainz T. Karg\thanksrefZeuthen A. Karle\thanksrefMadisonPAC U. Katz\thanksrefErlangen M. Kauer\thanksrefMadisonPAC A. Keivani\thanksrefPennPhys J. L. Kelley\thanksrefMadisonPAC A. Kheirandish\thanksrefMadisonPAC J. Kim\thanksrefSKKU T. Kintscher\thanksrefZeuthen J. Kiryluk\thanksrefStonyBrook T. Kittler\thanksrefErlangen S. R. Klein\thanksrefLBNL,Berkeley R. Koirala\thanksrefBartol H. Kolanoski\thanksrefBerlin L. Köpke\thanksrefMainz C. Kopper\thanksrefEdmonton S. Kopper\thanksrefAlabama D. J. Koskinen\thanksrefCopenhagen M. Kowalski\thanksrefBerlin,Zeuthen K. Krings\thanksrefMunich M. Kroll\thanksrefBochum G. Krückl\thanksrefMainz S. Kunwar\thanksrefZeuthen N. Kurahashi\thanksrefDrexel A. Kyriacou\thanksrefAdelaide M. Labare\thanksrefGent J. L. Lanfranchi\thanksrefPennPhys M. J. Larson\thanksrefCopenhagen F. Lauber\thanksrefWuppertal K. Leonard\thanksrefMadisonPAC M. Leuermann\thanksrefAachen Q. R. Liu\thanksrefMadisonPAC E. Lohfink\thanksrefMainz C. J. Lozano Mariscal\thanksrefMunster L. Lu\thanksrefChiba J. Lünemann\thanksrefBrusselsVrije W. Luszczak\thanksrefMadisonPAC J. Madsen\thanksrefRiverFalls G. Maggi\thanksrefBrusselsVrije K. B. M. Mahn\thanksrefMichigan Y. Makino\thanksrefChiba S. Mancina\thanksrefMadisonPAC I. C. Mariş\thanksrefBrusselsLibre R. Maruyama\thanksrefYale K. Mase\thanksrefChiba R. Maunu\thanksrefMaryland K. Meagher\thanksrefBrusselsLibre M. Medici\thanksrefCopenhagen M. Meier\thanksrefDortmund T. Menne\thanksrefDortmund G. Merino\thanksrefMadisonPAC T. Meures\thanksrefBrusselsLibre S. Miarecki\thanksrefLBNL,Berkeley J. Micallef\thanksrefMichigan G. Momenté\thanksrefMainz T. Montaruli\thanksrefGeneva R. W. Moore\thanksrefEdmonton M. Moulai\thanksrefMIT R. Nagai\thanksrefChiba R. Nahnhauer\thanksrefZeuthen P. Nakarmi\thanksrefAlabama U. Naumann\thanksrefWuppertal G. Neer\thanksrefMichigan H. Niederhausen\thanksrefStonyBrook S. C. Nowicki\thanksrefEdmonton D. R. Nygren\thanksrefLBNL A. Obertacke Pollmann\thanksrefWuppertal A. Olivas\thanksrefMaryland A. O’Murchadha\thanksrefBrusselsLibre E. O’Sullivan\thanksrefStockholmOKC T. Palczewski\thanksrefLBNL,Berkeley H. Pandya\thanksrefBartol D. V. Pankova\thanksrefPennPhys P. Peiffer\thanksrefMainz C. Pérez de los Heros\thanksrefUppsala D. Pieloth\thanksrefDortmund E. Pinat\thanksrefBrusselsLibre A. Pizzuto\thanksrefMadisonPAC M. Plum\thanksrefMarquette P. B. Price\thanksrefBerkeley G. T. Przybylski\thanksrefLBNL C. Raab\thanksrefBrusselsLibre M. Rameez\thanksrefCopenhagen L. Rauch\thanksrefZeuthen K. Rawlins\thanksrefAnchorage I. C. Rea\thanksrefMunich R. Reimann\thanksrefAachen B. Relethford\thanksrefDrexel G. Renzi\thanksrefBrusselsLibre E. Resconi\thanksrefMunich W. Rhode\thanksrefDortmund M. Richman\thanksrefDrexel S. Robertson\thanksrefLBNL M. Rongen\thanksrefAachen C. Rott\thanksrefSKKU T. Ruhe\thanksrefDortmund D. Ryckbosch\thanksrefGent D. Rysewyk\thanksrefMichigan I. Safa\thanksrefMadisonPAC S. E. Sanchez Herrera\thanksrefEdmonton A. Sandrock\thanksrefDortmund J. Sandroos\thanksrefMainz M. Santander\thanksrefAlabama S. Sarkar\thanksrefCopenhagen,Oxford S. Sarkar\thanksrefEdmonton K. Satalecka\thanksrefZeuthen M. Schaufel\thanksrefAachen P. Schlunder\thanksrefDortmund T. Schmidt\thanksrefMaryland A. Schneider\thanksrefMadisonPAC J. Schneider\thanksrefErlangen S. Schöneberg\thanksrefBochum L. Schumacher\thanksrefAachen S. Sclafani\thanksrefDrexel D. Seckel\thanksrefBartol S. Seunarine\thanksrefRiverFalls J. Soedingrekso\thanksrefDortmund D. Soldin\thanksrefBartol M. Song\thanksrefMaryland G. M. Spiczak\thanksrefRiverFalls C. Spiering\thanksrefZeuthen J. Stachurska\thanksrefZeuthen M. Stamatikos\thanksrefOhio T. Stanev\thanksrefBartol A. Stasik\thanksrefZeuthen R. Stein\thanksrefZeuthen J. Stettner\thanksrefAachen A. Steuer\thanksrefMainz T. Stezelberger\thanksrefLBNL R. G. Stokstad\thanksrefLBNL A. Stößl\thanksrefChiba N. L. Strotjohann\thanksrefZeuthen T. Stuttard\thanksrefCopenhagen G. W. Sullivan\thanksrefMaryland M. Sutherland\thanksrefOhio I. Taboada\thanksrefGeorgia F. Tenholt\thanksrefBochum S. Ter-Antonyan\thanksrefSouthern A. Terliuk\thanksrefZeuthen S. Tilav\thanksrefBartol M. N. Tobin\thanksrefMadisonPAC C. Tönnis\thanksrefSKKU S. Toscano\thanksrefBrusselsVrije D. Tosi\thanksrefMadisonPAC M. Tselengidou\thanksrefErlangen C. F. Tung\thanksrefGeorgia A. Turcati\thanksrefMunich R. Turcotte\thanksrefAachen C. F. Turley\thanksrefPennPhys B. Ty\thanksrefMadisonPAC E. Unger\thanksrefUppsala M. A. Unland Elorrieta\thanksrefMunster M. Usner\thanksrefZeuthen J. Vandenbroucke\thanksrefMadisonPAC W. Van Driessche\thanksrefGent D. van Eijk\thanksrefMadisonPAC N. van Eijndhoven\thanksrefBrusselsVrije S. Vanheule\thanksrefGent J. van Santen\thanksrefZeuthen M. Vraeghe\thanksrefGent C. Walck\thanksrefStockholmOKC A. Wallace\thanksrefAdelaide M. Wallraff\thanksrefAachen F. D. Wandler\thanksrefEdmonton N. Wandkowsky\thanksrefMadisonPAC T. B. Watson\thanksrefArlington C. Weaver\thanksrefEdmonton M. J. Weiss\thanksrefPennPhys C. Wendt\thanksrefMadisonPAC J. Werthebach\thanksrefMadisonPAC S. Westerhoff\thanksrefMadisonPAC B. J. Whelan\thanksrefAdelaide N. Whitehorn\thanksrefUCLA K. Wiebe\thanksrefMainz C. H. Wiebusch\thanksrefAachen L. Wille\thanksrefMadisonPAC D. R. Williams\thanksrefAlabama L. Wills\thanksrefDrexel M. Wolf\thanksrefMunich J. Wood\thanksrefMadisonPAC T. R. Wood\thanksrefEdmonton E. Woolsey\thanksrefEdmonton K. Woschnagg\thanksrefBerkeley G. Wrede\thanksrefErlangen D. L. Xu\thanksrefMadisonPAC X. W. Xu\thanksrefSouthern Y. Xu\thanksrefStonyBrook J. P. Yanez\thanksrefEdmonton G. Yodh\thanksrefIrvine S. Yoshida\thanksrefChiba T. Yuan\thanksrefMadisonPAC

(Received: date / Accepted: date)

Abstract

The IceCube Collaboration has observed a high-energy astrophysical neutrino flux and recently found evidence for neutrino emission from the blazar TXS 0506+056. These results open a new window into the high-energy universe. However, the source or sources of most of the observed flux of astrophysical neutrinos remains uncertain.

Here, a search for steady point-like neutrino sources is performed using an unbinned likelihood analysis. The method searches for a spatial accumulation of muon-neutrino events using the very high-statistics sample of about $497\,000$ neutrinos recorded by IceCube between 2009 and 2017. The median angular resolution is $\sim 1^{\circ}$ at 1 TeV and improves to $\sim 0.3^{\circ}$ for neutrinos with an energy of 1 PeV. Compared to previous analyses, this search is optimized for point-like neutrino emission with the same flux-characteristics as the observed astrophysical muon-neutrino flux and introduces an improved event-reconstruction and parametrization of the background. The result is an improvement in sensitivity to the muon-neutrino flux compared to the previous analysis of $\sim 35\%$ assuming an $E^{-2}$ spectrum. The sensitivity on the muon-neutrino flux is at a level of $E^{2}\mathrm{d}N/\mathrm{d}E=3\cdot 10^{-13}\,\mathrm{TeV}\,\mathrm{cm}^{-2}\,\mathrm{s}^{-1}$ .

No new evidence for neutrino sources is found in a full sky scan and in an a priori candidate source list that is motivated by gamma-ray observations. Furthermore, no significant excesses above background are found from populations of sub-threshold sources. The implications of the non-observation for potential source classes are discussed.

Keywords:

Neutrino IceCube point source

^†^†journal: Eur. Phys. J. C

1 Introduction

Astrophysical neutrinos are thought to be produced by hadronic interactions of cosmic-rays with matter or radiation fields in the vicinity of their acceleration sites Gaisser et al. (1995). Unlike cosmic-rays, neutrinos are not charged and are not deflected by magnetic fields and thus point back to their origin. Moreover, since neutrinos have a relatively small interaction cross section, they can escape from the sources and do not suffer absorption on their way to Earth. Hadronic interactions of high-energy cosmic rays may also result in high-energy or very-high-energy gamma-rays. Since gamma-rays can also arise from the interaction of relativistic leptons with low-energy photons, only neutrinos are directly linked to hadronic interactions. The most commonly assumed neutrino-flavor flux ratios in the sources result in equal or nearly equal flavor flux ratios at Earth Athar et al. (2006). Thus about 1/3 of the astrophysical neutrinos are expected to be muon neutrinos and muon anti-neutrinos.

In 2013, the IceCube Collaboration reported the observation of an unresolved, astrophysical, high-energy, all-flavor neutrino flux, consistent with isotropy, using a sample of events which begin inside the detector (‘starting events’) Aartsen et al. (2013a); Kopper (2017). This observation was confirmed by the measurement of an astrophysical high-energy muon-neutrino flux using the complementary detection channel of through-going muons, produced in neutrino interactions in the vicinity of the detector Aartsen et al. (2015a, 2016a); Haack (2017). Track-like events from through-going muons are ideal to search for neutrino sources because of their relatively good angular resolution. However, to date, the sources of this flux have not been identified.

In 2018, first evidence of neutrino emission from an individual source was observed for the blazar TXS 0506+056 Aartsen et al. (2018a, b). Multi-messenger observations following up a high-energy muon neutrino event on September 22, 2017 resulted in the detection of this blazar being in flaring state. Furthermore, evidence was found for an earlier neutrino burst from the same direction between September 2014 and March 2015. However, the total neutrino flux from this source is less than 1% of the total observed astrophysical flux. Furthermore, the stacking of the directions of known blazars has revealed no significant excess of astrophysical neutrinos at the locations of known blazars. This indicates that blazars from the 2nd Fermi-LAT AGN catalogue contribute less than about 30% to the total observed neutrino flux assuming an unbroken power-law spectrum with spectral index of $-2.5$ Aartsen et al. (2017a). The constraint weakens to about 40%-80% of the total observed neutrino flux assuming a spectral index of $-2$ Aartsen et al. (2018a). Note that these results are model dependent and an extrapolation beyond the catalog is uncertain. No other previous searches have revealed a significant source or source class of astrophysical neutrinos Aartsen et al. (2015b, c); Adrián-Martínez et al. (2016); Aartsen et al. (2016b, c, 2017b, 2017c, 2017d); Albert et al. (2017); Aartsen et al. (2017e, f).

Here, a search for point-like sources is presented that takes advantage of the improved event selection and reconstruction of a muon-neutrino sample developed in Aartsen et al. (2016a) and the increased livetime of eight years Haack (2017) between 2009 and 2017. The best description of the sample includes a high-energy astrophysical neutrino flux given by a single power-law with a spectral index of $2.19\pm 0.10$ and a flux normalization, at $100\,\mathrm{TeV}$ , of $\Phi_{100\,\mathrm{TeV}}=1.01^{+0.26}_{-0.23}\times 10^{-18}\,\mathrm{GeV}^{-1}\,\mathrm{cm}^{-2}\,\mathrm{s}^{-1}\,\mathrm{sr}^{-1}$ , resulting in 190 to 2145 astrophysical neutrinos in the event sample. Compared to the previous time-integrated point source publication by IceCube Abbasi et al. (2011); Aartsen et al. (2013b, 2014a, 2016b, 2017b), this analysis is optimized for sources that show similar energy spectra as the measured astrophysical muon-neutrino spectrum. Furthermore, a high-statistics Monte Carlo parametrization of the measured data, consisting of astrophysical and atmospherical neutrinos and including systematic uncertainties, is used to model the background expectation and thus increases the sensitivity.

Within this paper, the following tests are discussed: 1. a full sky scan for the most significant source in the Northern hemisphere, 2. a test for a population of sub-threshold sources based on the result of the full sky scan, 3. a search based on an a priori defined catalog of candidate objects motivated by gamma-ray observations Aartsen et al. (2017b), 4. a test for a population of sub-threshold sources based on the result of the a priori defined catalog search, and 5. a test of the recently observed blazar TXS 0506+056. The tests are described in Section 3.4 and their results are given in Section 4.

2 Data sample

Season	Start Date	Livetime / days	Events	Declination Range	$\log_{10}(E_{\nu}^{\mathrm{astro}}/\mathrm{GeV})$ Range	$\log_{10}(E_{\nu}^{\mathrm{atmos}}/\mathrm{GeV})$ Range
IC59	2009/05/20	353.39	21411	$0^{\circ}$ – $+90^{\circ}$	3.02 – 5.73	2.37 – 4.06
IC79	2010/06/01	310.59	36880	$-5^{\circ}$ – $+90^{\circ}$	2.96 – 5.82	2.36 – 4.04
IC2011	2011/05/13	359.97	71191	$-5^{\circ}$ – $+90^{\circ}$	2.89 – 5.76	2.29 – 3.98
IC2012	2012/05/15	331.35
IC2013	2013/05/02	360.45
IC2014	2014/05/06	367.96	367590	$-5^{\circ}$ – $+90^{\circ}$	2.91 – 5.77	2.29 - 3.91
IC2015	2015/05/18	356.18
IC2016	2016/05/25	340.95

Table 1: Data samples used in this analysis and some characteristics of these samples. For each sample start date, livetime, number of observed events, and energy and declination range of the event selections are given. The energy range, calculated using a spectrum of atmospheric neutrinos and astrophysical neutrinos, spans the central 90% of the simulated events. Astrophysical neutrinos were generated using the best-fit values listed in Section 1. Note that livetime values slightly deviate from Ref. Aartsen et al. (2016a); Haack (2017) as the livetime calculation has been corrected.

IceCube is a cubic-kilometer neutrino detector with 5160 digital optical modules installed on 86 cable strings in the clear ice at the geographic South Pole between depths of 1450 m and 2450 m Achterberg et al. (2006); Aartsen et al. (2017g). The neutrino energy and directional reconstruction relies on the optical detection of Cherenkov radiation emitted by secondary particles produced in neutrino interactions in the surrounding ice or the nearby bedrock. The produced Cherenkov light is detected by digital optical modules (DOMs) each consisting of a 10 inch photomultiplier tube Abbasi et al. (2010), on-board read-out electronics Abbasi et al. (2009) and a high-voltage board, all contained in a spherical glass pressure vessel. Light propagation within the ice can be parametrized by the scattering and absorption behavior of the antarctic ice at the South Pole Aartsen et al. (2013c). The detector construction finished in 2010. During construction, data was taken in partial detector configurations with 59 strings (IC59) from May 2009 to May 2010 and with 79 strings (IC79) from May 2010 to May 2011 before IceCube became fully operational.

For events arriving from the Southern hemisphere, the trigger rate in IceCube is dominated by atmospheric muons produced in cosmic-ray air showers. The event selection is restricted to the Northern hemisphere where these muons are shielded by the Earth. Additionally, events are considered down to $-5^{\circ}$ declination, where the effective overburden of ice is sufficient to strongly attenuate the flux of atmospheric muons. Even after requiring reconstructed tracks from the Northern hemisphere, the event rate is dominated by mis-reconstructed atmospheric muons. However, these mis-reconstructed events can be reduced to less than 0.3% of the background using a careful event selection Aartsen et al. (2016a); Haack (2017). As the data were taken with different partial configurations of IceCube, the details of the event selections are different for each season. At final selection level, the sample is dominated by atmospheric muon neutrinos from cosmic-ray air showers Aartsen et al. (2016a). These atmospheric neutrinos form an irreducible background to astrophysical neutrino searches and can be separated from astrophysical neutrinos on a statistical basis only.

In total, data with a livetime of 2780.85 days are analyzed containing about $497\,000$ events at the final selection level. A summary of the different sub-samples is shown in Tab. 1.

The performance of the event selection can be characterized by the effective area of muon-neutrino and anti-neutrino detection, the point spread function and the central 90% energy range of the resulting event sample. The performance is evaluated with a full detector Monte Carlo simulation Aartsen et al. (2017g).

The effective area $A_{\mathrm{eff}}^{\nu+\bar{\nu}}$ quantifies the relation between neutrino and anti-neutrino fluxes $\phi_{\nu+\bar{\nu}}$ with respect to the observed rate of events $\frac{\mathrm{d}N_{\nu+\bar{\nu}}}{\mathrm{d}t}$ :

\frac{\mathrm{d}N_{\nu+\bar{\nu}}}{\mathrm{d}t}=\int\mathrm{d}\Omega\int_{0}^{\infty}\mathrm{d}E_{\nu}\,A_{\mathrm{eff}}^{\nu+\bar{\nu}}(E_{\nu},\theta,\phi)\times\phi_{\nu+\bar{\nu}}(E_{\nu},\theta,\phi)\,,

(1)

where $\Omega$ is the solid angle, $\theta,\phi$ are the detector zenith and azimuth angle and $E_{\nu}$ is the neutrino energy. The effective area for muon neutrinos and muon anti-neutrinos averaged over the Northern hemisphere down to -5 degree declination is shown in Fig.1 (top).

At high energies, the muon direction is well correlated with the muon-neutrino direction ( $<0.1^{\circ}$ deviation above 10 TeV) and the muon is reconstructed with a median angular uncertainty $\Delta\Psi_{\nu}$ of about $0.6^{\circ}$ at $10\,\mathrm{TeV}$ . All events have been reconstructed with an improved reconstruction based on the techniques described in Ahrens et al. (2004); Schatto (2014). The median angular resolution is shown in Fig.1 (middle). The median angular resolution at a neutrino energy of $1\,\mathrm{TeV}$ is about $1^{\circ}$ and decreases for higher energies to about $0.3^{\circ}$ at $1\,\mathrm{PeV}$ .

The central 90% energy range is shown in Fig. 1 (bottom) as a function of $\sin\delta$ , with declination $\delta$ . Energy ranges are calculated using the precise best-fit parametrization of the experimental sample. The energy range stays mostly constant as function of declination but shifts to slightly higher energies near the horizon. The central 90% energy range of atmospheric neutrinos is about $200\,\mathrm{GeV}$ – $10\,\mathrm{TeV}$ .

Figure 1: Top: Muon neutrino and anti-neutrino effective area averaged over the Northern hemisphere as function of

\log_{10}

of neutrino energy. Middle: Median neutrino angular resolution as function of

\log_{10}

of neutrino energy. Bottom: Central 90% neutrino energy range for atmospheric (astrophysical) neutrinos as solid (dashed) line for each declination. Lines show the livetime weighted averaged of all sub-samples. Plots for individual seasons can be found in the supplemental material.

In Fig. 2, the ratio of effective area (top) and median angular resolution (bottom) of the sub-sample IC86 2012-2016 and the sample labeled 2012-2015 from previous time-integrated point source publication by IceCube is shown Aartsen et al. (2017b). The differences in effective area are declination dependent. When averaged over the full Northern hemisphere, the effective area produced by this event selection is smaller than that in Aartsen et al. (2017b) at low neutrino energies but is larger above $100\,\mathrm{TeV}$ . The median neutrino angular resolution $\Delta\Psi_{\nu}$ improves at $10\,\mathrm{TeV}$ by about 10% compared to the reconstruction used in Aartsen et al. (2017b) and improves up to 20% at higher energies. The event sample for the season from May 2011 to May 2012 has an overlap of about 80% with the selection presented in Ref. Aartsen et al. (2017b) using the same time range.

Figure 2: Ratio of effective area (top) and median angular resolution (bottom) of the sub-sample IC86 2012-2016 and the sample labeled 2012-2015 from previous publication of time-integrated point source searches by IceCube is shown Aartsen et al. (2017b).

3 Unbinned likelihood method

3.1 Likelihood & test statistics

The data sample is tested for a spatial clustering of events with an unbinned likelihood method described in Braun et al. (2008) and used in the previous time-integrated point source publications by IceCube Abbasi et al. (2011); Aartsen et al. (2013b, 2014a, 2016b, 2017b). At a given position $\vec{x}_{s}$ in the sky, the likelihood function for a source at this position, assuming a power law energy spectrum with spectral index $\gamma$ , is given by

\mathcal{L}=\prod_{i}^{\mathrm{events}}\left[\frac{n_{s}}{N}S_{i}(\vec{x}_{s},\gamma)+\left(1-\frac{n_{s}}{N}\right)B_{i}\right]\cdot P(\gamma)\,,

(2)

where $i$ is an index of the observed neutrino events, $N$ is the total number of events, $n_{s}$ is the number of signal events and $P(\gamma)$ is a prior term. $S_{i}$ and $B_{i}$ are the signal and background probability densities evaluated for event $i$ . The likelihood is maximized with respect to the source parameters $n_{s}\geq 0$ and $1\leq\gamma\leq 4$ at each tested source position in the sky given by its right ascension and declination $\vec{x}_{s}=(\alpha_{s},\delta_{s})$ .

The signal and background probability density functions (PDF) $S_{i}$ and $B_{i}$ factorize into a spatial and an energy factor

	$\displaystyle S_{i}(\vec{x}_{s},\gamma)$	$\displaystyle=$	$\displaystyle S^{\mathrm{spat}}(\vec{x}_{i},\sigma_{i}\|\vec{x}_{s})\cdot S^{\mathrm{ener}}(E_{i}\|\gamma)$		(3)
	$\displaystyle B_{i}$	$\displaystyle=$	$\displaystyle B^{\mathrm{spat}}(\vec{x}_{i})B^{\mathrm{ener}}(E_{i})\,,$		(4)

where $\vec{x}_{i}=(\alpha_{i},\delta_{i})$ is the reconstructed right ascension $\alpha_{i}$ and declination $\delta_{i}$ , $E_{i}$ is the reconstructed energy Abbasi et al. (2013) and $\sigma_{i}$ is the event-by-event based estimated angular uncertainty of the reconstruction of event $i$ Neunhöffer (2006); Abbasi et al. (2011).

A likelihood ratio test is performed to compare the best-fit likelihood to the null hypothesis of no significant clustering $\mathcal{L}_{0}=\prod_{i}B_{i}$ . The likelihood ratio is given by

TS=2\cdot\log\left[\frac{\mathcal{L}(\vec{x}_{s},\hat{n}_{s},\hat{\gamma})}{\mathcal{L}_{0}}\right]\,,

(5)

with best-fit values $\hat{n}_{s}$ and $\hat{\gamma}$ , which is used as a test statistic.

The sensitivity of the analysis is defined as the median expected 90% CL upper limit on the flux normalization in case of pure background. In addition, the discovery potential is defined as the signal strength that leads to a $5\,\sigma$ deviation from background in 50% of all cases.

In previous point source publications by IceCube Abbasi et al. (2011); Aartsen et al. (2013b, 2014a, 2016b, 2017b), the spatial background PDF $B^{\mathrm{spat}}$ and the energy background PDF $B^{\mathrm{ener}}$ were estimated from the data. Given the best-fit parameters obtained from Aartsen et al. (2016a) and good data / Monte Carlo agreement, it is, however, possible to get a precise parametrization of the atmospheric and diffuse astrophysical components, including systematic uncertainties. By doing this, it is possible to take advantage of the high statistics of the full detector simulation data sets which can be used to generate smooth PDFs optimized for the sample used in this work. Thus this parametrization of the experimental data allows us to obtain a better extrapolation to sparsely populated regions in the energy-declination plane than by using only the statistically limited experimental data. This comes with the drawback that the analysis can only be applied to the Northern hemisphere since no precise parametrization is available for the Southern hemisphere. Generating PDFs from full detector simulations has already been done in previous publications for the energy signal PDF $S^{\mathrm{ener}}$ , as it is not possible to estimate this PDF from data itself. The spatial signal PDF $S^{\mathrm{spat}}$ is still assumed to be Gaussian with an event individual uncertainty of $\sigma_{i}$ .

It is known from the best-fit parametrization of the sample that the data contain astrophysical events. The astrophysical component has been parametrized by an unbroken power-law with best-fit spectral index of $2.19\pm 0.10$ Haack (2017). In contrast to the previous publication of time-integrated point source searches by IceCube Aartsen et al. (2017b), which uses a flat prior on the spectral index in the range $1\leq\gamma\leq 4$ , this analysis focuses on those sources that produce the observed spectrum of astrophysical events by adding a Gaussian prior $P(\gamma)$ on the spectral index in Eq. 2 with mean 2.19 and width 0.10. As the individual source spectra are not strongly constrained by the few events that contribute to a source, the prior dominates the fit of $\gamma$ and thus the spectral index is effectively fixed allowing only for small variations. Due to the prior, the likelihood has reduced effective degrees of freedom to model fluctuations. As a result, the distribution of the test statistic in the case of only background becomes steeper which results in an improvement of the discovery potential assuming an $E^{-2}$ source spectrum.

However, due to the reduced freedom of the likelihood by the prior on the spectral index about 80% of background trials yield $\hat{n}_{s}=0$ and thus $TS=0$ . This pile-up leads to an over-estimation of the median 90% upper limit as the median is degenerate and the flux sensitivity is artificially over-estimated. Thus a different definition for the $TS$ is introduced for $\hat{n}_{s}=0$ . Allowing for negative $\hat{n}_{s}$ can lead to convergence problems due to the second free parameter of $\gamma$ . Assuming $\hat{n}_{s}=0$ is already close to the minimum of $\log\mathcal{L}$ , $\log\mathcal{L}$ can be approximated as a parabola. The likelihood is extended in a Taylor series up to second order around $n_{s}=0$ . The Taylor series gives a parabola for which the value of the extremum can be calculated from the first and second order derivative of the likelihood at $n_{s}=0$ . This value is used as test statistic

TS=-2\cdot\frac{\left(\left.\log\mathcal{L}^{\prime}\right|_{0}\right)^{2}}{2\left.\log\mathcal{L}^{\prime\prime}\right|_{0}}\,,\qquad\hat{n}_{s}=0\,,

(6)

for likelihood fits that yield $\hat{n}_{s}=0$ . With this definition, the pile-up of $\hat{n}_{s}$ is spread towards negative values of $TS$ and the median of the test statistic is no longer degenerate. Using this method, the sensitivity which had been overestimated due to the pile-up at $n_{s}=0$ can be recovered.

3.2 Pseudo-experiments

To calculate the performance of the analysis, pseudo-experiments containing only background and pseudo-experiments with injected signal have been generated.

In this search for astrophysical point sources, atmospheric neutrinos and astrophysical neutrinos from unresolved sources make up the background. Using the precise parametrization of the reconstructed declination and energy distribution¹¹1In Ref. Haack (2017), the reconstructed zenith-energy distribution has been parametrized, although, due to IceCube’s unique position at the geographic South Pole the zenith can be directly converted to declination. from Ref. Haack (2017), pseudo-experiments are generated using full detector simulation events. Due to IceCube’s position at the South Pole and the high duty cycle of >99% Aartsen et al. (2017g), the background PDF is uniform in right ascension.

As a cross check, background samples are generated by scrambling experimental data uniformly in right ascension. The declination and energy of the events are kept fixed. This results in a smaller sampled range of event energy and declination compared to the Monte Carlo-based pseudo-experiments. In the Monte Carlo-based pseudo-experiments, events are sampled from the simulated background distributions, and thus are not limited to the values of energy and declination present in the data when scrambling. P-values for tests presented in Section 4 are calculated using the Monte Carlo method and are compared to the data scrambling method for verification (values in brackets).

Signal is injected according to a full simulation of the detector. Events are generated at a simulated source position assuming a power law energy distribution. The number of injected signal events is calculated from the assumed flux and the effective area for a small declination band around the source position. In this analysis, the declination band was reduced compared to previous publication of time-integrated point source searches by IceCube Aartsen et al. (2017b), resulting in a more accurate modelling of the effective area. This change in signal modeling has a visible effect on the sensitivity and discovery potential, especially at the horizon and at the celestial pole. The effect can be seen in Fig. 3 by comparing the solid (small bandwidth) and dotted (large bandwidth) lines. The bandwidth is optimized by taking into account the effect of averaging over small declination bands and limited simulation statistics to calculate the effective area. As the bandwidth cannot be made too narrow, an uncertainty of about 8% on the flux limit calculation arises and is included in the systematic error.

3.3 Sensitivity & discovery potential

Figure 3: Sensitivity (dashed) and

5\sigma

discovery potential (solid) of the flux normalization for an

E^{-2}

source spectrum as function of the

\sin\delta

. For comparison, the lines from Aartsen et al. (2017b) are shown as well. The dotted line indicates the bandwidth effect discussed in Section 3.2.

The sensitivity and discovery potential for a single point source is calculated for an unbroken power law flux according to

\frac{\mathrm{d}N_{\nu_{\mu}+\bar{\nu}_{\mu}}}{\mathrm{d}E_{\nu}}=\phi_{100\,\mathrm{TeV}}^{\nu_{\mu}+\bar{\nu}_{\mu}}\left(\frac{E_{\nu}}{100\,\mathrm{TeV}}\right)^{-\gamma}\,.

(7)

In Fig. 3, the sensitivity and discovery potential as function of $\sin\delta$ is shown. Note that Fig. 3 shows $E_{\nu}^{2}\frac{\mathrm{d}N_{\nu_{\mu}+\bar{\nu}_{\mu}}}{\mathrm{d}E_{\nu}}=\phi_{0}E_{0}^{2}$ which is constant in neutrino energy for an $E^{-2}$ flux. The sensitivity corresponds to a 90% CL averaged upper limit and the discovery potential gives the median source flux for which a $5\sigma$ discovery would be expected. The flux is given as a muon neutrino plus muon anti-neutrino flux. For comparison, the sensitivity and discovery potential from the previous publication of time-integrated point source searches by IceCube Aartsen et al. (2017b) are shown. Despite only a moderate increase of livetime, this analysis outperforms the analysis in Aartsen et al. (2017b) by about 35% for multiple reasons: 1. the use of an improved angular reconstruction, 2. a slightly better optimized event selection near the horizon, 3. the use of background PDFs in the likelihood that are optimized on the parametrization from Haack (2017); Aartsen et al. (2016a) which improves sensitivity especially for higher energies, 4. the fact that due to the prior on the spectral index the number of source hypotheses is reduced which results in a steeper falling background $TS$ distribution, and 5. the use of negative $TS$ values which avoids overestimating the sensitivity, especially in the celestial pole region ( $\sin\delta\sim 1$ ), where the background changes rapidly in $\sin\delta$ . In Fig. 4, the differential discovery potentials for three different declination bands are shown.

Figure 4: Differential sensitivity (dashed) and

5\sigma

discovery potential (solid) flux for three different declinations. For high declinations and high energies, the effect of neutrino absorption within the Earth becomes visible. The flux is given as the sum of the muon neutrino and anti-neutrino flux.

3.4 Tested hypothesis

3.4.1 Full sky scan

A scan of the full Northern hemisphere from $90^{\circ}$ down to $-3^{\circ}$ declination has been performed. The edge at $-3^{\circ}$ has been chosen to avoid computational problems due to fast changing PDFs at the boundary of the sample at $-5^{\circ}$ . The scan is performed on a grid with a resolution of about $0.1^{\circ}$ . The grid was generated using the HEALPix pixelization scheme²²2Hierarchical Equal Area isoLatitude Pixelation of a sphere (HEALPix), http://healpix.sourceforge.net/ Górski et al. (2005). For each grid point, the pre-trial p-value is calculated. As the test statistic shows a slight declination dependence, the declination dependent $TS$ is used to calculate local p-values. $TS$ distributions have been generated for 100 declinations equally distributed in $\sin\delta$ . $10^{6}$ trials have been generated for each declination. Below a $TS$ value of 5, the p-value is determined directly from trials. Above $TS=5$ , an exponential function is fitted to the tail of the distribution which is used to calculate p-values above $TS=5$ . A Kolmogorov-Smirnov test Kolmogorov (1933); Smirnov (1939) and a $\chi^{2}$ test are used to verify the agreement of the fitted function and the distribution.

The most significant point on the sky produced by the scan is selected using the pre-trial p-value. Since many points are tested in this scan, a trial correction has to be applied. Therefore, the procedure is repeated with background pseudo-experiments as described in Section 3.2. By comparing the local p-values from the most significant points in the background sample to the experimental pre-trial p-value, the post-trial p-value is calculated. The final p-value is calculated directly from $\sim 3500$ trials³³3The background distribution of the local p-value $p_{\mathrm{local}}$ for the most significant point is described by $d\mathcal{P}=N(1-p_{\mathrm{local}})^{N-1}dp_{\mathrm{local}}$ , with an effective number of trials $N$ that is fitted to $241\,000\pm 9\,000$ . A rough approximation of this trial factor can be calculated by dividing the solid angle of the Northern hemisphere $\sim 2\pi$ by the squared median angular resolution. Considering that highest energy events dominate the sensitivity, we use $0.3^{\circ}$ for the median angular resolution. Thus we get $2\pi/(0.3^{\circ})^{2}\approx 229000$ effective trials, which is in the same order of magnitude as the determined value..

3.4.2 Population test in the full sky scan

Figure 5: Upper Panel: Number of local warm spots with p-values smaller that

p_{\mathrm{thres}}

as function of

p_{\mathrm{thres}}

. The observed number of local spots are shown as solid black line. The background expectation is shown as dashed line with

1\sigma

2\sigma

and

3\sigma

intervals corresponding to Poisson statistics. Lower Panel: Local Poisson p-value for given

p_{\mathrm{thres}}

. The most significant point is indicated by a dotted vertical line.

Due to the large number of trials, only very strong sources would be identified in a full sky scan, which attempts to quantify only the most significant source. However, the obtained $TS$ values can be tested also for a significant excess of events from multiple weaker sources without any bias towards source positions. This is done by counting p-values of local warm spots where the p-values are smaller than a preset threshold. An excess of counts with respect to the expectation from pure background sky maps can indicate the presence of multiple weak sources.

From the full sky scan, local spots with $p_{\mathrm{local}}<10^{-2}$ and a minimal separation of $1^{\circ}$ are selected. The number of expected local spots $\lambda$ with a p-value smaller than $p_{\mathrm{thres}}$ is estimated from background pseudo-experiments and shown in Fig. 5 as dashed line. The background expectation was found to be Poisson distributed. The threshold value is optimized to give the most significant excess above background expectation using the Poisson probability

p_{\mathrm{poisson}}=\exp(-\lambda)\sum_{m=n}^{\infty}\frac{\lambda^{m}}{m!}\,,

(8)

to find an excess of at least $n$ spots. Due to the optimization of the threshold in the range on $2<-\log_{10}p_{\mathrm{thres}}<\infty$ , the result has to be corrected for trials as well. To include this correction, the full sky scan population test is performed on background pseudo-experiments to calculate the post-trial p-value.

3.4.3 A priori source list

The detectability of sources suffers from the large number of trials within the full sky scan and thus individual significant source directions may become insignificant after the trial correction. However, gamma-ray data can help to preselect interesting neutrino source candidates. A standard IceCube and ANTARES a priori source list, containing 34 prominent candidate sources for high-energy neutrino emission on the Northern hemisphere has been tested Aartsen et al. (2017b), reducing the trial factor to about the number of sources in the catalog. The source catalog is summarized in Tab. 2. The sources were selected mainly based on observations in gamma rays and belong to various object classes. The sources from this list are tested individually with the unbinned likelihood from Eq. 2. For this test, p-values are calculated from $10^{6}$ background trials without using any extrapolation. Then the most significant source is selected and a trial-correction, derived from background pseudo-experiments, is applied. Note that some sources such as MGRO J1908+06, SS 433, and Geminga are spatially extended with an apparent angular size of up to several degrees, which is larger than IceCube’s point spread function. In such cases, the sensitivity of the analysis presented in this paper is reduced. E.g., for an extension of 1 degree, the sensitivity on the neutrino flux decreases by $\sim 20\%$ Aartsen et al. (2014a).

3.4.4 Population test in the a priori source list

Similar to the population test in the full sky scan, an excess of several sources with small but not significant p-values in the a priori source list can indicate a population of weak sources. Therefore, the $k$ most significant p-values of the source list are combined using a binomial distribution

P_{\mathrm{binom}}(k|p_{k},N)={{N}\choose{k}}p_{k}^{k}(1-p_{k})^{N-k}\,,

(9)

of p-values that are larger than a threshold $p_{k}$ . Here, $N=34$ is the total number of sources in the source list. The most significant combination is used as a test statistic and assessed against background using pseudo-experiments.

3.4.5 Monitored source list

IceCube and ANTARES have tested the a priori source list for several years with increasingly sensitive analyses Abbasi et al. (2011); Aartsen et al. (2013b, 2014a, 2017b). Changing the source list posterior may lead to a bias on the result. However, not reporting on recently seen, interesting sources would also ignore progress in the field. A definition of an unbiased p-value is not possible as these were added later. Therefore, a second list with sources is tested to report on an updated source catalog. In this work, this second catalog so far comprises only TXS 0506+056, for which evidence for neutrino emission has been observed.

3.5 Systematic uncertainties

The p-values for the tested hypotheses are determined with simulated pseudo-experiments assuming only background (see also Section 3.2). These experiments are generated using the full detector Monte Carlo simulation, weighted to the best-fit parametrization from Ref. Haack (2017). This parametrization includes the optimization of nuisance parameters accounting for systematic uncertainties resulting in very good agreement between experimental data and Monte Carlo. Because of this procedure, the p-values are less affected by statistical fluctuations that would occur when estimating p-values from scrambled experimental data as well as the effect of fixed event energies during scrambling. However, a good agreement of the parametrization with experimental data is a prerequisite of this method. As a cross check, p-values are also calculated using scrambled experimental data. These p-values are given for comparison in brackets in Section 4. We find that the two methods show very similar results confirming the absence of systematic biases.

The calculation of the absolute neutrino flux normalization based on Monte Carlo simulations is affected by systematic uncertainties. These uncertainties influence the reconstruction performance and the determination of the effective area. Here, the dominant uncertainties are found to be the absolute optical efficiency of the Cherenkov light production and detection in the DOMs Abbasi et al. (2010), the optical properties (absorption, scattering) of the South Pole ice Aartsen et al. (2014b), and the photo-nuclear interaction cross sections of high energy muons Bezrukov and Bugaev (1981); Bugaev and Shlepin (2003a, b); Bugaev et al. (2004); Abramowicz et al. (1991); Abramowicz and Levy ; Koehne et al. (2013).

The systematic uncertainties on the sensitivity flux normalization is evaluated by propagating changed input values on the optical efficiency, ice properties and cross section values through the entire likelihood analysis for a signal energy spectrum of $\mathrm{d}N/\mathrm{d}E_{\nu}\propto E_{\nu}^{-2}$ . Changing the optical efficiencies by $\pm 10\%$ results in a change of the flux normalization by $\pm 7.5\%$ . The ice properties have been varied by (+10%, 0%), (0%, +10%) and (-7.1%, -7.1%) in the values of absorption and scattering length. The resulting uncertainty of the flux normalization is $\pm 5.3\%$ . To study the effect of the photo-nuclear interactions of high energy muons, the models in Ref. Bezrukov and Bugaev (1981); Bugaev and Shlepin (2003a, b); Bugaev et al. (2004); Abramowicz et al. (1991); Abramowicz and Levy ; Koehne et al. (2013) have been used, which give a flux normalization variation of $\pm 5.1\%$ . Note, that these models are outdated and represent the extreme cases from common literature. Thus, the systematic uncertainty is estimated conservatively. The systematic uncertainties are assumed to be independent and are added in quadrature, yielding a total systematic uncertainty of $\pm 10.5\%$ for the $\nu_{\mu}+\bar{\nu}_{\mu}$ flux normalization. One should note that additionally, the modeling of point-like sources yields an uncertainty of about $\pm 8\%$ as discussed in Section 3.2.

Since the sample is assumed to be purely muon neutrino and muon anti-neutrino events, only $\nu_{\mu}+\bar{\nu}_{\mu}$ fluxes are considered. However, $\nu_{\tau}$ and $\bar{\nu}_{\tau}$ may also contribute to the observed astrophysical neutrinos in the data sample. Taking $\nu_{\tau}$ and $\bar{\nu}_{\tau}$ fluxes into account and assuming an equal flavor ratio at Earth, the sensitivity of the per-flavor flux normalization improves, depending on the declination, by 2.6% – 4.3%. The expected contamination from $\nu_{e}$ and $\bar{\nu}_{e}$ is negligible.

The relative systematic uncertainty is comparable with the systematic uncertainties quoted in previous publications of time integrated point source searches by IceCube Aartsen et al. (2017b). In addition, the systematic effect due to the chosen finite bandwidth is included in this analysis.

4 Results

No significant clustering was found in any of the hypotheses tests beyond the expectation from background. Both the full-sky scan of the Northern hemisphere and the p-values from the source list are compatible with pure background. The p-values given in this section are calculated by pseudo-experiments based on Monte Carlo simulation weighted to the best-fit parametrization of the sample (see Section 3.2). For verification, p-values calculated by pseudo-experiments from scrambled experimental data are given in brackets.

4.1 Sky scan

Refer to caption — Figure 6: Sky map of the local p-values from the sky scan in equatorial coordinates down to $-3^{\circ}$ declination. The local p-value is given as $-\log_{10}(p_{\mathrm{local}})$ . The position of the most significant spot is indicated by a black circle.

The pre-trial p-value map of the Northern hemisphere scan is shown in Fig. 6. The hottest spot in the scan is indicated by a black circle and is located at $\alpha=177.89^{\circ}$ and $\delta=23.23^{\circ}$ (J2000) with the Galactic coordinates $b_{\mathrm{gal}}=75.92^{\circ}$ , $l_{\mathrm{gal}}=-134.33^{\circ}$ . The best-fit signal strength is $\hat{n}_{s}=21.32$ ( $\Phi^{\nu_{\mu}+\bar{\nu}_{\mu}}_{100\,\mathrm{TeV}}=1.4\cdot 10^{-19}\,\mathrm{GeV}^{-1}\mathrm{cm}^{-2}\mathrm{s}^{-1}$ assuming $\hat{\gamma}=2.20$ ) with a fitted spectral index of $\hat{\gamma}=2.20$ close to the prior of $2.19$ . The $TS$ -value is $21.63$ which corresponds to $p_{\mathrm{local}}=10^{-5.97}$ . The post-trial corrected p-value is 26.5% (29.9%) and is thus compatible with background. A zoom into the local p-value landscape around the hottest spot position and the observed events is shown in Fig. 7. Events are shown as small circles where the area of the circle is proportional to the median $\log_{10}$ of neutrino energy assuming the diffuse best-fit spectrum. The closest gamma-ray source from the Fermi 3FGL and Fermi 3FHL catalogs Acero et al. (2015); Ajello et al. (2017) is 3FHL J1150.3+2418 which is about 1.1 degree away from the hottest spot. The chance probability to find a 3FGL or 3FHL source within 1.1 degree is 25%, which is estimated from all-sky pseudo-experiments. At the source location of 3FHL J1150.3+2418, the $TS$ value is 8.02 which is inconsistent with the best-fit point at the $3.6\,\sigma$ level, if assuming Wilks theorem with one degree of freedom Wilks (1938).

4.2 Population test in the sky scan

In Fig. 5, the number of spots with p-values below $p_{\mathrm{thres}}$ are shown together with the expectation from background. The most significant deviation was found for $p_{\mathrm{thres}}=0.5\%$ where 454.3 spots were expected and 492 were observed with a p-value of $p_{\mathrm{poisson}}=4.17\%$ . Correcting the result for trials gives a p-value of 42.0% (54.3%) and thus the result is compatible with background.

Figure 8: Single-flavor neutrino and anti-neutrino flux per source vs number of sources. An unbroken

E^{-2}

power law and equal fluxes of the sources at Earth are assumed. Solid lines show 90% CL upper limits and dashed lines indicate the sensitivity. Upper limits and sensitivity are calculated assuming that background consists of atmospheric neutrinos only and exclude an astrophysical component. Thus the limits are conservative, especially for small number of sources. For comparison, the results from Aartsen et al. (2017b); Glauch and Turcati (2017) are given. The dotted line gives the flux per source that saturates the diffuse flux from Ref. Haack (2017).

As no significant deviation from the background hypothesis has been observed, exclusion limits are calculated as 90% CL upper limits with Neyman’s method Neyman (1937) for the benchmark scenario of a fixed number of sources $N_{\mathrm{sources}}$ , all producing the same flux at Earth. Upper limits are calculated assuming that background consists of atmospheric neutrinos only, excluding an astrophysical component from background pseudo-experiment generation. Excluding the astrophysical component from background is necessary as the summed injected flux makes up a substantial part of the astrophysical flux in case of large $N_{\mathrm{sources}}$ . However, this will over-estimate the flux sensitivity for small $N_{\mathrm{sources}}$ . More realistic source scenarios are discussed in Section 5. This rather unrealistic scenario does not depend on astrophysical and cosmological assumptions about source populations and allows for a comparison between the analysis power of different analyses directly. The sensitivity and upper limits for $N_{\mathrm{source}}$ sources is shown in Fig. 8 together with the analyses from Aartsen et al. (2017b); Glauch and Turcati (2017)⁴⁴4The 90% CL upper limit from Ref. Aartsen et al. (2017b) has been recalculated to account for an incorrect treatment of signal acceptance in the original publication.. This analysis finds the most stringent exclusion limits for small number of sources to date. The gain in sensitivity compared to Ref. Aartsen et al. (2017b) is consistent with the gain in the sensitivity to a single point source.

4.3 A priori source list

Table 2: Results of the a priori defined source list search. Coordinates are given in equatorial coordinates (J2000). The fitted spectral index

\hat{\gamma}

is not given as it is effectively fixed by the introduced prior. As discussed in the text, negative

TS

values are assigned to sources with best-fit

\hat{n}_{s}=0

. Source types abbreviation: BL Lacertae object (BL Lac), Flat Spectrum Radio Quasar (FSRQ), Not Identified (NI), Pulsar Wind Nebula (PWN), Star Formation Region (SFR), Supernova Remnant (SNR), Starburst / Radio Galaxy (SRG), X-ray Binary and Micro-Quasar (XB/mqso).

Source	Type	$\alpha\,[\mathrm{deg}]$	$\delta\,[\mathrm{deg}]$	p-Value	$TS$	$\hat{n}_{s}$	$E^{2}\mathrm{d}N_{\nu_{\mu}+\bar{\nu}_{\mu}}/\mathrm{d}E\,[\mathrm{TeV}\,\mathrm{cm}^{-2}\,\mathrm{s}^{-1}]$
4C 38.41	FSRQ	248.81	38.13	0.0080	5.0893	7.69	1.27 $\cdot 10^{-12}$
MGRO J1908+06	NI	286.99	6.27	0.0088	4.7933	2.82	7.62 $\cdot 10^{-13}$
Cyg A	SRG	299.87	40.73	0.0101	4.7199	3.80	1.28 $\cdot 10^{-12}$
3C454.3	FSRQ	343.50	16.15	0.0258	2.9675	5.03	8.08 $\cdot 10^{-13}$
Cyg X-3	XB/mqso	308.11	40.96	0.1263	0.5695	4.33	8.20 $\cdot 10^{-13}$
Cyg OB2	SFR	308.09	41.23	0.1706	0.2554	2.82	7.64 $\cdot 10^{-13}$
LSI 303	XB/mqso	40.13	61.23	0.2056	0.1747	2.37	9.93 $\cdot 10^{-13}$
NGC 1275	SRG	49.95	41.51	0.2447	0.0230	0.50	6.96 $\cdot 10^{-13}$
1ES 1959+650	BL Lac	300.00	65.15	0.2573	0.0717	1.70	9.86 $\cdot 10^{-13}$
Crab Nebula	PWN	83.63	22.01	0.3213	-0.0197	0.00	4.74 $\cdot 10^{-13}$
Mrk 421	BL Lac	166.11	38.21	0.3460	-0.0205	0.00	5.79 $\cdot 10^{-13}$
Cas A	SNR	350.85	58.81	0.3808	-0.0169	0.00	7.01 $\cdot 10^{-13}$
TYCHO	SNR	6.36	64.18	0.3893	-0.0219	0.00	7.98 $\cdot 10^{-13}$
PKS 1502+106	FSRQ	226.10	10.52	0.3931	-0.1770	0.00	3.57 $\cdot 10^{-13}$
3C66A	BL Lac	35.67	43.04	0.4265	-0.1089	0.00	5.44 $\cdot 10^{-13}$
3C 273	FSRQ	187.28	2.05	0.4285	-0.3705	0.00	2.72 $\cdot 10^{-13}$
HESS J0632+057	XB/mqso	98.24	5.81	0.5017	-0.7603	0.00	2.82 $\cdot 10^{-13}$
BL Lac	BL Lac	330.68	42.28	0.5378	-0.4766	0.00	4.78 $\cdot 10^{-13}$
W Comae	BL Lac	185.38	28.23	0.5961	-1.0769	0.00	3.88 $\cdot 10^{-13}$
Cyg X-1	XB/mqso	299.59	35.20	0.6170	-1.0639	0.00	4.31 $\cdot 10^{-13}$
1ES 0229+200	BL Lac	38.20	20.29	0.6257	-1.6867	0.00	3.41 $\cdot 10^{-13}$
M87	SRG	187.71	12.39	0.7054	-2.9682	0.00	3.26 $\cdot 10^{-13}$
Mrk 501	BL Lac	253.47	39.76	0.7214	-1.9858	0.00	4.58 $\cdot 10^{-13}$
PKS 0235+164	BL Lac	39.66	16.62	0.7494	-3.5951	0.00	3.33 $\cdot 10^{-13}$
H 1426+428	BL Lac	217.14	42.67	0.7587	-2.5100	0.00	4.86 $\cdot 10^{-13}$
PKS 0528+134	FSRQ	82.73	13.53	0.7788	-4.4554	0.00	3.18 $\cdot 10^{-13}$
S5 0716+71	BL Lac	110.47	71.34	0.7802	-2.0711	0.00	8.02 $\cdot 10^{-13}$
Geminga	PWN	98.48	17.77	0.7950	-4.7785	0.00	3.41 $\cdot 10^{-13}$
SS433	XB/mqso	287.96	4.98	0.8455	-8.0055	0.00	2.71 $\cdot 10^{-13}$
M82	SRG	148.97	69.68	0.8456	-3.5574	0.00	8.04 $\cdot 10^{-13}$
3C 123.0	SRG	69.27	29.67	0.9056	-8.2916	0.00	4.11 $\cdot 10^{-13}$
1ES 2344+514	BL Lac	356.77	51.70	0.9518	-10.1395	0.00	5.28 $\cdot 10^{-13}$
IC443	SNR	94.18	22.53	0.9620	-16.4154	0.00	3.63 $\cdot 10^{-13}$
MGRO J2019+37	PWN	305.22	36.83	0.9784	-17.6070	0.00	4.54 $\cdot 10^{-13}$

Figure 9: Sensitivity (dashed) and

5\sigma

discovery potential (solid) of the flux normalization for an

E^{-2}

source spectrum as function of the

\sin\delta

. For comparison, the lines from Aartsen et al. (2017b) are shown as well. 90% CL Neyman upper limits on the flux normalization for sources in the a priori and monitored source list are shown as circles and squares, respectively.

The fit results of sources in the a priori source list are given in Tab. 2. The most significant source with a local p-value of 0.8% is 4C 38.41, which is a flat spectrum radio quasar (FSRQ) at a redshift of $z=1.8$ . Taking into account that 34 sources have been tested, a post-trial p-value of 23.7% (20.3%) is calculated from background pseudo-experiments which is compatible with background.

As no significant source has been found, 90% CL upper limits are calculated assuming an unbroken power law with spectral index of -2 using Neyman’s method Neyman (1937). The 90% CL upper limit flux is summarized in Tab. 2 and shown in Fig. 9. In case of under-fluctuations, the limit was set to the sensitivity level of the analysis. Note that 90% upper limits can exceed the discovery potential as long as the best-fit flux is below the discovery potential.

Interestingly, a total of three sources, 4C 38.41, MGRO J1908+06 and Cyg A, have a local p-value below or close to 1%. The p-value landscapes and observed events around these three sources are shown in Fig. 10 and Fig. 11.

4.4 Population test in the a piori source list

Figure 12: Local significance in Gaussian

\sigma

for binomial combinations of the

k

most significant sources in the a priori source list. Sources with

\hat{n}_{s}>0

and

\hat{n}_{s}=0

can be separated by the dashed vertical line.

The most significant combination of p-values from the a priori source list is given when combining the three most significant p-values, i.e. $k=3$ , with $2.59\sigma$ as shown in Fig. 12. The comparison with background pseudo-experiments yields a trial-corrected p-value of 6.6% (4.1%) which is not significant.

4.5 Monitored source list

Table 3: Results of the monitored source list search. The fitted spectral index

\hat{\gamma}

is not given as it is effectively fixed by the introduced prior. We use the abbreviation BL Lac for BL Lacertae objects.

Source	Type	$\alpha\,[\mathrm{deg}]$	$\delta\,[\mathrm{deg}]$	p-Value	$TS$	$\hat{n}_{s}$	$E^{2}\mathrm{d}N_{\nu_{\mu}+\bar{\nu}_{\mu}}/\mathrm{d}E\,[\mathrm{TeV}\,\mathrm{cm}^{-2}\,\mathrm{s}^{-1}]$
TXS 0506+056	BL Lac	77.38	5.69	0.0293	2.6475	7.87	6.19 $\cdot 10^{-13}$

The best-fit results for TXS 0506+056 in the monitored source list are given in Tab. 3. Note that the event selection ends in May 2017 and thus does not include the time of the alert ICECUBE-170922A IceCube Collaboration (2017) that led to follow-up observations and the discovery of $\gamma$ -ray emission from that blazar up to 400 GeV. The data, however, include the earlier time-period of the observed neutrino flare. The local p-value here is found to be 2.93%. This is less significant than the reported significance of the time-dependent flare in Aartsen et al. (2018a) but is consistent with the reported time-integrated significances in Aartsen et al. (2018a), when taking into account that this analysis has a prior on the spectral index of the source flux and does not cover the same time-range as in Aartsen et al. (2018a).

The local p-value landscape around TXS 0506+056 is shown in Fig. 11 together with the observed event directions of this sample.

5 Implications on source populations

Figure 13: 90% CL upper limits on the effective muon-neutrino luminosity within the energy range

10^{4}\,\mathrm{GeV}

–

10^{7}\,\mathrm{GeV}

at Earth and effective source density, derived from the hotspot population analysis and the sky scan.

The non-detection of a significant point-like source and the non-detection of a population of sources within the sky scan is used to put constrains on realistic source populations. In the following calculation, source populations are characterized by their effective $\nu_{\mu}+\bar{\nu}_{\mu}$ single-source luminosity $L_{\nu_{\mu}+\bar{\nu}_{\mu}}^{\mathrm{eff}}$ and their local source density $\rho_{0}^{\mathrm{eff}}$ . Using the software tool FIRESONG⁵⁵5FIRst Extragalactic Simulation Of Neutrinos and Gamma-rays (FIRESONG), https://github.com/ChrisCFTung/FIRESONG Taboada et al. (2017), the resulting source count distribution $\frac{\mathrm{d}N}{\mathrm{d}\Phi}$ as a function of the flux $\Phi$ for source populations are calculated for sources within $z<10$ and representations of this population are simulated. To calculate the source count distribution, FIRESONG takes the source density $\rho$ , luminosity distribution, source evolution, cosmological parameters, the energy range of the flux and the spectral index into account. Following Ref. Murase and Waxman (2016), sources are simulated with a log-normal distribution with median $L_{\nu_{\mu}+\bar{\nu}_{\mu}}^{\mathrm{eff}}$ and a width of 0.01 in $\log_{10}(L_{\nu_{\mu}+\bar{\nu}_{\mu}}^{\mathrm{eff}})$ which corresponds to a standard candle luminosity. The evolution of the sources was chosen to follow the parametrization of star formation rate from Hopkins and Beacom Hopkins and Beacom (2006) assuming a flat universe with $\Omega_{M,0}=0.308$ , $\Omega_{\lambda,0}=0.692$ and $h=0.678$ Ade et al. (2016). The energy range of the flux at Earth was chosen as $10^{4}\,\mathrm{GeV}$ – $10^{7}\,\mathrm{GeV}$ to calculate the effective muon neutrino luminosities of sources.

Generating pseudo-experiments with signal components corresponding to the flux distribution obtained from FIRESONG, 90% CL upper limits are calculated in the $\rho_{0}^{\mathrm{eff}}$ - $L_{\nu_{\mu}+\bar{\nu}_{\mu}}^{\mathrm{eff}}$ plane for various spectral indices assuming that background consists of atmospheric neutrinos only, as described in Section 8. The 90% CL upper limit is calculated based on the fact that the strongest source of a population does not give a p-value in the sky scan that is larger than the observed one. The 90% upper limits are shown as dashed lines in Fig. 13. In addition, 90% CL upper limits are calculated by comparing the largest excess measured with the population test in the sky scan. These 90% upper limits are shown as solid lines in Fig. 13. Populations that are compatible at the $1\sigma$ and $3\sigma$ level with the diffuse flux measured in Haack (2017) are shown as blue shaded band. 90% CL upper limits have been calculated assuming an $E^{-2}$ power-law flux. The same has been performed for an $E^{-2.19}$ power-law flux, which is the diffuse best-fit for this sample (this result can be found in the supplementary material). The computation of upper limits becomes very computing-intensive for large source densities. Therefore, the computation of the upper limits, resulting from the sky scan, are extrapolated to larger source densities (indicated by dotted line in Fig. 13). It can be seen that for large effective source densities and small effective luminosities, the limit resulting from the population analysis goes $\propto 1/L_{\nu_{\mu}+\bar{\nu}_{\mu}}^{\mathrm{eff}}$ which is the same scaling as one would expect from a diffuse flux. Indeed it is found that an excess of diffuse high-energy events, i.e. sources from which only one neutrino are detected, leads to a p-value excess in the population analysis. This is a result of taking the energy of the event into account in the likelihood. Limits from the hottest spot in the sky scan are a bit stronger for large effective luminosities while upper limits from the population test become stronger at about $L_{\nu_{\mu}+\bar{\nu}_{\mu}}^{\mathrm{eff}}\sim 10^{52}\frac{\mathrm{erg}}{\mathrm{yr}}$ .

6 Implications for individual source models

Table 4: Model rejection factors for source models in the source catalog. Given are source type, model reference, central energy range that contributes 90% to sensitivity, MRF sensitivity and MRF at 90% CL.

Type	Source Model	$\log_{10}(E/\mathrm{GeV})$	sensitivity	90% UL
Crab	Amato et. al Amato et al. (2003) $\Gamma=10^{4}$	1.5 - 9.0	23.38	31.47
	Amato et. al Amato et al. (2003) $\Gamma=10^{5}$	3.0 - 4.5	0.79	1.14
	Amato et. al Amato et al. (2003) $\Gamma=10^{6}$	4.0 - 5.5	0.16	0.21
	Amato et. al Amato et al. (2003) $\Gamma=10^{7}$	4.5 - 6.0	0.32	0.40
	Kappes et. al Kappes et al. (2007)	2.5 - 4.5	1.06	1.47
Blazar	3C273, Reimer Reimer (2015)	6.0 - 8.5	0.39	0.42
	3C454.3, Reimer Reimer (2015)	6.0 - 8.0	2.80	5.42
	Mrk421, Petropoulou et. al Petropoulou et al. (2015)	5.5 - 7.0	0.36	0.43
SNR	G40.5-0.5, Mandelartz et. al Mandelartz and Becker Tjus (2015)	3.5 - 5.5	1.45	4.57

In Section 4.3, constraints on source fluxes assuming $\mathrm{d}N/\mathrm{d}E_{\nu}\propto E_{\nu}^{-2}$ have been calculated. However, more specific neutrino flux models can be obtained using $\gamma$ -ray data. In pion decays, both neutrinos and $\gamma$ -rays are produced. Thus $\gamma$ -ray data can be used to construct models for neutrino emission under certain assumptions. Here, models for sources of the a priori source list are tested. For each model, the Model Rejection Factor (MRF) is calculated which is the ratio between the predicted flux and the 90% CL upper limit. In addition, the expected experimental result in the case of pure background is also calculated giving the MRF sensitivity. The energy range that contributes 90% to the sensitivity has been calculated by folding the differential discovery potential at the source position (similar to Fig. 4) with the flux prediction. Models for which the MRF sensitivity is larger than 10 are not discussed here.

The first source tested is the Crab Nebula, which is a Pulsar Wind Nebula (PWN) and the brightest source in TeV $\gamma$ -rays. Despite the common understanding that the emission from PWNe is of leptonic nature, see e.g. Aleksić et al. (2015a), neutrinos can be produced by subdominant hadronic emission. Predictions for neutrino fluxes from the Crab Nebula are proposed, e.g. by Amato et. al Amato et al. (2003) and Kappes et. al Kappes et al. (2007). The prediction by Amato et. al assumes pion production is dominated by p-p interactions and the target density is given by $n_{t}=10\mu M_{N_{\odot}}R_{\mathrm{pc}}^{-3}\,\mathrm{cm}^{-3}$ with $M_{N_{\odot}}$ the mass of the supernova ejecta in units of solar masses. Moreover, $R_{\mathrm{pc}}$ is the radius of the supernova in units of $\mathrm{pc}$ and $\mu$ is an unknown factor of the order of $1\leq\mu\leq 20$ that takes into account e.g. the intensity and structures of magnetic fields within the PWN. Here $\mu=20$ and a proton luminosity of 60% of the total PWN luminosity for Lorentz factors of $\Gamma=10^{4},\,10^{5},\,10^{6},\,10^{7}$ are used to provide a result that is model-independent and complementary to Amato et al. (2003). The model prediction by Kappes et al., assumes a dominant production of $\gamma$ -rays of the HESS $\gamma$ -ray spectrum Aharonian et al. (2006) by p-p interactions.

The model predictions, sensitivity and 90% CL upper limit are shown in Fig. 14 and are listed in Tab. 4. Sensitivity and upper limits are shown for the central energy range that contributes 90% to the sensitivity.

Figure 14: Differential source flux for the Crab Nebula. Solid lines show the model prediction, thick lines give the 90% CL upper limit and the dashed lines indicate the sensitivity flux. 90% CL upper limit and sensitivity are shown in the energy range that contributes 90% to the sensitivity.

For the model of Kappes et al., the sensitivity is very close to the model prediction while for Amato et al. with $\Gamma=10^{7}$ , the sensitivity is a factor of three lower than the prediction. The 90% CL upper limits are listed in Tab. 4. They are slightly higher but still constrain the models by Amato et al.

Another very interesting class of sources are active galactic nuclei (AGN). Here, the models being tested come from Ref. Petropoulou et al. (2015) for Mrk 421, a BL Lacertae object (BL Lac) that was found in spatial and energetic agreement with a high-energy starting event and from Ref. Reimer (2015) for 3C273 and 3C454.3 which are flat spectrum radio quasars (FSRQ). The models, sensitivities and 90% CL upper limits are shown in Fig. 15 and the MRF are listed in Tab. 4.

Figure 15: Differential source flux for 3C273, 3C454.3 and Mrk 421. Solid lines show the model prediction, thick lines give the 90% CL upper limit and dashed lines indicate the sensitivity flux. 90% CL upper limit and sensitivity are shown in the energy range that contributes 90% to the sensitivity.

The sensitivities for 3C273 and Mrk 421 are well below the model prediction and the 90% CL upper limits are at about 40% of the model flux. For 3C454.3, the sensitivity is a factor 2.8 above the model prediction. Since 3C454.3 is one of the few sources with a local p-value below $\sim 2.5\%$ , the 90% CL upper limit is much larger.

Another tested model was derived for the source G40.5-0.5 which is a galactic supernova remnant Mandelartz and Becker Tjus (2015). This supernova remnant can be associated with the TeV source MGRO J1908+06 which is the second most significant source in the a priori source catalog, although the association of G40.5-0.5 with MGRO J1908+06 is not distinct Aliu et al. (2006). In addition, the pulsar wind nebula powered by PSR J1907+0602 may contribute to the TeV emission of the MGRO J1908+06 region. However, here the tested model for the SNR G40.5-0.5 is adapted from Ref. Mandelartz and Becker Tjus (2015). The model, sensitivity and 90% CL upper limit are shown in Fig. 16 and are listed in Tab. 4.

Figure 16: Differential source flux for SNR G40.5-0.5. The solid line gives the model prediction, the thick line gives the 90% CL upper limit and the dashed line indicates the sensitivity flux. The 90% CL upper limit and sensitivity are shown in the energy range that contributes 90% to the sensitivity. G40.5-0.5 is associated with MGRO J1908+06.

The sensitivity of this analysis is a factor 1.4 above the model prediction and not yet sensitive to this model. As MGRO J1908+06 is the second most significant source in the catalog, with a local p-value of $<1\%$ , the upper limit lies nearly a factor of five above the model prediction.

7 Conclusions

Eight years of IceCube data have been analyzed for a time-independent clustering of through-going muon neutrinos using an unbinned likelihood method. The analysis includes a full sky search of the Northern hemisphere down to a declination of $-3^{\circ}$ for a significant hot spot as well as an analysis of a possible cumulative excess of a population of weak sources. Furthermore, source-candidates from an a priori catalog and a catalog of monitored sources are tested individually and again for a cumulative excess.

The analysis method has been optimized for the observed energy spectrum of high-energy astro-physical muon neutrinos Aartsen et al. (2016a) and a number of improvements with respect to the previously published search Aartsen et al. (2017b) have been incorporated. By implementing these improvements, a sensitivity increase of about 35% has been achieved.

No significant source was found in the full-sky scan of the Northern hemisphere and the search for significant neutrino emission from objects on a a priori source list results in a post-trial p-value of 23.7% (20.3%), compatible with background. Also the tests for populations of sub-threshold sources revealed no significant excess.

Three sources on the a priori source-list, 4C 38.41, MGRO J1908+06 and Cyg A, have pre-trial p-values of only about 1%. However, these excesses are not significant. The source TXS 0506+056 in the catalog of monitored sources has a p-value of 2.9 %. This is consistent with the time-integrated p-value in Aartsen et al. (2018a) for the assumed prior on the spectral index.

Based on these results, the most stringent limits on high-energy neutrino emission from point-like sources are obtained. In addition, models for neutrino emission from specific sources are tested. The model Amato et al. (2003) for the Crab Nebula is excluded for $\Gamma\geq 10^{6}$ as well as the predictions for 3C273 Reimer (2015) and Mrk 421 Petropoulou et al. (2015). In addition to these specific models, an exclusion of source populations as a function of local source density and single-source luminosity are derived by calculating the source count distribution for a realistic cosmological evolution model.

Acknowledgements.

The IceCube collaboration acknowledges the significant contributions to this manuscript from René Reimann. 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).

References

Gaisser et al. (1995) T. K. Gaisser, F. Halzen, and T. Stanev, Phys. Rept. 258, 173 (1995), [Erratum: Phys. Rept.271,355(1996)], arXiv:hep-ph/9410384 .
Athar et al. (2006) H. Athar, C. S. Kim, and J. Lee, Mod. Phys. Lett. A21, 1049 (2006), arXiv:hep-ph/0505017 .
Aartsen et al. (2013a) M. G. Aartsen et al. (IceCube Collaboration), Science 342, 1242856 (2013a), arXiv:1311.5238 .
Kopper (2017) C. Kopper (IceCube Collaboration), PoS ICRC2017, 981 (2017), arXiv:1710.01191 .
Aartsen et al. (2015a) M. G. Aartsen et al. (IceCube Collaboration), Phys. Rev. Lett. 115, 081102 (2015a), arXiv:1507.04005 .
Aartsen et al. (2016a) M. G. Aartsen et al. (IceCube Collaboration), Astrophys. J. 833, 3 (2016a), arXiv:1607.08006 .
Haack (2017) C. Haack (IceCube Collaboration), PoS ICRC2017, 1005 (2017), arXiv:1710.01191 .
Aartsen et al. (2018a) M. G. Aartsen et al. (IceCube Collaboration), Science 361, 147 (2018a), arXiv:1807.08794 .
Aartsen et al. (2018b) M. G. Aartsen et al. (Liverpool Telescope, MAGIC, H.E.S.S., AGILE, Kiso, VLA/17B-403, INTEGRAL, Kapteyn, Subaru, HAWC, Fermi-LAT, ASAS-SN, VERITAS, Kanata, IceCube, Swift NuSTAR), Science 361, eaat1378 (2018b), arXiv:1807.08816 .
Aartsen et al. (2017a) M. G. Aartsen et al. (IceCube Collaboration), Astrophys. J. 835, 45 (2017a), arXiv:1611.03874 .
Aartsen et al. (2015b) M. G. Aartsen et al. (IceCube Collaboration), Astrophys. J. 805, L5 (2015b), arXiv:1412.6510 .
Aartsen et al. (2015c) M. G. Aartsen et al. (IceCube Collaboration), Astrophys. J. 807, 46 (2015c), arXiv:1503.00598 .
Adrián-Martínez et al. (2016) S. Adrián-Martínez et al. (ANTARES Collaboration & IceCube Collaboration), Astrophys. J. 823, 65 (2016), arXiv:1511.02149 .
Aartsen et al. (2016b) M. G. Aartsen et al. (IceCube Collaboration), Astrophys. J. 824, L28 (2016b), arXiv:1605.00163 .
Aartsen et al. (2016c) M. G. Aartsen et al. (IceCube Collaboration), J. Inst. 11, P11009 (2016c), arXiv:1610.01814 .
Aartsen et al. (2017b) M. G. Aartsen et al. (IceCube Collaboration), Astrophys. J. 835, 151 (2017b), arXiv:1609.04981 .
Aartsen et al. (2017c) M. G. Aartsen et al. (IceCube Collaboration), Astrophys. J. 843, 112 (2017c), arXiv:1702.06868 .
Aartsen et al. (2017d) M. G. Aartsen et al. (IceCube Collaboration), A&A 607, A115 (2017d), arXiv:1702.06131 .
Albert et al. (2017) A. Albert et al. (ANTARES Collaboration & IceCube Collaboration & LIGO Scientific Collaboration & Virgo Collaboration), Phys. Rev. D 96, 022005 (2017), arXiv:1703.06298 .
Aartsen et al. (2017e) M. G. Aartsen et al. (IceCube Collaboration), Astrophys. J. 846, 136 (2017e), arXiv:1705.02383 .
Aartsen et al. (2017f) M. G. Aartsen et al. (IceCube Collaboration), Astrophys. J. 849, 67 (2017f), arXiv:1707.03416 .
Abbasi et al. (2011) R. Abbasi et al. (IceCube Collaboration), Astrophys. J. 732, 18 (2011), arXiv:1012.2137 .
Aartsen et al. (2013b) M. G. Aartsen et al. (IceCube Collaboration), Astrophys. J. 779, 132 (2013b), arXiv:1307.6669 .
Aartsen et al. (2014a) M. G. Aartsen et al. (IceCube Collaboration), Astrophys. J. 796, 109 (2014a), arXiv:1406.6757 .
Achterberg et al. (2006) A. Achterberg et al. (IceCube Collaboration), Astropart. Phys. 26, 155 (2006), arXiv:astro-ph/0604450 .
Aartsen et al. (2017g) M. G. Aartsen et al. (IceCube Collaboration), J. Inst. 12, P03012 (2017g), arXiv:1612.05093 .
Abbasi et al. (2010) R. Abbasi et al. (IceCube Collaboration), Nucl. Instrum. Meth. A 618, 139 (2010), arXiv:1002.2442 .
Abbasi et al. (2009) R. Abbasi et al. (IceCube Collaboration), Nucl. Instrum. Meth. A 601, 294 (2009), arXiv:0810.4930 .
Aartsen et al. (2013c) M. G. Aartsen et al. (IceCube Collaboration), Nucl. Instrum. Meth. A 711, 73 (2013c), arXiv:1301.5361 .
Ahrens et al. (2004) J. Ahrens et al. (AMANDA Collaboration), Nucl. Instrum. Meth. A 524, 169 (2004), arXiv:astro-ph/0407044 .
Schatto (2014) K. Schatto, Stacked searches for high-energy neutrinos from blazars with IceCube, Ph.D. thesis, Johannes Gutenberg-Universität (2014).
Braun et al. (2008) J. Braun et al., Astropart. Phys. 29, 299 (2008), arXiv:0801.1604 .
Abbasi et al. (2013) R. Abbasi et al. (IceCube Collaboration), Nucl. Instrum. Meth. A 703, 190 (2013), arXiv:1208.3430 .
Neunhöffer (2006) T. Neunhöffer, Astropart. Phys. 25, 220 (2006), arXiv:astro-ph/0403367 .
Górski et al. (2005) K. M. Górski et al., Astrophys. J. 622, 759 (2005), arXiv:astro-ph/0409513 .
Kolmogorov (1933) A. Kolmogorov, Inst. Ital. Attuari, Giorn. 4, 83 (1933).
Smirnov (1939) N. V. Smirnov, Bull. Math. Univ. Moscou 2, 3 (1939).
Aartsen et al. (2014b) M. G. Aartsen et al. (IceCube Collaboration), J. Inst. 9, P03009 (2014b), arXiv:1311.4767 .
Bezrukov and Bugaev (1981) L. B. Bezrukov and E. V. Bugaev, Sov. J. Nucl. Phys. 33, 1195 (1981).
Bugaev and Shlepin (2003a) E. V. Bugaev and Y. V. Shlepin, Nucl. Phys. B 122, 341 (2003a).
Bugaev and Shlepin (2003b) E. V. Bugaev and Y. V. Shlepin, Phys. Rev. D 67, 034027 (2003b), arXiv:hep-ph/0203096 .
Bugaev et al. (2004) E. Bugaev, T. Montaruli, Y. Shlepin, and I. Sokalski, Astropart. Phys. 21, 491 (2004), arXiv:hep-ph/0312295 .
Abramowicz et al. (1991) H. Abramowicz, E. M. Levin, A. Levy, and U. Maor, Phys. Lett. B 269, 465 (1991).
(44) H. Abramowicz and A. Levy, DESY 97-251, arXiv:hep-ph/9712415 .
Koehne et al. (2013) J.-H. Koehne et al., Comput. Phys. Commun. 184, 2070 (2013).
Acero et al. (2015) F. Acero et al. (Fermi-LAT Collaboration), Astrophys. J. Suppl. 218, 23 (2015), arXiv:1501.02003 .
Ajello et al. (2017) M. Ajello et al. (Fermi-LAT Collaboration), Astrophys. J. Suppl. 232, 18 (2017), arXiv:1702.00664 .
Wilks (1938) S. S. Wilks, Ann. Math. Statist. 9, 60 (1938).
Glauch and Turcati (2017) T. Glauch and A. Turcati (IceCube Collaboration), PoS ICRC2017, 1014 (2017), arXiv:1710.01179 .
Neyman (1937) J. Neyman, Phil. Trans. R. Soc. A 236, 333 (1937).
IceCube Collaboration (2017) IceCube Collaboration, GRB Coordinates Network, Circular Service 21916 (2017).
Taboada et al. (2017) I. Taboada, C. F. Tung, and J. Wood (HAWC Collaboration), PoS ICRC2017, 663 (2017), arXiv:1801.09545 .
Murase and Waxman (2016) K. Murase and E. Waxman, Phys. Rev. D 94, 103006 (2016), arXiv:1607.01601 .
Hopkins and Beacom (2006) A. M. Hopkins and J. F. Beacom, Astrophys. J. 651, 142 (2006), arXiv:astro-ph/0601463 .
Ade et al. (2016) P. A. R. Ade et al. (Planck Collaboration), A&A 594, A13 (2016), arXiv:1502.01589 .
Amato et al. (2003) E. Amato, D. Guetta, and P. Blasi, A&A 402, 827 (2003).
Kappes et al. (2007) A. Kappes, J. Hinton, C. Stegmann, and F. A. Aharonian, Astrophys. J. 661, 1348 (2007), arXiv:astro-ph/0607286 .
Reimer (2015) A. Reimer, PoS ICRC2015, 1123 (2015).
Petropoulou et al. (2015) M. Petropoulou, S. Dimitrakoudis, P. Padovani, A. Mastichiadis, and E. Resconi, Mon. Not. R. Astron. Soc. 448, 2412 (2015), arXiv:1501.07115 .
Mandelartz and Becker Tjus (2015) M. Mandelartz and J. Becker Tjus, Astropart. Phys. 65, 80 (2015), arXiv:1301.2437 .
Aleksić et al. (2015a) J. Aleksić et al. (MAGIC Collaboration), JHEAp 5-6, 30 (2015a), arXiv:1406.6892 .
Aharonian et al. (2006) F. Aharonian et al. (HESS Collaboration), A&A 457, 899 (2006), arXiv:astro-ph/0607333 .
Aliu et al. (2006) E. Aliu et al. (VERITAS Collaboration), Astrophys. J. 787, 166 (2014), arXiv:1404.7185 .

Appendix A Performance of individual sub-samples

The quality and statistical power of a sample, w.r.t. a search for point-like sources, can be characterized by the effective area of muon-neutrino and anti-neutrino detection, the point spread function and the central 90% energy range (see Section 2). As the data were taken with different partial configurations of IceCube, the details of the event selections are different for each season. In Fig. 1 the livetime average of all sub-samples is shown. In Fig. 17 the effective area, point spread function and central 90% energy range are shown for each sub-sample individually. The plot shows that - despite of different detector configurations and event selections - the characteristics of the event samples are similar.

Figure 17: Top: Muon neutrino and anti-neutrino effective area averaged over the Northern hemisphere as function of

\log_{10}

of neutrino energy. Middle: Median neutrino angular resolution as function of

\log_{10}

of neutrino energy. Bottom: Central 90% neutrino energy range for atmospheric (astrophysical) neutrinos as solid (dashed) line for each declination. Lines are labeled by there sub-season.

Appendix B Results for diffuse best-fit spectral index

An $E^{-2}$ power-law is often used as benchmark model and for a comparison between publications. However, the diffuse best-fit spectral index is $\gamma=2.19$ , which is, given the uncertainties is not consistent with $\gamma=2$ . Therefore, the sensitivity and discovery potential for single point sources are recalculated using this spectral index. In Fig. 18, the sensitivity and discovery potential for an $E^{-\gamma}$ spectum are shown with $\gamma=2.0$ and $\gamma=2.19$ . The flux normalization is evaluated at a pivot energy of $100\,\mathrm{TeV}$ . The sensitivity and discovery potential for the assumed spectral indices turn out to be very similar.

Figure 18: Sensitivity and discovery potential on the flux normalization at

100\,\mathrm{TeV}

for an

E^{-\gamma}

power-law spectrum. Lines are given for

\gamma=2.0

as in Fig. 3 and the diffuse best-fit spectral index of

\gamma=2.19

In addition also the 90% CL upper limit on source populations as described in Section 5 are recalculated for a spectral index of $\gamma=2.19$ . The upper limit are shown in Fig. 19. Comparing with Fig. 13, there is no strong indication of a dependence on the spectral index.

Figure 19: Same as Fig. 13 but for

\gamma=2.19