Bayesian Glitch Subtraction (BayesWave)

• Bayesian Glitch Subtraction (BayesWave) is a data analysis method that models both gravitational-wave signals and non-Gaussian glitches within a unified Bayesian framework.
• It employs wavelet-based glitch templates and trans-dimensional MCMC to dynamically adjust model complexity based on the data.
• This approach enhances signal recovery and reduces false alarms by jointly handling deterministic signals, glitches, and stochastic Gaussian noise.

Bayesian Glitch Subtraction (BayesWave) is a class of data analysis techniques for interferometric gravitational-wave (GW) detectors that simultaneously models astrophysical signals and non-Gaussian noise transients (“glitches”) within a unified Bayesian inference framework. BayesWave achieves robust glitch subtraction and signal recovery by decomposing the data into a deterministic signal model, a flexible non-Gaussian noise model using wavelets as glitch templates, and a stochastic (possibly non-stationary) Gaussian background. The number and properties of wavelet-based glitch templates are not fixed, but are determined dynamically using trans-dimensional Markov chain Monte Carlo (RJMCMC), enabling the model’s complexity to adapt to the data. This approach enables non-parametric mitigation of transient noise artifacts and supports high-fidelity recovery of gravitational-wave signals, especially in real, glitch-prone data from advanced detectors such as LIGO and Virgo.

1. Theoretical Framework and Motivation

The motivation for Bayesian glitch subtraction arises from limitations of standard methods in gravitational-wave data analysis. Matched filtering is statistically optimal when noise is stationary and Gaussian, but in practice, ground-based detectors exhibit a significant rate of non-Gaussian, non-stationary transients (glitches) that produce long tails in background distributions and can drive up detection thresholds or result in false alarms.

BayesWave introduces a fully Bayesian model where the observed data $d$ in each detector are represented as:

$d = s + g + n_G$

• $s$: The astrophysical gravitational-wave signal, coherently projected onto the network.
• $g$: A deterministic, non-Gaussian glitch component, modeled as a sum of wavelets.
• $n_G$: Colored Gaussian noise, possibly non-stationary, with a time–frequency-dependent power spectrum.

Unlike approaches that attempt to veto or exclude glitch-contaminated data, BayesWave explicitly models these features and fits both signal and noise (glitch + Gaussian) components simultaneously. This approach mitigates contamination and preserves maximal data sensitivity.

The likelihood function in this framework is constructed as:

$$p(d|\theta, M) = \exp\left[ -\frac{1}{2} (d - m(\theta, M) | d - m(\theta, M)) \right]$$

where $m(\theta, M)$ is the model of the data under hypothesis $M$ (signal+noise+glitch, etc.), and the inner product is defined with respect to the local noise power spectral density.

2. Wavelet-Based Noise and Glitch Decomposition

BayesWave employs wavelet transforms, specifically the discrete wavelet transform (DWT) or a continuous frame using Morlet-Gabor (sine–Gaussian) wavelets, to decompose data into localized time–frequency "pixels" or segments. Each basis function $\psi_{ij}$ has compact support in both time and frequency, making it well-suited to the localized structure of glitches.

For each time–frequency pixel, the noise model parameterizes the local noise amplitude and, if necessary, allows for non-stationary fluctuations. The glitch component is represented as a sum over selected wavelets:

$$g(t) = \sum_{n=1}^{n_G} A_n\, \psi_n(t; f_{0,n}, Q_n, t_{0,n}, \varphi_{0,n})$$

where each wavelet is defined by its amplitude $A_n$, central frequency $f_{0,n}$, quality factor $Q_n$, central time $t_{0,n}$, and phase $\varphi_{0,n}$. The number of active wavelets $n_G$ is unknown a priori and is itself a model parameter (see next section).

By working in a wavelet basis, BayesWave naturally captures the semi-coherent structure of transient noise excursions that appear as clusters in the time–frequency plane. The approach extends seamlessly to non-stationary noise by allowing pixelwise (or blockwise) variations in noise amplitude.

3. Explicit Glitch Modeling via Trans-Dimensional MCMC

A central innovation of BayesWave is the use of Reversible Jump Markov Chain Monte Carlo (RJMCMC) for trans-dimensional inference. The RJMCMC algorithm enables the sampler to move between models with different numbers of glitch wavelets; thus, the effective dimension of the glitch model is not fixed but determined by the data.

Consider models $\mathcal{M}_i$ with $n_G$ glitch wavelets and $\mathcal{M}_j$ with $n_G + 1$. The move acceptance probability is given by the Hastings ratio:

