There was an interesting paper by Seikel & Clarkson on astro ph recently (link) discussing the choice of covariance function and hyper parameter optimization methodology for Gaussian Process (GP) regression of w(z) and its derivatives (cosmology). The background here is that although we have a restricted class of “popular” cosmological models that are widely believed and which do a pretty good job at explaining our current observational datasets (WMAP+SDSS+2dFGRS+Union2+etc) none are particularly well-motivated with regard to their formulation of the mysterious dark energy component. Thus there is some interest in “non-parametric” reconstructions of the key cosmological distance/expansion functions, such as the equation of state (a univariate function), w(z), given only the available empirical data. The typical workhorse for this problem in modern statistics is, of course, the GP used as a regression tool. In this study Seikel & Clarkson have chosen to examine the methodological details of this procedure by way of numerical simulations with a particular focus on reconstructing mock forms of w(z) from mock datasets of 546 SNe.

The two key issues explored by the authors — the sensitivity of the posterior regression curve on the choice of covariance function (e.g. Matern, exponential, etc), and the effect on the same of deciding whether or not to marginalise over hyper parameters (as opposed to just optimizing *a la* empirical Bayes) — are too often ignored in astronomical studies, so it is good to see these highlighted here. In brief, both decisions can have a significant influence on the type of fit recovered: most notably, of course, in regions of sparse data (i.e. where we are harnessing the power of the GP for interpolation). However, although the authors do a nice job at exploring these issues and conduct detailed simulations to assess coverage under various scenarios I was left feeling that something was missing in terms of motivation and/or purpose; more precisely, I was left without a sense of what we need from these GP reconstructions, which in turn I might have guided the assessment stage more informatively. Do we need these non-parametric reconstructions of w(z) for planning future surveys robust against the assumption of a particular parametric model for w(z)? If so, do we need them to interpolate over the redshift desert, or to extrapolate beyond the range of our current datasets? Likewise, I would have also appreciated some design consideration: how does the expected distribution of observed SNe relate to our ability to constrain the hyper parameters of these covariance functions?

(And, since this is my pet interest, how would we fare if we went a step further and tried a model averaging of the posterior w(z)’s [via marginal likelihood] over all these covariance functions at once?)