Solving the PTA Data Analysis Problem with a Global Gibbs Scheme

The announcement in the summer of 2023 about the discovery of evidence for a gravitational wave background (GWB) using pulsar timing arrays (PTAs) has ignited both the PTA and the larger scientific community's interest in the experiment and the scientific implications of its findings. As a result, numerous scientific works have been published analyzing and further developing various aspects of the experiment, from performing tests of gravity to improving the efficiency of the current data analysis techniques. In this regard, we contribute to the recent advancements in the field of PTAs by presenting the most general, agnostic, per-frequency Bayesian search for a low-frequency (red) noise process in these data. Our new method involves the use of a conjugate Jeffrey's-like multivariate prior which allows one to model all unique parameters of the global PTA-level red noise covariance matrix as a separate model parameter for which a marginalized posterior-probability distribution can be found using Gibbs sampling. Even though perfecting the implementation of the Gibbs sampling and mitigating the numerical stability challenges require further development, we show the power of this new method by analyzing realistic and theoretical PTA simulated data sets. We show how our technique is consistent with the more restricted standard techniques in recovering both the auto and the cross-spectrum of pulsars' low-frequency (red) noise. Furthermore, we highlight ways to approximately characterize a GWB (both its auto- and cross-spectrum) using Fourier coefficient estimates from single-pulsar and so-called CURN (common uncorrelated red noise) analyses via analytic draws from a specific Inverse-Wishart distribution.