Tests of general relativity in the nonlinear regime: a parametrized plunge-merger-ringdown gravitational waveform model

The plunge-merger stage of the binary-black hole coalescence, when the bodies' velocities reach a large fraction of the speed of light and the gravitational-wave luminosity peaks, provides a unique opportunity to probe gravity in the dynamical and nonlinear regime. How much do the predictions of general relativity differ from the ones in other theories of gravity for this stage of the binary evolution? To address this question, we develop a parametrized waveform model, within the effective-one-body formalism, that allows for deviations from general relativity in the plunge-merger-ringdown stage. As first step, we focus on nonprecessing-spin, quasicircular black hole binaries. In comparison to previous works, for each gravitational wave mode, our model can modify, with respect to general-relativistic predictions, the instant at which the amplitude peaks, the instantaneous frequency at this time instant, and the value of the peak amplitude. We use this waveform model to explore several questions considering both synthetic-data injections and two gravitational wave signals. In particular, we find that deviations from the peak gravitational wave amplitude and instantaneous frequency can be constrained to about $20\%$ with GW150914. Alarmingly, we find that GW200129_065458 shows a strong violation of general relativity. We interpret this result as a false violation, either due to waveform systematics (mismodeling of spin precession) or due to data-quality issues depending on one's interpretation of this event. This illustrates the use of parametrized waveform models as tools to investigate systematic errors in plain general relativity. The results with GW200129_065458 also vividly demonstrate the importance of waveform systematics and of glitch mitigation procedures when interpreting tests of general relativity with current gravitational wave observations.

The plunge-merger stage of the binary-black-hole coalescence, when the bodies' velocities reach a large fraction of the speed of light and the gravitational-wave luminosity peaks, provides a unique opportunity to probe gravity in the dynamical and nonlinear regime. How much do the predictions of general relativity differ from the ones in other theories of gravity for this stage of the binary evolution? To address this question, we develop a parametrized waveform model, within the effective-one-body formalism, that allows for deviations from general relativity in the plunge-merger-ringdown stage. As first step, we focus on nonprecessing-spin, quasicircular black hole binaries. In comparison to previous works, for each gravitational wave mode, our model can modify, with respect to general-relativistic predictions, the instant at which the amplitude peaks, the instantaneous frequency at this time instant, and the value of the peak amplitude. We use this waveform model to explore several questions considering both synthetic-data injections and two gravitational wave signals. In particular, we find that deviations from the peak gravitational wave amplitude and instantaneous frequency can be constrained to about 20% with GW150914. Alarmingly, we find that GW200129 065458 shows a strong violation of general relativity. We interpret this result as a false violation, either due to waveform systematics (mismodeling of spin precession) or due to data-quality issues depending on one's interpretation of this event. This illustrates the use of parametrized waveform models as tools to investigate systematic errors in plain general relativity. The results with GW200129 065458 also vividly demonstrate the importance of waveform systematics and of glitch mitigation procedures when interpreting tests of general relativity with current gravitational wave observations.
Generally, tests of GR with GW observations have been developed following two strategies: theory independent and theory specific. The former assumes that the underlying GW signal is well-described by GR, and non-GR degrees of freedom (or parameters) are included to characterize any potential deviation. These tests use GW observations to check consistency with their nominal predictions in GR, and then constrain the non-GR parameters at a certain statistical level of confidence. Eventually, the non-GR parameters can be translated to the ones in specific modified theories of gravity, albeit there could be subtleties in doing it due to the choice of the priors and the actual parameters on which the measurements are done. By contrast, analyses that compare directly the data with proposed modified theories of gravity belong to the theory-specific framework of tests of GR.
Here, we focus on theory-independent tests of GR for BBHs. Historically, those tests have been proposed intro-ducing deviations in (or parametrizations of) the gravitational waveform, whether for the inspiral, the merger or the ringdown stages, in time or frequency domain. Those parametrizations are clearly not unique; neither they guarantee to fully represent the infinite space of modified gravity-theory waveforms. Furthermore, non-GR parameters may be degenerate with each other, limiting the study to a subset of them [9] or demanding the use of principal-component-analysis methods [25].
In this manuscript, we develop a parametrized time-arXiv:2212.09655v3 [gr-qc] 2 Aug 2023 domain IMR waveform model within the effective-onebody (EOB) formalism [52][53][54][55][56][57][58][59]. The EOB approach builds semianalytical IMR waveforms by combining analytical predictions for the inspiral [notably from post-Newtonian (PN), post-Minkowskian (PM), and gravitational self-force (GSF) approximations] and ringdown phases (from BH perturbation theory) with physically-motivated Ansätze for the plunge-merger stage. The EOB waveforms are then made highly accurate via a calibration to numerical relativity (NR) waveforms of BBHs. The EOB formalism relies on three key ingredients: the EOB conservative dynamics (i.e., a two-body Hamiltonian), the EOB radiation-reaction forces (i.e., the energy and angular momentum fluxes) and the EOB GW modes. Since the EOB waveforms are computed on the EOB dynamics by solving Hamilton's equations, in principle deviations from GR can be introduced in all the three building blocks, consistently. Here, for simplicity, following previous work [39,40] which focused on the ringdown stage, we introduce non-GR parameters in the plunge-merger-ringdown GW modes. We leave to future work the extension of the parametrization to the conservative and dissipative dynamics, notably by including in the EOB dynamics fractional deviations to the PN (as well as PM and GSF) terms, to NR-informed terms or specific new terms motivated by phenomena observed in modified gravity theories. We note that non-GR deviations in the EOB energy flux were implemented in Refs. [60,61], and the corresponding EOB waveforms were used in IMR consistency and other tests of gravity in Refs. [60][61][62].
Although the parametrized IMR model can in principle be constructed for precessing spinning BBHs, as first step, we consider nonprecessing BHs. There are two main EOB families, SEOBNR (e.g., see Refs. [63][64][65]) and TEOBResumS (e.g., see Refs. [66][67][68]). We consider here the former, and in particular we focus on the SEOBNRHM model developed in Refs. [63,64], which contains GW modes beyond the dominant quadrupole. We denote the parametrized version pSEOBNRHM. In Fig 1, we contrast a GR SEOBNRHM waveform with parameters similar to the first GW observation, GW150914, with a pSEOBNRHM waveform where the fractional deviations from GR are of the order of a few tens of percent. We can see that differences from GR occur just before, during, and after the merger stage, which is when the gravitational strain peaks.
The paper is organized as follows. In Sec. II, we describe how we build the pSEOBNRHM model starting from the baseline model SEOBNRHM, and introduce the non-GR parameters that describe potential deviations from GR during the plunge-merger-ringdown stage. In Sec. III, we study in detail the morphology of the parametrized waveform, and understand which parts of the waveform change when the non-GR parameters are varied one at the time. After discussing the basics of Bayesian analysis in Sec. IV, we perform a synthetic-signal injection study in Sec. V, and then apply our parametrized IMR model to real data in Secs. VI and VII, analyzing two events, GW150914 and GW200129. Finally, we summarize our conclusions and future work in Sec. VIII.
Unless stated otherwise, we work in geometrical units in which G = 1 = c.

