Variability of signal to noise ratio and the network analysis of gravitational wave burst signals

0.4 MB·English

by  Soumya D. Mohanty

Variability of signal to noise ratio and the network analysis of gravitational wave burst
signals
Soumya D. Mohanty,1,∗ Malik Rakhmanov,2,† Sergei Klimenko,2 and Guenakh Mitselmakher2
1The University of Texas at Brownsville, 80 Fort Brown, Brownsville, Texas, 78520, USA
2University of Florida, P.O.Box 118440, Gainesville, Florida, 32611, USA
The detection and estimation of gravitational wave burst signals, with a priori unknown polar-
ization waveforms, requires the use of data from a network of detectors. For determining how the
data from such a network should becombined, approaches based on themaximum likelihood prin-
ciple have proven to be useful. The most straightforward among these uses the global maximum
of the likelihood over the space of all waveforms as both the detection statistic and signal estima-
tor. However, in the case of burst signals, a physically counterintuitive situation results: for two
aligned detectors the statistic includes the cross-correlation of the detector outputs, as expected,
but this term disappears even for an infinitesimal misalignment. This two detector paradox arises
6 fromtheinclusion ofimprobablewaveformsinthesolution spaceofmaximization. Suchwaveforms
0 producewidelydifferentresponsesindetectorsthatarecloselyaligned. Weshowthatbypenalizing
0 waveforms that exhibit large signal-to-noise ratio (snr) variability, as the corresponding source is
2 moved on the sky, a physically motivated restriction is obtained that (i) resolves the two detector
n paradox and (ii) leads to a better performing statistic than the global maximum of the likelihood.
a Waveforms with high snr variability turn out to be precisely the ones that are improbable in the
J sensementionedabove. Thecoherentnetworkanalysismethodthusobtainedcanbeappliedtoany
8 network,irrespective of thenumberor themutual alignment of detectors.
1