$$H_{\mathcal{M}_i \to \mathcal{M}_j} = \frac{p(d|\theta_j, \mathcal{M}_j)\,p(\theta_j|\mathcal{M}_j)\,q(\theta_i|\mathcal{M}_i)}{p(d|\theta_i, \mathcal{M}_i)\,p(\theta_i|\mathcal{M}_i)\,q(\theta_j|\mathcal{M}_j)}\,|J_{ij}|$$

where $q$ denotes the proposal distributions, and $|J_{ij}|$ is the Jacobian for dimensional changes. This construction enforces a Bayesian Occam penalty: models must justify additional complexity (more glitch components) through improved likelihood.

Priors for the number and placement of wavelets are typically uniform over all combinations for a given $n_G$, e.g.,

$p(n_G|G_1) = \frac{1}{\binom{N}{n_G}}$

with $N$ the total number of wavelet pixels. Amplitudes have a prior $A_n \sim \mathcal{N}[0, 100 S_{ij}]$ to focus on significant transient power.

The algorithm thus adaptively “fits” just enough glitch structure to describe the excess transient noise, without overfitting random fluctuations or absorbing actual signals.

4. Model Selection and Performance Gains

Traditional searches either veto glitchy data or assume a fixed Gaussian background, practices that can lower sensitivity or result in elevated false-alarm rates when non-Gaussian tails dominate. BayesWave directly compares several alternative composite hypotheses:

• Gaussian noise only ("$N_0$")
• Gaussian noise plus detector-specific glitches ("$G_1$")
• Gaussian noise plus coherent astrophysical signal (and possibly glitches, "$S+G_1$")

Bayesian evidence (marginal likelihood) is numerically computed for each, with model selection naturally penalizing complexity. Strong discrimination between glitches and signals is achieved by exploiting the relative coherence or incoherence across the detector network.

Empirically, the BayesWave framework permits reliable recovery of GW signals at SNRs as low as $\sim 8$ in data containing realistic glitch backgrounds. In contrast, analyses that do not explicitly model glitches may either suffer false positives or lose sensitivity to such signals (Littenberg et al., 2010).

Residuals after glitch subtraction are demonstrably more consistent with the Gaussian noise hypothesis, supporting both lower detection thresholds and unbiased astrophysical parameter estimation.

5. Scaling to Detector Networks and Multimodal Applications

BayesWave is designed for multi-detector data. The mechanism exploits the physical coherence of GW signals (which appear consistently across interferometer network, up to time and polarization projection) versus the detector-local nature of glitches (often uncorrelated).

For $M$ detectors, the data vector is modeled jointly as:

$$\vec{d}(t) = [R_1 * h(t), \ldots, R_M * h(t)] + [g_1(t), \ldots, g_M(t)] + [n_{G1}(t), \ldots, n_{GM}(t)]$$

Multi-detector implementation ensures that astrophysical signals are preferred by Bayesian evidence if and only if the observed excess power is coherent in the way expected from gravitational waves, greatly reducing the rate of spurious glitch-induced triggers.

The framework is extensible: by adapting wavelet or signal models (e.g., generic burst, chirp, or polarizations), it can be used not only for compact binary coalescences, but for unmodeled bursts, non-GR polarization content, and advanced detector configurations.

6. Implications and Future Prospects

Integration of explicit glitch modeling via RJMCMC and wavelet decomposition enables robust, non-parametric cleaning of real interferometric GW data, crucially reducing detection thresholds in the presence of realistic, glitch-prone backgrounds. This approach paves the way for:

• Greater sensitivity to lower-SNR signals as instrumental artifacts are marginalized, not simply discarded.
• Quantifiable signal-versus-glitch discrimination, reducing misclassification rates and affording high-confidence detections (Littenberg et al., 2015).
• Automation of analysis pipelines: the Bayesian scheme is fully specified once likelihoods and priors are defined, minimizing subjective tuning.
• Incorporation into pipelines for multi-detector, low-latency, and multi-messenger GW astronomy, with the theoretical flexibility to include future detector upgrades and new noise models.

Limitations remain with very low-SNR, strongly overlapping or highly structured glitches (e.g., scattered light) when the prior penalizes low-amplitude or highly complex wavelet combinations. Ongoing research focuses on improving sampling efficiency, adapting priors, and integrating machine learning–informed glitch models.

Bayesian glitch subtraction, as pioneered by BayesWave and its descendants, represents a foundational advancement in gravitational-wave data analysis and continues to inform the development of robust inference strategies in the presence of non-Gaussian, non-stationary noise.