II. THE PARAMETRIZED PLUNGE-MERGER-RINGDOWN WAVEFORM MODEL
In this section we first review the GR waveform model developed within the EOB formalism. In Sec. II B, we explain how we deform this baseline model by introducing deformations away from GR in the plunge-mergerringdown phase.
A. A brief review of the effective-one-body gravitational waveform model The GW signal produced by a spinning, nonprecessing, and quasicircular BBH with component masses m 1 and m 2 , and total mass M = m 1 + m 2 , is described in GR by a set of eleven parameters, ϑ GR , given by where χ i (i = 1, 2) are the constant-in-time projections of each BH's spin vectors S i in the direction of the unit vector perpendicular to the orbital planeL, i.e., χ i = S i ·L/m 2 i , where |χ i | ⩽ 1, (ι, ψ) describe the binary's orientation through the inclination and polarization angles, (α, δ) describe the sky location of the source in the detector frame, D L is the luminosity distance, and t c and ϕ c are the reference time and phase, respectively. It is convenient to define the chirp mass M = M ν 3/5 , where ν = m 1 m 2 /M 2 is the symmetric mass ratio, the asymmetric mass ratio q = m 2 /m 1 , and the effective spin χ eff = (χ 1 m 1 + χ 2 m 2 )/M . We adopt the convention that m 1 ⩾ m 2 and thus q ⩽ 1.
The GW polarizations can be written in the observer's frame as where φ 0 is the azimuthal direction of the observer, where, without loss of generality, we set φ 0 = ϕ c , and −2 Y ℓm are the −2 spin-weighted spherical harmonics [73], ℓ is the angular number and |m| ⩽ ℓ is the azimuthal number of each GW mode, h ℓm .
We follow Refs. [40,74] and use as our baseline model (i.e., the waveform model upon which the non-GR deviation parameters are added) the time-domain IMR waveform developed in Refs. [63,64,70] within the EOB  [64,69,70] (solid line) for a face-on, nonspinning and quasicircular binary with GW150914-like mass-ratio q = m2/m1 ≈ 0.867, and detector-frame total mass M = m1 + m2 = 71.9 M⊙. The pSEOBNRHM waveform is generated with non-GR parameters values δ∆t = −0.2, δω = −0.4, and δA = 0.5. These parameters change respectively, in comparison to GR, the instant at which the GW amplitude peaks, the orbital frequency at this time instant, and the value of the peak amplitude. Both waveforms are phase aligned and time shifted around 20 Hz using the prescription of Refs. [59,69,71,72]. The details of how the waveform model is developed are in given Sec. II, and additional details about its morphology are presented in Sec. III.
As explained in Refs. [63,64], the SEOBNRHM waveform is constructed by attaching the merger-ringdown waveform, h merger−RD ℓm (t), to the inspiral-plunge waveform, h insp−plunge ℓm (t), at a matching time t = t ℓm match , where Θ(t) is the Heaviside step function and the value of t ℓm match is defined as value. We impose that the amplitude and phase of h ℓm (t) at t = t ℓm match are C 1 (i.e., they are continuously and differentiable at this time instant). The time t 22 peak is defined as where t Ω peak is the time in which the EOB orbital frequency peaks [80]. Calculations performed in the test-particle limit using BH perturbation theory found that the amplitude and the orbital frequency peak at different times, especially when the central BH has large spins [81][82][83][84]. This motivates the introduction of the time-lag parameter ∆t 22 peak in Eq. (2.5), which can be fitted against NR waveforms as function of the symmetric mass ratio ν and the BH's spins χ 1,2 (see Sec. II B in Ref. [63] for details). We impose the condition ∆t 22 peak ⩽ 0 to ensure that the attachment of the merger-ringdown waveform happens before the peak of the orbital frequency, and thus before the end of the binary's dynamics. For later convenience, we define Because we are interested in adding non-GR terms to h merger−RD ℓm (t), we now briefly review how the mergerringdown waveform is constructed. Further details can be found in Sec. IV E of Ref. [64]. The merger-ringdown mode is written as where σ ℓm0 are the complex-valued frequencies of the least damped quasinormal mode (QNM) of the remnant BH [85][86][87]. We define σ R ℓm0 = Im(σ ℓm0 ) < 0 and σ I ℓm0 = −Re(σ ℓm0 ) < 0. The functionsÃ ℓm andφ ℓm are given by [63] A ℓm = c ℓm 1,c tanh c ℓm 1,f (t − t ℓm match ) + c ℓm 2,f + c ℓm 2,c , (2.8a) where ϕ ℓm match is the phase of the inspiral-plunge mode h insp−plunge ℓm at t = t ℓm match . We see that Eqs. (2.8) depend on the set of parameters c ℓm i and d ℓm i (i = 1, 2), which are either constrained by imposing thatÃ ℓm ,φ ℓm are C 1 at t = t ℓm match (we append the subscript "c") or free parameters to be determined by fitting against NR waveforms (we append the subscript "f ").
We now impose that h ℓm is C 1 at t = t ℓm match . This yields two equations that relate the constrained coefficients c ℓm 1,c and c ℓm 2,c to the free coefficients c ℓm 1,f , c ℓm 2,f , to σ R ℓm0 and to the mode amplitude of h insp−plunge ℓm and its first time derivative at the matching time, namely |h insp−plunge ℓm (t ℓm match )| and ∂ t |h insp−plunge ℓm (t ℓm match )|. The equations are: We also obtain one equation that relates the constrained parameter d ℓm 1,c to the free coefficients d ℓm 1,f , d ℓm 2,f , to σ I The values of at t = t ℓm match are fixed by the so-called nonquasicircular (NQC) terms, N ℓm (t). The NQC terms describe nonquasicircular corrections to the modes during the late inspiral and plunge. The NQCs are a parametrized time series that is multiplied with the factorized PN GR modes, h F ℓm , such that the resultant time series is calibrated against NR simulations. They are crucial in guaranteeing a very good agreement of the SEOBNRHM amplitude and phase (relative to NR) during the late inspiral and plunge.
The GW modes in the inspiral-plunge part of the EOB waveform are given as where we refer the reader to Sec. IV C in Ref. [64] for details on how h F ℓm and N ℓm are constructed. For our purposes, it is sufficient to say that |h insp−plunge ℓm |(t ℓm match ), ∂ t |h insp−plunge ℓm (t ℓm match )|, and ω insp−plunge ℓm (t ℓm match ) are the same as the NR values of |h NR ℓm |, ∂ t |h NR ℓm |, and ω NR ℓm , at t = t ℓm match . The values of these three quantities are obtained for each BBH, from the Simulating eXtreme Spacetimes (SXS) catalog of NR waveforms [88], after which a fitting formula that depends on the symmetric mass ratio ν and spins χ 1 and χ 2 is obtained to interpolate over the parameter space covered by the catalog. Their explicit forms can be found in Ref. [64], Appendix B. At this point, we are left with the free parameters c ℓm i,f and d ℓm i,f (i = 1, 2) to fix. This is accomplished through fits against NR and Teukolsky equation-based waveforms [82,83], written also as functions of ν, χ 1 and χ 2 . The explicit form of these fits can be found in Ref. [64], Appendix C.

