4. Cosmological Consistency Checks

Gravitational lensing time delays provide a geometric probe of TEP at cosmological scales. Under TEP, the time delay between multiply-imaged quasars should exhibit scale-dependent structure—"temporal shear"—that varies with the measurement timescale. Analysis of the COSMOGRAIL dataset reveals evidence for this effect in the lens system DESJ0408, where scale-dependent delays are detected consistent with the screened soliton framework established in Section 2. These results are presented as a complementary time-domain constraint alongside the pulsar signal.

4.1 Data and Methods

The analysis uses publicly available light curves from the COSMOGRAIL collaboration, which has monitored gravitationally lensed quasars for over a decade. The sample includes:

4.1.1 Multi-Scale Delay Estimation

For each image pair, time delays are estimated at multiple variability timescales $\tau_w \in \{20, 40, 80, 160\}$ days. The analysis pipeline makes three critical methodological choices to ensure robustness against common lensing artifacts:

The pipeline implementation follows these steps:

1. Microlensing Detrending: Gaussian smoothing (σ = 200 days) removes slow trends uncorrelated between images.
2. Broadband Anchor: Cross-correlation of detrended light curves establishes the reference delay ($\Delta t_0$).
3. Bandpass Filtering: Difference-of-Gaussians isolates variability at each specific timescale $\tau_w$.
4. Constrained Search: The delay $\Delta t(\tau_w)$ is found by maximizing the CCF within $\Delta t_0 \pm 50$ days.

4.1.2 Temporal Shear Fitting

TEP predicts that enhanced gravitational time dilation induces a scale-dependent time delay, varying with the quasar variability timescale used for cross-correlation:

The quantity Γ—here termed "temporal shear"—measures the change in apparent time delay per decade of variability timescale. Under GR, Γ = 0. Under TEP, Γ ≠ 0.

Significance is assessed via σ = |Γ|/σ_Γ, where σ_Γ is the residual-based standard error. This empirical error estimate accounts for the correlated nature of the delay measurements across timescales.

With the estimator and error model defined, this pipeline is now applied to the COSMOGRAIL sample to search for the predicted temporal shear signature.

4.2 Results: Temporal Shear Measurements

4.3 Validation: Broadband Delays Match Literature

A critical validation is that the broadband (unfiltered) delay estimates match published literature values, indicating that the cross-correlation methodology is sound:

All broadband delays agree with Courbin et al. (2017) within 1σ. This consistency suggests that the temporal shear signal represents a real scale-dependent variation around the correct mean delay, rather than a systematic bias in the delay estimation.

4.4 Methodological Validation

4.4.1 Estimator Independence: ICCF vs Interpolation

A critical test of robustness is whether results depend on the choice of delay estimator. Two distinct estimator families commonly used in time-delay cosmography were compared:

• ICCF (Interpolated Cross-Correlation Function): Evaluates correlation on observed data points only. Robust to gaps.
• Linear Interpolation (Python numpy.interp): Linearly interpolates light curves onto a regular grid before correlation. Prone to gap artifacts.

The estimators agree within 1.1σ at the meta-analysis level. The high correlation (r=0.67) and sign agreement (77%) indicate that the measured temporal shear is not an artifact of a specific algorithm. ICCF is prioritized for the primary analysis due to its robust handling of seasonal gaps (see §4.1.1).

4.4.2 Gap Handling

Seasonal observing gaps (typically 3-6 months) can introduce spurious correlations if not handled correctly. The pipeline implements NaN-aware Gaussian smoothing and forbids interpolation across gaps exceeding 30 days, ensuring that no artificial structure is fabricated in the gaps.

4.5 Detection Summary and Future Prospects

The current analysis provides constraints on temporal shear signals with complex sign structure:

• DESJ0408: Γ = +32 ± 15 (A-B) and Γ = +35 ± 16 (B-D) days/decade
• RXJ1131: Statistically significant at 4–5σ but with mixed signs across pairs (see §4.2 for interpretation)
• J1004+4112: Γ = +4.9 ± 1.7 (A-B) days/decade (2.9σ) — cluster lens, smaller than TEP scaling predicts

The consistent positive signs for DESJ0408 are TEP-consistent (images in deeper potential show slower time flow). RXJ1131 shows mixed signs that may reflect different image configurations. The LSST era, with its deep multi-band monitoring of thousands of lensed quasars, will provide the statistical power needed to characterize these patterns.

4.5.1 TEP Predictions for Future Surveys

