On the classification of flaring states of blazars
Key Words.:X-ray: observations, galaxies: quasars
Aims:The time evolution of the electromagnetic emission from blazars, in particular high frequency peaked sources (HBLs), displays irregular activity not yet understood. In this work we report a methodology capable of characterizing the time behavior of these variable objects.
Methods:The Maximum Likelihood Blocks (MLBs) is a model-independent estimator which sub-divides the light curve into time blocks, whose length and amplitude are compatible with states of constant emission rate of the observed source. The MLBs yields the statistical significance in the rate variations and strongly suppresses the noise fluctuations in the light curves.
Results:We apply the MLBs for the first time on the long term X-ray light curves (RXTE/ASM) of Mkn 421, Mkn 501, 1ES 1959+650 and 1ES 2155-304, which consist of more than 10 years of observational data (1996-2007). Using the MLBs interpretation of RXTE/ASM data, the integrated time flux distribution is determined for each single source considered. We identify in these distributions the characteristic level as well as the flaring states of the blazars All the distributions show a significant component at negative flux values, most probably caused by an uncertainty in the background subtraction and by intrinsic fluctuations of RXTE/ASM. This effect interests in particular short time observations. In order to quantify the probability that the intrinsic fluctuations give rise to a false identification of a flare, we study a population of very faint sources and their integrated time flux distribution. We determine duty cycle or fraction of time a source spent in the flaring state of the source Mkn 421, Mkn 501, 1ES 1959+650 and 1ES 2155-304. Moreover, we study the random coincidences between flares and generic sporadic events such as high energy neutrinos or flares in other wavelengths.
Blazars are defined as Active Galactic Nuclei (AGNs) dominated by a highly variable component of non-thermal
radiation produced in relativistic jets pointed close to the line of sight
((Begelman et al.1984); (Urry et al.1995)).
One of their main characteristics is the flux variability on different time-scales:
from fast flares lasting few minutes to high states of several months.
Blazars are considered to be sites of energetic particle production and potential sources of
cosmic rays up to energies of at least eV.
The standard blazar spectral energy distribution (SED) shows two prevalent components: a hump at low-energy which peaks in the frequency range between infrared and X-ray bands, and a second hump at higher energy, proportionally shifted in the range from MeV up to TeV -rays. Two potential scenarios, the so-called leptonic and hadronic ones, have been proposed in order to model the SED. In leptonic models (e.g (Jones et al.1974); (Ghisellini et al.1996); (Mastichiadis et al.1997)), synchrotron emission from relativistic electrons is responsible for the first hump. Electrons in the jet plasma up-scatter low-energy photons to high energies via Inverse Compton, producing the second hump. In this scheme, the same electron population produces both components. In hadronic models (e.g (Mannheim and Biermann 1989); (Mannheim 1993); (Aharonian 2000)) protons are accelerated in the jet together with electrons. The synchrotron radiation produced by primary and proton-induced electrons contribute to the low-energy component. High-energy radiation originates from photo-meson interactions and from proton and muon synchrotron radiation. A comprehensive description of a Monte Carlo simulation of a stationary synchrotron proton blazar model, including all relevant emission processes, can be found in (Mücke and Prothereo 2000). In hadronic models, -ray production by pion photo-production results in simultaneous neutrino production. The decay of charged pions is the main neutrino production channel as discussed in ((Mücke et al.2003)).
The detection of very high energy neutrinos coming from blazars would be an unambiguous proof of the existence of baryonic loaded outflows and would indicate that blazars accelerate high energy cosmic rays. Neutrino telescopes (e.g. (Ahrens et al.2004); (Ackermann et al.2007); (Antipin et al.2007)) until now have not detected any extraterrestrial source of neutrinos in the TeV-PeV energy region. As discussed in (Rachen et al.1998), “the transience of energetic emission could improve the association of detected neutrinos with their putative sources, because one could use both arrival direction and arrival time information, allowing statistically significant statements even for total fluxes below the background level”. This is true under the assumption that neutrino production in HBLs is subjected to the same mechanisms at the base of the electromagnetic activity. Consequently, neutrino production and electromagnetic activity should show the same time modulation. The observation of time coincidences between electromagnetic flares and rare events, like neutrinos, represents a natural test to the hadronic scenario.
The main requirement to this approach is a clear definition and classification of the states of activity of the observed source. In this paper we discuss a procedure able to identify characteristic and flare states in a light curve. The estimator that best fits our requirements is the Maximum Likelihood Blocks (MLBs), since it is model-independent, it has been designed to identify blocks of data with a constant rate in variable periods, and it provides a statistical significance for each block. To test our approach, we perform a complete and detailed analysis on RXTE/ASM X-ray light curves. In particular, we analyze data from the brightest High energy peaked BLLacs (HBLs) ((Giommi and Padovani 1994)): Mkn 421, Mkn 501, 1ES 1959+650 and 1ES 2155-304. In the first part of this paper we describe the MLBs and how to separate flares from the characteristic level. Moreover, we introduce a definition for the duty cycle of the source. In the second part, we discuss the application of the method on RXTE/ASM data.
A variety of methods are used in astrophysics in order to assess the variability of a source and to qualify the character of the variability (periodic, correlated etc). It is not our intention here to review these methods. Each method is designed for a specific purpose. Often, data are affected by large uncertainties or the data spacing is rather inhomogeneous. The driving factors for the selection of a method are the goals of the analysis and the quality of the data to be analyzed. In our case we need a method that addresses the variability issue on light curves which are unevenly spaced and have short and long breaks, takes into consideration the statistical errors and possible unknown instrumental fluctuations on the measurements and gives a representation of the light curve in term of periods in which the data points are compatible with a constant level. A method that could satisfy these requirements is the Maximum Likelihood Blocks. The entire data analysis reported here is performed in ROOT ((Brun and Rademakers 1997)), an object-oriented data analysis framework. The only exception is for the Maximum Likelihood Blocks algorithm which is currently an IDL based program.
2.1 Representation of the light curve: Maximum Likelihood Blocks
The methods used in the study of temporal variability depend strongly
on the nature of the available data and of the signal of interest. In
all cases, the most basic step is the classification of the
time-series as “constant” or “variable”. Suitable and widely used
statistical tests include the Kolmogorov-Smirnov test for time-tagged
data (e.g. the arrival times of X- and -ray photons) and the
test for binned data (e.g. binned photon arrival times or
optical magnitudes). The next step in the analysis of light curves is
the characterization of their “shape”. We will use a simple and
model-independent approach that aims at dividing the light curve into
time intervals in which the source emission is compatible with a
constant level. An algorithm based on Bayesian statistics that
performs such a segmentation for data of different natures was
presented by (S98, (Scargle 1998)). In its form for
time-tagged data, this algorithm was used for example to characterize
the X-ray light curves of young pre-main sequence stars observed by
the Chandra Orion Ultra-Deep Project (COUP, (Getman et al.2005)) in
the Orion Nebula Cluster (ONC). A modified version of the S98
algorithm, based on a Maximum Likelihood rather than a Bayesian
approach, was recently employed in other studies of stellar X-ray
light curves ((Wolk et al.2005); (Favata et al.2005); (Stelzer et al.2006)). We will refer to this
algorithm as the Maximum Likelihood Blocks (MLBs) and we introduce
here a variant that is suitable for the analysis of binned light
Our algorithm is derived from the one presented by S98. Like S98, we tackle the problem of finding the best piecewise representation of a binned light curve in an iterative (and approximate) way: we begin by testing the data against a constant-flux model. If the model does not represent the data adequately we split the light curve into two parts at the most likely “change point”. We then repeat these two steps recursively on the resulting segments until all segments are compatible with constant emission. The fundamental difference with the algorithms presented by S98 lies in the statistics used to test if a light curve is variable and to find the most likely change point: rather than “marginal likelihoods” and “Bayes factors” (e.g. Eq. 7 and 48 in S98) we employ simple Likelihood functions, i.e. the probability densities of obtaining the observed data set given a parametric model.
As mentioned above the algorithm was first applied to time-tagged data. It is here adapted to the different statistical properties of binned time series. Our lightcurves can be described as a series of independent flux measurements, , each normally distributed about their mean values with standard deviations . The likelihood of a parametrized model, , of the lightcurve is maximized by minimizing the ), where are the model-predicted fluxes. In our case the model is either the single-segment representation or one of the possible two-segment representations of the light curve. We will refer to these models, respectively, as “1”, specified by one parameter, the constant flux level, and “2(j)” specified by three parameters: the change point “(j)” (more specifically the index of the last point in the first of the two segments) and the two flux levels. In this notation our algorithm reduces to: ) splitting the light curve if the minimum is such that the probability of obtaining a larger value is lower than a user defined confidence threshold (e.g. 1% or 0.1%); ) choosing the change point, , as the one that minimizes
2.2 Interpretation of the light curve: flares versus characteristic level
The goals of this analysis are the identification of the various levels of activity of a source and the separation between
bursting events (flares) and steady state period(s) (characteristic level(s)).
Sometimes periods of no variable activity are defined in the literature as “quiescent”.
As discussed for example in (Wolk et al.2005),
the meaning of quiescent emission is ambiguous. An apparently quiescent level can be due to
a superposition of numerous unresolved flares. Quiescent, as defined as inactive,
is therefore not appropriate to describe the level of activity in which the source spends most of the time.
We define the characteristic level as and the spread around it .
In order to determine the value of
we construct the distribution of the amplitude
and the duration of the single block . We call this
integrated time(T)-flux(r) distribution based on the MLBs interpretation (B): .
This provides the distribution of the total amount of time the source passes in a particular activity state.
The threshold above which a flux state
is defined as flare is then defined as .
Depending on , the probability that a selected flare state is
caused by a fluctuation of the characteristic level, by
an instrumental fluctuation or by a real enhancement of the photon emission from the source
can be fully assessed.
On the base of this definition of flares we can determine as well the frequency of flare states or duty cycle as:
In Sect. 4, the application of this method to RXTE/ASM data are reported.
The All-Sky Monitor (ASM) on board of the Rossi X-ray Timing Explorer
(RXTE) has been monitoring the X-ray sky routinely since March 1996.
During each orbit up to 80% of the sky is surveyed to a depth of 20-100 mcrab.
A source is observed roughly 10 times a day.
A set of linear least-square fits over 90 seconds observation periods,
by one of the three Scanning Shadow Cameras, yields
the source intensities in four energy bands (1.5-3, 3-5, 5-12, and 1.5-12 keV).
The intensities are usually given in units of the count rates expected if the sources were at the center
of the field of view in one of the cameras. In 1.5-12 keV band, the Crab Nebula
flux corresponds to about 75 ASM counts per second.
A detailed description of ASM can be found in (Levine et al.1996). RXTE standard data products
are collected directly from the HEASARC database.
We concentrate our study on RXTE/ASM data because this provides the longest light curves in X-ray of Mkn 421, Mkn 501, 1ES 1959+650 and 1ES 2155-304. However, for these kind of sources, the RXTE/ASM sensitivity is limited and data are affected by large errors. Moreover, the resolution and the background level of ASM observations depend on the Sun contamination or back-scattered solar X-rays and on the detector stability along the 10 years of data taking, (Wen et al.2006).
The results of the application of the MLBs to the RXTE/ASM data for the four HBLs considered are reported in Fig. 1 and in Fig. 2. Each change point identified by the algorithm has a statistical significance of at least . This operation is still not enough in order to characterize the behavior of a source. We haven’t yet quantified the characteristic noise of RXTE/ASM, then we can not distinguish if the change points identified by the MLBs are of instrumental nature or of a physics nature. In order to distinguish between these two scenarios and study the behavior of sources, we construct the integrated time(T)-flux(r) distribution , as discussed in Sect. 2.2. We first construct the for a set of very faint sources in order to study the intrinsic fluctuations of the instrument and finally we determine the for Mkn 421, Mkn 501, 1ES 1959+650 and 1ES 2155-304. On the base of these distributions we identify the states that can be considered flares with good confidence.
4.1 RXTE/ASM intrinsic fluctuations
In Fig. 1, we notice that the MLBs identify not only significant change points at positive amplitudes but also at negative ones. These negative fluctuations can be caused by uncertainties in the background subtraction and by intrinsic fluctuations of RXTE/ASM. In order to characterize such a component and its effect on the definition of flares, we have analyzed RXTE/ASM light curves for a set of very faint sources, since these are expected to spend most of their time at a flux level well below the ASM sensitivity. The sources considered, reported in Tab. 1, have a low X-ray monochromatic average emission (less then 0.6 Jy at 1keV), are randomly distributed in the sky and are at various redshifts. The flux distributions of the faint sources are all normal distributions, as expected for a random instrumental noise. On average, the normal distributions peak at rate ASM c/s and have a standard deviation of ASM c/s. Since the studied faint sources show similar flux distributions, we will use in the next just one of them for comparison; the source used is PKS 0118-272 and represent the average faint source in our sample. In Fig. 3, the flux distribution of PKS 0118-272 is compared with the flux distribution of the HBLs considered in this work. The distributions are normalized using the areas under the negative flux tails. In this way, we estimate the fraction of the HBLs flux distributions caused by the intrinsic fluctuations of RXTE/ASM ().
4.2 RXTE/ASM flare states
The for the Mkn 421, Mkn 501, 1ES 1959+650
and 1ES 2155-304 are reported in Fig. 3.
All the distributions differ significantly from a normal distribution
indicating that RXTE/ASM is indeed sensitive to high activity states of the sources considered.
For all the HBLs considered we observe that the flux distribution
shows a peak above the pure background distribution.
We define the central peak value (that corresponds also to the mode of the flux distribution)
and its standard deviation as and .
More sophisticated fitting procedures have been applied but, given the quality of the data of RXTE/ASM,
they did not provide a more precise estimation of and
. We observe that is next to the detection threshold of RXTE/ASM.
As discussed in Sect. 2.2,
is the threshold above which a flux state
is considered a flare.
Using this definition of flares we determine the frequency of flare states or duty cycle as described in Sect. 2.2. In Tab. 2 we report for the HBLs considered and the cases of N=1 and N=3. For the specific case of RXTE/ASM, we calculate as well the intrinsic fluctuation that affects the duty cycle as:
where is the flux distribution of a faint source (in this work PKS 0118-272). Results are reported in Tab. 2.
4.3 Example of application: correlation study between electromagnetic flares and neutrinos
As anticipated in the introduction, the study of the physics of HBLs develops through different approaches. One of the most frequently used sees the study of flare correlation among different wavelengths, for example X-ray and TeV- rays (Maraschi et al.1999). A study of the correlation among different messengers such as photons and neutrinos is a more recent interest (Achterberg et al.2005) and has been the motivation of this work. The significance of such correlations can be assessed only when the frequency of the electromagnetic flare states is determined, for example following the procedure described in this paper. In order to illustrate such a case, we study the distribution of coincidences between RXTE/ASM flare states and a set of N neutrinos or equivalent sporadic events. The flares are selected following the procedure described above for and the N neutrino events are uniformly distributed in the entire time period considered in this paper, approximately 10 years. The distributions of coincidences between the RXTE/ASM flare periods for Mkn 421, Mkn 501 and 1ES 1959+650 and the neutrino events are reported in Fig. 4. Depending on the number N of sporadic events we are able to determine the number of random coincidences multiplying the probability of a flare, or duty cycle of the source, and N. Once the random coincidence among flares and neutrinos is determined, then the statistical meaning of experimentally observed coincidences between flares and neutrinos can be determined. This eventually can hint at the association of detected neutrinos with an astronomical source even if N is below or at the level of the expected background.
The X-ray time behavior of Mkn 421, 1ES 1959+650, Mkn 501 and 1ES 2155-304 has been
characterized using approximately 10 years of data from RXTE/ASM.
The characteristic level and flaring states have been defined and values are reported in
Mkn 421 is the source that flares most often amongst those studied: for of the RXTE/ASM observations the source was in an active state (above ) and for of these it was in a very high state (above ). The probability that a flare state is caused by a fluctuation of the instrumental noise is marginal. This confirms the well known fact that Mkn 421 is an extremely variable HBL and quantifies for the first time its duty cycle even if this is only valid for RXTE/ASM. Mkn 501 flares often at low rates () and less often at high rates (). Also in this case, the intrinsic fluctuations do not significantly affect the flare states. 1ES 1959+650 flares more rarely in particular at high rates (only 2.6%). The systematic component affects one tenth of flares. In the case of 1ES 2155-304 nearly one fourth of flare is caused by an intrinsic fluctuation of RXTE/ASM. This simply means that this source is at the threshold limit for RXTE/ASM.
The study of the significance of correlations of flares among different wavebands and different messengers is foreseen for a future work.
|Source Name||X-ray at 1keV (Jy)||Mean (ASM c/s)||Sigma (ASM c/s)||z|
|Source Name||(ASM c/s)||(ASM c/s)||(%)||(%)||(%)||(%)|
Acknowledgements.E.R. and A.G. are funded by the Deutsche Forschungsgemeinschaft (DFG) through an Emmy Noether grant (RE 2262/2-1). This research has made use of data obtained through the High Energy Astrophysics Science Archive Center Online Service, provided by the NASA/ Goddard Space Flight Center. We acknowledge the RXTE/ASM and RXTE/PCA team for the X-ray data. Thanks to A. Taylor, S. Movit and S. Odrowski for their proofreading and useful comments.
- The minimization with respect to the flux levels is trivial and reduces to choosing the mean flux in the given time interval.
- Ackermann, M. et al. (IceCube Collaboration) 2008, AJ, 675, 1014
- Achterberg, A. et al. (IceCube Collaboration) 2005, Contribution to 29th International Cosmic Ray Conference (ICRC 2005). e-Print: astro-ph/0509330
- Aharonian, F.A. 2000, New A, 5, 377
- Ahrens, J. et al. (IceCube Collaboration) 2004, Astropart. Phys. 20, 507
- Antipin, K. et al. (Baikal Collaboration) 2007, arXiv:0710.3064v1
- Begelman, M.C., Blandford, R.D., & Rees, M.J. 1984, Rev. Mod. Phys., 56, 255
- Blazejowski, M., et al. 2005, ApJ, 630, 130
- Brun, R & Rademakers, F. 1997, NIMA, 389, 81
- Favata, F., Flaccomio, & E., Reale, et al. 2005, ApJS, 160, 469
- Getman, K. V., et al.2005, ApJS, 160, 319
- Ghisellini, G., Maraschi, L., & Dondi, L. 1996, ApJS, 120, 153
- Giommi, P. & Padovani, P. 1994, ApJ, 444, 567
- Jones, T.W., O’Dell, S.L., & Stein, W.A. 1974, ApJ, 188, 353
- Krawczynski, H., Coppi, P.S., Maccarone, T., et al. 2000, A&A, 353, 97
- Krawczynski, H., et al. 2002, MNRAS, 336, 721
- Krawczynski, H., et al. 2004, ApJ, 601, 151
- Levine, A.M., et al. 1996 ApJ, 469, 33
- Mannheim, K., & Biermann, P.L. 1989, A&A, 221, 211
- Mannheim, K. 1993, A&A, 269, 67
- Maraschi, L., et al. 1999, Astrop. Phys., 11, 189
- Mastichiadis, A., & Kirk, J. G. 1997, A&A, 320, 19
- Mücke, A., Protheroe, R.J. 2001, Astropart. Phys, 15, 121
- Mücke, A., Protheroe, R. J., et al. 2003, Astropart. Phys., 18, 593
- Rachen, J.P., et al. 1998, Phys. Rev. D, 58, 123005
- Scargle, J. D. 1998, ApJ, 504, 405 (S98)
- Stelzer, B., et al. 2006, A&A, 448, 293
- Urry, C.M., & Padovani, P. 1995, PASP, 107, 803
- Wen, L., et al. 2006 ApJS, 163, 372
- Wolk, S. J., Harnden, F. R., Flaccomio, et al. 2005, ApJS, 160, 423