Extending the search for muon neutrinos coincident with gamma-ray bursts in IceCube data

TL;DR

This study conducts an all-sky search for muon neutrinos coincident with gamma-ray bursts using IceCube data, aiming to detect evidence of ultra-high energy cosmic ray production, but finds no significant correlation.

Contribution

It extends previous IceCube searches by including additional years of data and a new Southern Hemisphere analysis, improving constraints on GRB neutrino models.

Findings

No significant neutrino-GRB correlation detected

New constraints placed on GRB neutrino production models

Enhanced sensitivity to high-energy neutrinos achieved

Abstract

We present an all-sky search for muon neutrinos produced during the prompt $\gamma$-ray emission of 1172 gamma-ray bursts (GRBs) with the IceCube Neutrino Observatory. The detection of these neutrinos would constitute evidence for ultra-high energy cosmic ray (UHECR) production in GRBs, as interactions between accelerated protons and the prompt $\gamma$-ray field would yield charged pions, which decay to neutrinos. A previously reported search for muon neutrino tracks from Northern Hemisphere GRBs has been extended to include three additional years of IceCube data. A search for such tracks from Southern Hemisphere GRBs in five years of IceCube data has been introduced to enhance our sensitivity to the highest energy neutrinos. No significant correlation between neutrino events and observed GRBs is seen in the new data. Combining this result with previous muon neutrino track searches and a search for cascade signature events from all neutrino flavors, we obtain new constraints for single-zone fireball models of GRB neutrino and UHECR production.

Figures

Tables

Table 1. Table 1: Three-year Northern Hemisphere track analysis coincident events. The duration T 100 subscript 𝑇 100 T_{100} , angular uncertainty σ 𝜎 \sigma (Fermi-GBM statistical-only uncertainties indicated by ∗ * ), and total γ 𝛾 \gamma -ray fluence ( erg ​ cm − 2 erg superscript cm 2 \mathrm{erg\,cm^{-2}} ) of GRBs with coincident events are described. Coincident events are summarized in terms of their time relative to the GRB start time T 1 subscript 𝑇 1 T_{1} , their angular uncertainty σ 𝜎 \sigma , angular displacement from the GRB location Δ ​ Ψ Δ Ψ \Delta\Psi , and reconstructed muon energy proxy (TeV). Event significance is estimated by their signal-to-background PDF ratio value 𝒮 / ℬ 𝒮 ℬ \mathcal{S}/\mathcal{B} (only events with 𝒮 / ℬ > 10 − 3 𝒮 ℬ superscript 10 3 \mathcal{S}/\mathcal{B}>10^{-3} are listed), while final GRB coincidence significance is given as pre-trials (post-trials) p-values relative to background-only test statistic distributions.

	Time	$\sigma$	$\Delta \Psi$	Fluence/Energy	Significance
GRB120612B	$T_{100}=63.24s$	∗${7.1}^{\circ }$		$2.06\times {10}^{-6}$	$p=0.049$ $\left(1\right)$
Event 1	$T_1+47.71s$	${5.3}^{\circ }$	${29.0}^{\circ }$	$\gtrsim 0.54$	$𝒮/ℬ=1.4$
GRB120911A	$T_{100}=28.58s$	${0.0003}^{\circ }$		$2.34\times {10}^{-6}$	$p=0.0044$ $\left(1\right)$
Event 1	$T_1+120.94s$	${4.6}^{\circ }$	${2.9}^{\circ }$	$\gtrsim 0.98$	$𝒮/ℬ=3.1$
GRB130116A	$T_{100}=66.82s$	∗${29.9}^{\circ }$		$9.27\times {10}^{-7}$	$p=0.076$ $\left(1\right)$
Event 1	$T_1+69.25s$	${0.5}^{\circ }$	${67.7}^{\circ }$	$\gtrsim 2.1$	$𝒮/ℬ=1.5$
GRB130318A	$T_{100}=137.99s$	∗${9.9}^{\circ }$		$3.41\times {10}^{-6}$	$p=0.021$ $\left(1\right)$
Event 1	$T_1+29.83s$	${0.6}^{\circ }$	${18.5}^{\circ }$	$\gtrsim 0.46$	$𝒮/ℬ=6.4$
Event 2	$T_1+44.58s$	${2.5}^{\circ }$	${48.2}^{\circ }$	$\gtrsim 0.32$	$𝒮/ℬ=0.024$
GRB130925B	$T_{100}=265.47s$	∗${4.1}^{\circ }$		$1.49\times {10}^{-5}$	$p=0.032$ $\left(1\right)$
Event 1	$T_1+108.8s$	${3.4}^{\circ }$	${12.6}^{\circ }$	$\gtrsim 0.70$	$𝒮/ℬ=16.3$
GRB131029B	$T_{100}=50.95s$	∗${5.8}^{\circ }$		$4.49\times {10}^{-6}$	$p=0.053$ $\left(1\right)$
Event 1	$T_1+50.49s$	${2.4}^{\circ }$	${18.2}^{\circ }$	$\gtrsim 0.68$	$𝒮/ℬ=4.9$
GRB131202B	$T_{100}=86.02s$	∗${2.2}^{\circ }$		$1.24\times {10}^{-5}$	$p=0.0069$ $\left(0.988\right)$
Event 1	$T_1+85.18s$	${2.1}^{\circ }$	${7.5}^{\circ }$	$\gtrsim 1.7$	$𝒮/ℬ=122.1$
GRB140404B	$T_{100}=26.63s$	∗${2.2}^{\circ }$		$8.18\times {10}^{-6}$	$p=0.026$ $\left(1\right)$
Event 1	$T_1-38.49s$	${5.4}^{\circ }$	${13.1}^{\circ }$	$\gtrsim 1.1$	$𝒮/ℬ=11.0$
GRB140521B	$T_{100}=46.59s$	∗${10.1}^{\circ }$		$2.75\times {10}^{-6}$	$p=0.051$ $\left(1\right)$
Event 1	$T_1+98.37s$	${1.6}^{\circ }$	${11.5}^{\circ }$	$\gtrsim 0.79$	$𝒮/ℬ=7.3$
GRB140603A	$T_{100}=138.24s$	∗${2.1}^{\circ }$		$1.86\times {10}^{-5}$	$p=0.025$ $\left(1\right)$
Event 1	$T_1+41.35s$	${1.1}^{\circ }$	${14.9}^{\circ }$	$\gtrsim 1.5$	$𝒮/ℬ=10.1$
Event 2	$T_1-33.78s$	${1.3}^{\circ }$	${38.7}^{\circ }$	$\gtrsim 2.1$	$𝒮/ℬ=0.026$
GRB141029B	$T_{100}=202.44s$	∗${1.0}^{\circ }$		$3.8\times {10}^{-5}$	$p=0.034$ $\left(1\right)$
Event 1	$T_1-10.33s$	${1.6}^{\circ }$	${11.7}^{\circ }$	$\gtrsim 0.70$	$𝒮/ℬ=6.4$
Event 2	$T_1-80.99s$	${1.0}^{\circ }$	${30.2}^{\circ }$	$\gtrsim 0.45$	$𝒮/ℬ=0.003$
GRB150428B	$T_{100}=161.8s$	${0.0003}^{\circ }$		$3.7\times {10}^{-6}$	$p=0.0020$ $\left(0.930\right)$
Event 1	$T_1+71.35s$	${2.9}^{\circ }$	${6.0}^{\circ }$	$\gtrsim 3.2$	$𝒮/ℬ=131.9$
GRB150428D	$T_{100}=32.51s$	∗${6.1}^{\circ }$		$1.53\times {10}^{-6}$	$p=0.024$ $\left(1\right)$
Event 1	$T_1-43.69s$	${4.4}^{\circ }$	${15.2}^{\circ }$	$\gtrsim 0.54$	$𝒮/ℬ=9.4$
GRB150507A	$T_{100}=63.49s$	∗${1.4}^{\circ }$		$1.52\times {10}^{-5}$	$p=0.039$ $\left(1\right)$
Event 1	$T_1+58.24s$	${1.1}^{\circ }$	${20.4}^{\circ }$	$\gtrsim 1.8$	$𝒮/ℬ=2.4$
Event 2	$T_1-74.44s$	${2.1}^{\circ }$	${36.3}^{\circ }$	$\gtrsim 0.69$	$𝒮/ℬ=0.0023$

