# Measuring stochastic gravitational-wave energy beyond general relativity

###### Abstract

Gravity theories beyond general relativity (GR) can change the properties of gravitational waves: their polarizations, dispersion, speed, and, importantly, energy content are all heavily theory-dependent. All these corrections can potentially be probed by measuring the stochastic gravitational-wave background. However, most existing treatments of this background beyond GR overlook modifications to the energy carried by gravitational waves, or rely on GR assumptions that are invalid in other theories. This may lead to mistranslation between the observable cross-correlation of detector outputs and gravitational-wave energy density, and thus to errors when deriving observational constraints on theories. In this article, we lay out a generic formalism for stochastic gravitational-wave searches, applicable to a large family of theories beyond GR. We explicitly state the (often tacit) assumptions that go into these searches, evaluating their generic applicability, or lack thereof. Examples of problematic assumptions are: statistical independence of linear polarization amplitudes; which polarizations satisfy equipartition; and which polarizations have well-defined phase velocities. We also show how to correctly infer the value of the stochastic energy density in the context of any given theory. We demonstrate with specific theories in which some of the traditional assumptions break down: Chern-Simons gravity, scalar-tensor theory, and Fierz-Pauli massive gravity. In each theory, we show how to properly include the beyond-GR corrections, and how to interpret observational results.

^{†}

^{†}preprint: LIGO-P1700234

## I Introduction

Besides transient signals, like those detected so far Abbott et al. et al. (2016a, b, c, 2017a, 2017b, 2017c, 2017d) by the Advanced Laser Interferometer Gravitational-wave Observatory (aLIGO) Aasi et al. (2015) and Virgo Acernese et al. (2015), gravitational-wave (GW) detectors are also expected to be sensitive to a persistent stochastic background Grishchuk (1976); Michelson (1987); Christensen (1990, 1992); Flanagan (1993); Allen (1996); Allen and Romano (1999); Abbott et al. et al. (2017e). This background signal is expected from primordial cosmological processes Starobinsky (1979); Easther et al. (2007); Barnaby et al. (2012); Cook and Sorbo (2012); Turner (1997); Easther and Lim (2006); Kamionkowski et al. (1994); Gasperini and Veneziano (1993a, b); Gasperini (2016); Caprini and Figueroa (2018), or the incoherent addition of myriad individually-unresolvable astrophysical sources, like compact binary coalescences Zhu et al. (2013); Marassi et al. (2011); Lasky et al. (2013); Rosado (2012); Zhu et al. (2011); Marassi et al. (2009); Buonanno et al. (2005); Sandick et al. (2006) or exotic topological defects Kibble (1976); Damour and Vilenkin (2005); Siemens et al. (2007); Sarangi and Tye (2002). Among many other rich scientific goals (see Maggiore (2000) for a review) detection of a stochastic background would provide an invaluable opportunity to study the fundamental nature of gravitational waves as they propagate over cosmological distances.

In the past decade or so, the formalism underlying stochastic GW searches has been extended to theories of gravity beyond general relativity (GR), primarily to account for the potential presence of nontensorial polarizations. Generic metric theories of gravity allow for up to six polarizations, corresponding to scalar (helicity 0), vector (helicity ) and tensor (helicity ) metric perturbations Eardley et al. (1973a, b). The effect of these extra polarizations on the stochastic background has been studied in particular for theories with scalar modes Maggiore and Nicolis (2000); Gasperini and Ungarelli (2001); Nishizawa and Hayama (2013), and in general for all possible modes in a theory-agnostic way Nishizawa et al. (2009, 2010). The problem of detecting nontensorial modes in the background has been studied in the context pulsar timing Lee et al. (2008); Chamberlin and Siemens (2012); Gair et al. (2015); Cornish et al. (2017), and GW measurements using astrometry O’Beirne and Cornish (2018). Beyond these proposals, a comprehensive data analysis framework has been recently implemented to search LIGO and Virgo data for GWs of any polarization, tensorial or otherwise, and some first upper limits have been placed on their amplitudes Callister et al. (2017); Abbott et al. et al. (2018).

The goal of searches for stochastic backgrounds, within GR or beyond, is to measure the amount of energy that the Universe contains in the form of gravitational waves. Consequently, treatments of stochastic GW signals are predominantly parametrized in terms of their effective energy-density spectrum [, defined in Eq. (29) below]. Such parametrization is only possible thanks to a standard set of assumptions about the properties of gravitational waves, the detectors, and the statistics of the background itself. Although generally justified within GR, the fundamental structure of beyond-GR theories may not always warrant all (or any) of those standard assumptions—even without considering modifications to specific emission mechanisms, or expected source populations. One must therefore be careful in applying the usual premises to searches for stochastic waves that aim to be theory agnostic, and should be aware that adopting any of these assumptions may come with additional observational restrictions.

Perhaps the most important example of an assumption that has been dubiously applied beyond GR concerns the form of the effective stress-energy of GWs. Multiple studies of stochastic signals beyond GR assume that the fractional energy density spectrum in GWs is related to the wave amplitudes in the same way as it is in GR Nishizawa et al. (2009, 2010); Nishizawa and Hayama (2013); Chamberlin and Siemens (2012); Yunes and Siemens (2013); Cornish et al. (2017); O’Beirne and Cornish (2018). Yet, as pointed out in Stein and Yunes (2011), the expression for the effective GW stress-energy need not be the same in all theories of gravity. This means that it is inadvisable to parametrize putatively model-independent searches for beyond-GR backgrounds assuming the GW energy density has the same functional form as in GR: doing so will result in the use of a quantity that should not generally be interpreted as the energy density in GWs. This is not only misleading, but (most importantly) can lead to incorrect comparisons between observational limits and theoretical predictions.

