🔭 astrophysics

Generalized Non-linear Bayesian Pulsar Timing with Enterprise

This study employs generalized non-linear Bayesian timing methods within the Enterprise package to analyze four pulsars, yielding refined mass constraints and providing the first evidence of intrinsic red noise in PSR J2043+1711 while demonstrating how physical priors improve the modeling of noise-parameter interplay.

Original authors: Andrew R. Kaiser, Jeffrey S. Hazboun, Maura A. McLaughlin, H. Thankful Cromartie, Emmanuel Fonseca, Joseph Simon, Stephen R. Taylor, Michele Vallisneri, Sarah J. Vigeland, Zaven Arzoumanian, Paul T. B

Published 2026-08-19

📖 8 min read🧠 Deep dive

CC BY 4.0

Original authors: Andrew R. Kaiser, Jeffrey S. Hazboun, Maura A. McLaughlin, H. Thankful Cromartie, Emmanuel Fonseca, Joseph Simon, Stephen R. Taylor, Michele Vallisneri, Sarah J. Vigeland, Zaven Arzoumanian, Paul T. Baker, Harsha Blumer, Paul R. Brook, Ismael Cognard, Megan E. DeCesar, Paul B. Demorest, Timothy Dolch, F. Adam Dong, Justin A. Ellis, Robert D. Ferdman, Elizabeth C. Ferrara, William Fiore, Nate Garver-Daniels, Peter A. Gentile, Deborah C. Good, Lucas Guillemot, Ross J. Jennings, Megan L. Jones, David L. Kaplan, Victoria M. Kaspi, Matthew Kerr, Aida Yu. Kirichenko, Michael T. Lam, Duncan R. Lorimer, Jing Luo, Ryan S. Lynch, Alexander McEwen, James W. McKee, Natasha McMann, Bradley W. Meyers, Arun Naidu, Cherry Ng, David J. Nice, Aditya Parthasarathy, Timothy T. Pennucci, Benetge B. P. Perera, Nihan S. Pol, Henri A. Radovan, Scott M. Ransom, Paul S. Ray, Brent J. Shapiro-Albert, Renée Spiewak, Ingrid H. Stairs, Kevin Stovall, Joseph K. Swiggum, Chia Min Tan, Shriharsh P. Tendulkar, Haley M. Wahl, WeiWei Zhu

