Make whitened residuals available

[abhisrkckl]

This partially overlaps with the stale PR #1000.

whitened_resids = (resids - noise_resids) / scaled_toa_error

If the fitting has been done properly, the whitened residuals should be normally distributed with unit variance.

[codecov]

Codecov Report

All modified and coverable lines are covered by tests ✅

Comparison is base (f2da697) 68.47% compared to head (cbb1924) 68.47%.

❗ Current head cbb1924 differs from pull request most recent head 7c869fb. Consider uploading reports for the commit 7c869fb to get more accurate results

Additional details and impacted files

@@           Coverage Diff           @@
##           master    #1671   +/-   ##
=======================================
  Coverage   68.47%   68.47%           
=======================================
  Files         104      104           
  Lines       24289    24294    +5     
  Branches     4336     4336           
=======================================
+ Hits        16631    16636    +5     
  Misses       6572     6572           
  Partials     1086     1086           

Files	Coverage Δ
src/pint/residuals.py	76.48% <100.00%> (+0.28%)	⬆️

☔ View full report in Codecov by Sentry.
📢 Have feedback on the report? Share it here.

[abhisrkckl]

Code:

from pint.models import get_model
from pint.simulation import make_fake_toas_uniform
from pint.fitter import Fitter
from pint.logging import setup as setup_log

from matplotlib import pyplot as plt
from astropy import units as u
import numpy as np

setup_log(level="WARNING")

m = get_model(
   StringIO(
        """
        RAJ    		05:00:00
        DECJ   		20:00:00
        F0     		100     1 
        F1     		-1e-14  1
        PEPOCH 		58000
        DM     		15
        EFAC tel gbt 	1.3
        TNRedAmp 	-14
        TNRedGam 	4.5
        TNRedC 		30
        """	
    )
)
t = make_fake_toas_uniform(
    57000, 
    60000, 
    1000, 
    model=m, 
    error=1 * u.us, 
    freq=u.Quantity([1400*u.MHz, 2500*u.MHz]), 
    add_noise=True, 
    add_correlated_noise=m.has_correlated_errors
)


ftr = Fitter.auto(t, m)
ftr.fit_toas()

# ===========

plt.subplot(211)
plt.errorbar(t.get_mjds(), ftr.resids.time_resids, ftr.resids.get_data_error(), ls="", marker="+")
plt.ylabel("resids (s)")
plt.xlabel("mjd")

plt.subplot(223)
plt.scatter(t.get_mjds(), ftr.resids.calc_whitened_resids(), marker="+")
plt.ylabel("wh. resids (s)")
plt.xlabel("mjd")

plt.subplot(224)
plt.hist(ftr.resids.calc_whitened_resids(), density=True, orientation="horizontal")

plt.show()

[dlakaplan]

Does this work for WB data?

[dlakaplan]

Is it better to have a separate method, or calc_time_resids() with a flag?

[abhisrkckl]

I think it should be a separate method since whitened residuals are dimensionless, whereas time resids have the dimentions of time.

[dlakaplan]

Then probably add a "See Also" to the other residual methods

[abhisrkckl]

OK.

[abhisrkckl]

Updated documentation.

[abhisrkckl]

Does this work for WB data?

It should. But I haven't tried it yet.

[dlakaplan]

Looks good to me. ready to go?

[abhisrkckl]

Looks good to me. ready to go?

Yes.