Besides this, some of the simplifying assumptions about the properties of the stochastic background that are usually justified in GR are not acceptable in general, and should not be extended to model-independent analyses. This is the case even without considering changes to the potential sources of the background in beyond-GR theories, which may themselves break more of the assumed symmetries. For instance, it is not reasonable to always assume that the usual linear GW polarization amplitudes will be statistically independent, as this will not be true unless the chosen polarization basis diagonalizes the kinetic matrix of the underlying theory of gravity. Similar arguments can be made about the assumptions that the polarizations are equipartitioned, or even that they have well defined phase velocities—let alone that they propagate at the speed of light.

In view of this, our goal is to straighten out the framework underlying searches for stochastic gravitational backgrounds, to make it formally valid and easily applicable to a large family of theories beyond GR. In Sec. II, we lay out a generic formalism for such searches, review the most commonplace assumptions in standard analyses, and evaluate their degree of applicability to other frameworks; along the way, we also clarify some relevant differences in conventions used by the theory and data analysis literatures. In Sec. III, we provide a series of examples of theories that break the premises behind one or more of these assumptions, and show the impact this has on the analysis—we focus on differences in the predicted form of the effective GW stress energy, but also discuss other problematic points. In particular, we use these examples to show how to go from the action defining a theory to (1) a relation between the fractional GW energy density spectrum and the correlation of polarization amplitudes, and (2) to the cross-correlation of GW detector outputs—which is the relevant observable for ground-based detectors. We review the derivation for general relativity in Sec. III.1, and then move on to Chern-Simons gravity in Sec. III.2, scalar-tensor theories in Sec. III.3, and Fierz-Pauli massive gravity in Sec. III.4. Finally, we offer a summary and conclusions in Sec. IV.

## Ii Formalism

In this section, we provide the framework required to search for stochastic GW backgrounds without assuming GR is correct. In Sec. II.1, we review the four-dimensional Fourier transform of a generic GW, lay out its decomposition into polarizations, and provide some useful identities for later use in Sec. III. In Sec. II.2, we focus on the properties of stochastic backgrounds, carefully reviewing the assumptions made in traditional analyses to determine whether they hold in theories beyond GR. In Sec. II.3, we describe the measurement process, including complications that may arise in generic theories. Finally, in Sec. II.4, we sketch the calculations needed to relate the effective stochastic GW energy in any given theory to the polarization amplitudes measurable by a detector.

Here, and throughout this document, spatial three-vectors are identified by an arrow (e.g. ), or a circumflex accent if they have unit norm (e.g. ). Four-vectors and higher-rank tensors are denoted by boldface, or abstract index notation (e.g. or ). For tensor coordinate components, spacetime Greek indices (, , , …) take values in the range 0–3, while spatial Latin indices (, , , …) span 1–3. We use metric signature , using for generic background metrics and for the Minkowski metric. Our conventions for the Levi-Civita tensor follow Misner et al. (1973): where is the determinant of the metric, and is the Levi-Civita symbol, with ; similarly, where is the determinant of the spatial metric, and . We normalize (anti-)symmetrizations as idempotent projection operations, e.g., and .

### ii.1 Decomposition of the metric perturbation

In any metric theory of gravity, as long as the observation region is small compared to the curvature radius, an arbitrary GW metric perturbation at a spacetime point may be expressed as a plane-wave expansion by the compact expression:

(1) |

integrating over all directions of propagation, and over both positive and negative frequencies. Here is the complex-valued Fourier amplitude for the wave-vector ; we let be the angular frequency, and the spatial wave-vector, implicitly defining as the sky location of the source. To simplify our notation in Eq. (1), we have defined the four-dimensional integral over the measure

(2) |

where is the Dirac delta function, and the last equality assumes an implicit integration over the magnitude of .
In order to write this, we assume that there is just one dispersion relation, , that determines the modulus of and implicitly defines .
^{1}^{1}1Some theories violate this assumption; for example, bimetric gravity Hassan and Rosen (2012) has one massless and one massive gravitational wave mode—we will allow for this briefly in Sec. II.3 only.
The dispersion relation is specific to the theory of gravity: for example, and in GR.

With the integration measure defined as in Eq. (2), in a local Lorentz frame (so that ), Eq. (1) can be recast in a form most common in stochastic GW literature (see, e.g., Flanagan (1993); Allen and Romano (1999); Romano and Cornish (2017)):

(3) |

where is the (potentially frequency-dependent) phase velocity of the wave ( in GR). Finally, to guarantee that be real, we must necessarily have

(4) |

where the asterisk indicates complex conjugation. In Appendix A, we elucidate the equivalence between Eqs. (1) and (3), derive the second equality in Eq. (2) and discuss differences between our Fourier conventions and those from the field theory literature.

For any given frequency and direction of propagation, the Fourier amplitudes may be written as a linear combination of at most six tensors corresponding to the six polarizations supported by generic metric theories of gravity Eardley et al. (1973a, b), even if the wave speed is slightly different from the speed of light Will (1993). Therefore, the most generic gravitational wave in this large category of theories may be written as a function of six independent amplitudes. We may, therefore, define six orthogonal polarization tensors, , such that