Table 2. Table 2: Southern Hemisphere track analysis ontime events, following the conventions of Table  1 .

	Time	$\sigma$	$\Delta \Psi$	Fluence/Energy	Significance
GRB110105A	$T_{100}=123.39s$	∗${2.0}^{\circ }$		$2.09\times {10}^{-5}$	$p=0.037$ $\left(1\right)$
Event 1	$T_1+102.0s$	${0.3}^{\circ }$	${13.1}^{\circ }$	$\gtrsim 15$	$𝒮/ℬ=2.2$
GRB110207A	$T_{100}=109.32s$	${0.0132}^{\circ }$		$4.4\times {10}^{-6}$	$p=0.00035$ $\left(0.540\right)$
Event 1	$T_1+87.4s$	${0.3}^{\circ }$	${0.9}^{\circ }$	$\gtrsim 12$	$𝒮/ℬ=271.6$
GRB111205A	$T_{100}=80.38s$	${0.1}^{\circ }$		$1.7\times {10}^{-4}$	$p=0.0023$ $\left(1\right)$
Event 1	$T_1+150.9s$	${18.7}^{\circ }$	${17.3}^{\circ }$	$\gtrsim 482$	$𝒮/ℬ=9.5$
GRB121127A	$T_{100}=3.51s$	${0.08}^{\circ }$		$9.34\times {10}^{-7}$	$p=0.00043$ $\left(1\right)$
Event 1	$T_1+2.42s$	${60.1}^{\circ }$	${79.5}^{\circ }$	$\gtrsim 175$	$𝒮/ℬ=0.85$
GRB121231A	$T_{100}=32.77s$	∗${6.5}^{\circ }$		$2.94\times {10}^{-6}$	$p=0.035$ $\left(1\right)$
Event 1	$T_1+66.5s$	${0.5}^{\circ }$	${13.9}^{\circ }$	$\gtrsim 24$	$𝒮/ℬ=4.2$
GRB130909A	$T_{100}=33.79s$	∗${17.2}^{\circ }$		$1.98\times {10}^{-6}$	$p=0.010$ $\left(0.989\right)$
Event 1	$T_1+14.9s$	${0.2}^{\circ }$	${19.4}^{\circ }$	$\gtrsim 53$	$𝒮/ℬ=30.6$
GRB130924A	$T_{100}=37.1s$	∗${6.0}^{\circ }$		$3.73\times {10}^{-6}$	$p=0.033$ $\left(1\right)$
Event 1	$T_1+92.6s$	${27.1}^{\circ }$	${8.0}^{\circ }$	$\gtrsim 72$	$𝒮/ℬ=1.3$
Event 2	$T_1+6.6s$	${0.4}^{\circ }$	${19.3}^{\circ }$	$\gtrsim 2.8$	$𝒮/ℬ=0.84$
GRB131119A	$T_{100}=34.8s$	∗${7.3}^{\circ }$		$1.85\times {10}^{-6}$	$p=0.025$ $\left(1\right)$
Event 1	$T_1-23.1s$	${0.4}^{\circ }$	${22.9}^{\circ }$	$\gtrsim 16$	$𝒮/ℬ=8.2$
GRB141012A	$T_{100}=37.64s$	∗${3.1}^{\circ }$		$6.64\times {10}^{-6}$	$p=0.014$ $\left(1\right)$
Event 1	$T_1+100.54s$	${11.5}^{\circ }$	${22.4}^{\circ }$	$\gtrsim 114$	$𝒮/ℬ=2.5$
GRB141013A	$T_{100}=82.43s$	∗${3.8}^{\circ }$		$8.81\times {10}^{-6}$	$p=0.017$ $\left(1\right)$
Event 1	$T_1+34.4s$	${17.6}^{\circ }$	${48.0}^{\circ }$	$\gtrsim 459$	$𝒮/ℬ=2.3$
GRB150222C	$T_{100}=74.75s$	∗${11.32}^{\circ }$		$3.84\times {10}^{-6}$	$p=0.020$ $\left(1\right)$
Event 1	$T_1+22.73s$	${0.3}^{\circ }$	${24.5}^{\circ }$	$\gtrsim 31$	$𝒮/ℬ=8.4$
Event 2	$T_1-61.2s$	${0.2}^{\circ }$	${52.3}^{\circ }$	$\gtrsim 50$	$𝒮/ℬ=0.0064$

Equations

$$p+\gamma \to {\Delta }^{+}\to n+{\pi }^{+}\to n+e^{+}+{\nu }_e+{\bar{\nu }}_{\mu }+{\nu }_{\mu }.$$

$${\dot{N}}_{\mathrm{signal}}={\int }_{\Omega }d{\Omega }^{\prime }\int d{E}_{\nu }{A}_{\mathrm{eff}}(E_{\nu },{\Omega }^{\prime })\times {\Phi }_{\nu }(E_{\nu },{\Omega }^{\prime }),$$

$$\mathcal{L}\left(n_s|n_b,\{{\mathbf{x}}_i\}\right)=P_N\prod _{i=1}^N\left[p_s\mathcal{S}\left({\mathbf{x}}_i\right)+p_b\mathcal{B}\left({\mathbf{x}}_i\right)\right],$$

$$P_N=\frac{{\left(n_s+n_b\right)}^N{e}^{-\left(n_s+n_b\right)}}{N!}.$$

$$\mathcal{T}=\ln \left[\frac{\mathcal{L}\left({\hat{n}}_s\right)}{\mathcal{L}\left(n_s=0\right)}\right]=-{\hat{n}}_s+\sum _{i=1}^N\ln \left[\frac{{\hat{n}}_s\mathcal{S}\left({\mathbf{x}}_i\right)}{⟨n_b⟩\mathcal{B}\left({\mathbf{x}}_i\right)}+1\right].$$

$${\mathcal{S}}_{\mathrm{space}}({\vec{x}}_i)=\frac{\kappa }{4\pi \sinh (\kappa )}{e}^{\kappa \cos (\Delta {\Psi }_{i,\mathrm{GRB}})}$$

$$\mathcal{T}=\sum _{c,s}{\mathcal{T}}_{c,s}.$$