B. Construction of the parametrized model
With this framework established, our strategy to develop a parametrized SEOBNRHM model (hereafter pSEOBNRHM) is the following. We will introduce fractional deviations to the NR-informed formulas for the mode amplitudes and angular frequencies at t = t ℓm match , i.e., and we will also allow for changes to t ℓm match by modifying the time-lag parameter ∆t GR ℓm [defined in Eq. (2.6)] as, where we constrain δ∆t ℓm > −1 to ensure that t ℓm match remains less than t Ω peak , and thus before the end of the dynamics, as originally required [64,69]. Equations (2.12) and (2.13) modify the constrained parameters c ℓm i,c and d ℓm i,c through Eqs. (2.9)-(2.10), and consequentlyÃ ℓm and ϕ ℓm that appear in the merger-ringdown waveform (2.7) and are given by Eqs. (2.8). It is important to emphasize that Eqs. (2.12) and (2.13) also modify the NQC coefficients which enter the inspiral-plunge waveform in Eq. (2.11). This is because both |h NR ℓm | and ω NR ℓm are used to fix some parameters in the explicit form of N ℓm . We refer the reader to Refs. [63,89] and in particular to  [40], while the merger deformation parameters δA ℓm , δω ℓm , and δ∆t ℓm are introduced here for the first time. As explained in Sec. II B, although we call these merger parameters, they do also affect the late inspiral-plunge part of the waveform. We quote under the column labeled "bound" the constraints on the parameter's values required by our waveform model.
We also introduce non-GR deformations to the QNMs, following the same strategy applied in Refs. [33,34,39,40,43,90]. It consists in modifying the QNM oscillation frequency and damping time, defined respectively for the zero overtone n = 0, as, according to the substitutions and we impose that δτ ℓm0 > −1 to ensure that the remnant BH is stable (i.e., it rings downs, instead of "ringingup" exponentially). Note that in Refs. [39,40], such deformations also concerned with the higher overtones, since the EOB model used for the merger-ringdown included higher overtones. Put it all together, we have the following set of plungemerger-ringdown parameters: intended to capture possible signatures of beyond-GR physics in the most dynamical and nonlinear stage of a BBH coalescence. We will casually refer to them as "non-GR" or as "deformation" (away from GR) parameters. In Table I, we summarize the ϑ nGR parameters, their meaning, and the constraints, if any, on their values. The GR limit is recovered when all parameters in ϑ nGR are set to zero.
The pSEOBNRHM model allows us to change the non-GR plunge-merger parameters ϑ merger nGR for each (ℓ, m) mode individually. Here, for a first study, we will assume that their values are the same across different modes, that is to say, δA ℓm = δA, δω ℓm = δω, and δ∆t ℓm = δ∆t , (2.17) for all the ℓ and m modes in the waveform model. This choice is motivated by the fact that in GW150914 there are no significant changes in the posterior distributions of the binary parameters when using all the modes and only the ℓ = m = 2 mode. As for the non-GR ringdown parameters ϑ RD nGR , we will assume that they are nonzero only for the least-damped (n = 0) (2, 2) mode. Under these assumptions, we have a 16-dimensional parameter space to work with, where the GR parameters ϑ GR are defined in Eq. (2.1). Some comments follow in order. First, the parametrized deformation of SEOBNRHM we have introduced is not unique. For instance, we could have added additional fractional changes to ∂ t |h NR ℓm |(t ℓm match ) or to the free parameters in the merger-ringdown waveform segment [see Eq. (2.8)]. We have found a compromise between the number of new parameters we can introduce and the physics we want to model; the optimal scenario being that of having the most flexible GW model that depends on the least number of deviation parameters. In our case, we find the parameters ϑ nGR defined in Eq. (2.16) to be sufficient for our purposes. Second, one may fear that by effectively "undoing" the NR calibration we would obtain nonphysical GWs. This is not the case, as shown in Fig. 1 and as we will see in Sec. III. Our model produces waveforms that are smooth deformations of the ones of GR and have sufficient flexibility to be applied in tests of GR (Secs. V and VI) and provide a diagnostic tool for the presence of systematic effects in GR GW models (Sec. VII).

