A $\gamma$-ray Quasi-Periodic modulation in the Blazar PKS 0301$-$243?

TL;DR

This paper reports a significant quasi-periodic gamma-ray signal in the blazar PKS 0301-243, suggesting possible evidence of a binary supermassive black hole based on analysis of Fermi-LAT data.

Contribution

First detection of a high-confidence gamma-ray quasi-periodic modulation in PKS 0301-243 using multiple analysis techniques.

Findings

Quasi-periodicity of 2.1±0.3 years detected at ~5σ significance.

The periodicity is observed in a data set roughly four times longer than the period.

Potential indication of a binary supermassive black hole in the blazar.

Abstract

We report a nominally high-confidence $\gamma$-ray quasi-periodic modulation in the blazar PKS 0301$-$243. For this target, we analyze its \emph{Fermi}-LAT Pass 8 data covering from 2008 August to 2017 May. Two techniques, i.e., the maximum likelihood optimization and the exposure-weighted aperture photometry, are used to build the $\gamma$-ray light curves. Then both the Lomb-Scargle Periodogram and the Weighted Wavelet Z-transform are applied to the light curves to search for period signals. A quasi-periodicity with a period of $2.1\pm0.3$ yr appears at the significance level of $\sim5\sigma$, although it should be noted that this putative quasi-period variability is seen in a data set barely four times longer. We speculate that this $\gamma$-ray quasi-periodic modulation might be evidence of a binary supermassive black hole.

Figures

Tables

Table 1. Table 1 AIC values for the ARIMA models fitting the 10-day-bin light curve.

ARIMA(p,d,q)		MA(q)
		d=0						d=1
		MA(0)	MA(1)	MA(2)	MA(3)	MA(4)	MA(5)	MA(0)	MA(1)	MA(2)	MA(3)	MA(4)	MA(5)
AR(p)	AR(0)	1475	1388	1370	1368	1366	1368	1425	1391	1365	1363	1365	1367
	AR(1)	1359	1361	1359	1359	1361	1363	1410	1363	1364	1365	1366	1368
	AR(2)	1361	1363	1360	1360	1362	1364	1397	1364	1363	1366	1360	1369
	AR(3)	1362	1360	1355	1362	1364	1364	1387	1365	1366	1368	1370	1371
	AR(4)	1362	1363	1362	1364	1365	1366	1383	1367	1368	1370	1360	1371
	AR(5)	1363	1362	1364	1356	1360	1368	1373	1368	1368	1370	1371	1373

Full text

1. INTRODUCTION

Blazars are a subclass of radio-loud active galactic nuclei (AGNs) whose relativistic jets almost point at observers (Urry & Padovani, 1995). It is generally thought that a central supermassive black hole (SMBH) provides the energy that powers the relativistic jet through BH spin or rotating accretion disk. The emissions from blazar are dominated by the nonthermal emission from the relativistic jet, extending from MHz radio frequencies to TeV $\gamma$-rays energies, also exhibiting variabilities at all energies on a wide range of timescales. The typical mutiwavelength spectral energy distribution (SED) is distinguished by two broad peaks: a synchrotron component peaking at infrared to X-ray bands and a Compton component peaking in MeV to GeV energies.

The periodic variabilities of blazars have been extensively investigated in optical band (e.g., Bai et al., 1998, 1999; Fan & Lin, 2000; Xie et al., 2008; Li et al., 2009; Urry, 2011; King et al., 2013; Zhang et al., 2014; Bhatta et al., 2016; Fan et al., 2016). An interesting case is OJ 287 with a $\sim$ 12-yr period cycle (Kidger et al., 1992; Valtonen et al., 2006). Searches for $\gamma$-ray quasi-periodic oscillations (QPOs) became possible after the launching of the Fermi Gamma-ray Space Telescope in 2008. So far, the Large Area Telescope (LAT; Abdo et al., 2009; Atwood et al., 2009) onboard * Fermi* has collected gamma rays over 8 yr. The possible quasi-periodic variabilities of blazars with period of $\sim$2-3 yr have been reported in $\gamma$-ray fluxes of several blazars (e.g., Sandrinelli et al., 2014; Ackermann et al., 2015; Sandrinelli et al., 2016a, b, 2017; Zhang et al., 2017a, b). In particular, PKS 2155$-$304 have been reported having significance of $\sim$ 4$\sigma$ (Zhang et al., 2017a). The quasi-periodic modulations in the blazars carry important information on the BH - jet system.

