Periodicity in the continua and broad line curves of a quasar E1821+643

1.1 MB·

by  Andjelka Kovacevic

Periodicity in the continua and broad line curves of a
quasar E1821+643
A. Kovaˇcevi´c1 • L. Cˇ Popovi´c1,2 •
A. I. Shapovalova 3 • D. Ili´c 1
7
1
0
2

n
a
J

6

] Abstract Here we present an in-depth analysis of the ing an insight into their structure even on the small-
A
periodicity of the continua and broad emission lines of est scales. Among the most interesting related issues
G a quasarE1821+643. We appliednon-parametriccom- is the quest for periodicity, because it could be a sig-
. posite models, the linear sum of stationary and non- nature of supermassive black hole binary (SMBHB,
h
p stationary Gaussian processes, and quantified contri- Lobanov & Roland2002). PeriodicityofAGNhasbeen
- bution of their periodic parts. We found important inthespotlightsincetheworkofSillanp˝aa˝ et al.(1988),
o
qualitativedifferencesamongthethreeperiodicsignals. where was noticed a regularity in the historical optical
r
t Periodsof∼2200and∼4500daysappearinbothcon- lightcurve ofthe blazarOJ 287. The largeflare events
s
a tinua5100˚A and4200˚A,aswellasintheHγ emission occurred12yearsinthelightcurveoftheblazar,during
[ line. Their integer ratio is nearly harmonic ∼ 12, sug- the period of almost one century of observations.
1 gestingthesamephysicalorigin. Wediscussthenature Up to a few years ago just some AGNs were con-
v of these periods, proposing that the system of two ob- sidered as long-term (or quasi1) periodic objects in
6 jects in dynamically interaction can be origin of two
optical and UV continuum (Hagen-Thorn et al. 2002;
6
largest periods.
5 Fan et al. 2002; Lainela et al. 1999; Valtonen et al.
1 2008; Graham et al. 2015a). Graham et al. (2015b)
Keywords (Galaxies:) quasars: supermassive black
0 made systematic search for the long-term variations
. holes; Methods: data analysis
1 of continuum in large sample of quasars and found
0 more than one hundred candidates. Thorough anal-
7
1 Introduction ysis of these results (see Vaughan et al. 2016) warned
1
: communityaboutcriticalinterpretationofderivedperi-
v Investigationofactivegalacticnuclei(AGN)variability odicities. Namely, AGN periodical flux variability may
i
X overpowers any other observational technique in giv- be explained by multiplex and violent physical mech-
r anisms. We list, as an examples, several possibilities
a
A.Kovaˇcevi´c (for the first three see Ackermann et al. 2015, and
references therein): pulsational accretion flow instabil-
Department of Astronomy, Faculty of Mathematics, University
ofBelgrade,Studentski trg16,11000Belgrade,Serbia, ities; jet precession, jet rotation, helical jet structure;
[email protected],tel[phone], fax[phone] anaccretion-outflowcouplingmechanism;orbitingdisk
L.Cˇ Popovi´c hotspots(Guilbert et al.1983). Also,thelasttwomen-
Astronomical Observatory Belgrade, Volgina7, 11060 Belgrade, tioned processes can be presented in SMBHB systems
Serbia,[email protected],
too, distorting and hiding periodic signals. Accord-
Department of Astronomy, Faculty of Mathematics, University
ingly, it becomes very important to be able to reveal
ofBelgrade,Serbia
the periodic signals in AGN time series from the noise
A.I.Shapovalova
that is hiding them.
Special Astrophysical Observatory of the Russian AS, Nizhnij
Arkhyz,Karachaevo-Cherkesia369167, Russia,[email protected]
D.Ili´c 1Some models of multiple black hole systems (e.g.
Department of astronomy, Faculty of mathematics, Univer- Sundeliusetal. (1997)) have shown that perturbations in
sity of Belgrade, Studentski trg 16, 11000 Belgrade, Serbia, suchsystemsmaynothavestrictlyperiodiccharacter,sowecan
[email protected] thinkofthem asquasi(pseudo)periodic.
2
UptonowinvestigationsofAGNledtoonlyfewnon- sparse and irregular sampled time series. GP are suit-
controversial explanations of SMBHB (e. g. Popovi´c abletoAGNlightcurves,sinceingeneraltheyarecom-
2012; Bon et al. 2012). Such pseudo-periodic signals plex and with no simple parametric form available. In
canbeusedinordertopredicttheepochsoffutureflux such situation our prior knowledge about right model
flareswithinknownconfidencelimits. Thiswouldallow can be limited to some general facts. For example, we
us to confirmor discardmultiple black-holecandidates mayjustknowthatourobservationsoriginatefromun-
andplanfuture observationsbasedonthe expectedac- derlying process which is discrete, non-smooth, that
tivity phase of the source. has variations over certain characteristic time scales
Among the binary black hole candidates is a quasar and/or has typical amplitudes. Surprisingly, it arises
E1821+643. Shapovalova et al. (2016), hereafter Pa- naturally to work with the infinite space of all func-
perI,presentedtheresultsofthefirstlong-term(1990- tions that have such general characteristics. As these
2014) optical spectroscopic monitoring of this object. functions arenotcharacterizedwith explicitsets ofpa-
Intheirstudy,applicationofthestandardLomb-Scar- rameters, these techniques are called non-parametric
gle method (LS, Lomb1976) suggesteda possible peri- modeling. Since the dominant tools for working with
odicity in the continuum and broad emission line light suchmodelisprobabilitytheory,theyarerecognizedas
curves, which may be caused by orbital motion. Here Bayesian non-parametric models, and particular mem-
we extended their periodicity analysis and investigate berofsuchmodelsareGPs. Moreover,AGNvariability
the signal/red noise distinction in the continua and amplitude is large and their luminosity is always posi-
emission lines time series of this object. We applied tive, allowing us to perform calculations using variable
composite non-parametric models, i.e. the linear com- l(t) = ln(Flux(t)), which can be assumed to be GP
bination of simple Gaussian processes, to the continua (Nipoti & Binney 2015).
at5100 ˚A,4200˚A,Hβ andHγ emissionlines. The aim The simplified version of GPs (namely continuous
of this paper is twofold: first, to investigate in depth time first-order autoregressive process (CAR(1))) has
periodic signals in the light curves of this object and been used to model AGN continuum light curves (see
second, to provide more clues on its broad line region Pancoastet al. 2014, and reference therein). These
(BLR) geometry and its possible SMBHB nature. worksusedOrnstein-Uhlenbeck(OU)covariancefunc-
A short description of used observations and insight tion:
into Gaussian processes are presented in the section
σˆ2
Data sets and Methods. In the section Results and K (∆t)= exp α0∆t
DiscussionwedisplayGaussianprocessesmodelsofour OU 2α0 −
lightcurves,theirinferredperiodsanddiscusstheirim-
where 1 is characteristic time scale, σˆ2 is variance of
plications on rough geometrical characteristics of BLR α0
driving noise, and ∆t is the time separating two obser-
of E1821+643.
vations (see Kelly et al. 2014). (Note here that this
kernel can be written also in the form
2 Data sets and Methods |x−x|
K (x,x)=σ2exp − ′ ,
OU ′ l
(cid:16) (cid:17)
We examine four diverse time series containing the in-
tegralfluxesofthebroadHβ,Hγ emissionlinesandthe wheretheparameterσ2isthevarianceandlisthechar-
fluxes for the continuum at the rest frame wavelength acteristiclengthscaleoftheprocess,i. e. howclosetwo
of 5100 ˚A and 4220 ˚A. As the temporal coverage and pointsxandx′ havetobe toinfluence eachother. The
sampling of data sets are important for an extensive OU process is a stationary GP. Since the information
search for a periodic behavior of an AGN, reader can abouttheunderlyingnatureofAGNlightcurvesiscru-
find this analysis in Paper I. cial for probing their variability, Kelly et al. (2011),
AGN time series are usually sparse and irregularly Andrae et al. (2013) and Zu et al. (2013) modeled the
sampled. Evenifeachtimeseriesisdefinedoveracom- light curves as a parameterizedstochastic process with
mon continuous time interval, they can be defined on compositecovariancefunctions. Followingthesamedi-
different collection of time points due to irregularity. rectionofinvestigations,we applied some ofstationary
Also, the distribution of intervals between time points andnon-stationarycompositecovariancefunctions. As
areusuallynotuniformandthenumberofobservations they are a linear combination of other simpler covari-
in different time series can be different. Learning any ancefunctions,weareabletoincorporateinsightabout
characteristic in such setting is difficult. periodicity, red noise and nonstationarity of underlay-
It has been shown that the Gaussian process (GP ing processes.
Rasmussen & Williams 2006) naturally conforms to
3
Inourresearchweemployedtwo kindofkernelfam- Many stochastic processes behave, at least for long
ilies: stationary - already mentioned OU precess, the ranges of time, like random walks with small but fre-
squared exponential (SE): quent jumps. The argument above implies that such
processes will look, at least approximately, on the
proper time scale, like Brownian motion. It is very
K (∆t)= σˆ2 exp (α0∆t)2, importantwhenprobingthe periodicityoflightcurves,
SE 2α −
0 totakeinto accountthe rednoise. Itisoftenpresented
which alternative representation is in the AGN continuum and broad emission lines. Such
noisecanproducespurioussignalsinapowerspectrum
K (x,x)=σ2exp − ||x−x′||2 , of periodogram. We point out that ’red noise’ vari-
SE ′ (cid:16) 2l2 (cid:17) ability of AGNs originates in physical processes of the
source and it is not result of measurement uncertain-
rational quadratic (RQ)
ties. The red noise simulations were produced by OU
σˆ2 α2|∆t|2 process having the variance of real light curves.
KRQ(∆t)= 2α (1+ 0α )−α, withα=2; Since all of these base kernels are positive semidefi-
0
nite kernels (i.e. the valid covariance functions can be
with an alternative form
defined) and closed under addition and multiplication,
|x−x|2 we can create richly structured and interpretable ker-
K (x,x)=σ2(1+ ′ ) αwithα=2;
RQ ′ αl2 − nels with their summation and/or multiplication. By
summing kernels, data can be modeled as a superpo-
non-stationary – the standard periodic kernel (StdPer)
sition of independent functions, representing different
structures. In the case of our time series models, sums
σˆ2 2sin2(π(∆t))
K (∆t)= exp − P of kernels can express superposition of different pro-
StdPer 2α0 (cid:16) (α10)2 (cid:17) cesses operating at different scales. Since we are ap-
plyingtheseoperationstothecovariancefunctions,the
with alternative version
composition of even a few base kernels can catch com-
2sin2(π(x−x′)) plexrelationshipsindatawhicharenotofsimplepara-
KStdPer(x,x′)=σ2exp − l2 P ; metric form. This enable us to search for the struc-
(cid:16) (cid:17)
tureinourdata,i.e. periodicity/rednoise,andcapture
and Brownian motion (Brw, red noise or Wiener pro-
them. We determine the number of probing and valid
cess)
periods by utilizing so called ’a two tier approach’ of
Vlachos et al. (2005). On the first tier the period can-
σˆ2
K (∆t)= min(∆t), didates are extracted from the periodogram analysis.
Brw 2α
0 On the second tier, these periods are verified by other
which alternative form is method (these authors used autocorrelation function).
In our case, on the first tier we extracted three can-
didate periods (see Paper I) by means of periodogram
KBrw(x,x′)=σ2min(x,x′). analysis. Here, on the second tier, we put these three
candidate periods into composite GP models. If these
The SE is infinitely differentiable allowing GP with
models provide better modeling result than their non
this covariance function to generate functions with no
periodic parts and if the optimized values of periods
sharp discontinuities. An alternative to the sum of SE
areclosetothe candidatevalues(withinafewhundred
kernels is RQ. It was shown that this kernel is equiva-
of days) then we can consider these values as a valid
lent to summing many SE kernelswith differentscales.
periods.
This produces variations with a range of time scales.
For our analysis, we used the following kernels:
Parameter α is roughness, it determines the relative
weighting of large and small scale variations. Under
condition α→∞, the RQ becomes identical to SE. In SE+ StdPer (1)
i
StdPer kernel hyper parameters P,l correspond to the iX=1,3
period and length scale of the periodic component of
RQ+ StdPer (2)
the variations. Inthis casel isrelatedtoP soithasno i
iX=1,3
dimension.
We applied Brownian motion because it is, in a cer- Brw+ StdPeri (3)
tain sense, a limit of rescaled simple random walks. iX=1,3
4
0 4.0 0 4.0 Forselectingamodelfromasetofcandidatemodels,
3.8 3.8
which is best supported by our data, the Bayes factor
20 3.6 20 3.6
3.4 3.4 (BF)(Jeffreys1961;Kass & Raftery1995)isused. The
40 3.2 40 3.2 BF of two models Model ,Model is just ratio of their
′X 3.0 ′X 3.0 i b
60 2.8 60 2.8 marginallikelihoodsML(Modeli),ML(Modelb)derived
2.6 2.6 using the same dataset, and it includes all informa-
80 80
2.4 2.4
tion about a particular Bayesian model. BF is pre-
2.2 2.2
0 20 40 60 80 0 20 40 60 80 ferredfor model selection criteria because they contain
X X
the information about a model prior. For interpreting
0 8.8 the values of BFs we used half-units of the logarith-
20 8.0 mic metric (LBF = log10(BF) = log10ML(Modeli)−
7.2
log ML(Model )), since we have already calculated
40 6.4 10 b
the marginal log likelihood of our models. The in-
′X 5.6
60 4.8 terpretation of LBF is that when 0 < LBF < 0.5,
4.0
80 3.2 there is an evidence against Modeli when compared
with Model , but it is only worth a bare mention.
2.4 b
0 20 40 60 80 When 0.5 ≤ LBF < 1 the evidence against Model
X i
is definite but not strong. For 1 ≤ LBF ≤ 2 the ev-
Fig.1 Thepriorcovariancematricesassociatedwithcom- idence is strong and for LBF > 2 it is decisive. In
posite kernelsgiven byequations1 (top left), 2(top right), our case we will compare models to the base model -
and 3 (bottom) are calculated from 100 artificial equidis-
red noise Model =M which is generated from OU
b OU
tant points. The lightest shade indicates the highest and
kernel. Results in this work have been calculated with
thedarkest shade lowest covariance.
1.0.7 versionof the Pythonic Gaussian process toolbox
GPy (http://github.com/SheffieldML/GPy).
whereStdPer,i=1,2,3areperiodickernelsforsearch-
i
ing for specific periods Per,i=1,2,3.
i
The structure of above defined kernels’ covariances
3 Results and Discussion
can be visualized by ’heatmap’ of the values in the co-
variancematrix. Fordemonstrationpurposeonly,arbi-
As our GP models have a number of hyperparameters
trary one dimensional input space has been discretized
(the covariance functions are composites, see Table 1),
equally, and the covariance matrices of Eq. 1, 2 and
we first assigned a prior values to some of hyper pa-
3 are constructed from this ordered list (see Figure 1).
rameters, obtained from our analysis of light curves in
The region of high SE and RQ covariance is depicted
PaperI.Inassigninginformativepriorstothreeperiods
asadiagonalwhiteband,reflectingthelocalstationary
inourGPmodels,weusedpriorsof4500,2000and1400
nature of these two kernels. Increasing the lengthscale
days (from the LS periodogram analysis), constraining
l increasesthe width ofthe diagonalband, since points
that the most important periods were on the order of
at larger distance from each other become correlated.
thousand days, rather than on larger or smaller scales
Since we summed periodic and non-periodic kernels
(e. g. seconds or millennia). So the range of priors to
it is important to quantify the periodicity of inferred
the periods were [0,1400],[0,2000]and [0,4500] days.
composite models. For this purpose, we applied the
From the heteroskedastic OU process and CARMA
method of Durrande et al. (2016) for the first time in
model (Kelly et al. 2014) random light curves were
astronomy. According to their prescription, let y and
p
sampledtoairregulartimeintervalsofreallightcurves,
y be the periodic and non-periodic components of
a
in order to to test if periodicity found by LS method
our model calculated onrandom variable of input data
is the product of random variations. In such a way,
points. To quantify the periodicity signal, we calculate
we also compare their effects as a base red noise mod-
periodicity ratio as:
els. The reason for trying this two types of models
is a radio-loud, X-ray nature of E1821+643, due to
Var(yp) which its light curves may not be well described by
S= (4)
Var(yp+ya) OU models. Moreover, Kozlowski (2016) has shown
that OU models can be a degenerate descriptors. A
Furthermore,Scannotbeviewedasapercentageofthe
CARMA models might serve better for this purpose,
periodicityofthesignalrigorouslysinceVar(y +y )6=
p a particularlyastherednoisetermsinthemcanproduce
Var(y )+Var(y ). Consequently, S can have values
p a quasi-periodicities which are missing from simple OU
greater than 1.
5
models. For the purposeofCARMA modeling we used 0
rednoiseOU
BrandonKellycarma packPythonpackageavailableon −50 Fcnt4200
rednoiseCARMA(4,1)
https://github.com/brandonckelly/carma pack. CARMA
−100
models lightcurveassumof(deterministic)autoregres-
sion plus (random) stochastic noise. The CARMA −150
modelorderinputisoptional. Weautomaticallychoose −200
the CARMA order (parameters (p,q)) by minimizing 0
rednoiseOU
the AIC (Aikake Information Criterion). Also, the %)) −50 Hγ
( rednoiseCARMA(5,0)
residualsfromtheone-step-aheadpredictionsstandard- y
ized by their standard deviation were uncorrelated. bilit −100
a
b
Random light curves were sampled with built in func- Pro −150
tions in GPy and farm pack packages. g(
o −200
We recalculated periods of each real and simulated L
0
curve(seeFigure2)bymeansofBayesianformalismfor rednoiseOU
−50 Fcnt5100
the generalized LS periodogram (BGLS Mortier et al. rednoiseCARMA(4,1)
−100
2015), which gives the probability of signal’s existence
−150
in the data. It can be seen that powers of BGLS peaks
−200
of OU and CARMA ’red noise’ curves are asymptoti-
−250
callyclosetoeachotheraswellastothe continua4200
−300
˚A 5100 ˚A and the Hγ line (Figure 2). 2.6 2.8 3.0 Log(P3e.2riods) 3.4 3.6 3.8
Due to this fact we will consider as base red noise
model OU. Also, in our calculations will be included
models given by Eq. 3 and Fig. 2 Bayesianperiodograms ofthesimulated’rednoise’
and real light curves on a logscale.
PM =K+ StdPer,K=OU,Brw.
K i To address the heteroskedasticity of data errors we
iX=1,3
have applied homoskedastic GP regression on loga-
Models PM ,PM consist of a set of signals rithmictransformeddependentvariableheteroscedastic
OU Brw
(StdPer) superimposed on the red noise (Brw,OU). GPy regression on non-transformed data. Both meth-
i
These models allow us to disentangle signals in the ods provide almost the same results, which is clearly
light curves from the ’red noise’ background in which seeninthecaseofthreeperiodiccompositemodels(Ta-
they are immersed. ble 1 and Table 2).
We initially set default hyperparameters of variance The one periodic component models (Table 3) are
to the variance of light curves and lengthscale between not favored since S are large while LBF are small and
0.5 and100 days. The variancedetermines the average vice versafor heteroscedasticregression(Shet aresmall
distance of our function away from its mean. Initially andLBFhet islarge). Despitetwoperiodiccomponents
weexpecttohavethesamevarianceasourlightcurves. models (Table 4) have large S and Shet they are not
We chose range of lengthscales between the minumum favoredsincetheir meanvaluesofLBFandLBFhet are
ofsamplingtimesofobservations(∼0.5days)andtwice smaller than corresponding values in Table 1 and Ta-
of mean of sampling times (∼ 100 days). The length- ble2respectivelly. Sothreeperiodiccomponentmodels
scaledeterminesthe lengthofthe ’wiggles’inourfunc- inferredfromhomoskedastic(Table1)andheterskedas-
tion, sowith largerlengthscale than100daysthe func- ticregressions(Table2)aredominant,providingalmost
tion would be smoother. the same results.
For each light curve, the values of these parameters Forsimplificationpurposeandfollowingprescription
areestimatedusingmaximumlikelihood. Theassigned of Nipoti & Binney (2015), we will refer on the results
valuesofthemodelsparametersmaynotberelevantfor from Table 1 in further discussion.
the currentdata(e. g. the confidenceintervalsofmod- In Table 1, we compare and contrast quantitatively
els can be too wide). So, we optimized them by maxi- how well each model captured the pattern of periodic-
mizing the marginal log likelihood of the observations, ity and random noise of each light curve. Before any
usingBroyden-Fletcher-Goldfarb-Shanno(BFGS)algo- quantification of the difference between models, from
rithm,whichistheoptionoftheGPytool. Itautomati- simple visual inspection of Figure 3 can be seen that
callyfindsacompromisebetweenmodelcomplexityand there is more uncertainty around the predictions of
data fit. MOU then other models. Actually, in regions of sparse
6
data, the predictivedistribution is unchangedfromthe equals the intervals between the first flare and third,
prior,whereasin data dense regions,the uncertanity is and second and last flare. The largest period is some-
smaller. This suggests worse performance of GP with- what greater than the interval between the first and
out periodic components. The predictions improve no- the last flare in R band curve since its monitoring pe-
tably as more structure is added to models (i. e. pe- riod is shorter than optical spectroscopic light curves.
riodic kernels). The similar characteristics, but not so However it equals the interval between the two local
pronounced, can be seen in the models of continuum maximuminopticalfluxesobservedaroundMJD51500
4200˚A (Figure 4). Also,Table1 showsthatevenfrom and 56000 (see Figure1-3,5 in Paper I). Association of
combinationofthreeperiodickernelsthelearnedmodel ∼ 4500 days period with ocurence of two local maxi-
will select two periods. This is not the case when we mumin opticallight curvesas wellas its comensurable
used M , M models for the continuum 5100 ˚A and relation with period of ∼ 2200 days indicate its physi-
1 2
PM1 model for the continuum 4200 ˚A. Perhaps ex- cal relevance, despite of observationalcoverage of 5700
OU
planation lies in the very nature of these signals (i.e. days(seeVaughan et al.2016,forstatisticalconditions
the shortest periodic signal is weakest). of reliable period detection).
Noneofourcompositemodelswereabletorepresent
the Hβ line. An example of such failure is depicted 3.1 The nature of periodicity in E1821+643
on Figure 6. We attribute it to the loss of periodic
signal. Furthermore, in Paper I was shown that the Herewewillanalyzepossiblephysicaloriginofinferred
time delay of Hβ is larger than the Hγ emission line, periods. It is interesting to see do we have enough in-
meaningthattheHβ lineoriginatesatlargerradiithan formation to discern whether the detected periods can
Hγ. Combining previous facts (failure of Hβ modeling arise from the orbital motion of distinct sources of ra-
and its larger distance of origin), it seems that at the diation(such as brightspots)within the accretiondisk
Hβ line distance some chaotic processes can exist (like or distinct physical objects in mutual dynamical inter-
outflows) destroying all periodic signals. We discuss action (which can be also binary black hole system).
reasons for this scenario in the Section 3.1. Through variety of numerical simulations it has been
Moreover our learned models did not select a single determined that configurations of SMBH binaries as-
period which is the least common multiple of three (or sume circumsecondary(”mini”) disk aroundsecondary
two larger) periods. This fact shows that such model component, beside circumbinary disk. Also, there are
withoneperiodcannotreproducethelightcurveswhile scenarioswhere mini-disk is seen as the main source of
the models with two periods can produce them very the emission. One of signature of mini-disk presence
well. Theresultsimplythatthe contributionsfromthe wouldbe ripples and/oroscillationsin the broadFeKα
two periods are significant, and the flare-lake events line as suggested by McKernan et al. (2013). How-
reported in Paper I are not random but have the same ever, these signatures will be possible to detect with
physical origin. future observational missions. The numerical study of
Thecontributionofperiodicpartofeachmodelises- Noble et al. (2012) shown that if the accretion rate to
timatedwiththeratioSgivenbyEq. 4. Forallmodels the inner regions of circumbinary disk were compara-
the periodic part is strong (S≥1). We emphasize that ble to that of ordinary AGN, the surface density, and
flare like events are associated with GP mean curve, therefore the luminosity, of such a circumbinary disk
whichensuresthatthe behaviorofthe GP modelscan- could approach AGN level. Most strikingly, the lumi-
not simply be interpreted as overfitting due to added nosityshouldbe modulatedperiodicallyatafrequency
periodic kernels. Also, Table 1 shows that the LBF ex- determined by the binary orbitalfrequency and the bi-
hibited preference for the models with periodic parts. nary mass ratio. Since E1821+643 is one of the most
Factor values larger than 2 on the log scale gives de- luminous quasars and we detected periodical signals in
10
cisive evidence in favor of composite models (with pe- optical light curves, we will put our discussion within
riodic parts). Thus, all models (except M ) expressed thescenarioofluminouscircumbinarydisk,beingaware
3
decisive evidence in favor of composite models. How- of other possibilites.
ever, M model of Hγ (Figure 5) and the continuum The reason to consider the first possibility is that
3
5100 ˚A still indicate that the data express substantial broad line profiles of this object have unusual shapes
to strong evidence in favor of composite models. The (thorough line shape analysis is given in Paper I). If
meanperiodsaveragedoverallmodelsaregiveninMP the accretion disk is non-axisymmetric, changes in the
row of Table 1. The shortest value of averagedperiods line shapes will be caused by the disk features orbiting
is comparablewith the consecutive intervals between 4 on dynamical time scales. The dynamical time scale
peaks in R band photometric light curve found in Pa- (orbitaltime) of the light emitting regionat the radius
per I (see their Fig.1). While medium period almost r of the disk can be approximated by
7
Modelposteriormean
Ln(Fcnt5100) 95%Cl
Kernel=RBF+StdPer1+StdPer2+StdPer3 observations
Kernel=RQ+StdPer1+StdPer2+StdPer3
Kernel=Brw+StdPer1+StdPer2
Kernel=OU+StdPer1+StdPer2
Kernel=OU
Fig. 3 Comparison of GP posteriors (dashed line - posterior mean function, shaded region - 95% confidence interval) of
thecontinuumat5100˚A (denotedoneachplotasmarkerswitherrorbars)fordifferentkernelcomposites. Weusednatural
logarithmic values of fluxes. Composite of Brownian and three periodic kernels provide a richer form of mean function,
which explains betterdata.
vents us to explain two largest detected periods as the
orbitalmotionwithinthediskoffeaturessuchasbright
2
tdyn ≃5×10−3(cid:18)10M8M (cid:19)(cid:18)rrg(cid:19)3 sthpeotHs.βFlianielucraenobfeGePxpmlaoindeedlsbtyoitesxtdreatcetriporearitoiodns.frTohme
⊙
degree of its disfiguration is easily diagnosed by com-
GM parison of the mean profiles of Hβ and Hγ (see Figure
r =
g c2 17 in Paper I). The mean Hβ profile exhibits more ex-
tendedredwingthenthemeanHγ. Itcanbecausedby
where r = GM is the gravitational radius around a
g c2 differentmechanisms. Anadditionalemissioninthefar
black hole of mass M, G is the gravitational constant
wing of the mean Hβ or even characteristics of obser-
and c the speed of light. By means of this relation we
vational instruments can result in such degradation of
can estimate orbital time of the Hβ and Hγ broad line
the Hβ line. Also, some processes like outflows can oc-
emitting region. Using their time lags (120 ld for Hβ
cur within the disk. We recall that Bogdanovi´cet al.
and 60ld for Hγ, see Table 10in PaperI), we estimate
(2007) proposed that a black hole ejection should be
thatcorrespondingvaluesoft totheirradiioforigin
dyn uncommonintheaftermathofgas-richmergers. Inthe
are ∼11.3 yr and ∼7.1 yr, respectively. These orbital
gas-rich galactic mergers, torques arising from gas ac-
times within the disk are close to two mean largestpe-
cretion align the spins of supermassive black holes and
riods detected by GP (Table 1). However, GP models
their orbital axis with a large-scale gas disk. On the
failed to detect any period in the Hβ line, which pre- otherhand,kickingoftheblackholemayhappenespe-
8
cially in merging galaxies that are relatively gas-poor.
Hβ Kernel=Brw+StdPer1+StdPer2+StdPer3 Modelposteriormean
Soifthekickspeedofthisobjectisabouttheescapeve- 95%Cl
observations
locityofthehostgalaxy(asreportedbyRobinson et al.
2010) then we should expect that disk is gas poor and
violent processes may contribute to it.
Close similarity between calculated orbital times
(t ∼11.3 yr and 7.1 yr, corresponding to the radii
dyn
of Hβ and Hγ origin within the disk), and two largest
periods detected by GP in both continua (at 4200 ˚A
and 5100 ˚A) and the Hγ line is not decisive factor fa-
voring bright spots within the disk. Because of this we
consider another possibility that dynamically related
physical objects (which can be also binary black hole Fig. 6 The same as in Figure 3 but for the H emission
β
system) can cause periodic signals. line and just one composite kernel. The shaded area is
Two larger periods (∼ 2200 and ∼ 4500 days, see enlarged around all points, proving that the GP fails to
Table 1) can originate from the real periodic process, describe daataset. It seems that chaotic processes destroy
periodic signals in the Hβ line, preventing GP models to
because they are related commensurably (their ratio is
∼ 1). Information about third signalof ∼1700 days is captureany such pattern in thedata.
2
weakly presented in both continua and is absent from
Hγ. We note that this period is incommensurate with following relation to estimate the rotational period of
respect to the two largest periods. It indicates that black hole (see Volvach et al. 2007):
this signal can be caused by some other process (like
flares) within the accretion disk. In such case, 1700
m +m 3T T cosθ
dayorbitalmotionwithinthediskwouldcorrespondto 1 2 = R Prec (5)
m 4 T2
disk’s radius of 2 ld. 2 O
In the simplest scenario the two largest periods de-
where m ,m are masses of two black holes in the sys-
1 2
tected by GP can be attributed to the periods of the
tem, T and T are rotational period of black hole
R Prec
orbital motion of any two objects in the close mutual
and precession period of its orbit, while T is orbital
O
dynamical interaction.
period. Theangleθ ishalf-angleoftheprecessioncone.
Robinson et al. (2010) reported that E1821+643
Assuming the binaries can be of equal mass m ∼m ,
1 2
might be an example of gravitational recoil, while re- and taking m + m = 2.6 × 109M from Paper I,
sultingSMBHismovingwithvelocityof∼2100kms−1. T = 12.2 y1r T =2 5.8yr then we o⊙btain T ∼ 7.4
Prec O R cosθ
Velocitiesofsuchorderofmagnitudehavebeenfoundin
yr. In close binary system cosθ ∼ 1 since θ converge
many theoretical works. For example, Gonzalez et al.
to small values. It is clear that the period of 4.8 yr
(2007)predictedkicksof∼2500kms−1,Campanelli et al. (seeMProwofTable1)is notrelatedto therotational
(2007) found kicks velocity up to ∼ 4000kms−1 by period. This is one more reason to believe that it can
means of empirical formula, and Lousto & Zlochower
originate from the random processes within the disk.
(2009) obtained ∼ 3300kms−1 for nearly maximal Using the third Kepler’s law, we can calculate the
spins of black holes in the system. It is believed that
radius binary black hole orbit a:
configuration,leading to such events (with high veloci-
tieskicks),comprisesofnearlyequalmasscircularbina- T 2 4π2a3
O
= (6)
ries black holes, with large equal and opposite spins in 1+z G(M +m)
(cid:16) (cid:17)
the orbital plane. Morevover, orbital angular momen-
Regaining the value T = 5.8 yr, we derived a ∼ 6×
tum and sum of two black holes’ spin-vectors precess O
1014m. Theupperboundofnaturalgravitationalwave
around total angular momentum during the in spiral
frequency is given by f < 1 c3 , where M is
(Barusse and Rezzola, 2009). If the conservation of GW 4√2πGM
total angular momentum is valid up to some approxi- the mass of the system (Hughes 2003). In our case
mation then the spin and orbital plane precess at the fGW ∼4.4×10−6Hz. Thereisnopossibilitytomeasure
same frequency. In this scenario, we can assume that gravitationalwavesofsuchlowfrequencyusingcurrent
the precessions of two spin-vectors and orbital angular ground-basedinstruments.
momentum is ∼ 12 yr , while the orbital period of bi-
nary system is ∼ 6 yr (see Table 1). We can use the
9
Table 3 Comparison oflight curvemodelswith oneperiodiccomponent,usingBayesfactor. Thenamesofcolumnshave
thesame meaning as in Table 1. Subscript stands for modeling with heteroscedastic errors.
het
Lc M NP Per ML S LBF NP Per ML S LBF
het het het het het
[day]
MOUa 3 - -125.76 - - 2 - 119.62 - -
00 PMOUb6 4734 -127.79 1.0 2.02 5 1693 117.950.29 -1.66
t51 M1c 6 2348-112.2890.84-13.48 5 623 110.770.68 -8.9
n
C M2d 7 4529 -126.60 0.63 0.83 6 4500 113 0.63 -6.62
M3e 5 1160 -126.39 0.62 6.8 4 3519 113.03 0.6 -6.58
0 MOU 3 - -63.41 - - 2 - 105.76 - -
420 PMOUb6 1841 -64.36 1.18 0.95 5 3532 102.630.46 -3.13
Cnt M1d 7 2321 -61.21 1.0 -2.2 6 2419 105.5 0.35 -0.26
M2e 5 4562 -64.04 0.95 0.63 4 4763 102.04 0.9 -3.36
M3e 5 4171 -63.79 1.13 0.38 4 4763 102.4 0.62 -3.36
MOU 3 - -56.75 - - 2 - 31.37 - -
PMOUb6 2010 -58.22 1.03 1.47 5 - - - -
Hγ M1d 7 1265 -49.08 4.23 -7.67 6 2196 30.71 0.38 -1.0
M2e 5 5062 -58.5 1.02 1.75 4 5024 27.62 0.48 -4.09
M3e 5 4595 -57.38 1.04 0.63 4 4371 29.34 0.54 -2.37
aOUmodelofLCareconsideredasbasedmodel
bPMOU assumesOU+StdPeri,i=1
cM1 assumesSE+StdPeri,i=1
dM2 assumesRQ+StdPeri,i=1
eM3 assumesBrw+StdPeri,i=1
Table 4 Comparison oflight curvemodels withtwoperiodiccomponent,usingBayesfactor. Thenamesofcolumnshave
thesame meaning as in Table 1. Subscript stands for modeling with heteroscedastic errors.
het
Lc M NPPer1 Per2 ML S LBFNPhetPer1hetPer2hetMLhetShetLBFhet
[day][day] [day] [day]
MOUa 3 - - -125.76 - - 2 - - 119.62 - -
00 PMOUb9 3784 1760-126.861.06 1.1 8 3573 1154 110.76 0.7 -8.85
t51 M1c 9 3505 1789-130.361.07 4.59 8 2396 3040 306.944.96187.324
n
C M2d 10 4744 1607-128.110.63 1.01 9 2903 1682 113.820.96 -5.79
M3e 8 4415 4699-126.770.99 1 7 4827 4722 110.590.94 -9.03
0 MOU 3 - - -63.41 - - 2 - - 105.76 - -
420 PMOUb9 3750 1163 -66.09 6.65 2.68 8 3693 1153 101.410.59 -4.35
Cnt M1c 10 4890 2373 -63.90 0.99 0.48 9 955 1144 103.230.87 -2.53
M2d 8 4848 2373 -69.72 1.37 6.31 7 4727 - 103.120.95 -2.64
M3e 9 4317 2298 -64.44 1.02 1.03 8 3534 - 102.570.53 -3.19
MOU 3 - - -56.75 - - 2 - - 31.37 - -
PMOUb9 2648 2081 -60.08 1.03 3.33 8 4509 1124 26.42 0.89 -4.94
Hγ M1c 10 4975 2156 -62.05 0.99 5.30 9 4891 2167 26.04 1.0 -5.33
M2d 8 4993 2163 -61.91 0.97 5.16 7 4982 2167 26.04 1.0 -5.33
M3e 9 4975 2156 -62.06 0.99 5.31 8 4982 2167 26.04 0.99 -5.33
aOUmodelofLCareconsideredasbasedmodel
bPMOU assumesOU+StdPeri,i=1,2
cM1 assumesSE+StdPeri,i=1,2
dM2 assumesRQ+StdPeri,i=1,2
eM3 assumesBrw+StdPeri,i=1,2
10
Ln(Fcnt4200)
Kernel=RQ+StdPer1+StdPer2
Kernel=Brw+StdPer1+StdPer2
Kernel=OU+StdPer1+StdPer2+StPer3
Kernel=OU
Fig. 4 ThesameasinFigure3butforthecontinuumat4200˚A.Notethatmodelingwaspossibleonlywithtwoperiodic
kernels.
4 Conclusions So we propose that two largest periods correspond
to the mutual dynamical interaction of two objects at
We have performed a set of GP modeling (Table 1) least. Also, we analyzed the subcenario if these two
of the continua at 4200, 5100 ˚A and the Hγ emis- objects are two kicked merging black holes. In such
sion line light curves of E1821+643. Following direc- situation,weassociatedperiodof∼12yrwiththe pre-
tion of mixed covariance stochastic modeling of AGN
cessionoforbitalangularmomentumofthesystemand
light curves, we used composite GP models consisting
the precessionofthe sumoftwo blackholes’spins dur-
of Ornstein - Uhlenbeck, the squared exponential, the
ing the infall. The orbital period of black hole within
rationalquadratic,Brownianandnon-stationary–peri-
this systemcan be associatedwith ∼6 yr. In suchset-
odickernels. Thisallowsustoquantifyperiodicsignals
ting,therotationalperiodofblackholewouldbenearly
in the light curves. 7.4 yr, where θ is half-angle of the precession cone.
We found three strong signals with periods of ∼ 12 cosθ
yr,∼6yr,and4.8yr(averagedvaluesoverallmodels).
Acknowledgements This work was supported by
The two largestperiods arein aharmonicrelationship,
the Ministry of Education and Science of Republic of
which indicates their common physical origin. Ana-
Serbia through the project AstrophysicalSpectroscopy
lyzing the scenario with orbital motion of bright spots
of Extragalactic Objects (176001) and RFBR (grants
withintheaccretiondisk,wefoundthatt scalescor-
dyn N12-02-01237a, 12-02-00857a, 12-02-01237a,N15- 02-
respondingtotheradiioforiginoftheHβ andHγ lines
02101).
arealmostequalto the two largestperiods detected by
GP models. However,the brightspot motionis not fa-
voredbecauseanyinformationaboutperiodicsignalsis
vanished from the Hβ line. Since the Hβ line arises at
larger radius than Hγ, it seems that here some chaotic
processes (like outflow) are presented that are destroy-
ing all periodic signals.
The third period of 4.8 yr is incommensurate with
respectto the twolargestperiods, suggestingits differ-
ent physical origin (such as flares at the radius of 2 ld
within the accretion disk of the system).