gwin.models.marginalized_gaussian_noise module

This module provides model classes that assume the noise is Gaussian and allows for the likelihood to be marginalized over phase and/or time and/or distance.

class gwin.models.marginalized_gaussian_noise.MarginalizedGaussianNoise(variable_params, data, waveform_generator, f_lower, psds=None, f_upper=None, norm=None, time_marginalization=False, distance_marginalization=False, phase_marginalization=False, marg_prior=None, **kwargs)

Bases: gwin.models.gaussian_noise.GaussianNoise

The likelihood is analytically marginalized over phase and/or time and/or distance.

For the case of marginalizing over phase, the signal, can be written as:

\[\tilde{h}(f; \Theta, \phi) = A(f; \Theta)e^{i\Psi(f; \Theta) + i \phi},\]

where \(\phi\) is an arbitrary phase constant. This phase constant can be analytically marginalized over with a uniform prior as follows: assuming the noise is stationary and Gaussian (see GaussianNoise for details), the posterior is:

\[\begin{split}p(\Theta,\phi|d) &\propto p(\Theta)p(\phi)p(d|\Theta,\phi) \\ &\propto p(\Theta)\frac{1}{2\pi}\exp\left[ -\frac{1}{2}\sum_{i}^{N_D} \left< h_i(\Theta,\phi) - d_i, h_i(\Theta,\phi) - d_i \right>\right].\end{split}\]

Here, the sum is over the number of detectors \(N_D\), \(d_i\) and \(h_i\) are the data and signal in detector \(i\), respectively, and we have assumed a uniform prior on \(phi \in [0, 2\pi)\). With the form of the signal model given above, the inner product in the exponent can be written as:

\[\begin{split}-\frac{1}{2}\left<h_i - d_i, h_i- d_i\right> &= \left<h_i, d_i\right> - \frac{1}{2}\left<h_i, h_i\right> - \frac{1}{2}\left<d_i, d_i\right> \\ &= \Re\left\{O(h^0_i, d_i)e^{-i\phi}\right\} - \frac{1}{2}\left<h^0_i, h^0_i\right> - \frac{1}{2}\left<d_i, d_i\right>,\end{split}\]

where:

\[\begin{split}h_i^0 &\equiv \tilde{h}_i(f; \Theta, \phi=0); \\ O(h^0_i, d_i) &\equiv 4 \int_0^\infty \frac{\tilde{h}_i^*(f; \Theta,0)\tilde{d}_i(f)}{S_n(f)}\mathrm{d}f.\end{split}\]

Gathering all of the terms that are not dependent on \(\phi\) together:

\[\alpha(\Theta, d) \equiv \exp\left[-\frac{1}{2}\sum_i \left<h^0_i, h^0_i\right> + \left<d_i, d_i\right>\right],\]

we can marginalize the posterior over \(\phi\):

\[\begin{split}p(\Theta|d) &\propto p(\Theta)\alpha(\Theta,d)\frac{1}{2\pi} \int_{0}^{2\pi}\exp\left[\Re \left\{ e^{-i\phi} \sum_i O(h^0_i, d_i) \right\}\right]\mathrm{d}\phi \\ &\propto p(\Theta)\alpha(\Theta, d)\frac{1}{2\pi} \int_{0}^{2\pi}\exp\left[ x(\Theta,d)\cos(\phi) + y(\Theta, d)\sin(\phi) \right]\mathrm{d}\phi.\end{split}\]

The integral in the last line is equal to \(2\pi I_0(\sqrt{x^2+y^2})\), where \(I_0\) is the modified Bessel function of the first kind. Thus the marginalized log posterior is:

\[\log p(\Theta|d) \propto \log p(\Theta) + I_0\left(\left|\sum_i O(h^0_i, d_i)\right|\right) - \frac{1}{2}\sum_i\left[ \left<h^0_i, h^0_i\right> - \left<d_i, d_i\right> \right]\]

For the case of marginalizing over distance, the signal can be written as,

\[\tilde{h}_{j} = \frac{1}{D} \tilde{h}_{j}^{0}\]

The distance can be analytically marginalized over with a uniform prior as follows: assuming the noise is stationary and Gaussian (see GaussianNoise for details), the likelihood is:

\[log L = -\frac{1}{2}\left<d-h|d-h\right>\]

We see that :math: <h|h> is inversely proportional to distance squared and :math: <h|d> is inversely proportional to distance. The log likelihood is therefore

\[\log L = -\frac{1}{2}\left<d|d\right> - \frac{1}{2D^{2}} \left<h|h\right> + \frac{1}{D}\left<h|d\right>\]

Consequently, the likelihood marginalised over distance is simply

\[\log L = \log\left(\int_{0}^{D}{L p(D) dD}\right)\]

If we assume a flat prior

\[\log L = \log\left(\int_{0}^{D}{\exp{\log L} dD}\right)\]

For the case of marginalizing over time, the signal can be written as,

\[\tilde{h}_{j} = \tilde{h}_{j}^{0} \exp\left(-2\pi ij\Delta ft\right)\]

The time can be analytically marginalized over with a uniform prior as follows: assuming the noise is stationary and Gaussian (see GaussianNoise for details), the likelihood is:

\[\log L = -\frac{1}{2}\left<d-h|d-h\right>\]

We note that :math: <h|h> and :math: <d|d> are time independent while :math: <d|h> is dependent of time

\[\left<d|h\right>(t) = 4\Delta f\sum_{j=0}^{N/2} \frac{\tilde{d}_{j}^{*} \tilde{h}_{j}^{0}}{S_{j}} \exp\left(-2\pi ij \Delta f t\right)\]

For integer timesteps :math: t=kDelta t

\[ \begin{align}\begin{aligned}\left<d|h\right>(k\Delta t) = 4\Delta f\sum_{j=0}^{N/2} \frac{ \tilde{d}_{j}^{*}\tilde{h}_{j}^{0}} {S_{j}}\exp(-2\pi \frac{ijk}{N}\\\left<d|h\right>(k\Delta t) = 2\Delta f\sum_{j=0}^{N} \frac{ \tilde{d}_{j}^{*}\tilde{h}_{j}^{0}} {S_{j}} \exp(-2\pi \frac{ijk}{N}\end{aligned}\end{align} \]

Using a FFT, this expression can be evaluated efficiently for all :math: k

\[\left<d|h\right>(k\Delta t) = 2\Delta f FFT_{k} (\frac{dh}{S})\]

since :math:: left<h|dright> = left<d|hright>^{*},

\[\left<d|h\right> + \left<h|d\right> = 4\Delta f FFT_{k} (\frac{dh}{S})\]

and so the likelihood marginalised over time is simply

\[\log{L} = \log\left(\int_{0}^{T} np.exp(\np.log(L) p(t))\right)\]

where p(t) is the prior. If we assume a flat prior then,

\[\log{L} = \log\left(\int_{0}^{T} np.exp(\np.log(L))\right)\]

Parameters:

• time_marginalization (bool, optional) – A Boolean operator which determines if the likelihood is marginalized over time
• phase_marginalization (bool, optional) – A Boolean operator which determines if the likelihood is marginalized over phase
• distance_marginalization (bool, optional) – A Boolean operator which determines if the likelihood is marginalized over distance
• marg_prior (list, optional) – An instance of pycbc.distributions which returns a list of prior distributions to be used when marginalizing the likelihood
• **kwargs – All other keyword arguments are passed to GaussianNoise.

classmethod from_config(cp, data, delta_f=None, delta_t=None, gates=None, recalibration=None, **kwargs)

Initializes an instance of this class from the given config file.

Parameters:

• cp (WorkflowConfigParser) – Config file parser to read.
• data (dict) – A dictionary of data, in which the keys are the detector names and the values are the data. This is not retrieved from the config file, and so must be provided.
• delta_f (float) – The frequency spacing of the data; needed for waveform generation.
• delta_t (float) – The time spacing of the data; needed for time-domain waveform generators.
• recalibration (dict of pycbc.calibration.Recalibrate, optional) – Dictionary of detectors -> recalibration class instances for recalibrating data.
• gates (dict of tuples, optional) – Dictionary of detectors -> tuples of specifying gate times. The sort of thing returned by pycbc.gate.gates_from_cli.
• **kwargs – All additional keyword arguments are passed to the class. Any provided keyword will over ride what is in the config file.

name = 'marginalized_gaussian_noise'