(5) |

where the sum is over six polarizations indexed by , and the ’s are the Fourier transforms of the six scalar fields, , encoding the amplitude of each mode, as defined by means of Eq. (1)

In order to study interactions between waves and detectors, it is usually convenient to pick a “synchronous” gauge^{2}^{2}2In a diffeomorphism invariant theory, one may always gauge transform into synchronous gauge by solving an initial value problem. If the theory is not diff-invariant, the Stückelberg trick can be used to restore the symmetry and then gauge transform. We provide an example of this in Sec. III.4.
such that the perturbation is purely spatial in the frame of interest (), and correspondingly so are the polarization tensors.
For instance, in an orthogonal frame in which the -axis is aligned with the direction of propagation (so that in that frame), we may write the six degrees of freedom as

(6) |

in terms of the linear tensor polarizations (, ), linear vector polarizations (x, y), and scalar breathing (b) and longitudinal (l) modes.

For the purpose of analyzing the output of multiple GW detectors, it is often convenient to write the polarization tensors in terms of unit vectors tangent and normal to the celestial sphere at each sky location. A standard linear polarization basis is given by

(7a) | ||||

(7b) | ||||

(7c) | ||||

(7d) | ||||

(7e) | ||||

(7f) |

where and are respectively the celestial polar and azimuthal coordinate vectors for a given source sky location determined by ; by design, these vectors satisfy . Other frame choices are possible, and multiple conventions abound in the literature.

As an example of the polarization decomposition of Eq. (5), consider theories in which gravitational perturbations carry spin-weight 2, like GR. In that case, we may choose to work with the two transverse-traceless linear polarization tensors corresponding to the plus () and cross () amplitudes shown in Eq. (6), and Eq. (5) becomes simply:

(8) |

Because the linear polarization tensors are real-valued by definition [cf. Eq. (7)], the reality condition for the amplitudes, Eq. (4), implies

(9) |

Alternatively, instead of the linear modes of Eq. (6), we could choose to work with eigenmodes of the helicity operator, i.e. the right- and left-handed circular polarization tensors (denoted “R” and “L” respectively). These modes satisfy an eigenvalue equation

(10) |

for (not summed on the RHS), where we have defined the factor , with the plus (minus) sign corresponding to the R (L) mode. Then, the circular polarization tensors can be written in terms of the ones for plus and cross as

(11) |

Using the circular tensors as a basis, we would write, instead of Eq. (8),

(12) |

and the reality condition, Eq. (4), would now imply (note the “L/R” subscript on the right hand side)

(13) |

instead of Eq. (9). The circular polarization modes can be similarly defined for the vector polarizations to obtain eigenmodes of helicity . On the other hand, the scalar modes have helicity 0, so in a sense are already circular.

For future reference, note that the spin-weight 2, spin-weight 1 and the transverse spin-weight 0 linear polarization tensors are normalized as usual such that, for a given direction of propagation,

(14) |

for , and the Kronecker delta; on the other hand, the longitudinal tensor satisfies . Similarly, the spin-weight 2 circular polarization tensors of Eq. (11) satisfy

(15) |

The basis tensors for the circularly-polarized vector modes also satisfy Eq. (15).

Although the two linear and circular bases discussed above are probably the most common in the GW literature (modulo normalizations), we are of course free to pick any other. For instance, in the analysis of differential-arm instruments, it is generally convenient to instead work with the traceless linear combination of and , since that is what such detectors can measure. Similarly, different theories may also define their own preferred polarization bases, given by the choice that diagonalizes their kinetic matrices.

### ii.2 Stochastic signals

In the case of stochastic signals, the Fourier amplitudes, , are, by definition, random variables and, as such, can be fully characterized by the moments of some (multivariate) probability distribution. Most standard searches for a stochastic GW background make the following assumptions about the random process that produced these amplitudes (see, e.g., Romano and Cornish (2017) for a review): the random process is (II.2) Gaussian, (II.2) ergodic, and (II.2) stationary, with no correlation between amplitudes from different (II.2) sky locations or (II.2) polarizations, and with (II.2) equipartition of power across polarizations; furthermore, the process is commonly (although not universally) assumed to be (II.2) isotropic. We break down these assumptions below, and introduce some important definitions along the way.

Stochastic backgrounds are expected to arise from primordial cosmological processes Starobinsky (1979); Easther et al. (2007); Barnaby et al. (2012); Cook and Sorbo (2012); Turner (1997); Easther and Lim (2006); Kamionkowski et al. (1994); Gasperini and Veneziano (1993a, b); Gasperini (2016); Caprini and Figueroa (2018), or by the incoherent superposition of a great number of signals from contemporary astrophysical events Zhu et al. (2013); Marassi et al. (2011); Lasky et al. (2013); Rosado (2012); Zhu et al. (2011); Marassi et al. (2009); Buonanno et al. (2005); Sandick et al. (2006); Kibble (1976); Damour and Vilenkin (2005); Siemens et al. (2007); Sarangi and Tye (2002); Maggiore (2000). The assumption (i) that the astrophysical background is produced by a Gaussian random process is motivated by the central limit theorem—this guarantees that the properties of any large number of incoherently-added GW signals will be normally distributed, regardless of the specific characteristics of any given source. A similar argument can be applied to primordial signals by considering the independent evolution of waves from causally-disconnected regions Allen (1996). Although waves from inflation will technically have non-Gaussianities, they will be small as long as inflation satisfied the slow-roll approximation Maldacena (2003); Caprini and Figueroa (2018); Allen (1996).

