Simulating a Gaussian stochastic gravitational wave background signal in pulsar timing arrays

TL;DR

The paper develops a unified frequency- and Fourier-domain framework for modeling a Gaussian SGWB in pulsar timing arrays by introducing transfer functions that map the SGWB power spectrum and HD spatial correlation to the observed Fourier coefficients of pulsar timing residuals. It derives explicit transfer functions for both unpolarized and circularly polarized signals, showing that the residual covariance is a convolution of the spectrum with these transfer functions, and validates the theory against standard point-source simulations and a covariance-based Gaussian simulation. The circular polarization analysis introduces Stokes parameters, anisotropy via a kinematic dipole, and cross-bin sine–cosine correlations as a distinctive signature, though practical detection remains challenging for realistic dipole speeds. Overall, the work provides a rigorous, robust framework for PTA signal modeling that highlights temporal correlations and guides future simulation-driven inference, including potential frame-independent polarization analyses.

Abstract

We revisit the theoretical modeling and simulation of a Gaussian stochastic gravitational wave background (SGWB) signal in a pulsar timing array (PTA). We show that the correlation between Fourier components of pulsar timing residuals can be expressed using transfer functions; that are indicative of characteristic temporal correlations in a SGWB signal observed in a finite time window. These transfer functions, when convolved with the SGWB power spectrum and spatial correlation (Hellings \& Downs curve), describe the variances and correlations of the pulsar timing residuals' Fourier coefficients. The convolutions are the exact frequency- and Fourier-domain representations of the time-domain covariance function. We derive explicit forms for the transfer functions for unpolarized and circularly polarized SGWB signals. We validate our results by comparing Gaussian theoretical expectation values with standard simulations based on point sources and our own covariance-matrix-based approach. The unified frequency- and Fourier-domain formalism provides a robust foundation for future PTA precision analyses and highlights the importance of temporal correlations in interpreting GW signals.

Figures (9)

• Figure 1: Representative transfer functions (\ref{['eq:alpha_filter']}-\ref{['eq:beta0_filter']}) for an SGWB signal in PTA observation with $T=15$ years, $f_1 \sim 2.1$ nanohertz, and $f_k \sim k f_1$.
• Figure 2: Simulated pulsar timing residuals due to a Gaussian SGWB signal with a power-law spectrum \ref{['eq:power_spectral_density_power_law']}; signal parameters: $A_{\rm gw}=10^{-15}$, $\gamma_{\rm gw}=13/3$, and $f_{\rm gw} \in (1, 100)$ nHz. Transfer functions (\ref{['eq:alpha_filter']}-\ref{['eq:beta_filter']}) are implemented in PTAfast. The STD shown are the average standard deviations of the timing residuals of the pulsar sample $\{a,b,c\}$; STD$\,=\{ \sqrt{ E[r(t)^2]-E[r(t)]^2 } \}\propto$ Amplitude of the time series, where $E$ is a time average operator and $\{ \cdots \}$ is the sample average over $a, b, c$. Note that $E[r(t)] = 0$ or that $\alpha_0=0$ in all cases.
• Figure 3: Fourier components of the timing residual correlation: (a) sample variances, (b) two-point correlation function, (c) cross-Fourier statistic $\{ \alpha_{k} \beta_{k} \}/\sqrt{ \{ \alpha_k^2 \} \{ \beta_k^2 \} }$, and (d) cross-Fourier spatial correlation $\{ \alpha_{a1} \beta_{b1} \}$. Green and red points refer to TEMPO2 (67 pulsars at NG positions NANOGrav:2023gor) and PTAfast (100 pulsars scattered anisotropically) simulations, respectively, and the error bars in subplots (a-c) represent the standard error across ${\cal O}\left( 10^3 \right)$ realizations. For subplot (d) the error bars shown are the error on the mean. Black points and dashed curves show the Gaussian theoretical expectation values based on Eqs. (\ref{['eq:xy_correlation']}-\ref{['eq:F_functional']}) and Eqs. (\ref{['eq:alpha_filter']}-\ref{['eq:beta0_filter']}) such as the HD curve in (b).
• Figure 4: (a) Temporal correlation in the timing residual correlation induced by a circularly polarized SGWB with a flat power spectrum, $I_{lm}(f)\sim$ Constant, $V_{lm}(f)\sim$ Constant, $T = 15$ years; $\overline{C}(t_a, t_b)$ (red solid) and $\overline{S}(t_a, t_b)$ (blue dashed) are temporal correlations associated with intensity and circular polarization, respectively, of the SGWB signal; (b) ORFs, $\gamma^I_{00}(\zeta)=\gamma^{\rm HD}(\zeta)$, $\gamma^I_{10}(\zeta)$, $\gamma^I_{11}(\zeta)$, and $\gamma^V_{11}(\zeta)$, in the computational frame for a circularly polarized SGWB Kato:2015byeSato-Polito:2021efuLiu:2022skj.
• Figure 5: Representative transfer functions for circular polarization \ref{['eq:cp_filter']} in a PTA residuals' Fourier components with $T=15$ years, $f_1 \sim 2.1$ nanohertz, and $f_k \sim k f_1$.