Original paper licensed under CC BY 4.0 (http://creativecommons.org/licenses/by/4.0/). ✨ This is an AI-generated explanation of the paper below. It is not written or endorsed by the authors. For technical accuracy, refer to the original paper. Read full disclaimer

Deep in the quiet of the cosmos, spinning like lighthouses in the dark, are neutron stars so dense that a single teaspoon of their material would weigh a billion tons on Earth. Some of these stars, known as pulsars, spin with a precision that rivals the atomic clocks on our best laboratory benches, emitting beams of radio waves that sweep past Earth with clockwork regularity. When astronomers catch these signals, they can measure the exact moment a pulse arrives, down to a fraction of a microsecond. This incredible precision turns these distant stars into cosmic laboratories, allowing scientists to test the fundamental laws of physics, including how gravity behaves in the most extreme environments imaginable. However, to use these stars as tools, scientists must first build a perfect model of their motion, accounting for everything from their spin to the gravity of any companion stars they might have. If the model is slightly off, or if the data contains hidden noise, the measurements of the star's properties can drift, leading to incorrect conclusions about the nature of matter itself.

A team of researchers has now developed a new way to build these models, moving beyond the traditional methods that have been used for decades. Instead of relying on a simplified, straight-line approximation of how a pulsar moves, they used a powerful, flexible approach that treats the timing data and the noise within it as a single, interconnected system. By applying this method to four specific pulsars, they found that their new technique often provides tighter, more accurate constraints on the stars' masses and orbits than previous studies. In some cases, they discovered that the stars were spinning slightly differently than previously thought, and for the first time, they found clear evidence of a specific type of internal noise in one of the pulsars that had been hiding in the data. This work suggests that by letting the computer explore the full complexity of the data rather than forcing it into a simpler shape, we can learn more about the heaviest objects in the universe.

The researchers focused their efforts on four binary pulsars, which are neutron stars orbiting a companion object, often a white dwarf. Because these stars are locked in a gravitational dance, their timing signals carry the imprint of their orbital motion. By analyzing how the pulses arrive, scientists can measure the mass of the pulsar and its companion. This is crucial because the mass of a neutron star tells us about the state of matter at densities that cannot be recreated on Earth. If a neutron star is too heavy, it might collapse into a black hole; if it is light enough, it reveals how nuclear forces hold up under such pressure. For years, scientists have used a method called generalized least-squares fitting to determine these masses. This technique finds the best set of parameters that minimizes the difference between the observed pulse times and the predicted times. However, this method assumes that any remaining errors, or residuals, are small and linear, meaning they can be treated as simple corrections to the main model.

The new study, led by Andrew Kaiser and colleagues, challenges this assumption by using a Bayesian approach within a software package called Enterprise. In this framework, the researchers do not just find the single best fit; they explore the entire landscape of possible solutions, weighing how likely each one is given the data. They tested three different ways of handling the timing model. The first method, which mirrors the traditional approach, treats the timing parameters as linear corrections that are mathematically removed from the data before analyzing the noise. The second method, which is the most computationally demanding, treats the timing model as a fully general, non-linear signal, recalculating the entire model at every step of the analysis. The third method is a hybrid, combining the two. By comparing these approaches, the team wanted to see if the traditional linear assumption was hiding important details or distorting the results.

When they applied these methods to the pulsar PSR J2043+1711, they found something unexpected. Previous studies had not found significant evidence of intrinsic red noise, which is a type of slow, random fluctuation in the pulsar's spin. However, using their fully general Bayesian method, the team detected this noise for the first time. They found that the pulsar's spin was not just perfectly steady but contained a subtle, low-frequency jitter. This discovery changed the picture of the pulsar's timing parameters, shifting the estimated values for its rotation speed and its companion's mass. The study showed that when this noise is ignored or handled with a simplified model, the timing parameters can shift to compensate, leading to slightly different conclusions about the star's physical properties.

For the pulsar PSR J1600–3053, the team found that their new method provided tighter constraints on the mass of the star. Previous measurements had suggested a mass of about 2.0 solar masses, but with a large margin of error. The new analysis narrowed this range significantly, providing a more precise estimate. Similarly, for PSR J0740+6620, a pulsar known for having one of the highest masses ever measured, the team's results were consistent with previous findings but offered even sharper constraints. They confirmed that this star has a mass of about 2.06 times that of our Sun, with a very small margin of error. This precision is vital because it helps rule out certain theories about what happens to matter inside a neutron star. If the mass is too high, it would break the laws of physics as we currently understand them, suggesting that our theories about the strong nuclear force need to be revised.

The most dramatic results came from the pulsar PSR J1640+2224. In previous studies, this pulsar had yielded a mass estimate that was surprisingly high, sometimes suggesting a value greater than three times the mass of the Sun, which would be physically impossible for a neutron star. The researchers suspected that the high mass was an artifact of the modeling method rather than a real property of the star. To test this, they imposed a physical limit on the mass in their model, restricting it to be less than three solar masses. When they did this, the estimated mass dropped to a more reasonable value of about 2.2 solar masses. This result highlights a critical point: if the mathematical model is not flexible enough to account for the complex interplay between the timing parameters and the noise, it can produce physically impossible results. By allowing the model to explore the full range of possibilities and then applying physical constraints, the team was able to recover a more realistic picture of the star.

The study also addressed the issue of how to handle the vast amount of data collected by radio telescopes. Pulsar timing arrays collect data from multiple telescopes over many years, and this data is often contaminated by noise from the interstellar medium, the space between the stars. The researchers found that their fully general method was better at distinguishing between the true signal of the pulsar's motion and the noise from the environment. In many cases, the traditional linear method absorbed some of the noise into the timing parameters, making the star's motion look slightly different than it actually was. By treating the timing model and the noise as a coupled system, the new method could separate them more effectively, leading to more accurate measurements.

One of the key takeaways from this work is that the choice of statistical method matters. The researchers found that while the traditional linear approach works well for many pulsars, it can fail when the data is complex or when the noise is significant. The fully general Bayesian method, while more computationally expensive, provides a more robust way to analyze the data. It allows scientists to see the full range of possible solutions and to understand how different parameters are related to one another. This is particularly important for future searches for gravitational waves, which rely on the precise timing of many pulsars. If the timing models are not accurate, the subtle ripples in spacetime caused by gravitational waves could be missed or misinterpreted.

The team also looked at how different data sets affected the results. They analyzed data from the NANOGrav collaboration, which includes observations spanning over a decade. They found that as more data is collected, the constraints on the pulsar masses become tighter, but the central values can sometimes shift. This suggests that the timing models are still evolving and that the best-fit values from earlier studies may not be the final answer. The researchers hope that with even more data, the models will stabilize, and the measurements will become consistent across different data sets.

In the end, this paper is about refining our tools for listening to the universe. By moving away from simplified assumptions and embracing the full complexity of the data, the researchers have shown that we can learn more about the most extreme objects in the cosmos. Their work provides a new standard for how pulsar timing should be done, ensuring that the measurements of neutron star masses and orbits are as accurate as possible. This is not just a technical improvement; it is a step toward understanding the fundamental nature of matter and gravity. As we continue to listen to these cosmic clocks, the methods developed in this study will help us hear the faintest whispers of the universe, from the spin of a single neutron star to the collision of black holes billions of light-years away.

Technical Summary: Generalized Non-linear Bayesian Pulsar Timing with Enterprise

Problem Statement
Traditional pulsar timing analysis relies on Generalized Least Squares (GLS) fitting, where timing parameters are estimated by minimizing $\chi^2$ between a model and observed Time of Arrivals (TOAs). In Pulsar Timing Array (PTA) searches for gravitational waves (GWs), the standard approach assumes that any remaining timing residuals can be represented by linear corrections to the original timing solutions. This allows for the analytical marginalization of "nuisance" timing parameters (such as spin frequency, position, and binary parameters) to reduce the dimensionality of the parameter space. However, this linearization assumes that deviations from the best-fit solution remain within a linear regime. This assumption may fail for parameters that are highly non-linear or when noise is incorrectly absorbed into timing parameters, potentially hiding combinations of parameters that do not obey the linear approximation. Furthermore, the interplay between intrinsic red noise (RN) and timing model parameters is often not fully explored in standard analyses.

Methodology
The authors implement a fully general, non-linear Bayesian timing framework using the Enterprise software package and its extensions. They analyze four binary millisecond pulsars: PSR J1600–3053, PSR J2043+1711, PSR J0740+6620, and PSR J1640+2224, using data from the NANOGrav 5-yr, 9-yr, and 12.5-yr data releases, as well as specific datasets from Fonseca et al. (2021) and Cromartie et al. (2020).

The study compares three distinct modeling approaches:

Model A: The standard approach where all timing model parameters are analytically marginalized and linearized around the GLS best-fit values.
Model B: A fully general approach where the timing model is treated as a deterministic signal. Timing residuals are recalculated at every Markov Chain Monte Carlo (MCMC) step using libstempo (a wrapper for TEMPO2) and PINT, allowing for the direct numerical sampling of all timing parameters, including non-linear binary parameters.
Model C: A hybrid approach where non-linear parameters (e.g., binary parameters) are sampled directly, while linear "nuisance" parameters (e.g., DMX, frequency-dependent terms, JUMPs) are analytically marginalized.

The authors employ uniform priors for most parameters but impose physical boundaries on eccentricity, inclination, and companion mass. For PSR J1640+2224, they specifically test the effect of restricting the pulsar mass prior to $m_p < 3 M_\odot$ . They also investigate the inclusion of intrinsic red noise (RN) characterized by amplitude ( $A$ ) and spectral index ( $\gamma$ ) to assess its impact on timing parameter posteriors.

Key Contributions and Results

Validation of Generalized Timing: The authors find generally good agreement between the maximum a posteriori (MAP) values of their fully general Bayesian methods (Models B and C) and the best-fit values from GLS and previous Bayesian studies. In many cases, the generalized methods provide more constrained parameter estimations.
Pulsar Mass Constraints: Using the most up-to-date data releases, the authors report tighter constraints on pulsar masses compared to published results:
- PSR J2043+1711: $m_p = 1.6(1) M_\odot$ (NANOGrav 12.5-yr).
- PSR J1600–3053: $m_p = 2.3^{+0.9}_{-0.7} M_\odot$ (NANOGrav 12.5-yr).
- PSR J0740+6620: $m_p = 2.06(6) M_\odot$ (Fonseca et al. 2021 data), which improves upon previous constraints.
Intrinsic Red Noise Detection: For the first time, the study finds evidence for intrinsic red noise in PSR J2043+1711. The inclusion of RN in the models shifts the timing parameter posteriors (particularly spin frequency $\nu$ and its derivative $\dot{\nu}$ ) away from the GLS best-fit values, demonstrating the covariance between noise and timing parameters.
Impact of Physical Priors on PSR J1640+2224: This pulsar has historically yielded mass estimates $>3 M_\odot$ . When the authors restrict the prior to $m_p < 3 M_\odot$ , the derived mass shifts to $m_p = 2.2(5) M_\odot$ (NANOGrav 12.5-yr). Without this restriction, the mass estimates remain high ( $m_p \approx 5 M_\odot$ ), highlighting the sensitivity of the solution to prior choices and the non-linear nature of the parameter space.
Model Comparison: The study demonstrates that while Model C (hybrid) often converges faster than Model B (fully general), Model B can provide tighter constraints on timing parameters. The authors note that for PSR J1640+2224, Model B struggled to converge fully within computational limits, suggesting that high-dimensional, non-linear spaces can be challenging for MCMC sampling without sufficient resources.

Significance
The paper argues that fully general Bayesian timing is a robust method for modeling the interplay between intrinsic noise and timing parameters, particularly in cases where linear approximations may be insufficient. By simultaneously sampling timing and noise parameters, the method can better account for noise absorption into timing models, leading to more accurate mass estimates and noise characterization. The authors conclude that these techniques are essential for future PTA analyses, as they allow for a more rigorous treatment of covariances between timing parameters and potential GW signals, ultimately improving the sensitivity of gravitational wave detection. The work emphasizes that timing models are evolving and that solutions in highly covariant spaces require careful handling of priors and non-linear effects to ensure stability and accuracy.

More like this