For Gaussian processes, all properties of the probability distribution are determined by its first two moments (correlation functions)—namely, the mean and power spectrum (respectively, the one- and two-point correlation functions).
The first moment of the distribution, the mean , will not appear explicitly in any of the expressions below, so we ignore it.^{3}^{3}3Some authors explicitly set this value to zero because the contribution from a nonvanishing mean would take the form of a coherent offset in the Fourier amplitudes as a function of frequency, which not only would be hard to justify physically, but would also hardly classify as “stochastic” (see, e.g., Romano and Cornish (2017); Allen and Romano (1999)).
Here and below, the expectation value, denoted by angle brackets , corresponds to ensemble averages, as well as space/time-averages by assumption (ii) of ergodicity.
The expectation of ergodicity itself comes from the assumption that the Universe is homogeneous (for more discussion on this topic, see Caprini and Figueroa (2018)).

The second moment of the distribution will end up being an important observable. In order to write down an expression for it, we can make use of assumptions (II.2) and (II.2). First, stationarity (iii) is motivated by the fact that observation times (order of months to years) are extremely small relative to the dynamical timescales intrinsic to the cosmological processes that could change the properties of the background (order of billions of years); therefore, any changes in the stochastic background would be unnoticeable to us. Formally, stationarity means that the first moment is constant, while the second moment depends only on time differences (see, e.g., Jaranowski and Królak (2009)). As a consequence, the Fourier transform of a stationary random variable can be shown to be such that amplitudes at different frequencies will be statistically independent and, therefore, uncorrelated (Appendix B).

Next, the assumption (iv) that amplitudes from different sky locations will be uncorrelated is justified for primordial waves because signals from different points in the sky are only coming into causal contact now at Earth, under ordinary topological assumptions. One could potentially search for nonstandard spatial topologies in a sufficiently “small” universe through angular correlations in gravitational waves Grishchuk et al. (1975), in much the same way as in the cosmic microwave background (CMB) Lachieze-Rey and Luminet (1995); Cornish et al. (1998). A small universe with nonstandard spatial topology would induce circles of excess correlation in both the CMB and gravitational-wave background. As there has been no evidence of this phenomenon in the CMB, in this article we consider primordial signals from different sky directions to be uncorrelated.

For contemporary (“astrophysical”) backgrounds, (II.2) comes from the assumption that the contributions from multiple sources throughout the sky (say, binary systems) are added “incoherently”—that is, sources are not perfectly aligned and timed as would be needed for signals from different directions to reach us with matching phase and amplitude evolution. Even though such astrophysical sources were in causal contact at some point in the past, they are embedded in chaotic astrophysical environments (with e.g. turbulent magnetohydrodynamics) with Lyapunov times sufficiently short that in practice, they can be treated as uncorrelated. In principle, strong gravitational lensing may introduce correlations between sky bins into the stochastic background, whether primordial or contemporary, but we can expect this effect to be negligible in practice Smith et al. (2018).

With assumptions (II.2) and (II.2) in place, we may write the second moment of the amplitude distribution in the form (Appendix B):

(16) |

This equation defines the one-sided cross-power spectral density, , for two signals, , sharing a wave-vector but with potentially different polarizations and . For linear polarizations, this quantity satisfies , because of the reality condition of Eq. (4). Since we are usually interested in the total measured power at a given frequency, regardless of sky direction, we also define the integral of over the sky,

(17) |

which carries units of . For , this is nothing more than the one-sided power spectral density (PSD) in polarization , which we denote . In general, for any real-valued random variable , the PSD can be approximated as twice the square of the band-limited Fourier transform Jaranowski and Królak (2009),

(18) |

in practice always computed for some long but finite integration time, , on the order of months to years for observations of the stochastic background. As usual, the factor of 2 in Eq. (18) accounts for the fact that this is the one-sided PSD, .

Assumption (v) that the different polarizations are statistically independent may be used to discard off-diagonal terms in the cross-power spectrum, so that . However, one must be careful with this simplification: the assumption is valid if and only if one works in a polarization basis that diagonalizes the kinetic matrix of the theory. Importantly, as we will show with specific examples, such a basis need not be the linear polarization basis used in most GR analyses. Even when working within GR, it is generally better, from a theoretical standpoint, to work in terms of the circular modes, as they are eigenstates of the helicity operator, and they might be produced with different intensities in the early universe Alexander et al. (2006); Seto and Taruya (2007, 2008).

Besides assuming that the polarizations are uncorrelated, it is also common to assume that there is equipartition of power between them—assumption (vi) in our list above. Under this presumption, the background is said to be unpolarized and the polarization PSDs may be written in terms of the total GW spectral density, , such that , where is the number of polarizations allowed to propagate in a given theory. In general, this assumption is only justified if the polarizations both diagonalize the kinetic matrix and interact similarly with matter, so that they are sourced in equal amounts. This not always the case: for example, in both massive gravity de Rham (2014); Hinterbichler (2012) and dynamical Chern-Simons gravity Jackiw and Pi (2003); Alexander and Yunes (2009); Sopuerta and Yunes (2009), different polarizations couple to sources with different strengths.