We present the results of searching for QPO in the $\gamma$-ray light curves of blazar PKS 0301$-$243. A clear quasi-periodic signal with a period cycle of $\sim 2.1$-yr is found at the significance level of $\sim$5$\sigma$, though given that the full data set is only 8.78 years long this signal can easily have arisen randomly (e.g. Press, 1978). The paper is organized as follows: the detailed LAT data analysis and the main results are reported in Section 2. In Section 3 we summary the results and present a brief discussion on the findings.

2. Observations and Analysis

PKS 0301$-$243 is a high-synchrotron-peaked blazar (HSP) with its synchrotron peak frequency \approx 10^{15}\Hz (Abramowski et al., 2013), and its redshift is 0.266 (Pita et al., 2012). The High Energy Stereoscopic System (H.E.S.S.) has detected TeV photons from this source (Abramowski et al., 2013).

The events are collected between 2008 August 4 and 2017 May 19 (Modified Julian Date, MJD: 54,682.66 – 57,892.66) in the energy range from 100 MeV to 500 GeV, and in a square region of interest (ROI) of $20^{\circ}\times 20^{\circ}$ centered at the position of PKS 0301$-$243. The position of the target is located at right ascension (R.A.) = $\rm 03^{h}~{}03.442^{m}$, declination (decl.) = $\rm-24^{h}~{}07.192^{m}$ (J2000; $l=214.621,~{}b=-60.177$). The analysis is performed with the Fermi Science Tools version v10r0p5 package which is provided in the Fermi Science Support Center (FSSC).111https://fermi.gsfc.nasa.gov/ssc/data/analysis/software/ The Pass 8 LAT data (Atwood et al., 2013) are used with keeping only the SOURCE class photon-like events (with options evclass = 128 and evtype = 3 in the Tool $gtselect$). To minimize the contamination due to the gamma-rays bright Earth limb, we exclude the events with zenith angles $\geqslant~{}90^{\circ}$. By running Tool $gtmktime$, we obtain the good time intervals with high-quality photons. All the data reductions follow the data analysis thread provided by FSSC 222https://fermi.gsfc.nasa.gov/ssc/data/analysis/scitools/. We adopt the instrumental response function (IRF) ‘P8R2_SOURCE_V6’ in the analysis. Two diffuse model files333https://fermi.gsfc.nasa.gov/ssc/data/access/lat/BackgroundModels.html, namely gll_iem_v06.fit and iso_P8R2_SOURCE_V6_v06.txt, are used to model the Galactic and extragalactic diffuse $\gamma$-rays . A binned maximum likelihood is adopted to fit the events in the whole time range with the model file generated with the script make3FGLxml.py. This file contains the information on the spectral parameters of all known 3FGL sources (Acero et al., 2015) in the ROI. The $\gamma$-ray spectrum of the target is power-law in the Fermi 3FGL. The best-fitting results are derived with Fermi Tool gtlike, and are saved as a new model file. We also use the spectra in Fermi 3FGL model file to fit the events in the square ROI. The integrated photon flux of the best-fitting results above 100 MeV is $\rm F_{0.1-500~{}GeV}=(4.2\pm 0.1)\times 10^{-8}~{}photons\ cm^{-2}\ s^{-1}$, and the index of power law is $1.90\pm 0.01$ with the TS value of 9391.7 (the results in this paper with statistical errors only). We construct the light curves based on this new model file.

