Blazar origin of some IceCube events

Abstract

Recently the ANTARES collaboration presented a time dependent analysis of a selected number of flaring blazars to look for upward going muon events produced from the charge current interaction of the muon neutrinos. We use the same list of flaring blazars to look for a possible positional correlation with the IceCube neutrino events. In the context of the photohadronic model we propose that the neutrinos are produced within the nuclear region of the blazar where Fermi accelerated high energy protons interact with the background synchrotron/SSC photons. Although we found that some objects from the ANTARES list are within the error circles of a few IceCube events, the statistical analysis shows that none of these sources have a significant correlation.

1 Introduction

Interactions of ultra high energy cosmic rays (UHECRs) with the background medium photons and protons produce high energy $\gamma$-rays and neutrinos. On their way to Earth the UHECRs can be deflected in the magnetic field and the high energy $\gamma$-rays can be absorbed. So both of these heavenly messengers will lose their directionality. On the other hand neutrinos will be directly pointing to the source, which is why neutrinos are considered as ideal cosmic messengers.

The IceCube detector located at South Pole in Antarctic ice is precisely built to look for high energy neutrinos (above few TeV) by measuring the Cherenkov radiation of the secondary particles created in each neutrino event. The energy deposited by each event and the secondary particles’ direction and topology can be calculated from the trail of the observed Cherenkov light. In 2012 the IceCube collaboration published two years of data (2010–2012) in which 28 neutrino events with energies between 30 and 1200 TeV were observed [1]. Twenty-one of these events are shower-like and the rest are muon tracks. In this analysis two events were PeV neutrino shower events. Adding a third year of analysis in a total 988-days data revealed a total of 37 events, of which nine are track events and the rest are shower events [2]. The shower events have larger angular errors (an average of $15^{\circ }$) than the track events (about $1^{\circ }$). These events have flavors, directions, and energies inconsistent with those expected from the atmospheric muon and neutrino backgrounds. So the study of arrival directions is helpful to find sources of high energy neutrinos and the relevant acceleration mechanism acting within the source.

The isotropic distribution of these IceCube neutrino events suggests a contribution from at least some extragalactic sources. There exist different types of potential astrophysical sources to produce UHECRs and hence high energy neutrinos and $\gamma$-rays. The list includes: $\gamma$-ray bursts (GRBs) [3], core of active galactic nuclei (AGN) [4], high energy peaked blazars (HBLs) [5–7], starburst galaxies [8], and sources from Galactic center [9]. In Ref. [5] many positional correlations of BL Lac objects and galactic pulsar wind nebulae with the IceCube events are shown. There are also nonstandard physics interpretations of these IceCube events from the decay of superheavy dark matter particles, the leptoquark interaction, and the decay of exotic neutrinos [10] (see [11] for a recent review).

Recently the ANTARES collaboration presented a time dependent analysis [12] to look for upward going muon tracks by charge current interaction of $\nu _{\mu }$ from flaring blazars selected from the Fermi-LAT and TeV $\gamma$-ray observed by ground based telescopes H.E.S.S, MAGIC, and VERITAS, respectively. In this analysis the most significant correlation was found with a GeV flaring blazar from the Fermi-LAT catalog. However, the post-trial probability estimate shows that the event was compatible with background fluctuations. In this work we would like to analyze the above list of Fermi-LAT flaring blazars to see if there is any correlation with the IceCube neutrino events. We use the unbinned maximum likelihood method (MLM) with two different values of the spectral index for our analysis of the positional correlation of these objects.

2 Candidates

Blazars are believed to be the most likely candidates to produce UHECRs and neutrinos [5–7]. These are extragalactic objects characterized by relativistic jets with a small viewing angle with respect to the line of sight and are powered by a supermassive black hole in the center of their respective galaxy. These objects are also efficient accelerators of particles through shock or diffusive Fermi acceleration processes with a power-law spectrum given as ${\mathrm{d}N}/{\mathrm{d}E} \propto E^{-\kappa }$, with the power index $\kappa \ge 2$ [13]. Protons can reach an ultra high energy through the above acceleration mechanisms. Fractions of these particles escaping from the source can constitute the UHECRs arriving on Earth. These objects also produce high energy $\gamma$-rays and neutrinos through pp and/or $p\gamma$ interactions [14]. The classification for these sources is according to the properties of their emission lines: if a strong broad emission line in the optical spectrum is present, it is classified as Flat Spectrum Radio Quasar (FSRQ), otherwise it is a BL Lacerate (BL Lac) object. Depending on the frequency of the first peak, the BL Lac objects are further classified into low (LBL), intermediate (IBL), and high energy (HBL) peaked objects.