If TEP is real, the geometric scaling predicts that high-redshift systems should show larger temporal shear. Current data cannot test this prediction due to insufficient precision, but future surveys targeting z_source > 3 systems with θ_E > 1.5″ would provide the most sensitive probe.

4.6 Systematic Considerations

Having established the presence of marginal signals, it is necessary to rigorously test whether they could arise from known astrophysical or instrumental systematics. The null result obtained in controls demonstrates that the pipeline does not generate false positives, but specific physical confounds like microlensing require detailed exclusion.

4.6.1 The Geometric Fingerprint: Excluding Local Systematics

A robust correlation is observed between the magnitude of the shear (|Γ|) and the source redshift ($r = 0.504$, $p = 0.014$). This supports a cosmological origin for the signal.

4.6.2 Quantitative Microlensing Exclusion

Beyond the redshift argument, the microlensing parameters required to mimic the observed signal can be quantified:

Standard microlensing does not naturally reproduce the observed temporal shear magnitude. The injection simulations show that typical microlensing trends produce |Γ| < 5 days/decade, approximately 10× weaker than the DESJ0408 signal. Reproducing the signal solely via microlensing would require stellar populations or velocities in strong tension with observational constraints.

4.6.3 Chromaticity Falsifier

The second important discriminator is multi-band chromaticity. TEP predicts achromatic temporal shear (Γ_blue ≈ Γ_red), whereas stellar microlensing generically induces chromatic distortions due to wavelength-dependent source sizes (Kochanek 2004).

An explicit multi-band test was performed using Q2237+0305 (Goicoechea et al. 2020; VizieR J/A+A/637/A89), which provides light curves in g, r, V, and I. For the A–B pair, the inferred chromatic differences are consistent with zero within the current uncertainties: ΔΓ(g−r) = 57 ± 169 days/decade and ΔΓ(V−I) = −2.6 ± 81 days/decade. A second multi-band check was performed using HE0435-1223 V and R light curves (Sorgenfrei et al. 2025; VizieR J/A+A/703/A250), finding ΔΓ(V−R) = 5.5 ± 59.0 days/decade for the A–B pair. A third multi-band check is possible in HE1104-1805 using BVRIJ light curves (Blackburne et al. 2015; VizieR J/ApJ/798/95). Using R as a reference band for the A–B pair, the results show ΔΓ(B−R) = 2.3 ± 161.0 days/decade, ΔΓ(I−R) = 143.9 ± 107.4 days/decade, and ΔΓ(J−R) = 202.1 ± 256.3 days/decade. A fourth independent dataset is available for Q2237+0305 from Vakulik et al. (2004; VizieR J/A+A/420/447), providing VRI photometry. This system has near-zero physical delays (~hours), so it serves as a null control for chromaticity: ΔΓ(V−R) = 0.0 ± 16.3 days/decade and ΔΓ(I−R) = 0.0 ± 17.4 days/decade. These constitute direct achromaticity checks in four systems, but they do not yet replace the critical test on the primary high-significance detections (DESJ0408 and PG1115), for which only R-band light curves are publicly available in the COSMOGRAIL releases (VizieR J/A+A/609/A71 and J/A+A/616/A183; Courbin et al. 2018; Bonvin et al. 2018). The public TDCOSMO data repository was also inspected (TDCOSMO/TD_data_public), which contains notebooks and mock imaging artifacts, and found no band-tagged monitoring light-curve time series for DESJ0408, PG1115, or J1206. Future multi-band monitoring of DESJ0408 and PG1115 remains a high priority for a robust chromaticity test.

4.7 Internal Consistency: Residual Cross-Correlation

If TEP is real, residuals from the Γ fit should show coherent structure across image pairs within the same system (shared path through the gravitational potential). If the signal were noise, residuals should be independent.

Cross-correlations of fit residuals between all pair combinations within each system are computed:

The significant positive mean correlation (p = 0.02) supports the interpretation that residuals share coherent structure, consistent with a real physical effect rather than independent noise.

4.8 High-Redshift Predictions

TEP predicts that |Γ| should scale with the geometric factor (1+z_S)/(1+z_L). Extrapolating from current detections yields the following predictions:

The current highest-redshift signal (DESJ0408, z_S = 2.375) shows Γ ≈ +36 days/decade (2.0σ), consistent with the linear extrapolation. RXJ1131 (z_S = 0.658) shows smaller but consistently positive values (Γ ~ +2 days/decade), suggesting the scaling includes both geometric and system-specific factors (lens mass, Einstein radius, image geometry).

5. Discussion