$${\Phi }_{\nu }(E_{\nu })={\Phi }_0\times \{\begin{array}{ll}{\epsilon }_b^{-1}{E}_{\nu }^{-1},& E_{\nu }\le {\epsilon }_b\\ E_{\nu }^{-2},& {\epsilon }_b<E_{\nu }\le 10{\epsilon }_b\\ E_{\nu }^{-4}{\left(10{\epsilon }_b\right)}^2,& 10{\epsilon }_b<E_{\nu },\end{array}$$

Full text

astroparticle physics, neutrinos, acceleration of particles, gamma-ray burst: general

††slugcomment: Submitted to Astrophysical Journal

1 Introduction

The sources of ultra-high energy cosmic rays (UHECRs) with energies above $10^{18}\,\mathrm{eV}$ remain unknown as intergalactic magnetic fields deflect UHECRs while they propagate through the universe. In their source, interactions of accelerated hadrons with matter or radiation are expected to produce high energy neutral particles, namely photons and neutrinos. As these particles are chargeless, they propagate through the universe undeflected, meaning they can be associated with astrophysical sources to elucidate the origin of UHECRs (Beatty & Westerhoff, 2009; Kotera & Olinto, 2011). Observations of $\gamma$-rays alone are insufficient to locate hadronic accelerators as purely electromagnetic processes may also produce $\gamma$-rays at these sources. Furthermore, $\gamma$-rays are not ideal messengers at the highest energies because their propagation is hindered by interactions with interstellar media or radiation, reducing their observable distance to the Local Group. On the other hand, neutrinos only interact through the weak force and gravity, allowing them to propagate from their source to Earth unimpeded. The detection of high energy neutrinos from an astrophysical object would constitute unambiguous evidence for hadronic acceleration, revealing UHECR sources (Learned & Mannheim, 2000; Halzen & Hooper, 2002; Anchordoqui & Montaruli, 2010; Anchordoqui et al., 2014).

One possible class of sources for UHECRs are gamma-ray bursts (GRBs), which release immense quantities of $\gamma$-ray radiation on time scales of $10^{-3}-10^{3}\,\mathrm{s}$. The predominant model for GRB phenomenology involves the release of a relativistic fireball (Piran, 2005; Mészáros, 2006; Fox & Mészáros, 2006)—a plasma of electrons, photons, and hadrons—that is triggered by a cataclysmic stellar collapse or binary system merger. Shock waves present within the fireball are capable of accelerating protons and electrons to very high energy through first-order Fermi-acceleration (Krymskii, 1977; Fermi, 1949). As the relativistic electrons are accelerated, they will radiate $\gamma$-rays that contribute to the observable prompt $\gamma$-ray flux once radiative pressure expands the fireball sufficiently so that the plasma becomes optically thin. Should protons in the fireball be accelerated with comparable efficiency and abundance to electrons in the fireball, the cosmological energy density of GRBs is sufficient to explain the measured UHECR flux (Waxman, 1995; Vietri, 1995). Additionally, these protons would interact with the ambient photon field to produce high energy neutrinos primarily through the $\Delta^{+}$ resonance:

[TABLE]

These neutrinos, called prompt neutrinos, would be observable in both temporal and spatial coincidence with the prompt $\gamma$-ray emission of GRBs.

The IceCube Neutrino Observatory (Achterberg et al., 2006; Aartsen et al., 2016b) is currently the most sensitive detector to astrophysical neutrinos. An astrophysical neutrino flux was discovered in neutrino interactions occurring within the detector volume (Aartsen et al., 2013, 2014c), while observations of $\nu_{\mu}+\bar{\nu}_{\mu}$ events from the Northern Hemisphere ($\delta>-5^{\circ})$ later confirmed the discovery (Aartsen et al., 2015c, 2016d). IceCube has not yet observed a neutrino signal associated with GRBs (Abbasi et al., 2010b, 2012a; Aartsen et al., 2015d, 2016a). These results are consistent with the non-detection in multiple years of analysis by AMANDA (Achterberg et al., 2007, 2008) and ANTARES (Adrián-Martínez et al., 2013a, b).

This paper presents a continued search for prompt $\nu_{\mu}+\bar{\nu}_{\mu}$ neutrino track events from GRBs with IceCube (Aartsen et al., 2015d) in three additional years of data, but introduces two additional components to the analysis: an analysis of each observed GRB individually, and an extension of the Northern $\nu_{\mu}+\bar{\nu}_{\mu}$ track event search to the Southern Hemisphere ($\delta\leq-5^{\circ})$, where IceCube is most sensitive to the highest energy neutrinos as Earth absorption attenuates this neutrino signal from the Northern Hemisphere (Connolly et al., 2011). Section 2 describes the prompt neutrino models tested in this analysis. Section 3 then reviews the IceCube detector and data acquisition. The neutrino candidate event characterization and selections performed for the separate Northern Hemisphere track and Southern Hemisphere track analyses are summarized in Section 4. The unbinned maximized likelihood method for discovery of a prompt neutrino signal from GRBs in both the stacked and per-GRB contexts is then outlined in Section 5. Finally, Section 6 presents the results of our all-sky track analysis, with Section 7 providing conclusions and an outlook for future neutrino searches from GRBs with IceCube.

2 GRB Prompt Neutrino Production

The search for neutrinos associated with GRBs in this paper considers only a flux of neutrinos observable during the prompt stage of $\gamma$-ray emission, and does not explicitly test precursor (Razzaque et al., 2003) or afterglow (Waxman & Bahcall, 2000; Murase & Nagataki, 2006) models. In the absence of an observed neutrino flux, we chose to place limits on two classes of GRB models. The first class normalizes the expected neutrino flux to that of the observed UHECR flux. The second class is more detailed and derives an expected neutrino flux from the details of $\gamma$-ray emission for each of the GRBs entering into the analysis.

The principal model describing the phenomenology of GRBs involves a beamed, relativistic fireball of electrons, photons, and hadrons released from a black hole central engine (Piran, 2005; Mészáros, 2006; Fox & Mészáros, 2006). In the standard internal shock fireball model, particle acceleration is achieved through first-order Fermi-acceleration at shock fronts created in collisions of shells of plasmas moving at different speeds within the fireball. The observed prompt $\gamma$-ray emission of GRBs can be produced through synchrotron radiation and inverse-Compton scattering of accelerated electrons or the decays of neutral pions to very energetic photons, while accelerated hadrons within the fireball could possibly escape as UHECRs.

If one assumes that the highest energy cosmic rays are produced through the escape of accelerated hadrons from GRBs, then a subsequent flux of neutrinos can be calculated (Waxman & Bahcall, 1997; Ahlers et al., 2011). The bulk Lorentz factor $\Gamma$ of the fireball is largely unknown and thought to be in the range $100\lesssim\Gamma\lesssim 1000$, though recent multi-wavelength observations of several long GRBs have found values for $\Gamma$ as low as $\gtrsim 10$ (Laskar et al., 2015). This value affects both the normalization and spectral break energy of neutrinos produced in GRB fireballs. A benchmark value of $\Gamma=300$ is taken in the literature (Waxman & Bahcall, 1997) to calculate average neutrino spectra assuming average $\gamma$-ray emission parameters, leading to neutrino spectra that are double broken power laws peaking around $100\,\mathrm{TeV}$. We present flux limits within a range of expected neutrino spectral break energies in Section 6.

Alternatively, one may assume that the predicted neutrino flux is related to the observed $\gamma$-ray emission of individual GRBs (Hümmer et al., 2012; Zhang & Kumar, 2013), allowing direct limits to be placed on emission-generating model parameters. Three representative models are tested, the internal shock fireball model (Hümmer et al., 2012; Zhang & Kumar, 2013), the photospheric fireball model (Rees & Mészáros, 2005; Murase, 2008; Zhang & Yan, 2011), and the Internal Collision-induced Magnetic Reconnection and Turbulence (ICMART) model (Zhang & Yan, 2011; Zhang & Kumar, 2013). Phenomenologically, they primarily differ in neutrino production radius from the GRB central engine, which scales the energy of the primary hadrons, and the density and number of neutrino-producing interactions (Zhang & Kumar, 2013). The neutrino spectra are calculated numerically using the Monte Carlo particle interaction generator SOPHIA (Mücke et al., 2000): protons are propagated in an ambient $\gamma$-ray field with spectrum derived from that measured at Earth, parameterized as a broken power law, and accounting for the redshift of the GRB and the bulk Lorentz factor $\Gamma$ of the fireball. Neutrinos are produced in a full simulation of possible $p\gamma$ interactions, accounting for synchrotron losses of interaction products. We do not, however, consider Fermi acceleration of these products—especially pions and muons—that might significantly enhance neutrino production in GRBs (Klein et al., 2013; Winter et al., 2014). The neutrino fluence of a given GRB at Earth is then determined accounting for the cosmological distance of the source and neutrino oscillations. All models considered here assume that proton acceleration occurs at a single location where $\gamma$-rays are also produced and emitted. In these cases, the predicted prompt neutrino fluence will scale linearly with the proton content of the fireball. When this acceleration location constraint is relaxed and a dynamic GRB outflow is considered, the predicted prompt neutrino fluence is significantly reduced, and well below the sensitivity of IceCube (Bustamante et al., 2015; Globus et al., 2015).

Information about each GRB is gathered from the Gamma-ray Coordinates Network (GCN)111http://gcn.gsfc.nasa.gov and the Fermi GBM database (von Kienlin et al., 2014; Gruber et al., 2014), and is compiled on a publicly accessible website222http://grbweb.icecube.wisc.edu (Aguilar, 2011). The temporal search window $T_{100}$ is defined by the interval between the earliest reported start time $T_{1}$ and the latest reported stop time $T_{2}$ among all observing satellites ($T_{100}=T_{2}-T_{1}$), while the burst localization is chosen from the most precise measurement reported. Similarly, the $\gamma$-ray fluence, break energy, and observed redshift are used as inputs to the neutrino emission calculation. In some cases these values are not measured, and in such cases we adopt conventions previously used by earlier analyses (Abbasi et al., 2012a; Aartsen et al., 2015d, 2016a). We distinguish short GRBs, $T_{100}$ $\leq$ $2\,\mathrm{s}$, from long GRBs, $T_{100}$ $>$ $2\,\mathrm{s}$. In both cases, if a measured photon fluence is unavailable, an average value of $10^{-5}\,\mathrm{erg\,cm^{-2}}$ is used. If the $\gamma$-ray break energy is unmeasured, we assume a value of $200\,\mathrm{keV}$ for long GRBs and $1000\,\mathrm{keV}$ for short GRBs. A redshift measurement is not available for all GRBs; for these GRBs, we use values of 2.15 for long bursts and 0.5 for short bursts. If the redshift was measured, the isotropic luminosity can be approximated from the redshift, photon fluence, and $T_{100}$ (Hümmer et al., 2012); otherwise, an average value of $10^{52}\,\mathrm{erg\,s^{-1}}$ for long bursts and $10^{51}\,\mathrm{erg\,s^{-1}}$ for short bursts is used. The variability time scale is generally unknown, so we use $0.01\,\mathrm{s}$ for long GRBs and $0.001\,\mathrm{s}$ for short GRBs, values that are consistent with current assumptions in the literature (Baerwald et al., 2011; Hümmer et al., 2012; Zhang & Kumar, 2013).

Using benchmark model parameters of $\Gamma$ $=$ $300$ and a baryonic loading—the ratio of fireball energy in protons to electrons—of $f_{p}$ $=$ $10$, the expected model fluxes are shown in Figure 1. We present the neutrino fluence calculations for the analyzed sample of GRBs as a quasi-diffuse flux, assuming an average of 667 potentially observable GRBs per year distributed over the full sky, following previous IceCube publications (Abbasi et al., 2012a; Aartsen et al., 2015d, 2016a). Similar spectra can be calculated for arbitrary values of $f_{p}$ and $\Gamma$.

3 IceCube

IceCube consists of 5160 digital optical modules (DOMs) placed at depths from 1450 m to 2450 m in the Antarctic ice shelf below the South Pole (Achterberg et al., 2006; Abbasi et al., 2012b). Each DOM consists of a photomultiplier tube (PMT) within a glass pressure sphere that detects the Cherenkov radiation of neutrino interaction products (Abbasi et al., 2010a). The photon signal measured by the PMT is digitized in the DOM and relayed to the surface computing laboratory (Abbasi et al., 2009). The detector is arranged in a hexagonal array of 86 vertical cables, each connected to 60 DOMs, called strings. The DOMs are placed on the strings at 17 m intervals, while the inter-string spacing is $\mathord{\sim}125\,\mathrm{m}$. Above the in-ice IceCube detector, 81 ice-tank pairs (each tank containing two DOMs) compose the IceTop surface array (Abbasi et al., 2013). In this paper, the IceTop array is used to veto likely atmospheric muon events from the sky above the detector whose air shower deposits Cherenkov radiation in the surface tanks.

Muons produced by $\nu_{\mu}$ or $\bar{\nu}_{\mu}$ charged current (CC) interactions—appearing as long tracks of light within the IceCube detector—are an especially convenient channel in searches for high energy astrophysical neutrino sources. This is due to their long propagation length in the ice, which increases the effective interaction volume of the detector and yields sub-degree resolution in reconstructed muon direction from the long lever-arm in the detector. Track events from astrophysical sources, however, are difficult to disentangle from IceCube’s primary backgrounds: atmospheric muons and $\nu_{\mu}+\bar{\nu}_{\mu}$ produced in cosmic ray extensive air showers. Atmospheric muons produced in the sky above IceCube trigger the detector at $\gtrsim 2\,\mathrm{kHz}$, while atmospheric $\nu_{\mu}+\bar{\nu}_{\mu}$ trigger the detector at $\mathord{\sim}20\,\mathrm{mHz}$ over the full sky. In recent searches for track events coincident with GRBs, only the Northern Hemisphere sky was analyzed to effectively remove the atmospheric muon background, with the irreducible atmospheric neutrino background remaining. In those analyses, declination $\delta>-5^{\circ}$ is chosen to define the Northern Hemisphere, as above these declinations no atmospheric muons can reach the detector. Such searches have a diminished sensitivity to the highest energy neutrino events compared to the Southern Hemisphere sky ($\delta\leq-5^{\circ}$), as Earth absorption becomes relevant for neutrinos with energy $\gtrsim 1\,\mathrm{PeV}$, though this is partially ameliorated by an increased sensitivity to signal near the analysis horizon. The continued non-detection of a neutrino signal from GRBs compels the extension of these track searches to the Southern Hemisphere.

IceCube data are searched for track events consistent with $\nu_{\mu}$ and $\bar{\nu}_{\mu}$ CC interactions. The previously published Northern Hemisphere GRB track analysis included data from three years of the partially completed IceCube detector and one year of the completed detector (Aartsen et al., 2015d) and included 506 GRBs. This has been extended to include three additional years of IceCube data between May 2012 and May 2015, during which 508 GRBs occurred in the Northern Hemisphere during good detector operation. The newly introduced Southern Hemisphere GRB track analysis has been applied to IceCube data between May 2010 and May 2015, during which IceCube operated for a year in a 79-string configuration and four years of the full IceCube detector. During good detector operation, 664 GRBs occurred in the Southern Hemisphere sky. In total, we searched for neutrino emission from 663 new GRBs, while 509 GRBs that were included in the all-sky cascade GRB analysis (Aartsen et al., 2016a) were analyzed for the first time in the $\nu_{\mu}+\bar{\nu}_{\mu}$ track channel.

4 Event Reconstruction and Selection

High quality track events were selected in IceCube data by the topology of the particle light deposition. Likelihood-based reconstructions fit a track hypothesis to the timing and position of PMT photoelectron pulses of a given event, accounting for photon scattering and absorption in the ice, to obtain the muon’s direction (Ahrens et al., 2004; Aartsen et al., 2014d). Expected light yield probability distribution functions (PDFs) are either defined analytically or from spline fits to simulated light yields from cascades or minimally ionizing muons (Whitehorn et al., 2013; Aartsen et al., 2014d). By parameterizing the behavior of the likelihood space as a function of muon direction near the reconstructed best-fit direction, the angular uncertainty in these reconstructions can be estimated, allowing further selection of high-quality track events. For use in the significance calculation in this analysis, the reconstructed angular uncertainty is determined using the Cramer–Rao lower bound (Cramer, 1945; Rao, 1945) on the covariance of angular direction measures from the inverse of the likelihood Fisher information matrix. Finally, the muon energy as it reaches the IceCube detector, as well as its individual stochastic losses, can also be reconstructed through similar likelihood-based fits to the measured light deposition. These algorithms yield reconstructed muon energies with a resolution of $\mathord{\sim}30\%$, and a resolution of the total deposited energy in the detector along the muon track of $10-15\%$ (Aartsen et al., 2014a).

The track samples in the Northern Hemisphere and Southern Hemisphere are obtained separately as the primary background in each is fundamentally different. In the Northern Hemisphere event sample, most events in low-level data are actually atmospheric muon events from the Southern Hemisphere that are misreconstructed to have an origin from the Northern Hemisphere. The event selection primarily focuses on removing these poorly reconstructed muons. In the Southern Hemisphere event sample, nearly all the events are well-reconstructed atmospheric muons, which must be separated from the muon signal produced in $\nu_{\mu}+\bar{\nu}_{\mu}$ CC interactions. In both cases, the background is characterized using events that are more than $\pm 2\,\mathrm{hr}$ (termed off-time data) from the prompt $\gamma$-ray emission of any analyzed GRB, which avoids contamination of a possible GRB neutrino signal. The event selection was optimized to maximize the retention of a Monte Carlo simulation of interactions of diffuse astrophysical $\nu_{\mu}+\bar{\nu}_{\mu}$ neutrinos with an $E^{-2}$ spectrum, giving the selection sensitivity to the wide range of neutrino production models of Section 2.

The Northern Hemisphere track sample was obtained following the same selection as the Northern Hemisphere GRB track analysis presented by Aartsen et al. (2015d). The background in this portion of the sky is dominated by atmospheric muons with misreconstructed direction. A selection optimized to well-reconstructed $\nu_{\mu}+\bar{\nu}_{\mu}$ signal can efficiently remove most of this background. The parameters that effectively distinguish these events have been described in previously published IceCube point source (Abbasi et al., 2011) and GRB (Aartsen et al., 2015d) searches. These parameters are used in a boosted decision tree (BDT) forest (Freund & Schapire, 1997), a multivariate machine learning algorithm that scores the effective signal-ness of a candidate event, to robustly separate background off-time data from signal simulation. By removing events below a certain BDT value, the most background-like events are eliminated (Aartsen et al., 2014d, 2016a), arriving at the final Northern Hemisphere event sample with a rate of $\mathord{\sim}6\,\mathrm{mHz}$. This data sample is dominated ($\mathord{\sim}80\%$) by the irreducible atmospheric $\nu_{\mu}+\bar{\nu}_{\mu}$ background from the Northern Hemisphere, with the remainder of the background being composed of misreconstructed atmospheric muons. The final BDT score cut was chosen such that the discovery and limit-setting potential of the stacked unbinned likelihood analysis (Section 5) was approximately optimized for a signal with an $E^{-2}$ spectrum produced by all analyzed GRBs, following Aartsen et al. (2016a). The cut may be further adjusted to coincide with the optimal discovery potential of the newly-introduced per-GRB analysis for an $E^{-2}$ spectrum produced by a single, randomly selected GRB in the sample, simulating a detectable neutrino fluence produced in a small subset of the analyzed GRBs. Though an $E^{-2}$ spectrum was used in the optimization, the GRB sensitivity was robust against variation of signal spectra in both Northern and Southern Hemisphere analyses. The $E^{-2}$ spectrum was therefore chosen for generality and consistency with previous IceCube GRB analyses.

The Southern Hemisphere track selection was modeled on recent IceCube point source analyses (Aartsen et al., 2014d, 2017) in this portion of the sky. This selection first removes the bulk of the low-energy atmospheric muon background through cuts on single parameters. Machine learning is then employed to reduce the atmospheric muon background further, especially the high-multiplicity bundles of muons produced concurrently in high-energy cosmic ray air showers and traversing the IceCube detector together. These muon bundles deposit large amounts of light in the detector, and are difficult to distinguish from single high-energy muons by reconstructions of event energy. A BDT forest is trained to distinguish atmospheric muons and muon bundles from the well-reconstructed muons of simulated interactions of $\nu_{\mu}+\bar{\nu}_{\mu}$ with an $E^{-2}$ spectrum. Parameters supplied to the BDT forest are those used in the point source analyses (Aartsen et al., 2014d, 2017), as well as a number of new parameters: 1) an azimuthal measure of the event to regularize artificially preferred directions in background events due to the IceCube detector geometry, 2) the distance of the estimated neutrino interaction vertex from the detector edge measured along the reconstructed track direction to select lower energy neutrino starting tracks (Aartsen et al., 2016c), and 3) the reconstructed muon energy and zenith. These additional parameters result in a more efficient signal selection than the Southern Hemisphere point source selection at all neutrino energies. A cut on the per-event BDT score yields the final event sample (optimized under the same procedure as the Northern Hemisphere event selection), with a background data rate of $2-3\,\mathrm{mHz}$ that is still dominated by atmospheric muons ($\mathord{\sim}4\%$ atmospheric neutrinos).

