I noticed the recent paper by Vallisneri & Haasteren on the arXiv entitled “Taming outliers in pulsar-timing datasets with hierarchical likelihoods and Hamiltonian sampling” in which the authors point out that Bayesian hierarchical models can, and perhaps should, be used to account for outliers when modelling pulsar-timing datasets with Gaussian processes. I wasn’t convinced that this point justifies an entire paper, since the idea of using such a model with a mixture of likelihoods for good and bad data points is well known already in astronomy as per Hogg, Bovy & Lang (as pointed out by V&H) and others (e.g. I would point to Taylor et al.); the only novelty here being that it’s described and applied to a single mock pulsar dataset. I was also not convinced by the authors’ arguments to show that the principled approach would surely beat the ad hoc alternative of 3-clipping, which was argued on the grounds that since the pulsar model includes estimation of the error term there is no single 3 threshold to work with. One would think that an simple marginalised estimator of model ‘discrepancy’ (marginalising over the unknown and underlying GP) for each data point would provide a tractable equivalent to threshold possible outliers. It would then be meaningful to compare the performance of this approach against the hierarchical model in terms of bias and runtime under a variety of bad data scenarios.
Having said that I was very interested to see that random Fourier features are being embraced by the pulsar timing community, as described in greater detail in an earlier paper by the same authors and another by Lentati et al. I was recently talking to Seth Flaxman about how odd it is that random Fourier bases aren’t more popular, with one possibility being that people find they need a huge number of components to get a close approximation to their targeted GP. But of course this is going to depend on the application, or whether you even care: like with the INLA/SPDE approximation we use a lot at MAP, in practical applications it may be sufficient to say ‘we adopt the approximate model as our model’ and then just make sure you do a damn good job testing its performance on mock data and out-of-sample cross-validation (etc.).