ROUGH DRAFT authorea.com/98745
Main Data History
Export
Show Index Toggle 1 comments
  •  Quick Edit
  • Detecting a stochastic gravitational wave signal

    This note discusses trying to detect a generic gravitational wave with an unknown waveform emitted from a particular sky position in data from two separate gravitational wave detectors. We define two slightly different approaches to this problem.

    The signal

    First we will define the gravitational wave signal at one timestamp, \(i\), observed in one detector, \(L\). We envisage two possible methods for this analysis with slightly different model definitions. The first uses

    \begin{equation} \label{eq:signal1}h_{i}^{L}=h_{0}\left({A_{+}}_{i}{F_{+}}_{i}^{L}(\psi_{i})+{A_{\times}}_{i}{F_{\times}}_{i}^{L}(\psi_{i})\right),\\ \end{equation}

    where \({A_{+}}_{i}\) and \({A_{\times}}_{i}\) are the signal amplitudes scale factors (which could be positive or negative) in the plus and cross polarisations at timestamp \(i\) (which would be the same for different detectors), \({F_{+}}_{i}^{L}(\psi_{i})\) and \({F_{\times}}_{i}^{L}(\psi_{i})\) are the detector’s antenna response to the plus and cross polarisations for a given polarisation angle \(\psi_{i}\)11Note that \(\psi\) could change between data points, so this is also indexed for the current timestamp., and \(h_{0}\) is an overall underlying gravitational wave amplitude. The second uses

    \begin{equation} \label{eq:signal2}h_{i}^{L}={A_{+}}_{i}{F_{+}}_{i}^{L}(\psi_{i})+{A_{\times}}_{i}{F_{\times}}_{i}^{L}(\psi_{i}),\\ \end{equation}

    where, in this case, the \({A_{+}}_{i}\) and \({A_{\times}}_{i}\) are the actual signal amplitudes in the plus and cross polarisations at timestamp \(i\).

    First Method

    \label{sec:method1}

    Here we will examine the details of the first method, which uses the signal model defined in Equation \ref{eq:signal1}.

    Now, if we had one data point for detector \(X\), \(d_{i}^{X}\), and assuming the noise in the detector is Gaussian with zero mean and standard deviation of \(\sigma_{i}^{X}\), then the likelihood for the data given the model is

    \begin{equation} \label{eq:singlelikelihood}p(d_{i}^{X}|h_{0},{A_{+}}_{i},{A_{\times}}_{i},\psi,I)=\frac{1}{\sqrt{2\pi}\sigma_{i}^{X}}\exp{\left(-\frac{(d_{i}^{X}-h_{i}^{X})^{2}}{2{\sigma_{i}^{X}}^{2}}\right)}.\\ \end{equation}

    We now add another detector, \(Y\), with data point \(d_{i}^{Y}\), where the \(i\) timestamp index in detector \(Y\) is actually indexing a time that is shifted with respect to that in detector \(X\) based on the time delay between detectors for the known signal position. So, e.g. \(t_{i}^{Y}=t_{i}^{X}+\Delta t_{i}(\alpha,\delta)\). This gives a joint likelihood of the data for the two detectors of

    \begin{equation} \label{eq:jointlikelihood}p(\{d_{i}^{X},d_{i}^{Y}\}|h_{0},{A_{+}}_{i},{A_{\times}}_{i},\psi_{i},I)=\frac{1}{2\pi\sigma_{i}^{X}\sigma_{i}^{Y}}\exp{\left(-\sum_{L=X,Y}\frac{(d_{i}^{L}-h_{i}^{L})^{2}}{2{\sigma_{i}^{L}}^{2}}\right)}.\\ \end{equation}

    We would like to get a posterior probability distribution on \(h_{0}\) alone (and in fact we also want the evidence marginalised over \(h_{0}\) too). If we assume that the \(A\) scale factors and \(\psi\) change on the timescale of individual data points then we want to marginalise over them for each point, e.g.

    \begin{align} p(\{d_{i}^{X},d_{i}^{Y}\}|h_{0},I)= & \int_{{A_{+}}_{i}}\int_{{A_{\times}}_{i}}\int_{0}^{\pi}p(\{d_{i}^{X},d_{i}^{Y}\}|h_{0},{A_{+}}_{i},{A_{\times}}_{i},\psi_{i},I)\times\notag \\ & \label{eq:h0likelihood} p({A_{+}}_{i}|I)p({A_{\times}}_{i}|I)p(\psi_{i}|I){\rm d}{A_{+}}_{i}{\rm d}{A_{\times}}_{i}{\rm d}\psi_{i},\\ \end{align}

    where \(p({A_{+}}_{i}|I)\) and \(p({A_{\times}}_{i}|I)\) are the priors on the scale factors and \(p(\psi_{i}|I)\) is the prior on \(\psi_{i}\).

    To get the joint likelihood over all the data we can just use the product of Equation \ref{eq:h0likelihood}, such that we have

    \begin{equation} \label{eq:fulljoint}p(\{\vec{d}^{X},\vec{d}^{Y}\}|h_{0},I)=\prod_{i=1}^{N}p(\{d_{i}^{X},d_{i}^{Y}\}|h_{0},I),\\ \end{equation}

    where \(N\) is the total number of data points. We can then get the posterior on \(h_{0}\) as

    \begin{equation} \label{eq:posterior1}p(h_{0}|\{\vec{d}^{X},\vec{d}^{Y}\},I)=\frac{p(\{\vec{d}^{X},\vec{d}^{Y}\}|h_{0},I)p(h_{0}|I)}{p(\{\vec{d}^{X},\vec{d}^{Y}\}|I)},\\ \end{equation}

    where \(p(h_{0}|I)\) is the prior on \(h_{0}\) and \(p(\{\vec{d}^{X},\vec{d}^{Y}\}|I)\) is the evidence for the data. The evidence, \(Z\), is given by

    \begin{equation} \label{eq:evidence1}Z=\int p(\{\vec{d}^{X},\vec{d}^{Y}\}|h_{0},I)p(h_{0}|I){\rm d}h_{0}.\\ \end{equation}

    Practical evaluation

    Depending on the choice of prior some of the marginalisations in Equation \ref{eq:h0likelihood} are analytical, but others will need to be evaluated numerically. To start with we can set a Gaussian priors on \({A_{+}}_{i}\) and \({A_{\times}}_{i}\), both with zero mean and standard deviations of \(\sigma_{A_{+,\times}}\)

    \begin{equation} \label{eq:priorAp}p({A_{+,\times}}_{i}|I)=\frac{1}{\sqrt{2\pi}\sigma_{A_{+,\times}}}\exp{\left(-\frac{{A_{+,\times}}_{i}^{2}}{2\sigma_{A_{+,\times}}^{2}}\right)}.\\ \end{equation}

    We can get a marginalised likelihood on \(h_{0}\) and \(\psi_{i}\) by multiplying Equation \ref{eq:h0likelihood} by these priors and integrating over \({A_{+}}_{i}\) and \({A_{\times}}_{i}\) between \(-\infty\) and \(\infty\)

    \begin{align} \label{eq:margApAc}p(\{d_{i}^{X},d_{i}^{Y}\}|h_{0},\psi_{i},I) & =\int_{-\infty}^{\infty}\int_{-\infty}^{\infty}p(\{d_{i}^{X},d_{i}^{Y}\}|h_{0},{A_{+}}_{i},{A_{\times}}_{i},\psi_{i},I)p({A_{+}}_{i}|I)p({A_{\times}}_{i}|I){\rm d}{A_{+}}_{i}{\rm d}{A_{\times}}_{i}, \\ & =\int_{-\infty}^{\infty}\int_{-\infty}^{\infty}\frac{\exp{\left(-\sum_{L=X,Y}\frac{(d_{i}^{L}-h_{i}^{L})^{2}}{2{\sigma_{i}^{L}}^{2}}\right)}\exp{\left(-\frac{1}{2}\left[\frac{{A_{+}}_{i}^{2}}{\sigma_{A_{+}}^{2}}+\frac{{A_{\times}}_{i}^{2}}{\sigma_{A_{\times}}^{2}}\right]\right)}}{4\pi^{2}\sigma_{i}^{X}\sigma_{i}^{Y}\sigma_{A_{+}}\sigma_{A_{\times}}}{\rm d}{A_{+}}_{i}{\rm d}{A_{\times}}_{i}\\ \end{align}

    Following, e.g. Appendix B of Haasteren et al. (2012) (and expanding out the model terms) we can put this equation into vector and matrix form to give

    \begin{equation} \label{eq:vh1}p(\mathbf{d}|h_{0},\psi_{i},I)=\int_{-\infty}^{\infty}\int_{-\infty}^{\infty}\frac{\exp{\left(-\frac{1}{2}\left(\mathbf{d}^{T}-\boldsymbol{\xi}\mathbf{M}\right)^{T}\mathbf{C}^{-1}\left(\mathbf{d}^{T}-\boldsymbol{\xi}\mathbf{M}\right)\right)}\exp{\left(-\frac{1}{2}\boldsymbol{\xi}^{T}\boldsymbol{\Sigma}_{0}^{-1}\boldsymbol{\xi}\right)}}{4\pi^{2}\sqrt{\left|\boldsymbol{\Sigma}_{0}\right|\left|\mathbf{C}\right|}}{\rm d}{A_{+}}_{i}{\rm d}{A_{\times}}_{i},\\ \end{equation}

    where we a data vector \(\mathbf{d}=\{d_{i}^{X},d_{i}^{Y}\}\), a design matrix

    \begin{equation} \mathbf{M}=h_{0}\left[\begin{array}[]{cc}{F_{+}}^{X}_{i}(\psi_{i})&{F_{+}}^{Y}_{i}(\psi_{i})\\ {F_{\times}}^{X}_{i}(\psi_{i})&{F_{\times}}_{i}^{Y}(\psi_{i})\end{array}\right]\\ \end{equation}

    and parameter vector \(\boldsymbol{\xi}=\{{A_{+}}_{i},{A_{\times}}_{i}\}\), and a data covariance matrix (and its inverse) of

    \begin{equation} \mathbf{C}=\left[\begin{array}[]{cc}\left(\sigma_{i}^{X}\right)^{2}&0\\ 0&\left(\sigma_{i}^{Y}\right)^{2}\end{array}\right],\mathrm{}\mathbf{C}^{-1}=\left[\begin{array}[]{cc}\left(\sigma_{i}^{X}\right)^{-2}&0\\ 0&\left(\sigma_{i}^{Y}\right)^{-2}\end{array}\right].\\ \end{equation}

    For the Gaussian priors on the amplitudes we use a covariance matrix (and its inverse) of

    \begin{equation} \boldsymbol{\Sigma}_{0}=\left[\begin{array}[]{cc}\sigma_{A_{+}}^{2}&0\\ 0&\sigma_{A_{\times}}^{2}\end{array}\right],\mathrm{}\boldsymbol{\Sigma}_{0}^{-1}=\left[\begin{array}[]{cc}\sigma_{A_{+}}^{-2}&0\\ 0&\sigma_{A_{\times}}^{-2}\end{array}\right].\\ \end{equation}

    This can be rearranged (to isolated the \(\boldsymbol{\xi}\) term) into the form given in equation (B2) of Haasteren et al. (2012) for which the solution is

    \begin{equation} \label{eq:margApAcmat}p(\mathbf{d}|h_{0},\psi_{i},I)=\frac{\sqrt{\left|\boldsymbol{\Sigma}\right|}}{4\pi^{2}\sqrt{\left|\boldsymbol{\Sigma}_{0}\right|\left|\mathbf{C}\right|}}\exp{\left(-\frac{1}{2}\left[\mathbf{d}^{T}\mathbf{C}^{-1}\mathbf{d}+\boldsymbol{\chi}^{T}\boldsymbol{\Sigma}^{-1}\boldsymbol{\chi}\right]\right)},\\ \end{equation}

    where

    \begin{equation} \label{eq:chi}\boldsymbol{\chi}=\left(\mathbf{M}^{T}\mathbf{C}^{-1}\mathbf{M}+\boldsymbol{\Sigma}_{0}^{-1}\right)^{-1}\left(\mathbf{M}^{T}\mathbf{C}^{-1}\mathbf{d}\right)\\ \end{equation}

    and

    \begin{equation} \boldsymbol{\Sigma}^{-1}=\mathbf{M}^{T}\mathbf{C}^{-1}\mathbf{M}+\boldsymbol{\Sigma}_{0}^{-1}.\\ \end{equation}

    For each data point, \(i\), we are not interested in the value of \(\psi_{i}\) so this can also be marginalised over. Again we require a prior, which we take as flat between 0 and \(\pi\), so

    \begin{equation} p(\psi_{i}|I)=\pi^{-1}.\\ \end{equation}

    From Equation \ref{eq:margApAcmat} we therefore get

    \begin{equation} \label{eq:margpsi}p(\mathbf{d}|h_{0},I)=\int_{0}^{\pi}p(\mathbf{d}|h_{0},\psi_{i},I)p(\psi_{i}|I){\rm d}\psi_{i}.\\ \end{equation}

    We cannot perform this integral analytically, so it must be performed numerically using the trapezium rule.

    As the values of \(h_{0}\) are common to the whole data set (i.e. all the pairs of detector data points) the values of Equation \ref{eq:margpsi} over a range of \(h_{0}\) values need to be produced. Next the product of all these likelihoods can be calculated for each of the \(h_{0}\) values and these can be used in Equations \ref{eq:fulljoint}, \ref{eq:posterior1} and \ref{eq:evidence1}.