Finally, the simplest searches for a stochastic background also adopt assumption (vii) of isotropy, in which case , by Eq. (17).
In GR, if one disregards the proper motion of the solar system, this assumption is expected to hold well for most foreseeable sources of a stochastic background detectable by existing ground-based observatories, since they are expected to originate from cosmological distances
Romano and Cornish (2017); Allen (1996); Allen and Romano (1999).
^{4}^{4}4The same will not necessarily be true for LISA, which will be sensitive to galactic stochastic sources, like the “confusion noise” from white-dwarf binaries Ruiter et al. (2010).
For cosmological sources, isotropy is likely also a good
assumption in many beyond-GR theories; however, isotropy should not be expected to hold
in theories with a preferred frame, which are intrinsically
anisotropic Kostelecky and Samuel (1989); Kostelecky and Potting (1991, 1995); Kostelecky (2004); Kostelecký and Mewes (2016); Jacobson and Mattingly (2001); Eling et al. (2004); Jacobson (2007).
For simplicity, the rest of this document will treat only the case of an isotropic background, but this does not affect the spirit of the results, which can be easily generalized to the anisotropic case.
For predictions of the angular power spectrum of astrophysical GR backgrounds see Cusin et al. (2018), and for corresponding observational limits that do not assume isotropy see Abbott et al. et al. (2017f).

Assuming both (II.2) an isotropic background and (II.2) uncorrelated polarizations, on top of (II.2) stationarity and (II.2) uncorrelated sky bins, Eq. (16) can be written directly in terms of the power spectral density,

(19) |

If one further assumed (II.2) equipartition, would be replaced with , as explained above. This is the form of the expression most common in recent literature about detection of stochastic gravitational-wave backgrounds (e.g., Romano and Cornish (2017)).

### ii.3 Detection

Because the output of ground-based GW detectors is largely dominated by stochastic instrumental and environmental noise Abbott et al. et al. (2016d); Covas et al. (2018), it is not possible to measure the power spectrum of the polarization amplitudes, , directly with a single detector at any level of interest. However, this quantity may be inferred by looking instead at the cross-correlation of the output of two or more instruments (see, e.g., Romano and Cornish (2017) for a comprehensive review of data analysis methods).

We assume that each GW detector has a purely linear response to gravitational waves. Therefore, in the Fourier domain, the response of detector to a plane wave must be expressible as

(20) |

for some tensor representing the detector’s frequency- and direction-dependent transfer function. This tensor encodes all relevant information about the detector and the physics of the measurement process Forward (1978); Linsay et al. (1983); Estabrook (1985); Schilling (1997); Rakhmanov (2005); Rakhmanov et al. (2008); Rakhmanov (2009); Essick et al. (2017) (for considerations specific to gravity beyond GR, see e.g. Eardley et al. (1973b); Tobar et al. (1999); Błaut (2012); Poisson and Will (2014); Isi et al. (2017); Gair et al. (2015)). The detector’s output (e.g. the calibrated current out of a photodiode) is a gauge-invariant observable. However, the metric perturbation is gauge-dependent; therefore, the detection tensor must also depend on the gauge choice, so that the overall gauge dependence on the RHS of Eq. (20) exactly cancels.

Assuming a basis of polarization states that have well defined phase velocities (i.e. they diagonalize the kinetic matrix of the theory), we may use Eq. (20) to write the Fourier transform of the signal at detector explicitly as a sum over polarizations and an integral over sky directions,

(21) |

defining the Fourier-domain response functions as the contraction between the detector and polarization tensors, , which must also be gauge-dependent.

The time-domain analogue of Eq. (20) is given by a convolution,

(22) |

with the location of detector , and its impulse response. Since is gauge-dependent, the same must be true for . For an ideal differential arm-length instrument, it is easiest to write down this detector tensor in a synchronous gauge ( in the detector frame), wherein the end test masses’ coordinate locations will not change Misner et al. (1973). In such a gauge, the resulting differential-arm detector tensor is the purely geometric factor

(23) |

with and spacelike unit vectors pointing along the detector arms. For real interferometric detectors, like LIGO and Virgo, Eq. (23) is valid only in the small-antenna limit (arm length GW wavelength) Schilling (1997); Rakhmanov (2005); Rakhmanov et al. (2008); Rakhmanov (2009); Essick et al. (2017).

For any realistic detector, the tensor of Eq. (23) will vary in time due to the motion of the instrument with respect to the inertial frame of the wave (e.g. due to Earth’s rotation, for ground-based observatories). However, for the cases we are interested in, we can take this variation to be slow with respect to the period of the waves, so that it can be ignored if Eq. (21) is implemented via short-time Fourier transforms. In this ideal “slow-detector” limit, we may then treat the response as time- and frequency-independent to write , and so Eq. (21) simplifies to

(24) |

with the frequency-independent antenna patterns defined in full analogy to our definition of above,

(25) |

For details on this simplification, and nuances applicable to anisotropic backgrounds, see Sect. IV in Allen and Ottewill (1997).

In the Fourier domain, the cross-correlation between the output of two detectors may then be written in terms of the second moment of the distribution of polarization amplitudes as

(26) | ||||

where, again, assumption (II.2) of ergodicity is tacitly implied. If we also assume, as we will throughout this paper, that the background is (II.2) stationary and (II.2) isotropic, and (II.2) that sky bins are uncorrelated, we may then use Eq. (19) to simplify this to

(27) |