III. WAVEFORM MORPHOLOGY
Having introduced our waveform model, we now discuss how each of the parameters ϑ merger nGR modify the GW signal in GR. An analogous exploration was done for ϑ RD nGR in Ref. [40], for this reason the present discussion is restricted to the merger parameters. In each of the following sections, we vary the parameters δA, δω, and δ∆t one at a time. We take the binary component masses and spins to be which are archetypal values of a GW150914-like event [91], the inclination to be ι = 0 and, for clarity, we show results only for h 22 . This is the dominant mode for such a quasicircular, nonspinning, and comparable-mass BBH. We end each section by showing how the waveform is modified when we apply the deformations, with the same values, simultaneously to all GW modes present in pSEOBNRHM. A. The amplitude parameter δA Let us start with δA, the amplitude parameter. In Fig. 2 we show the real part of h 22 (t), rescaled by the luminosity distance D L and total mass M , for two values of δA: 0.5 (top panel) and −0.5 (bottom panel). The dashed segment corresponds to t ⩽ t 22 match (i.e., the inspiral-plunge part of waveform), whereas the solid segment corresponds to t > t 22 match (i.e., the merger-ringdown part of the waveform). In both panels, the black curve corresponds to the GR signal (δA = 0) with the same binary parameters. Both the GR and non-GR waveforms have been shifted in time and aligned in phase around 20 Hz following the prescription of Refs. [59,69,71,72]. The amplitudes of the non-GR waveforms ±|h 22 | are shown by the dotted lines and form the envelope around Re(h 22 ).
Unsurprisingly, for positive values of δA, the amplitude |h 22 | increases relative to its GR value while keeping t 22 match ≈ 1704 M the same. The situation is more interesting for δA < 0. For the binary under consideration, we find that |h 22 | decreases for δA ≳ −0.31, but for δA ≲ −0.31, we see that δA pinches downwards the amplitude enough to result in a local minimum (which we will refer to as t 22 min ) and two maxima, located before and after t 22 min , with the global maximum happening at t 22 max < t 22 min . The values of both maxima are smaller than the GR peak amplitude. By construction, the matching time t 22 match is then shifted to earlier times relative to its GR value. We now consider δω, the frequency parameter. Figure 4 is analogous to Fig. 2, except that we now consider δω = 0.5 (top panel) and δω = −0.5 (bottom panel). We see that δω induces a time-dependent phase shift to the waveform, with its effects being most noticeable near the merger, and causing t match to happen later (earlier) relative to GR when δω > 0 (δω < 0), while keeping the peak amplitude unaffected.
In Fig. 5, we show an analogous version of Fig. 3, but now for δω. Once more, the top panel shows the real part of the strain, the middle panel the strain amplitude, and the bottom panel the instantaneous frequency. We focus on the region near the merger and we plot the GR curves (δω = 0) with black solid lines. In the top panel, we can see the phase differences between the non-GR and GR waveforms, which are the largest around the time of merger and ringdown. This is in part due to the δω itself, but also to the phase-shift and time-alignment procedure already mentioned, which we perform with respect to the GR waveform. The effect of the latter is small, as can be seen in the middle panel for the amplitude, where all curves nearly overlap in time. In the bottom panel, we note sharp changes to f when |δω| ≈ 0.5. They originate from us not imposing the continuity of the time derivative of ω NR ℓm at t = t ℓm match [63,64]. In the top panel, we clearly see the phase difference between the non-GR and GR waveform near the merger. This is partially due to the δω itself, but also to the phase-shift and time-alignment done with respect to the GR waveform. The effect of the latter is small as can be seen in the middle panel, which shows the amplitude. The sharp changes to f in the bottom panel for |δω| ≈ 0.5 originate from us not imposing the continuity of the time derivative of ω NR ℓm at t = t ℓm match .
around the time of merger, obtained by finely covering the interval δ∆t ∈ [−0.5, 0.5]. The GR prediction is shown by the black solid line. We see that the changes to the strain (top panel), its amplitude (middle panel), and its frequency evolution (bottom panel) are small. Therefore, δ∆t introduces changes to the GR waveform which are in general subdominant relative to those due to δA and δω. We also remark that ∆t GR ℓm is not very sensitive to the EOB calibration against NR waveform. Hence, the fractional changes we are introducing on ∆t GR ℓm are comparable with the NR fitting errors. This explains why this parameter affects the GR waveforms so little.

IV. PARAMETER ESTIMATION
In the previous section, we have introduced our waveform model and discussed the properties of the waveform morphology. Here, we summarize the Bayesian inference formalism used for parameter estimation of GW signals and synthetic-data studies. We describe the prior choices and the criteria for the GW event selection.

A. Bayesian parameter estimation
Our hypothesis, H, is that in the detector data, d, an observed GW signal is described by the waveform model pSEOBNRHM. The model pSEOBNRHM has a set of GR and non-GR parameters, as in Eqs. (2.1) and (2.16), where ϑ nGR = {δA, δω, δ∆t, δf 220 , δτ 220 } . (4.1) As said, we assume that the merger modifications are the same for all (ℓ, m) modes present in the model pSEOBNRHM.
The posterior probability distribution on the parameters of the model, ϑ, given the hypothesis, H, is obtained using Bayes' theorem, where P (ϑ|H) is the prior probability distribution, P (d|ϑ, H) is the likelihood function, and P (d|H) is the evidence of the hypothesis H. For a detector with stationary, Gaussian noise and power spectral density S n (f ), the likelihood function can be written as where the noise-weighted inner product is defined as whereÃ(f ) is the Fourier transform of A(t), and the asterisk denotes the complex conjugation, and S n (f ) is the one-sided power spectral density of the detector. The integration limits f low and f high set the bandwidth of the detector's sensitivity. We follow the LVK analysis and set f low = 20 Hz, while f high is the Nyquist frequency [19]. The posterior distributions are computed by using LALInferenceMCMC [92,93], a Markovchain Monte Carlo that uses the Metropolis-Hastings algorithm to survey the likelihood surface and is implemented in LALInference [94], part of the LALSuite software suite [95].

B. Prior choices
The prior distributions on the GR parameters are assumed to be uniform in the component masses (m 1 , m 2 ), uniform and isotropic in the spin magnitudes (χ 1 , χ 2 ), isotropic on the binary orientation, and isotropically distributed on a sphere for the source location with P (D L ) ∝ D 2 L . For the non-GR parameters, as explained in Sec. II B, the internal consistency of the pSEOBNRHM model requires that both δ∆t and δτ 220 are larger than −1 (cf. Table I). We use this fact to fix a common lower limit on the uniform priors on all ϑ nGR . We set 1 to be an upper limit on the uniform priors on the non-GR parameters. This was sufficient in most of our analysis, but in a few cases we found that the marginalized posteriors distributions for one or more non-GR parameters had support at ϑ nGR ≈ 1. In such cases we extended the priors' domains to ϑ nGR ∈ (−1, +2]. Even at this wider range, we did not find anomalies in the waveform.

C. Event selection
The pSEOBNRHM ringdown analysis performed in Ref. [13] selected GW events from the GWTC-3 catalog [19] which had a signal-to-noise ratio (SNR) ⩾ 8 in the inspiral and post-inspiral regimes. The requirement on the inspiral regime allows one to break the strong degeneracy between the total mass of the binary and the ringdown deviation parameters [39,40]. Among the GW events that meet this criteria, two stand out in terms of their constraining power on ϑ RD nGR , namely GW150914 [14,15] and GW200129 065458 (hereafter GW200129) [19]. These two events, with a median total source-frame masses of 64.5 M ⊙ and 63.4 M ⊙ , respectively, are among the loudest BBH signals to date with a median total network SNR of 26.0 and 26.8, respectively [18,19]. GW150914 was detected by the two LIGO detectors at Hanford and Livingston, whereas GW200129 was detected by the threedetector network of LIGO Hanford, Livingston, and Virgo.
We guide ourselves by this result and use these two events to investigate what constraints we can place on the merger-ringdown parameters. We remark that this SNR selection criteria may be too strong if we are interested in ϑ merger nGR only. We leave the study of the optimal SNR to constrain only the merger parameters to a future work.