The final expected signal event rate in this event selection can be determined from Monte Carlo simulation of $\nu_{\mu}+\bar{\nu}_{\mu}$ interactions through

[TABLE]

where $A_{\mathrm{eff}}(E_{\nu},\Omega^{\prime})$ is the effective area of neutrino interaction for an event selection, $\Phi_{\nu}(E_{\nu},\Omega^{\prime})$ is the signal neutrino flux, and the integral is performed over the analysis solid angle $\Omega$ and neutrino energy $E_{\nu}$ range. The effective areas, scaled to all-sky, of the Northern and Southern Hemisphere track selections are shown in Figure 2, compared to the all-sky cascade selection of Aartsen et al. (2016a). The Northern Hemisphere selection is demonstrated to be most sensitive to neutrinos with energy $\lesssim 1\,\mathrm{PeV}$, while the effective area of the Southern Hemisphere selection displays the enhanced sensitivity of this channel to neutrinos above a few PeV. The resonant scattering of $\bar{\nu}_{e}$ with electrons in ice at $6.3\,\mathrm{PeV}$ (Glashow, 1960) is seen in the all-sky cascade effective area, and is yet to be observed by IceCube.

5 Unbinned Likelihood Analysis

Given an ensemble of neutrino events and a set of GRBs, a statistical test is required to distinguish an observation of prompt neutrinos from expected backgrounds. For a sample of $N$ events coincident with GRBs, we calculate the significance of the coincidences by an unbinned likelihood with observed number of signal events $n_{s}$ of the form:

[TABLE]

where $p_{s}=n_{s}/(n_{s}+n_{b})$, $p_{b}=n_{b}/(n_{s}+n_{b})$, and $P_{N}$ is the Poisson probability of the observed event count $N$ given expected signal and background event counts $n_{s}$ and $n_{b}$, respectively:

[TABLE]

The index $i$ runs over the neutrino candidate events, and $\mathcal{S}$ and $\mathcal{B}$ respectively represent the combined signal and background probability density functions (PDFs) for event characteristics $\bm{x}_{i}$. Each of the signal and background PDFs is defined with respect to the time and direction relative to the GRBs, and with respect to event energy. The final test statistic is the logarithm of the likelihood, maximized with respect to $n_{s}$ (maximized at $\hat{n}_{s}$) divided by the background-only likelihood ($n_{s}=0$), which simplifies to:

[TABLE]

The average expected number of background events can be determined from off-time data, denoted as $\langle n_{b}\rangle$.

The time component of the signal and background PDFs, shown as a signal-to-background PDF ratio in Figure 3, is defined by the $T_{100}$ of each burst. The signal time PDF is constant during $T_{100}$, with Gaussian tails before and after the GRB prompt phase. The functional form of the Gaussian tails is chosen to have a smooth transition on either side, and the Gaussian standard deviation $\sigma_{T}$ is chosen to be the same as $T_{100}$, but limited to minimum and maximum values of $2\,\mathrm{s}$ and $30\,\mathrm{s}$, respectively. For simplicity, the signal time PDF is truncated after $\pm 4\sigma$ in each of the Gaussian tails. The background time PDF is constant in this search time window.