2.1. $\gamma$-ray light-curve

We use the two techniques, the maximum likelihood optimization (ML) and the exposure-weighted aperture photometry (AP; Corbet et al., 2007; Kerr, 2011), to construct the $\gamma$-ray light curves of PKS 0301$-$243. The 30-day-bin ML light curve is generated by employing the unbinned maximum likelihood fitting technique. In this step, the Tool gtlike is employed for each time-bin, and the events are selected in a circle ROI of $15^{\circ}$ centered at the coordinates of the target. We use the same parameter value as that in the new model file for all the sources in the ROI and freeze the spectral parameters except for the target. The 30-day-bin ML light curve is shown in the left upper panel of Fig. 1. For testing whether the power peaks vary with different length of time-bin, we also produce the ML light curve with the 10-day-bin, which is shown in the left upper panel of Fig. 2. The light curve can also be produced by the method of the exposure-weighted aperture photometry. In this method, we calculate the probabilities for each photon with the Fermi-Tool gtsrcprob, and then sum the probabilities of each photon within $1^{\circ}$ radius centered on the position of target for each 2.5-day-bin, in which the counts are weighted by its relative exposure for each time-bin. The AP light curve is shown in the left upper panel of Fig. 3.

We notice that there is an isolated large flare around MJD 55320. To avoid its impact on searching for quasi-periodic variability, we remove this flare in the following quasi-periodicity analyses.

2.2. Searching for quasi-periodic variability

2.2.1 Analyses on $\gamma$-ray data

We apply the two widely used methods, Lomb-Scargle Periodogram (LSP; Lomb, 1976; Scargle, 1982) and Weighted Wavelet Z-transform (WWZ; Foster, 1996), to the $\gamma$-ray light curves. For the 30-day-bin ML light-curve, three power spectra, LSP power, WWZ power and time-averaged WWZ power, are shown in Fig. 1. A strong peak near a period cycle of $2.1\pm 0.3$ yr appears, in which the maximum power is $>$ 18.6 times of the mean power value. The probability (Prob) for obtaining a power larger than the maximum power from the noise is $<1.58\times 10^{-9}$ (corresponding to a $>6.0\sigma$ significance level). The Prob(P $>$ Pn) is assessed through the formula: $\rm Prob(P>Pn)=(1-2\times\frac{Pn}{N-1})^{(N-3)/2}$ with the normalization from Horne & Baliunas (1986), where N=105 is the number of time-bin in the month-bin light curve. We correct the probability in the range of $1/3000\rm\ day^{-1}$ - $1/60\rm\ day^{-1}$ with the “trial factor = 50” (the number of sampled independent frequencies) (Zechmeister & Kürster, 2009), and find that the false-alarm probability (FAP) is less than $7.8\times 10^{-8}$, corresponding to $>5.4\sigma$. By fitting the power peak with Gaussian-function, we derive the period cycle of $763.3\pm 114.9$ days. The uncertainty of the period is evaluated based on the half width at half maximum (HWHM) of the Gaussian fitting.

In order to evaluate the precise significance of the signal, we use the method in Emmanoulopoulos et al. (2013) (also see Ackermann et al., 2015; Bhatta et al., 2016) to simulate light curves $3\times 10^{6}$ times based on the obtained best-fitting result of power spectral density (with the form of $P(f)\sim 1/f^{\alpha}+c$, where $c$ represents the Poisson noise level) and the probability density function of observed variation. We then derive the significance curves of 5-$\sigma$ and 4-$\sigma$ based on the simulations, which are shown in the lower right panel of Fig. 1. The significance of the signal is $\simeq 5.4\ \sigma$. We also calculate the power spectra of the 10-day ML light curve and 2.5-day-bin AP light-curve, which are shown in Fig. 2 and Fig. 3, respectively. In these two power spectra, we also find strong signals at $\sim$2.1 yr.