5.1 The Ladder of Evidence

The Temporal Equivalence Principle has been tested using two time-domain astrophysical probes that directly measure the rate of proper time accumulation. The results form a coherent "Ladder of Evidence" for potential-dependent modifications to gravitational time flow.

The identifiability of the pulsar signal is established not just by the detection of a residual, but by the quantitative exclusion of Newtonian systematics via the "Systematics Exclusion Matrix" (Section 3.6.3). Specifically, the observation of suppressed density scaling (slope 0.35) and the binary inversion (-0.31 dex) directly contradicts the predictions of standard mass segregation (slope > 0.8, positive binary residual).

Critically, both anomalies point to the same enhancement factor. Lensing requires α ~ 10⁶–10⁷. The pulsar offset, if interpreted as modified time dilation, implies a similar magnitude. This unification addresses Occam's Razor objections. While standard physics requires distinct mechanisms for the pulsar and lensing anomalies (dynamical complexity vs statistical fluctuation), TEP provides a unified parameterization ($\rho_c$) that accommodates both.

5.2 Cross-Scale Consistency with ρ_c

The universal critical density ρ_c ≈ 20 g/cm³, independently calibrated from terrestrial clock correlations, defines the screening threshold across all scales. Since ρ_c far exceeds astrophysical densities, essentially all extended gravitational systems are in the unscreened regime:

The key test is not whether ρ_c predicts specific length scales, but whether the saturation behavior is observed: in unscreened systems, TEP effects should not scale indefinitely with density. The pulsar channel confirms this with a 4$\sigma$ rejection of $\rho^2$ scaling.

5.3 Suppressed Density Scaling

The suppressed density scaling result (Section 3.4–3.5) provides strong evidence against standard dynamical contamination. The observed slope (0.29 ± 0.11) is significantly flatter than the Newtonian expectation (0.82 fiducial; ~0.72 with exact structures + segregation)—a >4$\sigma$ rejection. The signal saturates rather than scaling with density, consistent with a screening threshold at $\rho_c$.

5.4 Connection to Other TEP Evidence

The ~10⁶–10⁷ enhancement factor is consistent with previous TEP findings:

The consistency across Earth-based (GNSS, SLR), galactic (pulsars, UCD), and cosmological scales is notable. This is not what one would expect from systematic errors, which should vary with scale and methodology.

5.5 Synthesis: The Hierarchy of Evidence

A "ladder of evidence" is constructed prioritizing results that are robust to systematics (Null Controls) and spatially resolved. Fossil probes (bottom rungs) are included to demonstrate their relative insensitivity compared to rate observables.

5.6 Validation against Local Systematics

A critical test of the "Temporal Shear" anomaly is its correlation with geometry. If the signal were due to local instrumental effects (e.g., telescope temperature, orbital phase) or data processing artifacts (e.g., spline fitting, detrending), it should be independent of the lens system's redshift.

The geometric correlation provides strong evidence against systematic artifacts. If the temporal shear were due to telescope optics, pipeline errors, or microlensing, it should be independent of the lens system's cosmological geometry. Instead, the redshift dependence predicted by TEP is observed: higher z_source → larger |Γ|.

5.7 Systematics and Discriminants

5.7.1 Microlensing (Lensing Channel)

Microlensing by stars in the lens galaxy could produce time-variable magnification. However:

• Microlensing is chromatic; TEP is achromatic
• Microlensing does not scale with source redshift
• The observed correlation with (1+z_S)/(1+z_L) disfavors simple microlensing as the primary explanation

5.7.2 The Mass Sheet Degeneracy

The Mass Sheet Degeneracy (MSD) limits the precision of $H_0$ measurements by rescaling time delays by a constant factor $\lambda$ (Birrer et al. 2020). While MSD affects the absolute normalization of the time delay, it does not introduce scale-dependent structure ($\Gamma$). A pure MSD transformation $\Delta t' = \lambda \Delta t$ would scale the temporal shear magnitude but cannot generate a non-zero $\Gamma$ from a zero baseline, nor can it explain the observed redshift correlation.

5.7.3 Selection Effects (Pulsar Channel)

Could low-Ṗ pulsars be preferentially detected in GCs? No:

• Pulsars are detected by period, not Ṗ
• High-Ṗ pulsars are easier to time
• No mechanism for this selection has been proposed

5.7.3 Population Differences (Pulsar Channel)

Are GC MSPs intrinsically different from field MSPs? No known mechanism:

• Both form through similar recycling channels
• Both are spun up by accretion
• No theoretical basis for intrinsic difference

A matched comparison of field MSPs (Section 4.7) shows no difference between binary and isolated systems (p = 0.70), whereas cluster binary MSPs show a significant offset (p = 0.011). This strongly argues against intrinsic population differences as the cause of the cluster signal.

5.7.4 Cluster Acceleration: A Question of Magnitude

Pulsars moving through globular cluster potentials experience line-of-sight acceleration that contributes to observed Ṗ. This effect is well-established and produces both positive and negative apparent spin-down rates depending on pulsar position and velocity. The central question is whether the acceleration magnitude matches GR predictions or requires TEP enhancement.

Under GR, the time dilation correction from cluster acceleration is negligible:

where a_∥ is the line-of-sight acceleration and R is the cluster scale. Under TEP with α_eff ~ 10⁶–10⁷, the same physical acceleration produces an enhanced effect:

This is a 1–10% effect. If TEP is correct, what pulsar astronomers measure as "cluster acceleration" is already a TEP-enhanced time dilation effect. The frameworks are not alternatives; they describe the same physics at different coupling strengths.

5.7.5 Fossil Probe Limitations

Fossil observables (integrated quantities such as supernova light curve shapes and stellar ages) are expected to be insensitive to TEP due to astrophysical systematics that dominate at orders of magnitude above the predicted signal. Type Ia supernovae (Test G) show the expected anti-correlation with host velocity dispersion due to progenitor age effects—a result consistent with the rate/fossil distinction established in Section 1.2. See Appendix A for detailed analysis. The insensitivity of fossil channels reinforces that only rate observables (pulsar timing, lensing time delays) are sensitive to TEP modifications.

5.7.6 Laboratory and Solar System Constraints

Modified gravity theories with screening mechanisms are tightly constrained by laboratory atom interferometry and Lunar Laser Ranging (LLR). Atom interferometry excludes a wide range of chameleon/symmetron parameters in vacuum (Burrage et al. 2018). However, TEP posits a screening transition at $\rho_c \approx 20 \text{ g/cm}^3$. Laboratory vacuum chambers are embedded within the Earth's density field, which is well above $\rho_c$, ensuring the local environment is screened. The predicted enhancement ($\alpha_{eff} \sim 10^6$) applies only to extended systems with density below $\rho_c$ (e.g., cluster outskirts, galactic halos), consistent with the observed null results in dense Solar System regimes.

5.7.7 Consistency with Pulsar Timing Arrays

Pulsar Timing Arrays (PTAs) such as NANOGrav, EPTA, and the Fermi-LAT PTA (Xia et al. 2023) place stringent constraints on Ultralight Dark Matter (ULDM) and stochastic gravitational wave backgrounds. Since the Galaxy is an "extended configuration" with density $\rho \ll \rho_c$, field pulsars used in PTAs might be expected to exhibit strong TEP signatures. The absence of such signals is consistent with TEP for three reasons:

• Static vs. Oscillatory Nature: PTA constraints on ULDM (e.g., Xia et al. 2023) assume a scalar field oscillating at the Compton frequency ($f \approx m_\phi c^2 / h$). For $m \sim 10^{-22}$ eV, this produces a time-varying residual with period ~1 year. TEP, by contrast, posits a static or slowly-varying scalar background (soliton). A static modification to the local potential produces a constant shift in the pulsar's spin frequency ($\nu$) and spin-down rate ($\dot{\nu}$). Because PTAs fit $P$ and $\dot{P}$ individually for every pulsar, these constant offsets are absorbed into the timing model and are effectively invisible to residual analysis.
• Screened Earth Term: PTA searches for correlated signals rely on the "Earth term"—the component of the signal common to all pulsars due to the detector's (Earth's) motion or potential. However, the Solar System density ($\rho \gg \rho_c$) ensures the Earth is locally screened. Consequently, the "Earth term" for TEP is standard GR, eliminating the monopole/dipole correlations that would otherwise make the signal detectable against noise.
• Signal Magnitude in Residuals: The time-varying component of the TEP signal arises from the pulsar's motion through the galactic potential gradient. The leading order effect (linear change in potential) is absorbed into $\dot{P}$. The first non-absorbed term is the "jerk" ($\ddot{\nu}$), driven by the curvature of the galactic potential.
  
  Explicit calculation for a pulsar moving at $v \sim 220$ km/s through the Galactic potential:
  $\Delta t_{\text{TEP}} \approx \frac{1}{6} \frac{\alpha \ddot{\Phi}}{c^2} T_{\text{obs}}^3 \sim 1 \mu\text{s} \quad (\text{over 10 years})$
  This drift (~1 $\mu$s) is comparable to or smaller than the intrinsic "red noise" often observed in millisecond pulsars over decadal baselines and is far below the deterministic shifts absorbed into $\dot{P}$. Thus, TEP does not violate current PTA constraints.