Signal neutrinos from GRBs are expected to be spatially associated with the observed GRB location. We define a PDF following the first-order non-elliptical component of the Kent distribution (Kent, 1982):

[TABLE]

where $\Delta\Psi_{i,\mathrm{GRB}}$ is the opening angle between the reconstructed event direction and GRB location, and the concentration term $\kappa$ is given by $\kappa=\left(\sigma_{i}^{2}+\sigma_{\textrm{GRB}}^{2}\right)^{-1}$ in units of radians. The Kent distribution is normalized on the unit sphere and is more appropriate than the typical two dimensional Gaussian representation, especially for events with large uncertainties in the reconstructed direction. The two dimensional Gaussian distribution is recovered for large concentration parameters ($1/\sqrt{\kappa}\lesssim 10^{\circ}$). Representative examples of the Kent distribution with varying directional uncertainties are shown in Figure 4. Data from the off-time sample are used to characterize the background space PDF. Due to the azimuthal symmetry of the IceCube detector, the background can be sufficiently described using only the zenith angle, with PDF normalized over the solid angle of each analysis. A spline is fit to a histogram of background data in $\cos\left(\theta_{\mathrm{zenith}}\right)$. The data histograms and spline fits for the Northern and Southern Hemisphere muon neutrino searches are shown in Figure 5.

One of the most powerful characteristics expected to distinguish GRB neutrinos from atmospheric neutrinos is their respective energy spectra: prompt neutrinos are expected to be produced at high energies where the steeply falling atmospheric spectrum is significantly diminished. The reconstructed energy of muon tracks is used as a proxy for the incoming neutrino energy. Although the muon energy proxy is only a lower bound, it scales with the neutrino energy and is therefore still useful for distinguishing signal from background. For generality, the signal energy PDF is calculated using the reconstructed muon energies of simulated neutrino events with an $E^{-2}$ spectrum. The background energy PDF is taken directly from the off-time data sample’s reconstructed energy spectrum. The Northern and Southern Hemisphere reconstructed energy PDFs are shown in Figure 6 along with the binned PDF ratio values as a function of reconstructed muon energy. The binned PDF ratio is fit with a spline to generalize the ratio to arbitrary reconstructed energies. At high and low energies where the distributions become sparsely populated, the PDF ratio is conservatively limited to the value of the nearest bin with sufficient statistics. In the Northern Hemisphere analysis, the background is largely made up of atmospheric neutrinos. As such, the background energy PDF can be artificially extended to very high energies by using simulated neutrino events with an atmospheric spectrum (Honda et al., 2007). The same technique is not valid for the Southern Hemisphere analysis because the background is composed of atmospheric muons. The simulation of cosmic ray air shower events is significantly more computationally intensive, meaning the statistics of the off-time data sample greatly exceeds our atmospheric muon simulations at final cut level. Furthermore, the simulated atmospheric muon events do not include a simulation of the IceTop detector, invalidating the comparison of the simulation to the off-time data set, which includes an IceTop veto selection. Thus, only off-time data is used to characterize the background energy PDF in the Southern Hemisphere analysis.

In previous IceCube searches for prompt neutrinos from GRBs, the search was performed by stacking all GRBs in each year and channel (i.e. Northern Hemisphere track, Southern Hemisphere track, all-sky cascade). This method, however, diminishes the significance of a concentrated neutrino signal from a single GRB within the stacked sample, as the test statistic treats such events equivalently to if they were distributed among all the GRBs in the sample. To increase sensitivity of the analysis to a signal concentrated in individual GRBs, we have adopted the new strategy of calculating a test statistic $\mathcal{T}_{g}$ for each GRB $g$. We then determine the GRB for which the maximal $\mathcal{T}_{g}$ value is obtained (called the $\max(\{\mathcal{T}_{g}\})$ method). This approach improves the discovery potential for a signal from a single neutrino-bright GRB and naturally moves into real-time style searches, as each GRB would be treated individually upon detection. The $\max(\{\mathcal{T}_{g}\})$ method is preferred to the selection of the most significant per-GRB coincidence (calculated relative to the background-only $\mathcal{T}_{g}$ distribution for GRB $g$) like that done in the IceCube point source searches, as it is less computationally intensive and yields comparable signal discovery potential. The trials-corrected significance of the final $\max(\{\mathcal{T}_{g}\})$ is determined through a comparison with the background expectation of this test statistic calculated over the entire analyzed Northern and Southern Hemisphere GRB sample.

The stacked test statistic was also still calculated for each year and channel to possibly discover a weak neutrino signal distributed over multiple GRBs. This allows the results presented in this paper to be combined with previous results (Abbasi et al., 2012a; Aartsen et al., 2015d, 2016a). A final stacked test statistic is calculated as a sum over individual channel $c$ and year $s$ test statistics:

[TABLE]

This combined test statistic is used to calculate limits on the GRB neutrino models of Section 2 as it is less sensitive to possible background fluctuations than the per-GRB method.

