Searching for High-Energy Neutrinos from Ultra-Luminous Infrared Galaxies with IceCube

This work presents an IceCube search for high-energy neutrinos from Ultra-Luminous Infrared Galaxies (ULIRGs). ULIRGs are the most luminous infrared objects in the sky, with infrared luminosities exceeding 1012 solar luminosities. They are mainly powered by starbursts that exhibit star-formation rates larger than 100 solar masses per year. In addition, an active galactic nucleus (AGN) can also contribute significantly to the ULIRG luminosity output. The acceleration of hadrons, and consequently the production of neutrinos, can occur both in starburst and AGN environments. As such, ULIRGs form a source population that could be responsible for a significant fraction of the diffuse neutrino flux observed by IceCube. In this study we perform a stacking analysis on a representative sample of 75 ULIRGs with redshift z ≤ 0.13 using 7.5 years of IceCube data. We find no evidence for astrophysical neutrinos correlated with our selection of ULIRGs. We therefore compute upper limits on the contribution of the ULIRG source population to the diffuse neutrino observations, and also use these limits to constrain model predictions.


Introduction
Since the discovery of a diffuse flux of astrophysical high-energy neutrinos by the IceCube collaboration in 2013 [1], the origin of these neutrinos remains largely unknown. Ultra-luminous infrared galaxies (ULIRGs) form a source class that could be responsible for a significant fraction of the diffuse neutrino observations [2,3,4]. ULIRGs possess infrared (IR) luminosities L IR ≥ 10 12 L between 8-1000 µm, making them the most luminous objects in the IR sky [5]. They are mainly powered by starburst nuclei with star-formation rates 100 M yr −1 , with a possible secondary contribution from an active galactic nucleus (AGN). Both the starburst and AGN components could be responsible for hadronic acceleration and neutrino production [6,7]. In this work we present an IceCube search for high-energy neutrinos from ULIRGs. Additional details about this study can be found in a recently submitted publication [8].

Stacking Analysis of ULIRGs
We first obtain an initial selection of 189 ULIRGs taken from three catalogs mainly based on data of the Infrared Astronomical Satellite (IRAS ) [9,10,11]. Subsequently, we obtain a representative sample of the local ULIRG population, which we define as a complete sample of ULIRGs up to a given redshift. We therefore perform a redshift cut on our initial ULIRG selection by finding the redshift up to which the least luminous ULIRGs (L IR = 10 12 L ) can be observed given a conservative IRAS sensitivity at 60 µm of 1 Jy. As such, we obtain a final representative selection of 75 ULIRGs within a redshift z ≤ 0.13. The IceCube Neutrino Observatory is a 1 km 3 detector that consists of 5160 digital optical modules (DOMs) distributed over 86 strings deep in the South Pole ice. These DOMs detect the Cherenkov light emitted by secondary charged particles produced in the interactions of neutrinos (and antineutrinos, but we will not make the distinction in the following) with the ice in or near IceCube. In particular, charged-current interactions of muon neutrinos produce muons which leave track-like signatures in the detector and have an angular resolution below 1 • for muon energies exceeding 1 TeV. Therefore, for our analysis we use the gamma-ray follow-up (GFU) data set [12], which consists of well-reconstructed tracks. The GFU sample contains 7.5 years of data collected with the full 86-string IceCube configuration between 2011-2018. However, the GFU sample consists mostly of background events, since tracks are also signatures of atmospheric muons (detection rate of several kHz) and atmospheric neutrinos (detection rate of several mHz) produced in cosmic-ray air showers. Quality cuts and neutrino event selections reduce the all-sky event rate of the GFU sample to the atmospheric neutrino level, namely 6.6 mHz, for a total of ∼1.5 million events.
To search for astrophysical neutrinos from our selection of M = 75 ULIRGs in the data that is dominated by atmospheric background, we perform a maximum-likelihood analysis [13]. In addition, we apply a source stacking method in order to enhance the sensitivity of the analysis [14]. The likelihood is given by Here, we fit for the number of astrophysical signal events, n s , and the power-law spectral index γ of the energy spectrum. The probability distribution functions (PDFs) of the background, B i , and signal of source k, S k i (γ), are evaluated for each event i ∈ {1, 2, ..., N } of the GFU sample. The background PDF is constructed from scrambled data, while the signal PDF is generated from simulations. For the spatial part of the signal PDF we assume a two-dimensional Gaussian centered around each source k, while for its energy component we assume an unbroken E −γ power-law spectrum. The stacking weight w k ∝ t k r k is proportional to the detector response of source k, r k , and a theoretical weight t k . We set the theoretical weight equal to the total IR flux of source k, F IR = L IR /4πd 2 L , with d L the luminosity distance of the source determined from its redshift measurements. We then construct a test statistic Here,n s andγ are the best-fit parameters that maximize the likelihood. We can distinguish a possible signal component from ULIRGs in the data from a background-only scenario, and also obtain a p-value of the analysis. We test the performance of our analysis by simulating signal events according to an unbroken power-law spectrum, normalized at E 0 = 10 TeV. We define the sensitivity at 90% confidence level (CL) as the amount of pseudo-signal required to obtain a p-value ≤ 0.5 in 90% of simulated pseudo-experiments. In addition, we define the 3σ (resp. 5σ) discovery potential as the amount of pseudo-signal required to obtain a p-value ≤ 2.70 × 10 −3 (resp. ≤ 5.73 × 10 −7 ) in 50% of simulated pseudo-experiments. The 90% sensitivity, 3σ, and 5σ discovery potentials of our ULIRG stacking analysis are shown in Fig. 1. In particular, the right panel of Fig. 1 indicates that the analysis is more sensitive towards harder spectra. This effect follows from the fact that harder spectra are easier to distinguish from the relatively soft E −3.7 spectrum of the atmospheric background.