The ANTARES collaboration searched for high energy cosmic muon neutrinos using the data taken during the period August 2008 to December 2012. The collaboration selected 41 very bright and variable Fermi-LAT blazars with significant time variability and having the flux ${>}10^{-9}$ photons ${\mathrm{cm}^{-2}\,\mathrm{s}^{-1}}$ for the $\gamma$-ray energy above 1 GeV. They have also selected seven TeV flaring objects reported by H.E.S.S., MAGIC, and VERITAS telescopes with the expectation that the TeV $\gamma$-rays may be correlated with the neutrino events. From the 41 Fermi blazar list, 33 are FSRQs, seven are BL Lacs and one is unknown. Similarly from the list of 7 TeV flaring blazars, one is FSRQ, and six are HBLs. It shows that both FSRQs and HBLs are probable sources of very high energy neutrinos and can be possible sources for some of the IceCube events. It is suggested that UHECRs are accelerated in the inner jet of the FSRQ and interact with the background from the broad-line region (BLR), synchrotron radiation or the photon from the accretion disk [15–17].

In a previous article [7] we proposed that photohadronic interactions of the Fermi accelerated high energy protons with the background photons in the nuclear region of the HBLs and AGN are responsible for some of the IceCube events. These objects were observed in multi-TeV $\gamma$-rays and some had also flaring. In this model it is assumed that the flaring of a blazar in high energy $\gamma$-ray occurs within a compact and confined region with a comoving radius $R'_f$ inside the blob of radius $R'_b$ [18] (henceforth $'$ implies a jet comoving frame). In the inner region, the photon density $n'_{\gamma ,f}$ is very high compared to the photon density $n'_{\gamma }$ in the outer region i.e. $n'_{\gamma ,f} \gg n'_{\gamma }$. Fermi accelerated high energy protons undergo photohadronic interaction with the seed photons in the inner region in the self-synchrotron Compton (SSC) regime through the intermediate $\Delta$-resonance. On the other hand, in a normal blazar jet, the photohadronic process is not an efficient mechanism to produce multi-TeV $\gamma$-rays and neutrinos, because $n'_{\gamma }$ is low, which makes the optical depth $\tau _{p\gamma }\ll 1$. But the assumption of a compact inner jet region overcomes this problem where the optical depth of the $\Delta$-resonance process is $\tau _{p\gamma }=n'_{\gamma ,f} \sigma _{\Delta } R'_f$ and $n'_{\gamma ,f}$ is unknown. We can estimate the photon density in this region by assuming that the Eddington luminosity is equally shared by the jet and the counter jet in the blazar. For a given comoving photon energy $\epsilon '_{\gamma }$ in the synchrotron/SSC regime we can get the upper limit on the photon density as $n'_{\gamma ,f} \ll L_{Edd}/(8\pi R'^2_f \epsilon '_{\gamma })$. Also by comparing the proton energy loss time scale $t'_{p\gamma }\simeq (0.5\,n'_{\gamma ,f}\sigma _{\Delta })^{-1}$ and the dynamical time scale $t'_{d}=R'_f$ we can estimate $n'_{\gamma ,f}$, so that the production of multi-TeV $\gamma$-rays and neutrinos takes place. In order not to have over production of neutrinos and $\gamma$-rays, we can assume a moderate efficiency (a few percents) by taking $\tau _{p\gamma } < 1$, which gives $n'_{\gamma ,f} < (\sigma _{\Delta } R'_f)^{-1}$. In this work we assume $1~\%$ energy loss of the UHE protons in the inner region on the dynamical time scale $t'_d$ corresponding to an optical depth of $\tau _{p\gamma } {\sim }0.01$ and $n'_{\gamma ,f}\sim 2\times 10^{10}\, R'^{-1}_{f,15}\, \mathrm{cm}^{-3}$. Here the inner blob radius $R'_f$ is expressed as $R'_f=10^{15} R'_{f,15}\, \mathrm{cm}$ and $R'_{f,15} {\sim }1$ [7].

In the photohadronic interaction, the intermediate $\Delta$-resonance produced will give both high energy neutrinos and $\gamma$-rays, and the relation between the seed photon and the neutrino energy is given by

$$\begin{aligned} E_{\nu }\epsilon _\gamma =0.016 \frac{\Gamma \delta }{(1+z)^2}\, \mathrm{GeV}^2, \end{aligned}$$

(1)

where $E_{\nu }$ and $\epsilon _{\gamma }$ are, respectively, the observed neutrino energy and the background photon energy. The source is located at a redshift z and the bulk Lorentz factor of the jet is $\Gamma$. The Doppler factor is given by $\delta$. But for FSRQ and BL Lac objects $\Gamma \simeq \delta$. So if z and $\Gamma$ of a blazar are known we can estimate the $\epsilon _{\gamma }$ from the given $E_{\nu }$. The neutrino flux is given as [19]

$$\begin{aligned} F_{\nu }=\sum _{\alpha } \int _{E_{\nu 1}(1+z)}^{E_{\nu 2}(1+z)}dE_{\nu } E_{\nu } J_{\nu _{\alpha }} (E_{\nu }), \end{aligned}$$

(2)

where for all neutrino flavors $\alpha$ (e, $\mu$, and $\tau$), a power-law spectrum of the form

$$\begin{aligned} J_{\nu _{\alpha }}(E_{\nu }) = A_{\nu _{\alpha } }\left( \frac{E_{\nu }}{100\, \mathrm{TeV}}\right) ^{-\kappa } \end{aligned}$$

(3)

is taken. The normalization constant $A_{\nu _\alpha }$ is given by

$$\begin{aligned} A_{\nu _\alpha }=\frac{1}{3} \frac{N_{\nu }}{T \Sigma _{\alpha } \int _{E_{\nu _1}}^{E_{\nu _2}} \mathrm{d}E_{\nu } A_{\mathrm{eff},\alpha } (E_{\nu } ) \left( \frac{E_{\nu }}{100\, \mathrm{TeV}}\right) ^{-\kappa } }, \end{aligned}$$

(4)

where $N_{\nu }$ is the number of neutrino events and $A_{\mathrm{eff},\alpha }$ is the effective area for different neutrino flavors. The energy integrals are done in the limit 25 TeV to 2.2 PeV. The time period $T=988$ days is used [1] for the calculation of normalization constant.

3 Unbinned maximum likelihood method

To identify the possible sources of IceCube events we employ the Unbinned Maximum Likelihood Method (MLM) [20] to find a spatial correlation between the blazar sample under consideration and the IceCube events. The signal and the background weights are not separable for an object and both contribute to the likelihood function, which is given by the product of the individual probability densities for the IceCube events as [21]

$$\begin{aligned} \mathcal{L} (n_{s}, \mathbf {x}_{s}) = \prod ^{N}_{i=1} \left[ \frac{n_s}{N} S_i (\mathbf {x}_{s}) + \left( 1-\frac{n_s}{N}\right) B_i \right] , \end{aligned}$$

(5)

where N is the number of IceCube events we take into account, $n_s/N$ is the weight associated with the signal probability density function (PDF); its values vary between 0 and 1. The background PDF depends on the neutrino energy and the declination, which is expressed as

$$\begin{aligned} B_i = \mathcal{B}(E_i,\delta _i). \end{aligned}$$

(6)

The background is constructed from the integrated effective areas of the IceCube 79 strings configuration [22]. The neutrino effective area depends on the detector geometry and the absorption of the neutrinos by the Earth. The background PDF takes into account the contribution from the atmospheric muon neutrinos. Above $\sim$100 TeV, neutrinos from the decay of charm hadrons $D^{\pm }, D^0$ contribute to the background neutrino flux known as prompt flux. Equal numbers of neutrinos and anti-neutrinos of electron and muon flavors are produced in this process. However, the prompt flux is poorly understood in the high energy limit. For the background calculation we also include the contribution from the prompt background [23, 24].

The signal PDF is defined as the product of a spatial term and the energy term as shown below:

$$\begin{aligned} S_i=\mathcal{S}_i(|\mathbf{x}_i- \mathbf{x}_s|,\sigma _i)\, \mathcal{E}_i (E_i,\delta _i,\kappa ), \end{aligned}$$

(7)

where we have defined

$$\begin{aligned} \mathcal{S}_i (\mathbf {x}_{s})= \frac{1}{2 \pi \sigma ^{2}_i} \mathrm{e}^{- \frac{| x_i - x_s |^2}{2 \sigma _i^2} }, \end{aligned}$$

(8)

which is a Gaussian function [25]. In Eq. (8), $| x_i - x_s |^2$ is the space angle difference between the source and the reconstructed event direction and $\sigma _i$ is the standard deviation of the ith IceCube angular error distribution. We also define

$$\begin{aligned} \delta \chi ^2=\frac{|x_i-x_s|^2}{\Delta }. \end{aligned}$$

(9)

The value of $\delta \chi ^2 \le 1$ signifies that the object is inside the median angular error $\Delta$ of the IceCube event. The signal energy PDF $\mathcal{E}_i$ depends on the event energy, the spectral index $\kappa$, and the declination. Here we use $\kappa =2$ and 2.5 for our analysis.

The ANTARES analysis takes into account both the temporal and the energy dependence of the flaring events whereas our analysis is independent of the time. The observed IceCube events can be modeled by taking into account two hypotheses: (1) the events could be produced by atmospheric muons and the muon neutrinos (background), or (2) from an astrophysical source which also includes the background contribution. A good test of compatibility is the ratio of these two hypotheses. We can take the ratio of the likelihood with the background of unique weight ($n_s=0$) and the maximized likelihood of the second hypothesis with the corresponding $n_s$ values defined as $n_s=n^*_s$. Now to evaluate each point source we use this Test Statistic (TS) taking minus twice the log of the likelihood ratio,

$$\begin{aligned} \mathrm{TS} = - 2\, \text {log} \left[ \frac{\mathcal{L} (n_s = 0)}{\mathcal{L} (n_s=n^*_s)} \right] . \end{aligned}$$

(10)

For this procedure we use a full-sky IceCube events. For our present analysis, we take into account 36 events out of reported 37 events (event 32 is excluded in the present analysis because its energy and direction are not reported). We calculate the significance of each source location, running 10,000 simulations in which the declination of each IceCube sample event is fixed but the right ascension is randomized. The p value is calculated as the number of simulations with $\mathrm{TS}_{(\mathrm{sim})} \ge \mathrm{TS}$ divided by the total number of simulations for a given source, where $\mathrm{TS}_{(\mathrm{sim})}$ is the TS value obtained from the simulation. Also, the posterior p value for each object is estimated as the fraction of the randomized simulations that yields an equal or higher TS value for at least one of the 41 ANTARES sources. The compatibility of the second hypothesis depends on the estimate of the posterior p value. If the posterior p value is close to unity then it is consistent with the background.

4 Results

In the context of recent IceCube results, we analyzed the 41 flaring blazars taken from the Fermi-LAT catalog which were previously studied by the ANTARES collaboration to look for possible temporal and spatial correlation [12]. We have also analyzed the 7 TeV flaring objects as discussed by ANTARES collaboration for the possible spatial correlation with the IceCube events. In fact all these seven objects are there in the TeVCat [26] which we had already analyzed in Ref. [7] and found that the only HBL, PG 1553+113 has the positional correlation with the IceCube event 17. So we do not discuss these seven flaring objects here any more. For our analysis of the possible correlation of IceCube events with the ANTARES sources we use the unbinned MLM and two different values of spectral index $\kappa =2$ and 2.5. We also do the separate analysis with and without the contribution from the prompt flux coming from the charm hadron decay. Our results are summarized in Table 1.

All the 28 shower events with their individual errors and the eight track events are shown in the sky map with equatorial coordinates in Fig. 1. The positions of ten FSRQs and two BL Lac objects are also shown in the sky map.

4.1 Spectral index $\kappa =2$

From the 41 Fermi blazars of ANTARES list, 32 objects have TS $> 0$ for the spectral index $\kappa =2$ without the prompt contribution to the background. However, this number reduces to 19 when we include the charm contribution.

Table 1 The objects which are in the error circles of the IceCube events (ID in the third column) are given in the first column. Below each object we also put their coordinates, right ascension, and declination ( R.A., Dec.) in degrees (this table is given in equatorial coordinates). The second column gives the type of object and below this we also give its redshift (z) and the bulk Lorentz factor ($\Gamma$). In the fourth column, the $\delta \chi ^2$ of the object is given. In the fifth and the sixth columns the deposited neutrino energy $E_{\nu }/\mathrm{TeV}$ and the corresponding seed photon energy $\epsilon _{\gamma }/\mathrm{keV}$ are given. In columns seven and eight the values of the $n^*_s$ and TS are given from the maximum likelihood method. In columns nine and ten the p value and the posteriori p value (post p value) are also shown. The last three objects are without (upper value) and with ( lower value) the prompt contribution to the background PDF

From the above 32 objects 12 are within the median angular error of at least one IceCube event having $\delta \chi ^2 \, < 1$. The FSRQ PKS 0235-618 is the only object associated with two IceCube events (7, 20). The FSRQs PKS 2326-502, PKS 0208-512, PKS 0235-618, and PKS 0244-470 are within the error circle of event 7, while the FSRQs 3C454.3, B22308+34, CTA102, and PKS 2227-08 are within the error circle of event 34. Another two FSRQs, PMNJ2345-1555 and PMNJ 2331-2148, are within the error circle of the IceCube event 21. The BL Lac objects OJ287 and PKS0805-07 are coincident with the events 26 and 27, respectively. All the relevant parameters of the above objects are shown in Table 1.

The posterior p values of all the above 12 objects are $\ge$99 %. This shows that our result (without the prompt contribution to the atmospheric background) is consistent with the background fluctuation.

By including the prompt contribution to the background we found that 19 objects have TS $> 0$ of which only three objects, two BL Lac objects (OJ287, PKS 0805-07) and one FSRQ (3C454.3) are within the median angular error of three IceCube events (26, 27, 34). These three objects are shown in the table.

We observed that the background photon energy $\epsilon _{\gamma }$ for most of the events are below <40 keV, which shows that the photon density $n'_{\gamma ,f}$ can be large in the inner region of the jet. By assuming a conservative $1~\%$ energy loss by the UHE protons we get a photon density in the inner region of $n'_{\gamma ,f}\sim 2\times 10^{10}\, \mathrm{cm}^{-3}$ which has a radius $R'_f\sim 10^{15}\, \mathrm{cm}$. The estimate of the $R'_f$ value depends on the outer blob radius $R'_b$, while the latter parameter is adjusted to fit the spectral energy distribution (SED) in the leptonic model of the objects. However, for most of the objects $R'_b > 10^{15}$ cm is taken to fit the SED [7]. So, here we take $R'_f\sim 10^{15}\, \mathrm{cm}$ for the estimation of $n'_{\gamma ,f}$. The simulation shows that $0< {\text{ T }S} \,< 1$ for all the objects.

The diffuse neutrino flux $F_{\nu }$ for all these objects is $2.31\times 10^{-9}\,$ $\mathrm{GeV}\,\mathrm{cm}^{-2}\,\mathrm{s}^{-1}\,\mathrm{sr}^{-1}$ because for $\kappa =2$, the integral in Eq. (2) is independent of the redshift, so the $F_{\nu }$ is the same for all the objects. On the other hand it varies slightly for $\kappa \ne 2$. The high posterior p value for all these objects shows that our result is consistent with the background fluctuation. We also repeated the simulation for $\kappa =2.5$. As $\kappa$ changes, the flux reduces, and we found that none of the 41 objects satisfy the condition TS $>0$.

5 Discussion

The ANTARES collaboration looked for possible temporal and spatial correlation of 41 flaring objects selected from the Fermi-LAT catalog. We analyzed the same objects for the possible spatial correlation with the IceCube events. For our analysis, we take into account the energy dependence of both the background and the signal constructed from the data of the 79 IceCube string configuration. We consider two different values of the spectral index, 2 and 2.5 and also analyzed our results with and without the prompt contribution to the atmospheric neutrino flux. We observed that, from the 41 flaring objects, for $\kappa =2$, the MLM gives 12 objects (without prompt flux contribution) and three objects (with prompt flux contribution) within the error circle of some IceCube events. For these objects we have also estimated the neutrino flux. However, for all these possible candidates, the TS value is very small, which leads to very high posterior p values $\ge$99 % and is consistent with the background fluctuation. It is possible that the high energy neutrino flux from these objects is much below the IceCube limit, or blazars may not have a sufficiently powerful central engine to produce very high energy cosmic rays. So most of the events in IceCube might be from some other type of sources. We have to wait for more data to look for a possible correlation of FSRQs and BL Lac objects with the IceCube events.

Fig. 1

The sky map is shown in equatorial coordinates with 37 IceCube events and their individual errors (only for shower events). Here $+$ corresponds to a shower event and the $\times$ sign corresponds to a track event with their corresponding event ID. We have also shown the positions and names of the blazars which are within the median angular error of the IceCube events and have a TS value >0. The objects in blue color are FSRQs and in red color are BL Lac