The background-only and background-plus-signal expectations of both stacked and per-GRB analyses are determined from Monte Carlo pseudo-experiments following the same methodology as described by Aartsen et al. (2016a). The sensitivity, both differential and integrated, of the stacked method to a per-flavor quasi-diffuse $E^{-2}$ neutrino spectrum is shown in Figure 7. This sensitivity is calculated for each individual search channel, as well as the final combined sensitivity. The Northern Hemisphere track analysis (combining the results of Aartsen et al. (2015d) with this paper’s extension to three additional years) is seen to be the most sensitive neutrino detection channel. The all-sky cascade and Southern Hemisphere track channels converge in sensitivity to the Northern Hemisphere track within a factor of a few at energies $\gtrsim 1\,\mathrm{PeV}$, while the Southern Hemisphere track analysis is the most sensitive GRB analysis to date for neutrinos $\gtrsim 10\,\mathrm{PeV}$. Each individual channel has sufficient sensitivity to detect a neutrino signal should the per-flavor quasi-diffuse GRB neutrino flux be comparable in magnitude to the measured IceCube astrophysical neutrino flux of $\mathord{\sim}10^{-8}\,\mathrm{GeV\,cm^{-2}\,sr^{-1}\,s^{-1}}$.

6 Results

The final event sample was searched in coincidence with the 508 GRBs of the three-year Northern Hemisphere sample and the 664 GRBs of the five-year Southern sample. Both per-GRB and stacked per-year and channel test statistics were calculated to discover a neutrino signal from GRBs. The results of the per-GRB analysis are presented for the Northern and Southern Hemisphere analyses in Tables 1 and 2, respectively. Here, basic information about the GRBs and coincident events are described, including their timing, angular uncertainty $\sigma$, angular separation $\Delta\Psi$, the measured $\gamma$-ray fluence of the GRB, and the estimated energy of the coincident event. The significance of the coincidences is summarized in two ways. Event signal-to-background PDF ratio values used in the test statistic calculation are provided to estimate relative event importance. The significance of the per-GRB test statistic is then given as a p-value calculated from that GRB’s expected background-only test statistic distribution, constituting that GRB’s pre-trials p-value. In parentheses, the post-trials p-value of this GRB coincidence is given, calculated relative to the combined three-year Northern Hemisphere track and five-year Southern track analysis $\max(\{\mathcal{T}_{g}\})$ test statistic distribution expected from background, respectively.

The most significant coincidence (in both pre-trials and post-trials p-value) was found in the Southern Hemisphere analysis coincident with GRB110207A, a Swift-localized long GRB ($T_{100}=109.32\,\mathrm{s}$) observed at a declination of $-10.8^{\circ}$. This event occurred during the $T_{100}$ of the GRB and had a reconstructed direction within $1^{\circ}$ of the GRB, with a moderate reconstructed muon energy of $E_{\mu}\gtrsim 12\,\mathrm{TeV}$, yielding a signal-to-background PDF ratio of $\mathcal{S}/\mathcal{B}=271.6$. The pre-trials significance is $p=3.5\times 10^{-4}$, making it the single most significant coincidence with a GRB to date in any IceCube GRB neutrino search. Although the event was within $1^{\circ}$ of the GRB location, the angular uncertainty of this event and GRB were $0.3^{\circ}$ and $0.01^{\circ}$, respectively. Combined, these lead to a $\mathord{\sim}3\sigma$ offset in the signal space PDF, reducing the significance of the coincidence. Monte Carlo simulations and reconstructions were performed of muons with similar energy, origin, and light deposition topology to the measured event, establishing that the reconstructed angular uncertainty of $0.3^{\circ}$ is consistent with the median angular resolution of the simulated muons of $0.24^{\circ}$. Furthermore, a full likelihood scan of a more detailed angular reconstruction, which accounts for muon stochastic losses, was performed on this event to verify the quality of the reconstructed direction (Aartsen et al., 2014a). It was found that the two reconstructions are consistent with each other, while the GRB110207A location is $>5\sigma$ from the advanced reconstructed direction, supporting that this event is inconsistent with the GRB location. Additionally, the post-trials significance of this event is $p=0.535$, making it consistent with the background-only hypothesis. Considered together, this event is concluded to be a background coincidence event.

Two additional coincident events were observed in the Northern Hemisphere track analysis that had event significances of $\mathcal{S}/\mathcal{B}\gtrsim 100$: one event in coincidence with GRB131202B, a Fermi-GBM localized long GRB ($T_{100}=86.02\,\mathrm{s}$) at a declination of $21.3^{\circ}$, and one in coincidence with GRB150428B, a Swift localized long GRB ($T_{100}=161.8\,\mathrm{s}$) at a declination of $4.1^{\circ}$. Both events occurred during the $T_{100}$ of the GRBs, and had reconstructed deposited energies above $1\,\mathrm{TeV}$. Due to the short tracks these events produced, each had a relatively large angular uncertainty between $2-3^{\circ}$. The opening angle between each GRB and event pair was a greater than $2\sigma$ deviation with respect to the signal space PDF. Though the pre-trials significances of these coincidences were $0.0069$ and $0.0020$ for GRB131202B and GRB150428B, respectively, correcting for trials these are $0.988$ and $0.930$. The remaining coincident events of Tables 1 and 2 are low significance coincidences, as measured by the event signal-to-background PDF ratios and post-trials p-values. In summary, the set of per-GRB coincidences observed is taken to be consistent with background.

The only coincidence that contributes significantly to a non-zero per-year and channel stacked test statistic is the one coincident with GRB110207A. The significance of the combined Northern and Southern Hemisphere track stacked analysis test statistic is $p=0.42$. Combined with the previously published four-year Northern Hemisphere track (Aartsen et al., 2015d) and three-year all-sky cascade (Aartsen et al., 2016a) analyses, the stacked analysis has a final significance of $p=0.60$, consistent with the background-only hypothesis. For the GRB sample analyzed in this paper, the benchmark standard fireball, photospheric fireball, and ICMART models were expected to yield 2.75, 4.66, and 0.10 neutrino events, respectively. When combined with previously published searches, these models are expected to yield 6.51, 11.02, and 0.25 neutrino events, respectively. Though a number of events have been found temporally coincident with GRBs, none has appeared to be a particularly compelling signal and they have occurred at a rate consistent with background.

Having found results consistent with background, limits can be placed on neutrino production models in GRBs. These amount to calculating the Neyman upper limit (Neyman, 1937) on the flux normalization of these models by determining the fraction of Monte Carlo pseudo-experiments in which such a model would yield a test statistic at least as extreme as that observed. For example, a model can be excluded at the 90% confidence level (CL) should it result in 90% of pseudo-experiments with $\mathcal{T}\geq\mathcal{T}_{\mathrm{obs}}$. Limits calculated account for systematic uncertainties in the ice model, DOM efficiency, and interaction cross sections, which translate to a 10%–20% uncertainty in model limits. The effect of these systematic uncertainties in calculated model limits is determined in a model-dependent way, as their effect is found to be much more pronounced at low energy than at high energy.

Constraints were first determined for a generic double broken power law neutrino flux of the form

[TABLE]

as a function of first break energy $\varepsilon_{b}$ and quasi-diffuse spectral normalization $\Phi_{0}$. These limits are presented in Figure 8 as excluded regions in this parameter space. Two models of neutrino production in GRBs where GRBs are assumed to be the sole origin of the measured UHECR flux are provided in this parameter space: the neutron escape model of Ahlers et al. (2011) and the proton escape model of Waxman & Bahcall (1997), which has been updated with recent measurements of the UHECR flux (Katz et al., 2009). Both models are excluded at over 90% confidence level (CL) with most of the model assumption phase space excluded at over the 99% CL. A thorough reconsideration of whether GRBs can be the sources of UHECRs from Baerwald et al. (2015) shows that the internal shock fireball model is still plausible if cosmic ray protons can efficiently escape the fireball with a low pion-production efficiency for a range of $f_{p}$ and $\Gamma$, which predict neutrino fluxes below the current limits.

