The global structure of the time delay likelihood

1. Problem Statement

The paper addresses a fundamental pathology in the likelihood surface for time delay inference, a critical observable in lensed-quasar cosmography (for measuring the Hubble constant, $H_0$) and reverberation mapping of active galactic nuclei.

• The Core Issue: Standard inference methods (both local optimizers/MCMC and global samplers like Nested Sampling) struggle to recover the true time delay. The likelihood function, when modeled using Gaussian Processes (GPs) for light curve variability, develops a global "W"-shaped structure.
• The Pathology:
  • There is a global maximum at the true time delay.
  • However, there are spurious modes (local maxima) at the boundaries of the observation window.
  • Crucially, the likelihood exhibits a gradual rise towards the edges of the time delay prior.
• Consequences:
  • Global Samplers (e.g., Nested Sampling): Unless configured with extremely high computational budgets (many "live points"), these algorithms are steered away from the true central peak toward the edge modes due to the gradual likelihood gradient.
  • Local Samplers/Optimizers: These tend to get trapped in local modes or, if they favor small delays implicitly, induce a systematic bias toward larger values of $H_0$.
  • Data Density Paradox: Counterintuitively, increasing the density of data points over a fixed time span sharpens the true peak but also sharpens the edge modes, making convergence failure more likely if sampling settings are not adjusted.

2. Methodology

The authors employ a combination of analytical derivation, toy model simulations, and numerical investigations using real data.

• Theoretical Framework:
  • They model intrinsic variability using stationary Gaussian Processes (GPs) with exponential kernels (damped random walk models).
  • They define the time-delayed function $h(t) = f(t + \Delta t)$ and analyze the marginal likelihood $p(y_1, y_2 | \Delta t)$.
  • Analytical Insight: They decompose the log-likelihood into a "goodness of fit" ($\chi^2$) term and a "complexity penalty" term. They demonstrate that the boundary rise is driven by the complexity penalty. As the time delay $\Delta t$ increases, the overlap between the two light curves shrinks. The GP effectively extrapolates the first light curve into the noise of the second. Since the extrapolated region has wide error bars (reverting to the prior), the data points in the non-overlapping region are statistically "ignored" (fit as noise), causing the likelihood to rise linearly with $\Delta t$.
• Toy Model Simulations:
  • Synthetic light curves were generated with known hyperparameters (amplitude, length scale $\ell$, noise $\sigma$, and true delay $\Delta t_{true}$).
  • They analyzed the data-averaged log-likelihood and its variance to confirm the "W" shape and the dominance of the edge modes in sampling.
• Numerical Investigation:
  • Nested Sampling (NS): Used PolyChord and BlackJAX implementations to sample the posterior. They systematically varied the number of live points and observation window sizes.
  • Convergence Metrics: Defined "unconverged" runs as those where the posterior mean deviated significantly from the true mode (verified via Laplace approximation).
  • Real Data Application: Applied the analysis to the lensed quasar HE 0435-1223, whitening the light curves to isolate the kernel structure.
  • Cosmological Priors: Investigated whether lens model and cosmological priors ($H_0$, $\Omega_m$) could naturally suppress the large delay modes.

3. Key Contributions

1. Identification of the "W"-Shaped Likelihood: The paper identifies that the boundary-driven rise in likelihood is an intrinsic feature of time delay inference with finite observation windows and stationary GPs, not a flaw in specific software or data. It arises because the problem is intrinsically extrapolative.
2. Mechanism Explanation: The authors clarify that the likelihood rise is not due to a better fit at the edges, but rather a reduction in the effective number of constraints. As $\Delta t$ increases, the GP treats more data points as "noise" (since they fall outside the correlation length of the shifted curve), reducing the complexity penalty and artificially inflating the likelihood.
3. Sampling Guidelines:
   • Live Points: The primary remedy is a drastic increase in the number of live points (or particles in SMC) in global samplers. The probability of initializing a live point in the narrow true mode decreases as data density increases; thus, a higher budget is required to ensure the true mode is discovered and retained.
   • Failure of Local Methods: Optimizers and local MCMC are shown to be prone to bias, often favoring small delays which translates to a positive bias in $H_0$.
4. Ineffectiveness of External Priors: The study demonstrates that standard cosmological and lens model priors (even when narrowing the $H_0$ range) are insufficient to suppress the spurious large-delay modes. The volume of the parameter space supporting large delays remains significant.
5. Regularized Likelihood Proposal: Theoretically, one could construct a "counterterm" prior to cancel the edge rise, but the authors argue this is impractical as it requires prior knowledge of the true delay and risks over-regularization. They advocate for better sampling strategies instead.

4. Key Results

• Likelihood Structure: The log-likelihood $L(\Delta t)$ is globally "W"-shaped. The true delay is the global maximum, but the likelihood rises towards the boundaries ($\pm T_{range}$).
• Data Density Effect: Increasing the number of data points ($N_{data}$) for a fixed time span sharpens the central peak but also increases the slope of the rise toward the edges. This makes the "W" shape more pronounced and convergence harder for standard settings.
• Convergence Failure:
  • With standard settings (e.g., 75 live points), Nested Sampling frequently fails to find the true mode, instead aggregating at the edges.
  • Increasing live points by a factor of 10 (e.g., to 250+) significantly reduces the fraction of unconverged runs.
  • Sequential Monte Carlo (SMC) and Parallel Tempering (PT) exhibit similar failure modes, driven by the same volume-dominated edge modes.
• Real Data Confirmation: Analysis of HE 0435-1223 confirms the "W" shape exists in real data, with edge likelihoods comparable to the central peak.
• Bias in $H_0$: Methods that implicitly favor small delays (or fail to resolve the true mode) result in a systematic bias toward larger $H_0$ values.

5. Significance and Implications

• Robust Cosmology: As surveys like LSST (Vera C. Rubin Observatory) increase the number of known lensed quasars by orders of magnitude, the demand for automated, reliable time delay inference will surge. This paper warns that standard pipelines will fail without specific adaptations.
• Algorithmic Requirements: The paper establishes that fully Bayesian inference for time delays requires rigorous convergence checks and significantly higher computational resources (more live points) than typically used. It challenges the assumption that global samplers automatically solve multimodal problems.
• Model Independence: While derived using GPs, the authors argue the pathology is model-agnostic, stemming from the finite observation window and the need for extrapolation.
• Future Work: The findings necessitate a shift in how time delay analyses are conducted, moving away from reliance on visual matching or simple optimization toward robust, high-budget global sampling with strict diagnostics.

In summary, Kroupa and Handley reveal that the difficulty in measuring time delays is not merely a data quality issue but a fundamental statistical challenge arising from the geometry of the likelihood surface. They provide the theoretical explanation and practical guidelines (specifically regarding live point budgets) necessary to achieve unbiased, robust cosmological measurements.