where we have defined the generalized overlap reduction function for polarizations , and detectors , ,

(28) |

in terms of the phase factor , which acquires a potential frequency dependence through the phase velocities. If there is one dispersion relation shared by all polarizations (true throughout the rest of this paper), the exponent in Eq. (28) can be written as , in terms of the separation between detectors . The overlap reduction functions encode all relevant information pertaining GW polarizations and speed, as well as detector geometry. The specific definition and normalization chosen here are intended to facilitate generalization of the analysis beyond GR, and are not necessarily standard (see, e.g., Sect. 5.3 of Romano and Cornish (2017) for a review of these functions and their properties).

Because the noise in different instruments will generally be statistically independent Abbott et al. et al. (2016d); Covas et al. (2018), by cross-correlating the output of a pair of detectors, one may directly measure the signal cross-correlation of Eq. (27), and hence infer the polarization power spectra (as proposed by Grishchuk (1976); Michelson (1987), and studied in multiple works since). In a theory that allows for independent polarizations, there will be up to different terms (only if the correlation matrix is diagonal), and at least as many detector pairs (“baselines”) will be needed to break all degeneracies between them.

### ii.4 Energy density

Searches for a stochastic gravitational-wave background attempt to measure the Universe’s total energy density in gravitational waves as a function of frequency. However, inferring this quantity from direct observables requires theoretical assumptions. Furthermore, the equivalence principle precludes being able to localize energy density in gravitational waves, so this is in fact an effective energy density. We elaborate on these important points below, and sketch the general procedure to link the effective GW energy density to observables at the detector in (almost) any given theory. Concrete examples of how to apply this are provided in Sec. III.

With an eye to cosmology, the quantity of interest in stochastic searches is usually chosen to be the log-fractional spectrum of the effective GW energy density Michelson (1987); Christensen (1990, 1992); Flanagan (1993); Allen and Romano (1999),

(29) |

with the effective GW energy density as a function of frequency, and the critical density required to close the universe,

(30) |

where is the present Hubble parameter Allen and Romano (1999). Presenting results of a stochastic background search in terms of this quantity facilitates their cosmological interpretation. More importantly, using an energy density (however parametrized) allows for direct comparison with theoretical models: in order to predict the properties of the GW background, one computes the typical GW power emitted by the system of interest (e.g., compact binaries, cosmic strings, or primordial fluctuations) and then obtains an energy spectrum by incoherently adding many such contributions (e.g., using the quadrupole formula with merger rates from population synthesis) Starobinsky (1979); Easther et al. (2007); Barnaby et al. (2012); Cook and Sorbo (2012); Turner (1997); Easther and Lim (2006); Kamionkowski et al. (1994); Gasperini and Veneziano (1993a, b); Gasperini (2016); Caprini and Figueroa (2018); Zhu et al. (2013); Marassi et al. (2011); Lasky et al. (2013); Rosado (2012); Zhu et al. (2011); Marassi et al. (2009); Buonanno et al. (2005); Sandick et al. (2006); Kibble (1976); Damour and Vilenkin (2005); Siemens et al. (2007); Sarangi and Tye (2002); Maggiore (2000).

However, GW detectors do not measure the effective physical GW energy density, but rather the amplitude of the waves at each instrument. In particular, searches for a stochastic background are sensitive to the (incoherent) strain amplitude power, Eq. (16). This will remain true for future detection methods, like space missions or pulsar timing. In the case of ground-based observatories, as outlined in Sec. II.3, the stochastic strain amplitudes are probed through the cross-correlation of detector outputs across a network, Eq. (26). Thus, whatever the detection method, we will need an object that relates gravitational-wave amplitudes to energies—a mapping that is theory-dependent.

The frequency-domain effective stress-energy tensor (ESET) for gravitational waves lets us translate between the more accessible two-point amplitude correlation function, Eq. (16), [or the two-detector-output cross-correlation, Eq. (26)] and the GW contribution to the energy density, Eq. (29). In GR the ESET is given by a simple expression first derived by Isaacson Isaacson (1968a, b) (see Sec. III.1 below), which enables stochastic searches to be parametrized directly in terms of Michelson (1987); Christensen (1990, 1992); Flanagan (1993); Allen and Romano (1999). Interestingly, the same relationship has been assumed to hold in most stochastic GW data analysis schemes that allow for departures from GR Nishizawa et al. (2009, 2010); Nishizawa and Hayama (2013); Chamberlin and Siemens (2012); Yunes and Siemens (2013); Cornish et al. (2017); O’Beirne and Cornish (2018), even though the Isaacson formula will not necessarily hold in arbitrary theories Stein and Yunes (2011). Using the Isaacson formula when inappropriate will lead to a mistranslation between detector cross-correlations and GW energy densities. This is not only misleading, but can also lead to errors when deriving constraints on theories from observations.

In the context of any specific theory of gravity, the ESET can be derived directly from the action. The ESET is given by a space-time average of the variation of the second-order perturbation of the action with respect to the background (inverse) metric Stein and Yunes (2011),

(31) |

where the double angular brackets indicate an averaging procedure over a spacetime region on the order of several wavelengths (e.g. Brill-Hartle averaging, though other procedures Zalaletdinov (1996) agree when there is a separation of length scales). We briefly summarize the approach here; we refer the interested reader to Stein and Yunes (2011) for more exposition.

The second-order Lagrangian is obtained from the action after perturbing the metric and other dynamical fields via