V. RESULTS: SYNTHETIC-SIGNAL INJECTION STUDIES
In this section, we use pSEOBNRHM to perform a number of synthetic-signal injection studies. As we saw in Sec. II, pSEOBNRHM is a smooth deformation of the GR waveform model SEOBNRHM, which is recovered when all ϑ nGR parameters are set to zero. This allows us to explore different scenarios that differ from one another on whether the GW signal and the GW model used to infer the parameters of this signal are described by GR (ϑ nGR = 0) or not (ϑ nGR ̸ = 0). We summarize these possibilities in Table III.
To prepare the GW signal we need to fix ϑ = ϑ GR ∪ ϑ nGR . In all cases, we use values of ϑ GR illustrative of a GW150914-like BBH as in Table II. We set all non-GR parameters to the same value, ϑ nGR = 0.1, whenever the injected signal is non-GR. By working exclusively with the pSEOBNRHM waveform model, we avoid introducing systematic errors due to waveform modeling in our analysis. We also employ an averaged (zero-noise) realization of the noise to avoid statistical errors due to noise. The resulting GW signal is then analyzed with the power spectral density S n (f ) of the LIGO Hanford and Livingston detectors both at design sensitivity [96]. In all cases, we set the distance to the binary to be such that the total network SNR is approximately 100.
In Sec. V A, we do a preliminary analysis where both injected and model waveforms are described by GR. This allows us to access the accuracy with which different binary parameters can be recovered from the data in the detector network. With these results as a benchmark, we can then proceed to inject a non-GR waveform and analyze it with a GR model. This allows us to study the systematic error introduced on the inferred binary parameters by assuming a priori that GR is true, while nature may not be so (the so-called fundamental bias). In Sec. V B, we inject a GR waveform and try to recover its parameters with a non-GR model. This allows us to answer how much the non-GR parameters can be constrained given an event consistent with GR. Finally, in Sec. V C, we use non-GR waveforms as both our injection and our model. This answers whether we can detect the presence of the non-GR parameters in our signal.

Model GR
non-GR Injection GR Sec. V A Sec. V B non-GR Sec. V A Sec. V C TABLE III. Summary of the synthetic-signal injection simulations performed in Sec. V. The label "GR " refers to the SEOBNRHM waveform model, whereas the label "non-GR" refers to the pSEOBNRHM waveform model, where all merger-ringdown parameters are set deviate in 10% deviations relative to their corresponding GR values.

A. Fundamental biases on binary parameters
We first explore the presence (or not) of biases in the inference of binary parameters when the template waveform model assumes GR, while the injected GW signal is non-GR [31,97]. For this purpose, we first inject a synthetic GR GW signal with SNR = 98 and recover the binary parameters with a GR model. By doing this exercise first, we gain an idea on the accuracy with which the parameters of the binary (cf. Table II) can be recovered in our set up. Next, we repeat the same analysis but now using as our synthetic GW signal the one obtained with pSEOBNRHM. The signal is prepared using the same binary parameters ϑ GR shown in Table II with SNR = 104, but now we let ϑ nGR = 0.1.
The results of our two analyses are shown in Fig. 8. We show the one-and two-dimensional posterior distributions of a subset of the intrinsic binary parameters, namely, the mass ratio q, the detector-frame chirp mass M and the effective spin χ eff . In all panels, the "true" (injection) values of these parameters are marked by the vertical and horizon lines. We see that in the case of a non-GR injection (solid curves), the posterior distributions of the parameters are shifted from the injected values and from the posterior distributions in the case of a GR injection (dashed curves). We attribute the differences in the 90% contours of the posterior distributions to the fact that in the non-GR injection a smaller value of the chirp mass M is inferred. This suggests that the GR waveform that best fits the data has a longer inspiral and this makes the inference of the other binary parameters more precise. The recovered SNR from the GR analysis of the non-GR signal is almost the same as the injected one, i.e., SNR = 104. Hence, if a GW signal with deviations from GR would be analyzed by current GR templates, the GW event would be interpreted as a BBH in GR with different values of the binary parameters.

B. Constraints on deviations to general relativity
We now inject a synthetic GW signal in GR using the parameters ϑ GR in Table II with SNR = 98. We analyze the signal using the pSEOBNRHM waveform model, allowing both ϑ GR in Eq. (2.1) and ϑ nGR in Eq. simulates a scenario where we have a GW event consistent with GR and we want to understand which constraints it places on the non-GR parameters in our waveform model. We summarize the results of the analysis in Fig. 9, where we show the one-and two-dimensional posterior probability distributions of the merger-ringdown parameters ϑ nGR . We find that the marginalized posterior distributions of the non-GR parameters are consistent with the corresponding injected values in GR, which are indicated by the markers. We can infer that a GW150914-like event with SNR = 98 would constrain the deformation parameters in the range between 5% (for δA and δf 220 ) and 20% (for δτ 220 ) at 90% credible level. In Appendix A, Fig. 15, we show the posterior distributions on the intrinsic binary parameters.
The best constrained parameter is the amplitude, δA, whereas the less constrained parameter is the time shift, δ∆t. For the latter, we obtain a posterior distribution that has support onto a wide range of the prior. This is perhaps unsurprising due to the small deviations caused by δ∆t in the waveform in comparison with δω (compare Figs. 4  and 6). We also observe a negative correlation between these two parameters and hence increasing precision on one is likely to increase uncertainty on the other; see the δ∆t-δω panel in Fig. 9. Together, these results suggest that considering δA and δω is sufficient, if one is interested in doing a test of GR only in the plunge-merger stage of the binary's coalescence.

C. Detecting deviations from general relativity
We now study whether we can detect the presence of the non-GR parameters. To do so, we inject a synthetic GW signal where the binary parameters are shown in Table II, SNR = 104, and we set the merger-ringdown parameters to be 10% larger than their corresponding GR values.
We summarize the outcome of our parameter estimation in Fig. 10, where we show the one-and two-dimensional posterior distributions for the ϑ nGR parameters. We see that all posteriors are consistent with the injected values, indicated by the markers. Moreover, the posteriors for ϑ nGR have support at their null, GR value. The exceptions are the amplitude δA and the QNM frequency δf 220 parameters, which have no support at their GR values at 90% credible level. This suggests that these two parameters are the most promising ones in signaling the presence of beyond-GR physics for GW150914-like binaries. In fact, we will see this suggestion taking place in our analysis of GW200129 in Sec. VII.