1 I. INTRODUCTION imum serves as an estimate of the true value of θ. The
v maximum,overθ, ofthe likelihoodratio(LR) functional
6 p(x|θ)/p(x|θ )yieldsadetectionstatistic. Thus,thelike-
Several interferometric gravitational wave (GW) de- 0
7 lihood based approach conveniently yields both a detec-
0 tectors are now operational around the world [1, 2, 3, 4]
tion statistic and source parameter estimates.
1 in addition to an already operating set of resonant mass
Forbinaryinspiralsignals,the polarizationwaveforms
0 detectors[5]. Sincegravitationalwavedetectorscanmea-
6 sureboththeamplitudeandphaseofincomingsignals,it can be computed in terms of a small set of parame-
0 ters describing the binary, such as the masses of the
ispossibletocombinethedatafromsuchanetworkofde-
/ two stars, the orbital inclination etc. Likelihood based
c tectors to not only detect GW sources better, compared
q to single detectors, but also estimate the waveforms of network analysis is fairly well understood in the case of
- such known waveforms [7, 8]. In this paper, we con-
r thetwoindependentpolarizationcomponentsoftheGW
g signal along with the direction to the source. sider likelihood based network analysis for burst signals
: that are short duration deterministic signals but have
v A straightforward approach to signal detection and
unreliable or simply unknown waveforms. Astrophysical
i estimation using data from a network of detectors is
X sources of burst signals include Supernovae [9], Gamma
through the likelihood [6] functional. The likelihood
r Ray Bursts [10] and the merger of compact objects such
a functionalofdata xis the jointprobabilitydensity func- as Neutron stars and Black Holes [11].
tionp(x|θ), whereθ isavectorofparameterswhose true
Althoughwedonothavepriorknowledgeofburstsig-
value for any given instance of x is not known a pri-
nal waveforms, it is possible to extend likelihood based
ori. In the case of GW detection, x is the collection
network analysis to burst signals in a formal sense. The
of outputs from all the detectors in a network and the
basicidea[12,13,14,15]istotreateachsampleofaburst
parameters are the sky position of the GW source and
signal as a parameter and maximize the likelihood over
the waveforms of the two GW polarization components,
all the parameters. We refer to the resulting detection
h (t) andh (t). The problemofsignaldetectionis con-
+ × method as standard likelihood (SL). Other approaches
cernedwithdecidingwhetherxisarealizationofp(x|θ )
0 to network analysis for burst signals have also been pro-
or p(x|θ 6= θ ), where θ denotes the signal absent case.
0 0 posedintheliterature[16,17]. Themethod[16]proposed
Theestimationproblemdealswithinferingthetruevalue
byGu¨rselandTintoisbasedonextremizingafunctional
of θ once a detection has been made. According to the
thatisconstructedpurelyoutoftheresponseofeachde-
maximum likelihood [6] approach to inference, the value
tectorinanetworkanddoesnotexplicitlydependonthe
θ at which the likelihood functional attains its max-
max polarizationwaveforms. A minimum of three misaligned
detectors are required in order to use the Gu¨rsel-Tinto
method whereas likelihood based methods do not have
this limitation.
∗Electronicaddress: [email protected]
†Currentlyat: DepartmentofPhysics,PennsylvaniaStateUniver- Althoughit is quite elegant,the SL methodhas a fun-
sity,UniversityPark,PA16802,USA damentalproblemthatwasfirstnotedin[14]andsystem-
2
aticallyanalyzedin[15](henceforth,KMRM).Calledthe portant qualitative differences in how the two methods
two detector paradox, the problem is most easily stated process data. For example, unlike constraint likelihood,
for a network of two detectors. When the detectors are the divergence from SL of the method introduced here
perfectly aligned, the SL detection statistic naturally in- depends on the actual data.
cludes the cross-correlation of detector outputs. Since Thepaperisorganizedasfollows. Wereviewthestan-
the responsetoa givenGW sourceis identicalinaligned dard likelihood method in Section II, establishing our
detectors, this term serves to increase the sensitivity of basic notation and conventions on the way. Section III
the overall detection statistic. However, it disappears contains a discussion of the two detector paradox. Sec-
when the SL approach is used for even infinitesimally tionIVconsiderstheSNRvariabilityofsignals,discusses
misaligned detectors. This is counterintuitive since we how it affects their detectability and defines a measure
expect the cross-correlation term to remain important for it. Section V contains the derivation of the detec-
for small misalignments. tion statistic and source parameter estimators that re-
Theoriginofthetwodetectorparadoxisthatallwave- sult from incorporating SNR variability. Results from
forms for h (t) and h (t) are allowed as solutions when simulations are presented in Section VI.
+ ×
maximizing the likelihood even though some of them
would be improbable in that they will produce widely
differentresponsesintwonearlyaligneddetectors. Itfol- II. STANDARD LIKELIHOOD ANALYSIS FOR
lows that this problem can be solved by removing such BURST SIGNALS
improbable waveforms from the space of allowed solu-
tions. In KMRM, this was achieved by imposing certain
We fix a Cartesian coordinate system with its origin
constraints on the space of waveforms yielding the class
at the center of the Earth, Z axis pointing at the North
of constraint likelihood methods. Modifying the likeli-
Pole and some arbitrary convention for the X and Y
hood functional with a waveform dependent weight fac-
axes. Thiswillbeourfiducialframeofreference. Letthe
tor, such as a Bayesianprior [18], is another way of sup-
sky position of an incoming GW signal in this frame be
pressing undesirable solutions. In the following we use
denotedbyθ,thepolarangle,andφ,theazimuthalangle.
the term “restriction” to encompass all approaches that
AnincomingGWsignalismostconvenientlyexpressedin
remove improbable solutions or reduce their influence.
the transverse traceless (TT) guage [21] associated with
Thenon-trivialaspectofsettinguprestrictionsisthat the source direction. The Z axis of the TT frame is
theperformanceoftheresultingstatisticmustnotbede- orientedalongthedirectiontotheGWsource. Weadopt
gradedand,ifpossible,shouldbe better thanthe SL ap- theconventionthattheX axisoftheTTgaugeliesinthe
proach. Additionally,the restrictionsmustbe applicable plane formed by the TT frame Z axis and the Z axis of
to all astrophysicalburst sources otherwise the resulting the fiducialframe. The GW signalis specified by its two
statistic willbe better foraparticulartype ofsourcebut independentpolarizationcomponentsh (t)andh (t)in
+ ×
may be worse for others. If sufficiently reliable knowl- the TT frame.
edge for a particular class is available, then it can be In this paper, we take into accountthe differentorien-
used[19,20]asanadditionalsetofrestrictionstofurther tations and geographic locations of currently operating
enhancesensitivity. TheconstraintsproposedinKMRM interferometric detectors but treat them as identical in
(seeSectionIII)satisfytheabovecondition,namelyinde- sensitivity. For instance, a GEO detector in the paper
pendencefromassumptionsaboutsourcedynamics. Sim- willmeanadetectorhavingthe samelocationandorien-
ulationsshowthatconstraintlikelihoodmethodsperform tation as GEO600 but with sensitivity similar to that of
significantly better than SL, especially in reconstructing LIGOdetectors. Thissimplifiesthe algebraconsiderably
the sky position of a burst source. However, the con- whithout altering the main result which is the relative
straintsusedwereafirstexampleandarenotnecessarily comparison of methods.
uniqueoroptimal. Itisimportanttoinvestigateindepen- Since real data is digital, we shift to the discrete
dent avenues of formulating generic restrictions in order time domain and from here on denote any time depen-
togainmoreinsightintothis promisingapproachtonet- dent quantity g(t) as g[j], where t = jδ, δ being the
work analysis. sampling interval. A finite sequence of time samples
We investigate the issue of generic restrictions from a {g[0],g[1],...,g[m−1]}willbe denotedbyg. Thescalar
physical point of view. It is found that the more im- response of a detector s[k] to an incoming GW is given
probable a set of waveforms h+ and h×, the faster is by,
the change in detectability of the corresponding source
as a function of sky position. By measuring the rate of s[k]=Wij[k]D , (1)
ij
changeofsignal-to-noiseratio(SNR)asafunctionofsky
position, therefore,we canrestrictsolutions that areun- where Wij is the symmetric trace free tensor that de-
desireable in the maximization of the LR. Monte Carlo scribestheGWinthetransverse-tracelessgaugeandD
ij
simulations show that the detectionstatistic that results is an equivalent tensor associated with the detector [22].
from the above restriction performs roughly the same as (The Einstein summation convention is in force in the
constraint likelihood methods. However, there exist im- aboveexpression.) Theresponses [k]oftheith detector,
i
3
i=1,...,N , is expressed as, The natural logarithm of the likelihood ratio (LR),
d
Λ(X|θ,φ,S), for an N detector network is given by
d
s [k]=F (θ,φ)h [kδ−δ ]+F (θ,φ)h [kδ−δ ], (2)
i +,i + i ×,i × i
1 NdN−1
where F+,i and F×,i depend on the rotations involvedin Λ(X|θ,φ,S) = 2 xi[j]2−(xi[j]−si[j])2 (9)
going from the TT frame associated with a detector to Xi=1Xj=0(cid:0) (cid:1)
the wave frame associatedwith the fiducial frame and δi N−1 1
is the travel time for the signal from the origin of the = X[j]·S[j]− kS[j]k2 , (10)
(cid:20) 2 (cid:21)
fiducial frame to the detector. For a given θ, φ, we can Xj=0
treat the detectors as co-located (δ = 0) at the origin
i
of the fiducial frame since the output from eachdetector where the detectors outputs are implicitly time shifted
can be shifted in time to compensate for δ . as explained above.
i
For a network of detectors, it is convenient to define The maximum likelihood ratio approach can be ex-
quantities that are analogous to those for a single de- tended to burst signalsby substituting Eq. 4 into Eq.10
tector. Thus, we can define a network response vector, and maximizing Λ(X|θ,φ,S) over each h+[j], h×[j] in-
dependently. Since the resulting equations involve only
quantitiesatthesamesamplej,thisinvolvesmaximizing
s [k] the summand in Eq. 10 independently for each j. From
1
S[k]= ..  . (3) now onwards we will focus on only one time instant and
.
dropthetimeindexwhenthereisnoscopeforconfusion.
 sNd[k] The quantity to be maximized is,