5.8 Key Discriminating Tests

A useful distinction is between tests that would strongly constrain the lensing-based interpretation and tests that primarily refine particular phenomenological scalings.

5.9 Limitations and Robustness

To aid critical evaluation, the primary limitations, parameter sensitivities, and failure modes of the analysis are explicitly identified.

5.10 Future Tests

6. Conclusions

This work presents time-domain astrophysical tests of the Temporal Equivalence Principle at intermediate gravitational scales (10⁵–10¹² M_☉). Analysis of 380 millisecond pulsars (182 GC, 198 field) with measured spin-down rates provides spatially-resolved evidence for environmental anomalies in pulsar spin-down rates, validated by independent controls and consistent with the universal critical density ρ_c ≈ 20 g/cm³ calibrated from terrestrial observations.

6.1 Summary of Findings

6.2 The Primary Detection: Pulsar Timing

Analysis of 380 MSPs with measured spin-down rates (Freire GCpsr + ATNF cross-match) reveals an environmental signal in globular cluster pulsars that satisfies three independent criteria consistent with TEP:

• Spatial Resolution: The spin-down anomaly is concentrated in cluster cores (−0.33 dex for inner binaries, p = 0.054) and absent in the outskirts (−0.09 dex, p = 0.63), directly tracking gravitational potential depth.
• Environmental Isolation: The Field Binary Control supports an environmental origin rather than an intrinsic one—the binary vs isolated difference vanishes in the galactic field (p = 0.70).
• Suppressed Density Scaling: While standard dynamics predicts residuals scaling strongly with density (slope ≈ 0.82 fiducial; ~0.72 with exact structures + segregation), the observed slope is only 0.35 ± 0.09—a 4.0σ rejection. All 15 clusters with sufficient statistics show positive controlled residuals (+0.02 to +0.33 dex), consistent with a universal environmental enhancement that saturates rather than scaling with density.

6.3 The Geometric Evidence: Gravitational Lensing

Analysis of the COSMOGRAIL dataset provides marginal geometric evidence for temporal shear at halo scales:

• Marginal Signal in DESJ0408: The B-D image pair shows temporal shear of $|\Gamma| \approx 36$ days/decade (2.0σ), with geometric sign inversion consistent with saddle-point topology.
• Geometric Verification: The sign inversion between A-B and B-D pairs is explained by the saddle-point topology of Image B, ruling out simple systematics.
• Injection-recovery tests demonstrate estimator linearity (mean bias 0.3 days).
• The geometric correlation with source redshift (r = 0.504) further supports a cosmological origin.

6.4 Cross-Scale Consistency

The convergence of time-domain evidence across scales is noteworthy:

The single parameter ρ_c defines a consistent screening threshold across all scales: systems with ρ ≪ ρ_c (all astrophysical environments) show saturation behavior, while Earth (ρ ~ ρ_c) shows a transition. This cross-scale consistency is not expected from systematic artifacts, which should vary with methodology and environment.

6.5 The Critical Path: Key Tests

The TEP hypothesis can be further constrained or strengthened by specific near-term observations. The following tests are prioritized:

Timeline: Both tests are achievable within 12–24 months: (1) chromaticity test on DESJ0408 using LSST commissioning data or targeted proposals; (2) detailed N-body simulations of Terzan 5 and 47 Tuc using existing computational resources. These measurements will significantly refine the TEP interpretation and constrain the parameter space.

6.6 Final Statement

This work investigated the hypothesis that intermediate-scale anomalies reflect modified temporal structure rather than dark sector physics. The data provide a split verdict: a robust detection in the pulsar channel and a stringent constraint in the lensing channel.

The Verdict: Pulsar timing reveals a spatially-resolved signal in globular cluster cores that deviates from standard Newtonian dynamics (4.0$\sigma$ suppression of density scaling) while tracking gravitational potential depth. Simultaneously, gravitational lensing places an upper limit of |Γ| ≲ 60 days/decade on temporal shear at halo scales, ruling out extreme modifications while remaining consistent with the screened parameter space implies by the pulsar signal.