In order to further check the reliability of the quasi-periodic signal, we fit $\gamma$-ray light curve with autoregressive integrated moving average (ARIMA) models (Box & Jenkins, 1976; Hamilton, 1994; Chatfield, 2003) to assess whether the signal is consistent with a stochastic origin of autoregressive noise. We use the Akaike Information Criterion (AIC; Akaike, 1973) to select the best-fit model. In Table 1, we show the AIC values for 72 ARIMA models fitting the 10-day-bin $\gamma$-ray light curve. One can see that the ARIMA (3,0,2) model, with the minimum AIC value of 1355444We note that the AIC values of several models [e.g., ARIMA(1,0,0), ARIMA(1,0,2), and ARIMA(1,0,3)] are very close to 1355. This indicates that more data are needed to clearly distinguish these models., is the best-fit one. In Fig. 4, we show the standard residuals and the auto-correlation function (ACF) of the residuals for the best-fit model. It can be seen that there is a spike at the lag of 660 days that exceeds the 95% confidence limit. This marginal evidence indicates that the $\gamma$-ray quasi-periodic variability may not be produced by such type of stochastic processes.

We fold the events within a square region of interest (ROI) of $20^{\circ}\times 20^{\circ}$ centered at the position of PKS 0301$-$243 into 15 uniform bins based on orbital phase with the phase zero corresponding to MJD 54,682.66. We then fit the data in each phase bin by using the above best-fitting model-file to obtain the phase-resolved likelihood results. In Fig. 5, one can see that this folded light curve varies with the phase, indicating substantial variability in the source brightness (see the upper panel of Fig. 5); but no variability appears in its spectral shape (see the lower panel of Fig. 5).

2.2.2 Analyses on optical and X-ray data

We also search for quasi-periodic signal in the optical and X-ray data from this source. The long-term optical data from the Catalina Sky Surveys covering from 2005 October to 2013 October and daily averaged X-ray data from Swift-BAT covering from 2005 February to 2017 January are shown in the upper panels of Fig. 6 and Fig. 7, respectively. The LSP powers of the optical data and X-ray data are respectively shown in the lower panels of Fig. 6 and Fig. 7. No obvious peak is found in the corresponding powers. It is noted that the X-ray data are weakly variable.

3. SUMMARY AND DISCUSSION

Possible $\gamma$-ray QPOs have been reported in several blazars (e.g., Sandrinelli et al., 2014; Ackermann et al., 2015; Sandrinelli et al., 2016a, b, 2017; Zhang et al., 2017a, b). However, the significance of the claimed QPOs is not very high. In this paper, we report the first detection of $\gamma$-ray quasi-periodic modulation at a nominal confidence level of $\sim 5\sigma$ in PKS 0301$-$243. No quasi-periodic modulation is found in its optical and X-ray data.

In PG 1553$+$113, the quasi-periodic variabilities in optical and $\gamma$-ray data have the same period cycle (Ackermann et al., 2015). In PKS 2155$-$304, the periods of optical and $\gamma$-ray quasi-periodic variabilities are different (Sandrinelli et al., 2014). In PKS 0426$-$380, no optical quasi-periodic variability is found (Zhang et al., 2017a). The lack of optical and X-ray quasi-periodic variabilities may be because of the optical and X-ray originating from the different region that does not contribute $\gamma$-rays. If the lack of optical and X-ray quasi-periodic variabilities is confirmed by futuer long-term monitoring, it would challenge the most popular one-zone blazar emission model in which optical, X-ray and $\gamma$-ray emissions are assumed to be produced in the same region (Zhang et al., 2017a).

The mechanism causing the $\gamma$-ray quasi-periodic modulation in blazars is poorly understood. Given that the $\gamma$-rays are produced in the jet, two possibilities may account for the $\gamma$-ray quasi-periodic variabilities in blazars (e.g., Ackermann et al., 2015): (i) pulsational accretion flow instabilities may induce a quasi-periodic injection of plasma into the jet, hence a quasi-periodic modulation appears in the $\gamma$-ray flux from the jet; and (ii) the Doppler magnification factor changes periodically caused by jet precession/rotation.