Vectors in RNd, the space of real N -tuples, will be de- 1
d λ[j]=X[j]·S[j]− kS[j]k2 , (11)
noted by boldface capital letters. By convention, ev- 2
ery vector A ∈ RNd will be a column vector. We de-
fine the Euclidean scalar product in this space, A·B = overthevectorSinRNd,thespaceofrealNd-tuples. Let
ATB. A sequence of vectors will be denoted by A = the maximum of λ occur at Smax. If the maximization
{A[0],...,A[m−1]}. allows all possible S, then it is obvious that Smax = X.
For given θ, φ, However, we are not allowed to choose any arbitrary S
in the evaluation of λ since from Eq. 4, S is a vector
S[k]=F(θ,φ)(cid:18) hh×+[[kk]](cid:19) , (4) tThhaits dimeppelinedssthoantoSnlmyutswtoliepainraamtewteorsdihm+e[jn]siaonndalhp×la[jn]e.
determined completely by the antenna pattern functions
where the matrix F is defined as, of the detectors and passing through the origin of RNd.
Though the maximization of λ is restricted to S lying
F (θ,φ) F (θ,φ) in the response plane, in the maximum likelihood ratio
+,1 ×,1
F(θ,φ)= .. ..  . (5) approach there are no further constraints on the choice
. . of S within this plane. We call the method of maximiz-
F+,Nd(θ,φ) F×,Nd(θ,φ)  ing λ overS with no further restrictions as the standard
likelihood (SL) method [11, 13, 14, 15, 19]. We denote
Given a network response vector S[k], the polarization the value of Λ (Eq. 10), obtained by summing the max-
components can be obtained as, imized values of λ corresponding to each time instant,
by Λ (θ,φ). This is the detection statistic of the SL
SL
h+[k] = H−1FT(θ,φ)S[k], (6) method.
(cid:18)h×[k](cid:19)
H(θ,φ) = FT(θ,φ)F(θ,φ), (7)
III. THE TWO DETECTOR PARADOX
providedH−1 exists. (NotethatHisasymmetric2-by-2
matrix and, hence, has two real eignevalues.) Inthetwodetectorcase,thenetworkdatavectorXlies
Letthenoiseintheith detectorben . Inthefollowing in the response plane. When the detectors are aligned,
i
we will assume that n is a zero mean, white, Gaussian we are constrained to use a network response vector S
i
process, i.e., E[n [j]] =0 and E[n [j]n [k]]=δ δ . In lying along the 45◦ line since the responses of the two
i i m im jk
analogy to the network response vector S[k], we define a detectorsmustbe identical,s1[k]=s2[k]=s[k]forallk.
network data vector X[k], Maximizing λ over the response s[k] yields,
1
x1[k] n1[k] s = (x1+x2), (12)
X[k]= .. = .. +S[k]. (8) 2
. . 1
λb ∝ (x +x )2 , (13)
xNd[k]  nNd[k] 2 1 2
b
4
where λ denotes the maximum value and λ = λ(s). It IV. VARIABILITY OF SNR OVER THE SKY
follows that,
b b b
Though our physical intuition suggests that an in-
N−1
finitesimal misalignment between two detectors should
Λ (θ,φ)∝ x [j]2+x [j]2+2x [j]x [j] . (14)
SL 1 2 1 2 not result in very different responses s ≃ s , this does
1 2
Xj=0 (cid:0) (cid:1) not forbid the existence of such sources. As discussed
above, mathematically, h (t) and h (t) can be com-
+ ×
Thus, the SL detection statistic contains a cross- pletelyarbitraryand,inparticular,canbesuchthats is
1
correlation term, x1[j]x2[j], which agrees with our very different from s2. However, we expect such sources
intuition: since thPe responses of the two detectors are to be improbable. Is their a quantitative basis for our
identical, the cross-correlation term will always have a expectation? In this section we show that, indeed, there
positive definite mean in the presence of a signal and, is a natural argument that “probable” GW sources, for
therefore, should be useful for detection. which s ≃ s , should occur more often than “improb-
1 2
However, the same procedure, when carried out for able” ones. Note that the categorization of a source as
eveninfinitesimallymisaligneddetectors,yieldsΛ that probableorimprobablealsodependsonthemutualalign-
SL
has no cross-correlationterm, mentofdetectors. Forexample,twodetectorsorientedat
45◦ to eachother canoften have very different responses
s =x ; s =x , (15) forthesamesource. Asshownbelow,ourargumenttakes
1 1 2 2
thisfactorintoaccount. Wefirstpresentthetwodetector
λ ∝ x2+x2 . (16)
b b1 2 case and then consider an arbitrary network.
b
This clearly violates what one expects on physical
groundssincethe cross-correlationtermshouldremaina
useful discriminant against noise when the detectors are A. Two detector case
only infinitesimally misaligned. However,for the class of
probable waveforms, the SL statistic becomes less sensi-
Consider a gravitational wave source located at (θ,φ)
tive than one in which a cross-correlationterm is added
with signal polarization components h [k] and h [k].
by hand. Thus, the SL approach fails to deliver a phys- + ×
Let the detector response vector be S[k]. Now consider
ically acceptable detection statistic in the two detector
an identical source located at (θ′,φ′). By this we mean
case. Following KMRM we call this problem the two
thatthefiducialframe,alongwiththeGWsource,bero-
detector paradox.
tatedrigidly leavingthe detectors in their originalplace.
Theoriginofthetwodetectorparadoxliesinthemax-
This leads to a network response S′[k] that is related to
imization of λ over all waveforms, h and h , including
+ × S[k] by,
thosethatcanproduceverydifferentresponsess ands
1 2
foreveninfinitesimally misaligneddetectors. Suchwave-
forms can lead to small or even zero cross-correlation S′ =F(θ′,φ′)F−1(θ,φ)S. (17)
between the responses. Physical intuition suggests that
such waveforms should be imrobable. Mathematically, Now,assumethatweusematchedfiltering[23]to detect
however, such waveforms are allowed since the matrix the source at the instant k. Thus, our detection statistic
F is invertible, even if it is nearly singular, and some will be X[k].S[k] and X[k].S′[k] for the two sources re-
h+ and h× can always be obtained for any given s1 and spectively. The detectability of a source under matched
s2. As a result, the SL approach“throws”out the cross- filtering is measured by the signal to noise ratio (SNR),
correlationtermsincethis termcannowcontribute pure ρ,definedastheratioofthemeanofthedetectionstatis-
noise. Geometrically, the two detector paradox comes tic to its standard deviation. For the two sources above,
about because the data vector X lies in the response the respective SNR will be ρ = kS[k]k and ρ′ = kS′[k]k.
plane andthe SL approachdoesnotconstrainthe choice Notethatweareconsideringamorerestrictedversionof
of S in this plane. Hence, the solution S = X that re- matchedfiltering thanthe conventionalonein whichthe
moves cross-correlation is an allowed one in this case. detection statistic would be X[k].S[k] and the SNR
This suggests that restricting the allowbed solutions in would be ρ2 = kS[k]k2. WePreturn to this point later
the response plane can resolve the two detector para- in this section.P
dox. The SL detection statistic recovers pair-wise cross- We denote the relative difference in SNR for the two
correlations in the case of three or more misaligned de- source locations as,
tectorssinceXnowliesoutsidetheresponseplane. How-
ever, the two detector paradox is just an extreme mani-
ρ′−ρ
festation of a problem that occurs for any arbitrary net- ǫ(θ,φ,θ′,φ′;S) = , (18)
ρ
work: depending on the geometry of the network, some
b
solutions are improbable and their inclusion in the max-
imization of the LR reduces its sensitivity. where S is the unit vector oriented along S. For small
b
5
displacements, θ′ =θ+δθ and φ′ =φ+δφ, we get, is assumed that the detectors are colocatedat the origin
of the fiducial frame and coplanar with angles between
∂ǫ ∂ǫ
ǫ(θ,φ,θ′,φ′;S) = δθ +δφ , (19) their respective X arms of α. The antenna patterns are,
b ∂θ′(cid:12)(cid:12)(cid:12)θφ′′==θφ, ∂φ′(cid:12)(cid:12)(cid:12)θφ′′==θφ, in this case,
∂ǫ = ∂ k(cid:12)F(θ′,φ′)F−1(θ(cid:12),φ)Sk2 , 1(1+cos2θ)cos2φ cosθ sin2φ
∂θ′(cid:12)(cid:12)(cid:12)θφ′′==θφ, ∂θ′q b (cid:12)(cid:12)(cid:12)θφ′′==θφ, F =12(1+cos2θ)cos2(φ+α) cosθ sin2(φ+α) .
(cid:12) ∂F(θ,φ) (cid:12) 2 
= S· F−1(θ,φ)S , (22)
(cid:18) ∂θ (cid:19)
Eachplotisoverallorientationsofthedetectorresponse
∂ǫ = Sb· ∂F(θ,φ)F−1(θ,φ)Sb . (20) vector S. First, we see that for almost coaligned detec-
∂φ′(cid:12)θ′=θ, (cid:18) ∂φ (cid:19) tors,α=0.1◦, Ifwelimitthe searchforthe maximumof
(cid:12)(cid:12)φ′=φ b b λ(c.f.,Ebq.11)overS suchthatG(θ,φ;S)≤10,say,then
(cid:12)
The norm of the gradient onlydetectorresponseswithorientationveryclosetothe
45◦ and −45◦ line wbill be allowed. Thibs means that the
2
G(θ,φ;S) = S· F F−1S + instantaneous response in the two detectors should be
φ
(cid:20)(cid:16) (cid:16) (cid:17)(cid:17) nearly equal. From this example, we see that waveforms
b b b 2 1/2 with low SNR variability are also the ones that lead to
S· FθF−1S , (21) nearlyidenticalresponsesinthetwodetectorsand,asdis-
(cid:16) (cid:16) (cid:17)(cid:17) (cid:21)
cussed above, are the ones which have relatively higher
b b
F = ∂F/∂θ and F = ∂F/∂φ, is a convenient measure chances of being detected. This bears out our physi-
θ φ
of the variability of the SNR of a source as a function of cal intuition that sources that produce very different re-
its location. sponses in nearly aligned detectors should be rare. As α
AhighvalueofG(θ,φ;S)impliesthatthesourcewould is increased, one sees that not only does G decrease but
have to occur at very special positions on the sky to be that it also approachesthe same value for all orientation
b
detectable. This is because if it has an SNR above some and, hence, all orientations of S are now allowed. The
detection threshold at one location, it may rapidly fall pattern opens up fully when the detectors are at 45◦ to
below it if it were to occur at a slightly different loca- eachother. Asαapproaches90◦b,theresponsesareagain
tion. Consider two classes of sources, class A with a increasingly restricted around the ±45◦ lines.
large associated G(θ,φ;S) and class B with low values ThoughFig. 1 shows that SNR variabilityis a natural
of the same. Assume that events from both classes oc- measureforrestrictingtheorientationofS,theparticular
b
cur with the same comovingrate and that they have the measure G(θ,φ;S) allows both ±45◦ lines as the locii of
same distribution on the sky. Then, due to the rapid restriction. Fornearlyaligneddetectors,bhowever,weex-
fluctuationofthe SNRforsourceclassA,the areaofthe pect detector resbponsesto be in-phase whichmeans that
sky over which its members will be detectable would be the responses should lie in the vicinity of the +45◦ line.
significantly less than that of class B. It follows that the
We introduce another function, Γ(θ,φ;S), that is also a
observed rate of the class A sources would be much re-
measure of SNR variability but has better properties in
duced in comparison to class B. Therefore, even though b
this respect.
the two sources may be equivalent in terms of, say, the
energy emitted in GWs, the high SNR variability source 2 2
Γ2(θ,φ;S) = F F−1S + F F−1S , (23)
would be detectable less often. Since a larger volume θ φ
(cid:13) (cid:13) (cid:13) (cid:13)
of space would required to detect such sources, it will b (cid:13)(cid:13) b(cid:13)(cid:13) (cid:13)(cid:13) b(cid:13)(cid:13)
naturally be harderto detect at a givenfalse alarmrate. In fact, from the Cauchy Schwartz inequality, it follows
Wethusarriveatthemainpointofthissection. Given that
a pair of misaligned detectors, we should include only
those responses S in our search for the maximum of the G(θ,φ;S)≤Γ(θ,φ;S). (24)
LRforwhichG(θ,φ;S)islowinsomewelldefinedsense.
b b
Hence, Γ(θ,φ;S) is an upper bound on SNR variabil-
In this way, the detectability of sources with an intrin-
b
sically low observable rate will be reduced but that of ity. Figure 2 shows polar plots of Γ(θ,φ;S) for fixed θ,
b
more frequently observable sources will be enhanced. It φ and different misalignment angles α. Now, for small
is important to note that the above argument does not misalignement (α = 0.1◦), one sees that obnly responses
require any prior assumptions about the astrophysics of along the +45◦ line are allowed while for misalignement
sources. AlsonotethatGdependsonlyontheorientation close to 90◦, only responses that are close to −45◦ are
ofS andnotits magnitude. Thatis, itis only concerned allowed. It makes sense, therefore, to use Γ(θ,φ;S) as a
with the shape of the waveforms of h+ and h× and not measure of SNR variability instead of G(θ,φ;S).
b
the distance to the source. Intheabovediscussion,weconsideredthedetectionof
Figure 1 shows polar plots of G(θ,φ;S) as a function asignalatasingleinstantusingamatchedfiltber,leading
of S for fixed θ, φ and different detector alignments. It tothedetectionstatisticX[k]·S[k]. Burstsignalscanlast
b
b
6
Misalignment=0.1deg; theta=60deg; phi=30 Misalignment=10deg; theta=60deg; phi=30
90 600 90 10
120 60 120 60
400
150 30 150 5 30
200
180 0 180 0
210 330 210 330
240 300 240 300
270 270
Misalignment=45deg; theta=60deg; phi=30 Misalignment=89.9deg; theta=60deg; phi=30
90 2 90 600
120 60 120 60
400
150 1 30 150 30
200
180 0 180 0
210 330 210 330
240 300 240 300
270 270
FIG. 1: Plots of G (Eq. 21) for different misalignment angles.
overseveralsamplingintervalsand,ingeneral,thedetec- (θ,φ) and (θ′,φ′) is,
tion statistic would in fact be a sum X[k]·S[k] overa
range of k. In this case, a high SNRPvariability response ǫ(θ,φ,θ′,φ′;S)=kF(θ′,φ′)H−1(θ,φ)FT(θ,φ)Sk−1.
S[k] can be balanced by a very low variability response (25)
S[k′]atsomeotherinstant,thusleadingtoalowvariabil- As before, forba small displacement of the sourcbe,
ity for the overallSNR =[ kSk2]1/2. This would allow
a broader range of signals tPo be used for maximizing the ∂ǫ = S· F H−1FTS , (26)
oinvsetraanlltalnikeeoluihsooSdNRΛ v(ac.rfi.a,bEilqit.y1i0s).simOpulrycahociocensoefrvuastiinvge ∂θ′(cid:12)(cid:12)(cid:12)θφ′′==θφ, b (cid:16) θ b(cid:17)
(cid:12)
one since signals with low instantaneous SNR variability ∂ǫ
= S· F H−1FTS , (27)
wMiollrenoevceers,saarsildyishcauvsseeldoweavrlaireira,bpihliytysicoafltihnetiuritoivoenrasullgSgNesRts. ∂φ′(cid:12)(cid:12)(cid:12)θφ′′==θφ, b (cid:16) φ b(cid:17)
(cid:12)
that the responses in closely aligned detectors should be
and our measure of SNR variability will be defined as,
nearly identical at each instant of time and not just on
the average. 2 2
Γ2(θ,φ;S) = F H−1FTS + F H−1FTS (.28)
θ φ
(cid:13) (cid:13) (cid:13) (cid:13)
b (cid:13)(cid:13) b(cid:13)(cid:13) (cid:13)(cid:13) b(cid:13)(cid:13)
Here,S isaunitvectorintheRNd spacebutlyinginthe
B. General network of detectors
response plane.
b
Foralgebraicsimplicity, itis convenientto expressthe
We will now generalize the SNR variability measure responsevectorSintermsofabasisintheresponseplane
Γ(θ,φ;S), introduced in Eq. 23, to the case of an arbi- ratherthanthefullRNd space. Inordertodoso,weneed
trarynetworkofdetectors. The relativedifference in the to fix two orthonormal vectors q and q on this plane.
0 1
b
observed SNR of the same source between two positions If e and e are the two orthonormal eigenvectors of H,
0 1
b b
b b
7
Misalignment=0.1deg; theta=60deg; phi=30 Misalignment=10deg; theta=60deg; phi=30
90 1500 90 15
120 60 120 60
1000 10
150 30 150 30
500 5
180 0 180 0
210 330 210 330
240 300 240 300
270 270
Misalignment=45deg; theta=60deg; phi=30 Misalignment=89.9deg; theta=60deg; phi=30
90 3 90 1500
120 60 120 60
2 1000
150 30 150 30
1 500
180 0 180 0
210 330 210 330
240 300 240 300
270 270
FIG. 2: Plots of Γ (Eq. 23) for different misalignment angles.
correspondingtoeigenvaluesν andν respectively,then where ρ = kSk. From the above it is obvious that the
0 1
it is easily seen that for, solution for S that maximizes λ is the projection of X
onto the response plane.
1
q = F(θ,φ)e (29)
0 0
|ν |
0
b p1 b
q = F(θ,φ)e (30)
1 1
|ν |
1
q ·qb = epTHe =0. b (31)
0 1 0 1
kq0k = kq1k=1. (32) V. PENALIZED MAXIMIZATION OF
b b b b
LIKELIHOOD RATIO
Since, q and q hbave the sbame form as the detector re-
0 1
sponses (c.f. Eq. 4), they are confined to the response
plane. bLet γ bbe the angle between S and q0. Then, We have shown that SNR variability, as measured by
Γ(θ,φ;γ) defined in Eq. 28, correlates well with how im-
S = kSk(q0cosγ+q1bsinγ)b, (33) probable a solution of the LR maximization problem is.
Improbablewaveformsh [k]andh [k]thatleadtovery
The SNR variability mebasure Γ(θ,bφ;S) then becomes a + ×
differentresponsesinnearlyaligneddetectorsalsocorre-
function Γ(θ,φ;γ) of γ.
b spondtoafastvariationinthedetectabilityofthesource
Recall Eq. 11 for the definition of the quantity λ that
as it is displaced on the sky. Such sources naturally be-
wewanttomaximizeovertheresponsevectorS. Wecan
comerarecomparedtoprobableones. InthisSection,we
now rewrite it as,
constructadetectionstatisticthatincorporatesΓ(θ,φ;γ)
1 to restrict the influence of improbable solutions on the
λ=(X·q0)ρcosγ+(X·q1)ρsinγ− ρ2 , (34) maximization of the LR.
2
b b
8
A. Algebraic simplifications m1 respectively. From Eq. 34 and 42, it follows that,
We begin by making some algebraic simplifications. b X·q 1
λ = ρ cosγ sinγ 0 − ρ2 ,
Substituting Eq. 33 into 28, we get, (cid:18)X·q1 (cid:19) 2
(cid:0) (cid:1) b
1
cosγ = ρ(α cosψ+β sinψ)b− ρ2 ,
Γ2(θ,φ,γ) = cosγ sinγ M , (35) x x 2
(cid:18) sinγ (cid:19)
(cid:0) (cid:1) = ρl cos(ψ−φ )− 1ρ2 , (43)
A C x x
M = , (36) 2
(cid:18)C B (cid:19) l = α2 +β2 , (44)
x x x
1
A = |ν0|(cid:0)kFθe0k2+kFφe0k2(cid:1) , (37) φx = ptan−1(cid:18)αβxx(cid:19) . (45)
1 b b
B = kF e k2+kF e k2 , (38)
θ 1 φ 1
|ν |
1 (cid:0) (cid:1) The maximization of λ has to be performed over ρ and
1 b b ψ. Thus, the SL solution is ρ = l , the length of the
C = [(F e )·(F e )+ x
|ν0ν1| θ 0 θ 1 projection of X onto the response plane, and ψ = φx
p(F e )·(F eb)] . b (39) which is the angle between the projection and the m0
φ 0 φ 1
vector.
We now perform anotherbchange obf basis in the response B. Penalized maximization b
plane. Let the eigenvalues and eigenvectors of M be µ ,
0
µ1 ≤ µ0 and m0 and m1 respectively. Since m0 and m1 Now we introduce one possible way to incorporate
form an orthonormalbasis, Γ(θ,φ,ψ) into the maximization of λ. We call this
b b b b
methodpenalized maximum likelihood ratio(PMLR).We
cosγ =m cosψ+m sinψ , (40) propose to maximize λ′,
(cid:18) sinγ (cid:19) 0 1
for some angle ψ betweenbthe responbse vector S andm . λ′ =λ−αΓ2(θ,φ,ψ), (46)
0
Then,
b where α ≥ 0 is a parameter under our control. The
Γ2(θ,φ;γ)=µ0cos2ψ+µ1sin2ψ . (41) extra term αΓ2(θ,φ,ψ) acts as a penalty that prevents
thesolutionfromreadilyconvergingtothedatavectorX.
Henceforth, we denote Γ2(θ,φ;γ) as Γ2(θ,φ;ψ). Thisresolvesthe twodetectorparadox(c.f., SectionIII).
Using the orthonormal basis vectors m and m , we
0 1 Once the solution ρ, ψ that maximizes λ′ in Eq. 46 is
can state
b b found, we find the value of λ at this solution (not λ′).
X·q Our final data functibonabl will be the sum of the values
(cid:18)X·q10 (cid:19)=αxm0+βxm1 , (42) of λ, found as above, over all samples of the data X.
b Reintroducing the time index and denoting the resulting
b b
where α and β debnote the components along m and functional as Λ (θ,φ), we get
x x 0 PL
b
1
Λ (θ,φ)= max ρl [j]cos(ψ[j]−φ [j])− ρ[j]2−α(µ cos2ψ[j]+µ sin2ψ[j]) . (47)
PL x x 0 1
Xj ρ,ψ (cid:20) 2 (cid:21)
The global maximum of Λ (θ,φ), over θ and φ, will data X projected on the 2D response plane. Also, the
PL
serve as the final detection statistic. importance of the penalty term diminishes as l−2. If a
x
The solution ρ, ψ that maximizes λ′ is easily found. signal is present, the latter increases with an increase
in the instantaneous signal amplitude. Hence, whenever
ρb=bl cos(ψ−φ ), (48) the signal amplitude is strong enough, the effect of the
x x
penaltytermisreducedatthatinstantandthewaveform
sin(2φ )
tan(2ψb) = b x . (49) estimators ψ and ρ converge to the SL estimates.
cos(2φ )−2α(µ −µ )/l2
x 0 1 x
b The termbµ0−µb1, which is always positive by defini-
Notethatwhenα=0,werecoverthestandardlikelihood tion, is a measure of how misaligned the detectors in a
ratio statistic Λ (θ,φ) wherein the solution for S is the network are. This difference is reduced with an increase
SL
9
inmisalignment,hence reducingtheeffectofthe penalty VI. RESULTS
term. Fig. 3 shows µ −µ as a function of sky position
0 1
for two networks: (1) three LIGO detectors: one at Liv-
ingston,andtwo atHanford,and (2)one LIGO detector
The PMLR method was implemented in Matlab [24]
at Livingston, one GEO600 and one VIRGO detector.
anditsperformancewasquantifiedwiththeMonteCarlo
Fig 3 also shows the ratio of the differences in eigenval-
simulation code previously used in [15]. We give a brief
uesforthetwonetworks. Itisclearthatthedifferencein
recapitulation of the simulation scheme: Binary Black
eigenvaluesforthefirstnetwork,whichconsistsofnearly
Hole merger waveforms from [25] are injected into short
aligned detectors, is larger than that for the second net-
stretches of white noise, one waveform per stretch. The
work over a large fraction of the sky.
signalscorrespondtoface-onbinariesscattereduniformly
The statistic Λ (θ,φ) solves the two detector para-
PL on the sky with randomly selected polarization angles.
dox by continuously bridging the alignedandmisaligned
ThesignalshaveavariablematchedfilteringSNRdueto
detectors cases. From Eq. 49, it is clear that as µ −µ
0 1 their sky position and polarization, but an overallfactor
increases,i.e., as the detectorsbecome more aligned,the
Gisusedtoscaletheamplitudesofthetwopolarizations
denominatorincreasesinthenegativedirectionwhilethe
h (t) and h (t). The sky-averaged SNR of the injected
+ ×
numerator, sin2φ , stays the same. Hence, 2ψ → π
x signals is approximately related to G as SNR ≃ 2.3G.
which implies ψ → π/2. From Eq. 41, the penalty func-
Further details can be found in [15].
tionΓ2hasitslowestvalueforψ =π/2. Thisisthedirec-
tion in the response plane along which all the responses
Consider, first, some case studies to give a qualitative
in a network get oriented as the detectors are aligned.
feel for the performance of PMLR compared to the SL
Therefore, with progressive alignment of the detectors
and constraint likelihood methods. Fig 4 compares SL,
in a network, the solution ψ is forced towards equal re- constraint likelihood (hard constraint) and PMLR (with
sponses in all detectors, thus leading smoothly into the α = 1500) for the same signal using the three detector
b
completely aligned case. In contrast, for the standard network: Hanford, Livingston and GEO. The true sky
likelihood statistic ΛSL(θ,φ) (=ΛPL when α=0), there position of the source was θ = 50◦, φ = 280◦. The best
exists a fundamental discontinuity betweenthe solutions estimates for sky position (in degrees) were (1) θ =76◦,
for aligned and misaligned cases: For the aligned case φ=32◦ for SL and(2) θ =48◦,φ=280◦ for PMLRand
ψ =π/2, while ψ =φx for the misaligned case. constraintlikelihood. Notethat,comparedtoSL,thesky
Finally, the penalty term, Γ2, depends only the direc- maps for constraint likelihood and PMLR are less noisy
tbion of the instabntaneous network response vector S. It andshow more contrastbetweenthe regionin which the
does not restrict its length. Hence, the magnitude of S signal is located and the rest of the sky.
at two distinct time instants can be arbitrarily different.
Thus, even white noise response waveforms are allowed The simulation code finds out the fraction of all in-
assolutions. Theonlyrequirementisonthemutualcon- jectedbinariesdetectedbythePMLR(orother)method.
sistency of detector responses at each instant of time. The detection statistic for all methods is the maximum
value over the respective skymap (e.g. Fig 4). Fig 5
shows the receiver operating characteristics (ROC) for
C. Connection with Bayesian analysis SL, constraintlikelihood(both hardandsoft constraint)
andPMLR(withα=150and1500). Themeasuredfrac-
The penalized likelihood approach, Eq. 46, has a tion ofdetected signals is plotted againstthe false alarm
straightforward connection to Bayesian [18] analysis. In rate corresponding to the detection threshold used. It is
the latter approach, the posterior degree of belief over clear from the figure that PMLR can vary continuously
the space of waveforms, or equivalently ρ and ψ (c.f., in performance, ranging from SL for low values of α to
Eq. ref), is obtained using Bayes’ law, the hard constraint method for high values.
p(ρ,ψ)p(X|ρ,ψ) The simulation code also estimates the error in sky
p(ρ,ψ|X)= , (50)
p(X) position reconstruction. Fig. 6 shows the fraction of the
injected signals that were detected within 8 degrees of
wherep(X|ρ,ψ)issimplythelikelihoodλinEq.11. The theirtrueposition,forarangeofthe scalefactorG. The
maximum a posteriori Bayesian estimates for the true ρ PMLR method (α=1500)has practically the same per-
and ψ would be the ones for which p(ρ,ψ|X) is maxi- formance as the hard constraint method and both per-
mized. Thisisequivalenttomaximizingthelogarithmof form somewhat better than the soft constraint method.
p(ρ,ψ|X)and,hence,thesumofλandlogp(ρ,ψ). Com- The performance of SL, however, is significantly worse
parisonwithEq.46showsthat1/Γ2 actsasaprior. The than all the other methods. This can also be seen quali-
difference is that this prior is not normalized in Eq. 46 tatively in the sky maps of the various methods (Fig. 4)
but is instead multiplied by α. Of course, α can be cho- wherethe SL skymapalwayshasless “contrast”around
sen to be the integral of 1/Γ2 over ψ thus reducing the the true signal location compared to the constraint like-
penalized maximization to a true Bayesiananalysis. lihood and PMLR methods.
10
FIG. 3: The difference in eigenvalues µ0 and µ1 plotted as a function of sky-position angles θ and φ. (All panels display the
logarithmtobase10ofµ0−µ1.) ThetopleftpanelisforanetworkconsistingofthreeLIGOdetectors: twoatHanfordandone
at Livingston. The top right panel is for a network consisting of one LIGO detector at Livingston, one GEO and one VIRGO
detector. The first network (top left), which is more aligned, shows much larger values of µ0−µ1 than the second (top right)
which has more misalignement between thedetectors. In thebottom left panel, thepart of thesky whereµ0−µ1 for thefirst
network is ≤ 10 times that for the second network has been suppressed. It is clear that this is a smaller fraction of the sky
compared tothepart whereµ0−µ1 for thefirst networkis ≥10 timesthat ofthesecond network. Inthebottom right panel,
only thepart of thesky where µ0−µ1 for thesecond network is greater than or equalto that for thefirst network is shown.
VII. CONCLUSION sky position reconstruction accuracy for normal signals.
PMLR contains a free parameter, the penalty factor
There exists a serious problem, called the two detec- α. For α = 0 PMLR simply reduces to the SL method.
tor paradox, in the formal application of the standard However,method is notvery sensitive to the exactvalue
likelihood method to the detection of burst GW sig- ofα,andtheskymapofthelikelihoodconvergesrapidly
nals. We resolve this problem by imposing a restriction asαismadelarge. Theroleofαneedstobeinvestigated
on the detector responses allowed as candidate solutions further.
for maximizingthe likelihoodfunctional. The restriction Resultsfromnumericalsimulationsshowthatthe per-
is imposed through a penalizing function that measures formance of the PMLR method can be tuned, using α,
the variability of the SNR of a candidate solution when to vary over the range covered by the standard likeli-
the correspondingsourceis displacedslightly onthe sky. hood and the constraint likelihood methods introduced
HighSNRvariabilitysourcesleadtoverydifferentdetec- inKMRM.Forthe type ofsourcepopulationusedinthe
torresponsesevenwhenthedetectorsarecloselyaligned. simulations, a high value of α yields the same perfor-
The standard likelihood method allows such sources as mance as the hard constraint method. The ROC curves
perfectly valid candidate solutions. However, such ab- are mildly better compared to standard likelihood but
normalsourceshave an inherently reduced observability. the accuracyof source direction reconstructionis signifi-
The proposed method (PMLR) sacrifices the detectabil- cantly better.
ity of suuch sources while boosting the sensitivity and PMLRis alsorelatedtothe Bayesianapproach. Work