VI. ANALYSIS OF GW150914: CONSTRAINTS ON THE PLUNGE-MERGER-RINGDOWN PARAMETERS
Having gained some intuition on the role of the mergerringdown parameters in the synthetic-signal injections presented in Sec. V, we now apply the pSEOBNRHM model to the analysis of real GW events. Our analysis, here and in Sec. VII, uses the power spectral density of the detectors from the Gravitational Wave Open Science Center (GWOSC) [98], and calibration envelopes as used for the analyses in Ref. [13]. We will start with GW1501914, the first GW event observed by the LIGO-Virgo Collaboration [14].
We will focus our analysis to two subsets of mergerringdown parameters due to the smaller SNR of this event (and of GW200129) in comparison to the SNR ≈ 100 scenarios studied in the previous section. First, we have seen that the time-shift parameter δ∆t is the hardest parameter to constrain, and that it has wide posteriors even at such large SNRs. This motivates us to consider, among the merger parameters, only to perform a "merger test of GR". Second, we performed a parameter estimation of GW150914, using all ϑ nGR parameters in Eq. (4.1). We found correlations between the frequency parameter δω and the QNM deformations parameters δf 220 and δτ 220 . Moreover, we also did a series of synthetic-signal injection studies using the binary parameters listed in Table II, with SNR= 26, and in Gaussian noise. In some of these cases, we also found correlations between δω and δf 220 and δτ 220 . In summary, these correlations arise either when the GW event has low SNR or due to noise. This suggest using ϑ nGR = {δA, δf 220 , δτ 220 } , (6.2) to perform a "merger-ringdown test of GR". In Fig. 11 we show the results of our merger test of GR. The corner plot shows the one-and two-dimensional posterior probability distributions of δA and δω. The posterior distributions are consistent with the null value predicted in GR. We obtain from GW150914, δA = −0.01 +0. 27 −0.19 , and δω = 0.00 +0.17 −0.12 , (6.3) at 90% credible level. This shows that we already constrain deviations from GR around the merger time of BBH coalescences to about 20% with present GW events. Figure 12 is a similar plot, but for the merger-ringdown test of GR. Once more, we find that the inferred values of the non-GR parameters are consistent with GR, δA = 0.03 +0. 29 −0.20 , δf 220 = 0.041 +0.151 −0.084 , δτ 220 = 0.04 +0.27 −0.29 , (6.4) at 90% credible level. The bound on the amplitude parameter is similar to the one obtained in the merger test, shown in Eq. (6.3). Also, the bounds on the ringdown parameters are similar to those obtain in Ref. [40] (δf 220 = 0.05 +0.11 −0.07 and δτ 220 = −0.07 +0.26 −0.23 ), which had only these two quantities as its non-GR parameters. In Appendix A, Fig. 16, we show the posterior distributions on the intrinsic binary parameters for both tests of GR.
When interpreting our inferences on these parameters, it is important to note that the statistical error in our analysis (≈ 20%) is larger than the systematic error due to fitting |h NR ℓm | and ω NR ℓm against NR data, which is at most around 4% with current models [63,64], depending on where one is in the η-χ eff parameter space. In fact, we see that the median values of δA and δω fall within this fitting error. In conclusion, we can claim to have placed a constraint on these non-GR parameters with GW150914.