Note that in our case the $\gamma$-ray photon index does not vary with the phase (see the lower panel of Fig. 5). The gamma-ray photon index is mainly determined by the high-energy electrons distribution. For HSP, the electron cooling is inefficient (e.g., Ghisellini & Tavecchio, 2008; Yan et al., 2014), and the electron distribution is mainly governed by the acceleration mechanism in the jet. This result indicates that the process yielding the QPO in the $\gamma$-ray flux would not have an impact on the acceleration process. The first possible origin for the QPO outlined above would have an impact on the energy outflow efficiency which is relative to the acceleration process in the jet (e.g., Ackermann et al., 2015). Therefore, our results may prefer to the second origin, i.e., jet precession. The jet precession could be the result of a helical jet (e.g, Rieger, 2004; Komossa & Zensus, 2016). Furthermore, a binary SMBH system would be involved in the formation of a helical jet (e.g, Komossa & Zensus, 2016; Sobacchi et al., 2017). Within such a scenario, the observed 2.1 yr period is the orbital time, and the equivalent intrinsic orbital time $P_{\rm int}$=$P_{\rm obs}/(1+z)$. The central SMBH of PKS 0301$-$243 is $\sim 8\times 10^{8}\ M_{\odot}$ (Ghisellini et al., 2010). Assuming the total mass of the binary SMBH of $10^{9}\ M_{\odot}$, the binary system size would be \sim 0.006\pc. At this stage, gravitational wave emission may be non-negligible in carrying away the energy.

In the jet precession model, the issue of the lack of optical and X-ray quasi-periodic variabilities could be resolved if the optical and X-ray radiations originate from a large region where the Doppler boosting is weak. Systematic sample study on QPOs at different electromagnetic frequencies in blazars may reveal deep physics of the jet (e.g., Sandrinelli et al., 2016a).

The $\gamma$-ray QPO in PKS 0301$-$243 is the first detection of such a kind of signal in blazars at a confidence level of $\sim 5\sigma$. Since there were barely four nominal quasi-periods in the currently Fermi-LAT data, this result certainly requires confirmation. Fortunately, our claim for a QPO should be tested rather soon, as the next flux maximum would be expected in 2018.

Finally, we would like to stress these claimed $\gamma$-ray QPOs in blazars are different from the X-ray QPO in BH X-ray binaries and narrow-line Seyfert 1 galaxy (Zhang et al., 2017b). For the X-ray QPOs, there is an inverse linear relation between QPO frequency and BH mass (e.g., Abramowicz et al., 2004; Török, 2005; Remillard & McClintock, 2006; Pan et al., 2016). This relation spans from stellar-mass to SMBH. No such relation is found in $\gamma$-ray QPO in blazars (Fig. 8). It seems that the intrinsic period of $\gamma$-ray QPO in blazars is independent on the SMBH mass. Moreover, the relation of the $\gamma$-ray QPO frequency-BH mass significantly deviates from the inverse relation found in the X-ray QPOs. The X-ray and $\gamma$-ray QPOs provide us different insights into the BH - jet system.

We thank the anonymous referees for useful and constructivecomments. Part of this work is based on archival data, software or online services provided by the ASI Science Data Center (ASDC). We acknowledge the financial support from the 973 Program of China under grant 2013CB837000, the National Natural Science Foundation of China (NSFC-11433004, NSFC-11525313, NSFC-11573060, NSFC-11573026 and NSFC-11661161010), and the Key Laboratory of Astroparticle Physics of Yunnan Province (No. 2016DG006). D.-H. Yan is very grateful to Prof. Yefei Yuan (USTC) for helpful discussions. The work of D.-H. Yan is supported by the CAS “Light of West China” Program.