(32) | ||||

(33) |

and collecting terms in the action order-by-order in the small parameter . This gives the expansion

(34) | ||||

where means both and are present. At order , the action generates the ordinary nonlinear background equations of motion for and . At order , the action is purely a “tadpole” term which vanishes when are on shell, and therefore does not contribute to any equations of motion. The same is true for the second-order perturbations , which appear linearly in . However, appear quadratically in : the quadratic action then generates the linear equations of motion for when varied with respect to ; at the same time, the variation with respect to will be a quadratic functional of , and results in the ESET.

From now on we drop the order-counting superscript, letting , since we will not encounter . In a local Lorentz frame whose time direction is aligned with the Hubble flow, we can define the position-space effective GW energy density as

(35) |

where the double argument is just to remind us that is a quadratic functional of . To use Eq. (29), we want in momentum space, so we need to make use of a plane-wave expansion like Eq. (1). The result will always be a momentum-space integral of the form

(36) |

where the (gauge-dependent) tensor encodes information about the kinetic matrix of the theory in momentum space, and we have used Eq. (4) to write . Notice that here we have replaced the spacetime averaging of Eq. (31) with ensemble averaging, based on assumption (II.2) ergodicity. When the two-point function is of the form of Eq. (16), the double integral will collapse to a single integral, and the physical energy density will be related to the power spectral density , with some potentially nontrivial frequency dependence arising from (we will see several examples below).

When this double integral collapses to a single integral, we can then define the fractional energy density per frequency bin via

(37) |

With this definition of , and the relationship between the energy density Eq. (36) and a two-point function like Eq. (16), it will be possible to relate the power spectral density to the cosmological fractional energy density , Eq. (29). The relationships between all these key quantities are illustrated in Fig. 1.

Once we have this, we may work directly with ; in particular, data analysis searches usually assume a power-law model like

(38) |

for some spectral index , and the characteristic amplitude at some arbitrary reference frequency . This is how LIGO generally parametrizes its searches, e.g. Abbott et al. et al. (2017e); for a discussion of the validity of this parametrization, see Callister et al. (2016).

## Iii Example theories

In this section, we show how different gravitational theories imply different functional relations between the effective fractional energy density spectrum, in Eq. (29), the strain cross-power spectrum, in Eq. (16), and, consequently, the cross-correlation between detector outputs, Eq. (26) As discussed in Sec. II.3, this last quantity is the relevant observable for ground-based instruments, on which we focus. The relationships between all the key quantities are illustrated in Fig. 1. Along the way, we also discuss the expected statistical properties of the polarization amplitudes under each framework, as required purely by the basic structure of the theory (that is, not considering specific source models).

We first demonstrate the procedure by rederiving the standard GR expressions from the Einstein-Hilbert action (Sec. III.1), and then offer a series of beyond-GR examples for which the analogous result is different: we consider the case of Chern-Simons gravity, a theory which is not parity-symmetric (Sec. III.2); this is followed by Brans-Dicke gravity, the prototypical example of a scalar-tensor theory (Sec. III.3); finally, we study Fierz-Pauli gravity (Sec. III.4), in which the graviton is endowed with a mass. The last two examples support nontensorial modes of the metric perturbation (see Sec. II.1).

For all the examples we consider, we find it reasonable to simplify our equations by assuming the stochastic background is (II.2) Gaussian, (II.2) ergodic, (II.2) stationary and (II.2) isotropic, with (II.2) no correlation between different sky locations. In all cases, then, we find that we can write the cross-correlation between the output of two ideal differential-arm detectors and in the form

(39) | ||||

where the sum is over some polarization basis that diagonalizes the kinetic matrix of the theory. Here the ’s are the generalized overlap reduction functions of Eq. (28), is the effective fractional energy spectrum in polarization defined by analogy to Eq. (29), and is a model-dependent factor encoding deviations from GR. In Einstein’s theory, for tensor polarizations and vanishes otherwise, as we show below.

Many of the results in this section are derived on a flat background, and will therefore be erroneous in a cosmological setting. However, because of the vast separation of scales between the gravitational wavelength and the Hubble parameter today , the error between the flat space results and the cosmologically-correct results will be of fractional order . This correction has been explicitly computed in GR Bonga and Hazboun (2017), and while we are not aware of the same computation in beyond-GR theories, it should remain true as long as the theory of gravity respects the separation of scales.

### iii.1 General relativity

The vacuum Einstein field equations can be derived from the Einstein-Hilbert (EH) action,

(40) |

where , is the determinant of the metric , and is the Ricci scalar Misner et al. (1973). We may now expand the metric around some background, , as in Eq. (32). The source-free linearized equations of motion, on a flat background (so that Riemann vanishes), and in the transverse-traceless gauge ( and ), take the simple form

(41) |

where is the d’Alembertian with respect to the background metric. Equation (41) leads to the standard geometric optics approximation to GW propagation, from which it follows that GWs show no birefringence and always propagate at the speed of light.

Focus now on the second-order perturbation of Eq. (40). On a flat background, the second-order Lagrangian density is given by Stein and Yunes (2011); MacCallum and Taub (1973)

(42) |

where all derivatives are taken with respect to the background metric , and is the trace-reversed metric perturbation. This piece of the Lagrangian density corresponds to in Eq. (34).