These findings do not constitute proof of TEP. They do, however, present a coherent "Ladder of Evidence" in which independent time-domain probes converge on a consistent picture. Importantly, the identifiability of the pulsar signal against "incomplete dynamical modeling" is established by specific falsification criteria: standard mass segregation predicts steeper density scaling ($\Gamma > 0.8$) and higher acceleration for binaries, while the data show suppressed scaling ($\Gamma \approx 0.35$, 4$\sigma$) and a binary inversion (-0.31 dex). This pattern specifically excludes the class of standard dynamical heating models.

The critical path forward requires pushing lensing precision to match the sensitivity already achieved in the pulsar channel, and performing full N-body verification of the suppressed density scaling result to rigorously test standard dynamical explanations.

6.7 Data and Code Availability

The complete data tables (including the full GC pulsar compilation and COSMOGRAIL light curves) and the Python analysis pipeline used to generate all figures and statistics in this work are available in the GitHub repository: https://github.com/matthewsmawfield/TEP-COS.

The repository includes a comprehensive reproduction guide (see README.md) to facilitate independent verification of the results. The analysis is fully containerized and reproducible, allowing researchers to verify the "Suppressed Density Scaling" and "Temporal Shear" results directly from the raw catalogs.

References

External References

Data Availability

Analysis code: https://github.com/matthewsmawfield/TEP-COS

MaNGA data: https://www.sdss.org/dr17/manga/

COSMOGRAIL light curves: Available from VizieR/CDS (see individual paper DOIs above)

GC Pulsar Catalog: Paulo Freire's GC Pulsar Catalog (MPIfR)

ATNF Pulsar Catalogue: https://www.atnf.csiro.au/research/pulsar/psrcat/

Acknowledgments

This work uses data from the Sloan Digital Sky Survey IV (SDSS-IV). Funding for SDSS-IV has been provided by the Alfred P. Sloan Foundation, the U.S. Department of Energy Office of Science, and the Participating Institutions.

The COSMOGRAIL light curves were obtained through the dedicated monitoring programs of the COSMOGRAIL collaboration. The author acknowledges the use of data from the Swiss 1.2m Euler telescope at La Silla Observatory, the SMARTS 1.3m telescope at Cerro Tololo, and the Mercator telescope at La Palma.

Pulsar timing data are from the ATNF Pulsar Catalogue (Manchester et al. 2005) and the comprehensive globular cluster pulsar catalog maintained by P. Freire (MPIfR). The author thanks the pulsar timing community for making these data publicly available.

Appendix A: Astrophysical Systematics in Fossil Probes

This appendix documents tests of fossil observables (integrated quantities rather than rates) that are expected to be insensitive to TEP due to astrophysical systematics. These results are included for completeness and to establish the distinction between rate-sensitive and fossil channels.

A.1 The Supernova Stretch Tension (Test G)

The most significant challenge to naive TEP predictions emerges from the Type Ia Supernova channel. TEP predicts that if time flows slower in deep potentials, standard clocks should appear to tick slower. Therefore, SNe in high-σ hosts should exhibit broader light curves (higher stretch $x_1$).

A.2 Why Fossil Probes Are Insensitive

The distinction between rate and fossil observables is fundamental to TEP phenomenology:

TEP modifications at the ~10⁻⁵ level are swamped by astrophysical scatter at the ~10⁻¹ level in fossil observables. The Supernova result is therefore expected under TEP—it is not a contradiction but a validation of the framework's predictions about observable sensitivity.

Appendix B: Data and Code Availability

To facilitate reproduction and independent verification of these results, the exact data snapshots, selection queries, and analysis scripts used in this work are provided via the public repository: https://github.com/matthewsmawfield/TEP-COS.

B.1 Catalog Snapshots

B.2 Selection Queries

# Standard Millisecond Pulsar (MSP) Definition
P_spin < 30 ms
P_dot_intrinsic > 0 (where available)
Not in binary with massive companion (> 10 M_sun)

# Cluster Association
Use Freire catalog "Cluster" field.
Filter out foreground contaminants identified in literature.
    

B.3 Analysis Code

All analysis steps are encapsulated in Python scripts available in the scripts/ directory. Key reproduction scripts include:

• scripts/steps/step_5_10_pulsar_population_controls.py: Implements the exact matching procedure for pulsar controls.
• scripts/steps/step_5_33_hierarchical_density_scaling.py: Runs the hierarchical mixed-effects model for density scaling.
• scripts/steps/step_3_0_cosmograil_temporal_shear.py: Performs the lensing time delay and temporal shear analysis.