VII. THE CASE OF GW200129: THE IMPORTANCE OF WAVEFORM SYSTEMATICS AND DATA-QUALITY IN TESTS OF GENERAL RELATIVITY
We now turn our attention to GW200129 and, following what we have learned in the previous section, we first consider pSEOBNRHM with only δA and δω as non-GR parameters. We show the one-and two-dimensional marginalized posteriors of these parameters with the black solid curves in the left panel of Fig. 13. We see that while our inferred value of δω (δω = −0.002 +0.097 −0.082 at the 90% credible level) is consistent with GR, our inferred value of δA (δA = 0.44 +0. 38 −0.28 at the 90% credible level) exhibits a gross violation of GR.
Have we found a strong evidence of violation of GR in GW200129? Assuming that this is not the case, the apparent violation of GR could be either due to statistical errors or to systematic errors. To explore the first possibility, we perform a series of synthetic-data injection studies. As our first step, we do a parameter-estimation study in zero noise, where the injected GW signal is generated with SEOBNRHM and we use the binary parameters corresponding to the maximum likelihood point from the GWTC-3 data release by the LVK [99] analysis of GW200129. The LVK analysis was done separately with two quasicircular and spin-precessing waveform models, SEOBNRv4PHM [65] and IMRPhenomXPHM [100], employing different parameter estimation libraries, RIFT [101][102][103] and Bilby [104,105], respectively. Here, as a reference, we use the maximum likelihood point of the analysis that employed the IMRPhenomXPHM model, and we expect the results to be qualitatively similar had we used SEOBNRv4PHM. More specifically, because the SEOBNRHM model we are using is nonprecessing, we use only the masses and luminosity distance from the maximum-likelihood point.
ingly, consistent with GR. We also repeat this analysis for ten Gaussian noise realizations, using the same synthetic GW signal (yellow solid curves in the left panel of Fig. 13). Consistent with the expectations, two noise realizations yield marginalized posteriors on δω and δA which are not consistent with GR at 90% credible level (shown by the thicker yellow solid curves). It is worth observing how the Gaussian noise curves have qualitatively the same shapes (spreads), with the two outliers being shifted away from (δA, δω) = (0, 0). This is an expected behavior consistent with the stationary, Gaussian assumption of statistical noise. These results, hence, disfavor the possibility that the violations of GR we are observing are due to Gaussian noise or due to the particular binary parameters inferred for this event. The latter alternative would have been quite unlikely in the first place, because both GW200129 and GW150914 have similar binary parameters and SNRs, and we have already found that GW1501914 is consistent with GR in Sec. VI (cf. Fig. 11). As our next step, we perform two additional parameter estimation runs, in zero noise, but now generating our synthetic GW signal with the SEOBNRv4PHM [65] and the NRSur7dq4 [106,107] waveform models. Both models allow for spin precession, unlike our pSEOBNRHM. Hence, we can study if the GR deviations we are finding are due to systematic errors in the GW modeling. Once again, the maximum-likelihood point of the LVK analysis of GW200129 using IMRPhenomXPHM was used, but this time with the binary in-plane spin components included. We show our results in the right panel of Fig. 13. The one-and two-dimensional posterior distributions of δω and δA are shown in dash-dotted curves for the NRSur7dq4 injection and with dotted curves for the SEOBNRv4PHM injection. For reference, we also include the posterior distribution associated to the SEOBNRHM injection (dashed curves) and to the data from GW200129 (solid curves). We see that these two spin-precessing GW signals, when analyzed in zero noise, are also in disagreement with GR, when analyzed with our nonprecessing non-GR model. We also see that our results using NRSur7dq4 (which compares the best against NR simulations in its regime of validity) are in good agreement with what we obtain by analyzing the GW200129 data. These results, compared with those obtained from the SEOBNRHM injection, suggest that the presence of spin precession in the GR signal, biases us to find a false evidence for beyond-GR effects when we use a nonprecessing non-GR model.
Is this the full story? In Ref. [108], Payne et al. revisited the evidence of spin precession in GW200129 [109]. They concluded that the evidence for spin-precession originates from the LIGO Livingston data, in the 20-50 Hz frequency range, alone. This range coincides with the frequency range that displays data quality issues, due to a glitch in the detector that overlapped in time with the signal [19]. By reanalyzing the GW200129 data with f low > 50 Hz (while leaving LIGO Hanford data intact and not using Virgo data), they showed that the evidence in favor of spin precession in this event disappears. See Ref. [108] for a detailed discussion. Moreover, a reanalysis of the LIGO Livingston glitch mitigation showed that the difference between the spin-precessing and nonprecessing interpretations of this event is subdominant relative to uncertainties in the glitch subtraction [108]. Since we have used the glitch-subtracted data in our parameter estimation, we are then led to the second conclusion of our study of this event, namely that: issues with data quality can introduce biases in non-GR parameters, to an extent that one can find significant false violations of GR in GW events detected with present GW observatories. See Ref. [110] for a recent study of this issue.
Furthermore, we repeat here the analysis we have performed for GW150914 where we considered ϑ nGR = {δA, δf 220 , δτ 220 } as our non-GR parameters. For the discussion that follows, we assume that GW200129 is an unmistakable genuine spin-precessing BBH. We show our results in Fig. 14. We see that while our inferred values of δf 220 and δτ 220 are consistent with GR at 90% confidence level, our inference of the amplitude parameter, δA = 0.50 +0. 46 −0.30 at 90% credible level, remains inconsistent with GR. Moreover, this value hardly changes from our {δA, δω}-study, i.e., δA = 0.44 +0. 38 −0.28 , at the same credible level.
This result is interesting for two reasons. First, it indicates that the systematic error caused by spin-precession mismodeling is robust to the inclusion of deformations For the latter, we used the maximum-likelihood point of LVK's original analysis of GW200129 which employed the IMRPhenomXPHM model to generate the synthetic GW signal. We performed the parameter estimation of these injections in zero noise (dashed curves) and in ten Gaussian noise realizations (yellow solid curves). Right panel: Similar, but having generated two additional GW200129-like synthetic signals with NRSur7dq4 (dot-dashed curves) and with SEOBNRv4PHM (dotted curves). Both models include spin precession effects. Observe how the posteriors distributions are in tension with GR [marker at (δA, δω) = (0, 0)] when we include spin-precession effects in the synthetic data and we recover with a nonprecessing and non-GR waveform model. to the ringdown QNM frequencies, at least for this event. Second, there is a commonality between our finding for GW150914 (see Fig. 12) and GW200129 (see Fig. 14) namely, that in both cases the posterior distributions of δf 220 and δτ 220 are consistent with GR, despite the larger parameter space due to the inclusion of δA. In the case in which one considers only δf 220 and δτ 220 as non-GR parameter, the consistency with GR had already been established in Ref. [40], and in particular in Ref. [13]; see Sec. VIII, Fig. 14 there. 1 Our analysis of these two GW events with the new pSEOBNRHM waveform model suggests the following: the model would be able to detect deviations from nonprecessing quasicircular GW signals in the plunge-merger-ringdown which otherwise would not be seen when having deformations to the ringdown only.
We close our discussion of GW200129 with two remarks. First, data-quality issues aside, we can think of our spin-precessing injection studies as illustrative of what could happen in upcoming LVK observation runs. 1 The LVK Collaboration also does an independent analysis of the ringdown using pyRing. This analysis lead to an odds ratio log 10 O nGR GR = −0.09 for GW200129, the largest among all events studied [11]. A positive value would quantify the level of disagreement with GR.
By doing so, we have then demonstrated the existence of a systematic error on the non-GR parameters caused by spin-precession mismodeling. 2 Second, although we have proposed pSEOBNRHM as a means of constraining (or detecting) potential non-GR physics in BBH coalescences, we can also interpret the merger parameters as indicators of our ignorance in GR waveform modeling. 3 More concretely, in a hypothetical scenario where GW modelers did not know that BBH can spin precess, an analysis of GW200129 with pSEOBNRHM would suggest that their model of the peak GW-mode amplitudes is insufficient to describe this event and hence be an indicative of new, nonmodeled binary dynamics that was absent in their waveform model. They would not be able to say that spin precession is the missing dynamics, but they would at least realize that something is missing. FIG. 14. The one-and two-dimensional posterior distribution functions for δA, δf220, and δτ220 for GW200129. All contours indicate 90% credible regions. We see that while our inferred values for δf220 and δτ220 are consistent with GR, δA is not.