Similar constraints were calculated for simple power law spectra consistent with IceCube’s observed astrophysical neutrino flux (Aartsen et al., 2014c, 2015c, 2015a, 2016d), concluding that $\lesssim 0.4\%$ of the astrophysical neutrino flux can be the result of a GRB prompt, quasi-diffuse flux assuming no spectral breaks. This constraint is weakened to a $\lesssim 1\%$ contribution should there be a low-energy spectral break in the astrophysical neutrino flux below $100\,\mathrm{TeV}$.

We also calculated limits for the numerical models of neutrino production in GRBs, where the expected measurable neutrino fluence is determined from the per-GRB $\gamma$-ray spectrum parameters. First, upper limits (90% CL) are calculated for the internal shock fireball, photospheric fireball, and ICMART models using benchmark parameters of the fireball baryonic loading $f_{p}=10$ and bulk Lorentz factor $\Gamma=300$. These are presented in Figure 9, scaling the model fluences to a per-flavor quasi-diffuse flux. Both the internal shock and photospheric fireball models are strongly constrained. The ICMART model significantly reduces the expected neutrino production in GRBs and remains beyond the sensitivity of the combined analysis.

These limits are extended to arbitrary values for $f_{b}$ and $\Gamma$ in the numerical models. Assuming all GRBs in the analyzed sample have identical values for $f_{p}$ and $\Gamma$, limits are presented in Figure 10 as exclusion regions in a scan of $f_{p}$ and $\Gamma$ parameter space. Here, the internal shock and photospheric fireball models are shown to be excluded at the 99% CL for benchmark model parameters. The 90% CL upper limits of all models are improved by about a factor of two compared to those presented in the all-sky cascade analysis (Aartsen et al., 2016a) with the inclusion of this new three year Northern Hemisphere and five year Southern sky $\nu_{\mu}+\bar{\nu}_{\mu}$ analysis. The primary regions in these models that still cannot be constrained require small baryonic loading and large bulk Lorentz factors. The ICMART model is limited in a much smaller interval of possible bulk Lorentz factors ($100<\Gamma<400$) as this model is much less well constrained; only regions of large baryonic loading and small bulk Lorentz factors can be meaningfully excluded.

7 Conclusions

We have performed a search for muon neutrinos and anti-neutrinos in coincidence with 1172 GRBs in IceCube data. This analysis consisted of an extension of previous Northern Hemisphere track analyses to three more years of data, and aa additional search for $\nu_{\mu}+\bar{\nu}_{\mu}$ induced track events in the Southern Hemisphere in five years of IceCube data, which improves the sensitivity of the analysis to neutrinos with energy above a few PeV. Taken together, these searches greatly improve IceCube’s sensitivity to neutrinos produced in GRBs when combined with previous analyses. A number of events were found temporally coincident with these GRBs, but were consistent with background both individually and when stacked together. New limits were therefore placed on prompt neutrino production models in GRBs, which represent the strongest constraints yet on the proposal that GRBs are the primary source of UHECRs during their prompt phase. General models of neutrino emission were first constrained as a function of spectral break energy and flux normalization, excluding much of the current model phase space where GRBs during their prompt emission are assumed to be the sole source of UHECRs in the universe at the 99% CL. Furthermore, models deriving an expected prompt neutrino flux from individual GRB $\gamma$-ray spectral properties were constrained as a function of GRB outflow hadronic content and Lorentz factor $\Gamma$. Models of prompt neutrino production that have not yet been excluded require GRBs to have much lower neutrino production efficiency, either through reduced hadronic content in the outflow, increased $\Gamma$-factor, or acceleration regions much farther from the central engine than the standard internal shock fireball model predicts. This analysis also does not meaningfully address the possible GRB production of neutrinos during their precursor or afterglow phases.

The continuing exclusion of the internal shock fireball model, as well as its photospheric extension, is not altogether surprising, as the radiative efficiency of these models has long been suggested to be insufficient to yield the observed $\gamma$-ray spectra (Fan & Piran, 2006; Zhang et al., 2007; Kumar & Zhang, 2015), unless the distribution of shell $\Gamma$-factors within the fireball is unrealistically broad (Beloborodov, 2000; Kumar & Zhang, 2015). Furthermore, Baerwald et al. (2015) have self-consistently constrained the hypothesis that GRBs are the source of UHECRs under the single-zone internal shock fireball model from the measured UHECR spectrum, $\gamma$-ray measurements, and IceCube limits to GRB neutrino production and cosmogenic neutrinos. These constraints in addition to null results in new IceCube prompt neutrino and cosmogenic neutrino searches increasingly require unphysically large baryonic loading factors or fireball bulk Lorentz factors that may be in tension with multi-wavelength measurements of $\Gamma\gtrsim 10$ for some GRBs (Laskar et al., 2015). Multi-zone internal shock fireball models of neutrino production remain beyond the sensitivity of this work (Bustamante et al., 2015; Globus et al., 2015), and thus are unconstrained.

This paper has introduced a new method for analyzing GRBs on an individual basis, which is adaptable to near real-time analyses for neutrino production in detected GRBs. Though the analysis in this paper improved constraints on neutrino production in GRBs, such constraints operate under the assumption of roughly uniform production across GRBs. Should a rare subclass of GRBs produce a significant neutrino signal, it may still be discoverable with fast follow-up by IceCube and multi-wavelength observations. The all-sky $\nu_{\mu}+\bar{\nu}_{\mu}$ CC interaction channel investigated in this paper is especially promising for this purpose. In addition, the proposed IceCube extension, IceCube-Gen2 (Aartsen et al., 2014b; Ahlers & Halzen, 2014), would increase the detector’s sensitivity to transient astrophysical neutrino sources and possibly reveal GRB neutrino production below IceCube’s current sensitivity. The continued non-detection of a prompt neutrino signal, however, will increasingly disfavor GRBs as a source of UHECRs.

We acknowledge the support from the following agencies: U.S. National Science Foundation-Office of Polar Programs, U.S. National Science Foundation-Physics Division, University of Wisconsin Alumni Research Foundation, the Grid Laboratory Of Wisconsin (GLOW) grid infrastructure at the University of Wisconsin - Madison, the Open Science Grid (OSG) grid infrastructure; U.S. Department of Energy, and National Energy Research Scientific Computing Center, the Louisiana Optical Network Initiative (LONI) grid computing resources; Natural Sciences and Engineering Research Council of Canada, WestGrid and Compute/Calcul Canada; Swedish Research Council, Swedish Polar Research Secretariat, Swedish National Infrastructure for Computing (SNIC), and Knut and Alice Wallenberg Foundation, Sweden; German Ministry for Education and Research (BMBF), Deutsche Forschungsgemeinschaft (DFG), Helmholtz Alliance for Astroparticle Physics (HAP), Research Department of Plasmas with Complex Interactions (Bochum), Germany; Fund for Scientific Research (FNRS-FWO), FWO Odysseus programme, Flanders Institute to encourage scientific and technological research in industry (IWT), Belgian Federal Science Policy Office (Belspo); University of Oxford, United Kingdom; Marsden Fund, New Zealand; Australian Research Council; Japan Society for Promotion of Science (JSPS); the Swiss National Science Foundation (SNSF), Switzerland; National Research Foundation of Korea (NRF); Villum Fonden, Danish National Research Foundation (DNRF), Denmark