Results and Interpretation
The stacking analysis yields no evidence for astrophysical neutrinos originating from ULIRGs. We fit a total number of astrophysical neutrino eventsn s = 0, such thatγ remains undetermined.  This result corresponds to a p-value = 1.0. We therefore set upper limits at 90% CL, corresponding with the 90% sensitivity in Fig. 1, on the stacked flux of our 75 selected ULIRGs within z ≤ 0.13, Φ z≤0. 13 νµ+νµ . These limits are extrapolated to upper limits on the diffuse flux of the full ULIRG population op to a redshift z max , Here, we assume that our sample of 75 ULIRGs is representative for the complete ULIRG source population. The factor ξ z effectively integrates over the redshift evolution of the sources [15]. In this work, we adopt the parameterization of [4] for the ULIRG source evolution, H(z) = (1+z) m , with m = 4 for z ≤ 1 and m = 0 for z > 1. The left panel of Fig. 2 shows the integral limits for unbroken E −2.0 , E −2.5 , and E −3.0 powerlaw spectra, which are compared to the diffuse neutrino observations of [16,17]. The E −2.0 and E −2.5 limits constrain the contribution of the ULIRG source population up to energies of ∼3 PeV and ∼600 TeV, respectively. Analogously, the right panel of Fig. 2 shows the quasidifferential limits, found by computing the E −2.0 limit in each bin of energy decade. We find that the contribution of ULIRGs to the diffuse observations is constrained between 10-100 TeV and 0.1-1 PeV.
We also compare our ULIRG upper limits with model predictions of starburst reservoir models. First, the comparison of our E −2.0 limit with the model by He et al. [2], who consider a neutrino flux powered by hypernovae in ULIRGs, is shown in the left panel of Fig. 3. We find that their prediction lies at the level of our upper limit. More data is required in order to validate or constrain this model in a follow-up analysis. Second, we compare our E −2.12 limit with the model of Palladino et al. [3] in the right panel of Fig. 3. In this case, we parameterize the source evolution of the sources according to the star-formation rate [18], H(z) = (1 + z) m , with m = 3.4 for z ≤ 1 and m = −0.3 for 1 < z ≤ 4. In this model, a generic population of hadronically-powered gamma-ray galaxies (HAGS), such as ULIRGs, is proposed as a source class responsible for the bulk of the diffuse neutrino observations. Our results exclude ULIRGs as the sole population of HAGS responsible for the diffuse neutrino flux.
Last, we compare our results to the beam-dump model by Vereecken & de Vries [4], who consider the AGN component of ULIRGs as the main source of high-energy neutrinos.  Their model strongly relies on the electron-to-proton luminosity ratio f e , which carries a large uncertainty. By fitting the model to our E −2.0 upper limit, we obtain an order-of-magnitude estimation of a lower limit on this f e . We find a lower limit for the electron-to-proton luminosity ratio f e 10 −3 . This value is comparable with previously obtained lower limits on f e in the context of obscured AGN analyzed with IceCube [19,20].

Conclusion
We presented the results of an IceCube stacking search for high-energy neutrinos from a representative local sample of 75 ULIRGs within z ≤ 0.13 using 7.5 years of data. No evidence for astrophysical neutrinos originating from these ULIRGs was found in the data. Upper limits were set on the stacked flux of this representative local ULIRG sample, which were subsequently extrapolated to upper limits on the diffuse neutrino flux of the full ULIRG source population. The integral limits for unbroken E −2.0 and E−2.5 constrained the contribution of ULIRGs to the diffuse neutrino observations up to ∼3 PeV and ∼600 TeV, respectively. The quasi-differential limits constrained the ULIRG contribution between 10-100 TeV and 0.1-1 PeV.
In addition, we compared our results with model predictions. First, we found that the prediction by He et al. [2] is at the level of our upper limit. A follow-up study with more years of data is therefore required to validate or exclude this model. Second, in the context of the work by Palladino et al. [3], we excluded ULIRGs as the sole population of HAGS that could be responsible for the diffuse neutrino observations. Last, we estimated a lower limit on the electron-to-proton luminosity ratio f e , which is the most uncertain parameter in the AGN beam-dump model by Vereecken & de Vries [4]. Our lower limit, f e 10 −3 , is consistent with previously reported values.