VIII. DISCUSSIONS AND FINAL REMARKS
We presented a time-domain IMR waveform model that accommodates parametrized deviations from GR in the plunge-merger-ringdown stage of nonprecessing and quasicircular BBHs. This model generalizes the previous iterations of the pSEOBNRHM model [38][39][40], which included deviations from GR in the inspiral phase or modified the QNM frequencies only, by introducing deformations parameters ϑ merger nGR that, for each GW mode, can change the time at which the GW mode peaks, the mode frequency at this instant, and the peak mode amplitude. This new version of pSEOBNRHM reduces to the state-ofthe-art SEOBNRHM model [63,64,70] for nonprecessing and quasicircular BBHs in the limit in which all deformations parameters are set to zero.
We used pSEOBNRHM to perform a series of injections studies for GW150914-like events exploring (i) the constraints that one could place on these non-GR parameters, (ii) the biases introduced on the intrinsic binary parameters in case nature is not described by GR and we model the signal with a GR template, and, finally, (iii) we studied the measurability of these non-GR parameters.
We also used pSEOBNRHM in a reanalysis of GW150914 and GW200129. For GW150914, we found that the deviations from the GR peak amplitude and the instantaneous GW frequency can already be constrained to about 20% at 90% credible level. For GW200129, we found an interesting interplay between spin precession and false violations of GR that manifests as a ∼ 2σ deviation from GR in the peak amplitude parameter. By interpreting the evidence for spin precession in this event as due to data-quality issues in the LIGO Livingston detector [19,108], we found a further a connection between data-quality issues and false violations of GR [110].
These results warrant further studies on the systematic bias due to spin precession in tests of GR. In the context of plunge-merger-ringdown test, this could be achieved by extending the SEOBNRv4PHM waveform model [65] to include the same set of non-GR parameters ϑ nGR used here. It is also natural to explore which systematic effects higher GW modes [111] and binary eccentricity can introduce in tests of GR. For the latter, see Ref. [112] for work in this direction for IMR consistency tests [60,113] and Ref. [114] in the context of deviations in the post-Newtonian (PN) GW phasing [31,36,115]. It would also be interesting to investigate these issues in the context of the ringdown test within the EOB framework employed by LVK Collaboration [13] and which relies on pSEOBNRHM [39,40]. This could be done by adding non-GR deformations to the SEOBNRv4EHM waveform model of Ref. [116]. It would also be important to investigate whether pSEOBNRHM can be used to detect signatures of non-GR physics, as predicted by the rapidly growing field of NR in modified gravity theories (see e.g., Refs. [117][118][119][120][121][122][123][124][125][126]); some of which predict nonperturbative departures from GR only in lateinspiral and merger ringdown [127][128][129][130]. One could also study what the theory-agnostic bounds we obtained with GW150914 on the amplitude and GW frequency imply to the free parameters of various modified gravity theories.
The deformations parameters ϑ merger nGR in our pSEOBNRHM model should have an approximate correspondence to the phenomenological deviation parameters (from NR calibrated values) in the "intermediate region" of the IMRPhenom waveform model used in the TIGER pipeline [35][36][37] of the LVK Collaboration [9][10][11][12]. Such a mapping could be derived through synthetic injection studies. This work only introduced non-GR parameters in the EOB GW modes and only during the plunge-mergerringdown. Importantly, and more consistently, in the near future we will extend the parametrization to the EOB conservative and dissipative dynamics.
The interplay between GW waveform systematics, characterization and subtraction of nontransient Gaussian noises in GW detectors, and non-GR physics will become increasingly important in the future. Planned groundbased [131,132] and space-borne GW observatories [133] will detect GW transients with SNRs that may reach the thousands depending on the source. Having all these aspects under control is a daunting task that will need to be faced if one wants to confidently answer the question "Is Einstein still right?" [134] in the stage of BBH coalescences where his theory unveils its most outlandish aspects.

ACKNOWLEDGMENTS
We thank Héctor Estellés, Ajit Kumar Mehta, Deyan Mihaylov, Serguei Ossokine, Harald Pfeiffer, Lorenzo Pompili, Antoni Ramos-Buades, and Helvi Witek for discussions. We also thank Nathan Johnson-McDaniel, Juan Calderón Bustillo, and Gregorio Carullo for comments on this work. We acknowledge funding from the Deutsche Forschungsgemeinschaft (DFG) -Project No. 386119226. We also acknowledge the computational resources provided by the Max Planck Institute for Gravitational Physics (Albert Einstein Institute), Potsdam, in particular, the Hypatia cluster. The material presented in this paper is based upon work supported by National Science Foundation's (NSF) LIGO Laboratory, which is a major facility fully funded by the NSF. This research has made use of data or software obtained from the Gravitational Wave Open Science Center (gwosc.org), a service of LIGO Laboratory, the LIGO Scientific Collaboration, the Virgo Collaboration, and KAGRA. LIGO Laboratory and Advanced LIGO are funded by the United States National Science Foundation (NSF) as well as the Science and Technology Facilities Council (STFC) of the United Kingdom, the Max-Planck-Society (MPS), and the State of Niedersachsen/Germany for support of the construction of Advanced LIGO and construction and operation of the GEO600 detector. Additional support for Advanced LIGO was provided by the Australian Research Council. Virgo is funded, through the European Gravitational Observatory (EGO), by the French Centre National de Recherche Scientifique (CNRS), the Italian Istituto Nazionale di Fisica Nucleare (INFN) and the Dutch Nikhef, with contributions by institutions from Belgium, Germany, Greece, Hungary, Ireland, Japan, Monaco, Poland, Portugal, Spain. KAGRA is supported by Ministry of Education, Culture, Sports, Science and Technology (MEXT), Japan Society for the Promotion of Science (JSPS) in Japan; National Research Foundation (NRF) and Ministry of Science and ICT (MSIT) in Korea; Academia Sinica (AS) and National Science and Technology Council (NSTC) in Taiwan [98].

Appendix A: Estimation of the intrinsic binary parameters
In this Appendix, we compare the posterior distributions on the intrinsic binary parameters obtained using the SEOBNRHM and pSEOBNRHM waveform models. This complements the results shown in Sec. V B (Fig. 15) and Sec. VI (Fig. 16). For simplicity, we focus on the total mass M , the mass ratio q, the effective spin χ eff and the luminosity distance D L . Figure 15 shows the posterior distributions on the intrinsic binary parameters for a GR signal with the properties shown in Table II FIG. 15. The one-and two-dimensional posterior distributions on the intrinsic binary parameters of the total mas M , the mass ratio q, the effective spin χ eff and the luminosity distance DL for a GR injection with the parameters in Table II. The parameter estimation is performed using SEOBNRHM (solid curves) and pSEOBNRHM (dashed curves) waveform models. All contours indicate 90% credible regions and the vertical lines mark the inferred median values for each parameter.
with the parametrized plunge-merger-ringdown waveform model pSEOBNRHM (cf. Sec. V B). In both cases, the SNR is 98. We see that the 90% confidence intervals of the posterior distributions in the two analyses overlap in the parameter space. The most important difference is that the 90% credible intervals are wider in the pSEOBNRHM analysis. There are also changes to the median values of the binary parameters, as can be seen through the vertical lines in the plot.
To be more precise, the posteriors of the pSEOBNRHM model have a tail, most evidently in the mass ratio q. For the mass ratio, at 90% confidence interval, we find q = 0.87 +0.12 −0.10 (for the SEOBNRHM recovery) and q = 0.84 +0.15 −0.20 (for the pSEOBNRHM recovery). The broader posteriors, and tails, are due to the fact that the pSEOBNRHM model has five additional parameters with respect to SEOBNRHM. Qualitatively, by increasing the number of parameters in the model we increase the number of possible waveforms that match, to some extent, the injected signal. This will be most evidently seen in Fig. 16 which we discuss below. Figure 16 shows the posterior distributions on the intrinsic binary parameters for our analyses of GW150914 (cf. Sec. VI). The solid curves are obtained when we use SEOBNRHM for parameter estimation, whereas the dashed and dotted curves are obtained when we use pSEOBNRHM with non-GR parameters {δA, δω} ("merger test of GR") and {δA, δf 220 , δτ 220 } ("merger-ringdown test of GR"), respectively. The figure is similar to Fig. 15, discussed above. Again, we see that the 90% confidence intervals of the posterior distributions in the three analyses overlap in the parameter space. However, here we can see more explicitly how the increase of extra non-GR parameters in the waveform model