We now apply the transverse-traceless gauge conditions and evaluate the perturbations on-shell (that is, we enforce the first-order equations of motion). Then, varying with respect to the inverse background metric, , as in Eq. (31), we obtain, far away from sources,

(43) |

This is the well-known expression for the effective stress-energy carried by a gravitational wave, first derived by Isaacson (and, consequently, known as the Isaacson formula) Isaacson (1968a, b).

We may now use Eq. (43) to relate and , as outlined in Sec. II.4. The Isaacson expression implies that, in a local Lorentz frame,

(44) |

where we have used the fact that the transverse-traceless metric perturbation will be purely spatial. Plugging in the plane-wave expansion of Eq. (1), using the reality condition of Eq. (4), and invoking (II.2) ergodicity, we may rewrite this as

(45) |

This means that, in GR and in our gauge, takes the form of Eq. (36) with

(46) |

It is convenient at this point to expand the Fourier amplitudes into polarizations. Because GR is parity-symmetric, in this theory all modes are generated and propagate equally, so one is free to choose between linear and circular polarizations; however, working with the former is slightly simpler because the corresponding polarization tensors, Eqs. (7a) and (7b), are real-valued. Then, summing over , ,

(47) |

We now use the fact that the Fourier amplitudes are given by a random process to simplify our expression for via Eq. (16). Following common practice and for the sake of simplicity, we will assume that the stationary Gaussian background is also (II.2) isotropic and (II.2) unpolarized, with equal contributions from the linear polarizations. Letting the total PSD in tensor polarizations be with , this means

(48) |

Again, the assumption (II.2) of equipartition is justified because GR conserves parity. The correlation of the Fourier polarization amplitudes, Eq. (19), then becomes

(49) |

With this in place, and noting that Eq. (14) implies when summing over , the effective energy density of Eq. (III.1) simplifies to

(50) |

Comparing with Eq. (37), we can immediately read off , and then, from the definition of , Eq. (29), we conclude

(51) |

As discussed in Sec. II.3, the actual observable for stochastic-background searches in data from ground-based observatories is the cross-correlation between the outputs of pairs of detectors. For an isotropic background, this is given by Eq. (27), which can be written in terms of the fractional energy density by means of Eq. (51):

(52) |

where we have defined the total tensor overlap-reduction function as . This is the desired expression relating the observable strain cross-correlation to the fractional effective-energy density spectrum, that will be predicted by theory. Eq. (52) is used in most LIGO and Virgo searches for a stochastic background, via parametrizations like the power-law of Eq. (38). Comparing to Eq. (39), and recalling with , we see that in GR, for tensor polarizations, and vanishes otherwise, as expected.

### iii.2 Chern-Simons gravity

Chern-Simons (CS) theory is an extension of GR with motivations ranging from anomaly-cancellation in curved spacetime, low-energy limits of both string theory and loop quantum gravity, effective field theory of inflation, and more Delbourgo and Salam (1972); Eguchi and Freund (1976); Alvarez-Gaume and Witten (1984); Green and Schwarz (1984); Ashtekar et al. (1989); Campbell et al. (1991, 1993); Jackiw and Pi (2003); Weinberg (2008); Taveras and Yunes (2008); Alexander et al. (2008); Mercuri and Taveras (2009). The theory is characterized by the presence of a parity-odd, axion-like scalar field, which couples to curvature through a parity-odd interaction (see Alexander and Yunes (2009) for a review). The ESET in this theory was derived in Stein and Yunes (2011), in an asymptotically-flat spacetime and approaching future null infinity (). As noted before, by promoting flat-space results to a cosmological setting, we are making an extremely small error of fractional order . Below, we provide a sketch of this derivation and show what the result implies for the stochastic background.

As a consequence of its lack of parity symmetry, CS gravity generally predicts birefringent propagation and generation of the metric perturbations, so that one of the circular tensor polarizations is amplified at the expense of the other Alexander and Martin (2005). Consequently, as is true for any theory lacking parity symmetry, it is not appropriate to assume that the stochastic background is unpolarized Contaldi et al. (2008). Furthermore, as we will see, the nondynamical version of the theory predicts an expression for the effective GW stress-energy different from the Isaacson formula of Eq. (43), and consequently differs from GR via a factor of in Eq. (39).

In the absence of matter, CS gravity is given by the Einstein-Hilbert action of Eq. (40), plus terms describing the axion-curvature coupling (), and dynamics () of the scalar field Jackiw and Pi (2003); Alexander and Yunes (2009),

(53) | ||||

(54) | ||||

(55) |

In the above, is the constant determining the coupling of the CS field to the gravitational sector, while controls the kinetic energy of the scalar; is the Pontryagin density, which is defined in terms of the Riemann tensor, , by

(56) |

with the Levi-Civita tensor. This term is parity-odd, and gives CS gravity much of its richness.

Studying the dynamics of the theory, one may show that gravitational waves in CS gravity will present the same tensor (spin-weight 2) propagating degrees of freedom as in GR Jackiw and Pi (2003); Sopuerta and Yunes (2009). On a flat background and in Lorenz gauge (), metric perturbations follow the first-order equations of motion Stein and Yunes (2011); Alexander and Yunes (2009),

(57) | ||||

where we split into a smooth background piece and a perturbation , and is the stress energy sourced quadratically by ,

(58) |

Again on a flat background, CS gravity admits an approximately traceless gauge Stein and Yunes (2011), so that can be replaced by in these equations, as was done for GR in Eq. (41).