Apparent cosmic acceleration from Type Ia supernovae

Abstract

Parameters that quantify the acceleration of cosmic expansion are conventionally determined within the standard Friedmann–Lemaître–Robertson–Walker (FLRW) model, which fixes spatial curvature to be homogeneous. Generic averages of Einstein's equations in inhomogeneous cosmology lead to models with non-rigidly evolving average spatial curvature, and different parametrizations of apparent cosmic acceleration. The timescape cosmology is a viable example of such a model without dark energy. Using the largest available supernova data set, the JLA catalogue, we find that the timescape model fits the luminosity distance-redshift data with a likelihood that is statistically indistinguishable from the standard spatially flat Λ cold dark matter cosmology by Bayesian comparison. In the timescape case cosmic acceleration is non-zero but has a marginal amplitude, with best-fitting apparent deceleration parameter, |$q_{\scriptstyle 0}=-0.043^{+0.004}_{-0.000}$|⁠. Systematic issues regarding standardization of supernova light curves are analysed. Cuts of data at the statistical homogeneity scale affect light-curve parameter fits independent of cosmology. A cosmological model dependence of empirical changes to the mean colour parameter is also found. Irrespective of which model ultimately fits better, we argue that as a competitive model with a non-FLRW expansion history, the timescape model may prove a useful diagnostic tool for disentangling selection effects and astrophysical systematics from the underlying expansion history.

cosmological parameters, dark energy, cosmology: observations – cosmology: theory

1 INHOMOGENEOUS ALTERNATIVES TO DARK ENERGY

One of the most important discoveries in cosmology (Riess et al. 1998; Perlmutter et al. 1999) is the observation that the luminosity distances and redshifts of Type Ia supernovae (SneIa) are well matched to the expansion history of a spatial homogeneous and isotropic Friedmann–Lemaître–Robertson–Walker (FLRW) model only if the Universe began an epoch of accelerated expansion late in cosmic history. Since gravity with matter obeying the Strong Energy Condition is universally attractive, this demands a cosmological constant, Λ, or some other unknown source of spatially homogeneous dark energy with an equation of state, P = wρc², that subsequent tests find to be close to the cosmological constant case, w = −1. Independently of SneIa, since the late 1990s other data sets including the cosmic microwave background (CMB) and galaxy clustering statistics, have been found to independently require late epoch cosmic acceleration in the FLRW model.

Despite this success, the nature of dark energy remains a mystery for fundamental physics. Furthermore, a number of puzzles remain (Buchert et al. 2016; Bull et al. 2016) which range from significant anomalies (Cyburt, Fields & Olive 2008) to lesser tensions (Delubac et al. 2015). While some puzzles may inevitably boil down to statistical sampling biases, systematic problems in data reduction or potentially unknown astrophysics, they may also in part be due to incorrect cosmological model assumptions. It is therefore important to consider alternatives to the standard cosmology, particularly any alternatives which do not change known gravitational physics but which reconsider assumptions that are not demanded theoretically in light of observations. The assumption of spatial homogeneity and isotropy falls in this category. While the remarkable isotropy of the CMB points to an initial state with a very high degree of smoothness, the late epoch Universe encompasses a complex cosmic web of structures (van de Weygaert et al. 2016). It is dominated in volume by voids that are threaded and surrounded by clusters of galaxies distributed in sheets, knots and filaments.

Spatial homogeneity is first encountered in some statistical sense in averages on scales ≳70–120 h⁻¹ Mpc based on the 2-point galaxy correlation function (Hogg et al. 2005; Scrimgeour et al. 2012). This indicates that we can expect cosmology to be modelled by some averaged description on larger scales, but does not guarantee that such a description is given at all times by the FLRW geometry, however. Demanding FLRW geometry on arbitrarily large scales in space and time presupposes one particular answer to the fitting problem (Ellis 1984; Ellis & Stoeger 1987; Wiltshire 2011). The fitting problem is an unanswered fundamental question in general relativity: on what scale(s) are matter and geometry dynamically coupled via Einstein's equations? Since general relativity is a causal theory, the speed of light yields a natural upper limit. However, most of the energy density of the present epoch Universe in the energy–momentum tensor on the right hand side of Einstein's equations is in non-relativistic matter subject to speeds much smaller than that of light. Consequently, the effective causal scale for density fluctuations is the matter horizon (Ellis & Stoeger 2009), the comoving region which has significantly contributed matter to the local physical environment of a given observer.

The matter horizon gives a natural effective causal upper limit for direct applicability of the Einstein equations. Regardless of what the scale is, if Einstein's equations apply only on small scales then generic large-scale averages do not lead to the Einstein equations for a prescribed averaged matter source. Instead one obtains averaged cosmological models with backreaction (Buchert 2000, 2001; Buchert & Räsänen 2012), whose evolution generically can differ significantly from that of the FLRW geometry. Since the cosmic expansion history is phenomenologically very well described by a FLRW model with a small number of parameters, then any physically relevant model of backreaction must be subject to some simplifying principles that restrict many possible mathematical choices for defining the average expansion history and its relation to local observables (Buchert & Carfora 2002, 2003).

The timescape scenario (Wiltshire 2007a,b, 2009; Duley, Nazer & Wiltshire 2013) provides such principles. Within Buchert's scalar averaging formalism (Buchert 2000, 2001), the relationship between statistical volume averages and local observables is restricted by applying the Cosmological Equivalence Principle (Wiltshire 2008), a generalization of the Strong Equivalence Principle. The clocks and rulers which best describe average cosmic evolution (the ‘bare’ or ‘volume–average’ quantities) in Buchert's statistical formalism are reinterpreted, and differ systematically from local clocks and rulers of observers in environments where the matter density differs significantly from the volume average. Given the smooth initial state of the Universe, this only becomes significant late in cosmic history, when spatial curvature gradients grow between bound structures (which by definition are greater than critical density) and the voids which dominate by volume.

The timescape model has passed a number of independent observational tests (Leith, Ng & Wiltshire 2008; Smale 2011; Smale & Wiltshire 2011; Duley et al. 2013; Nazer & Wiltshire 2015). Its distance–redshift relation is very close to that of particular Λ cold dark matter (ΛCDM) models over small redshift ranges, but effectively interpolates (Wiltshire 2009, 2014) between spatially flat ΛCDM cosmologies with different values of |$\Omega _{\scriptstyle M0}$| and |$\Omega _{\scriptstyle \Lambda 0}$| over larger redshift ranges. In particular, when the timescape model is fit to the angular diameter distance of the sound horizon in the CMB then the spatially flat ΛCDM model with the closest comoving distance at redshifts z ≳ 100 has a value of |$\Omega _{\scriptstyle M0}$| 15–27 per cent lower than that of the spatially flat ΛCDM model with the closest comoving distance at redshifts z < 1.2 (Wiltshire 2009, 2014).

Geometric tests of the timescape expansion history are most developed (Wiltshire 2009), and give rise to measures (Clarkson, Bassett & Lu 2008; Räsänen, Bolejko & Finoguenov 2015) which will definitely distinguish both the timescape model and other alternatives (Larena et al. 2009; Lavinto, Räsänen & Szybka 2013) from the standard FLRW model using Euclid satellite data (Sapone, Majerotto & Nesseris 2014). On the other hand, tests of the CMB anisotropy spectrum in the timescape model are at present limited by systematic uncertainties of 8–13 per cent in parameters which relate to the matter content (Nazer & Wiltshire 2015). This is a consequence of backreaction schemes having not yet been applied to the primordial plasma.

2 SUPERNOVA REDSHIFT-DISTANCE ANALYSIS

In the case of the redshift range probed by SneIa the difference between the timescape and ΛCDM cosmologies is comparable to the systematic uncertainties that arise between different methods for fitting the light curves of SneIa to obtain ‘standard candles’. In particular, in the last full analysis of the timescape model (using data available in 2010), Smale & Wiltshire (2011) found significant differences between data reduced by the MLCS2k2 (Multi Light-Curve Shape) fitter (Jha, Riess & Kirshner 2007) and the SALT/SALT2 (Spectral Adaptive Light-Curve Template) fitters (Guy et al. 2005, 2007). While the relative Bayesian evidence was sometimes ‘positive’ (but not very strong), the conclusion as to which cosmological model fitted better depended on the light-curve fitting method. Consequently, the empirical nature of light-curve fitting may mask effects due to the underlying expansion history if this deviates from the FLRW geometry.

The significantly larger Joint Light-curve Analysis (JLA) SneIa catalogue (Betoule et al. 2014) now makes possible a renewed comparison of the timescape and ΛCDM models, as well as further investigation of the systematics of light-curve fitting. Recently, Nielsen, Guffanti & Sarkar (2016, NGS16) have used the JLA catalogue to reinvestigate systematic issues associated with SneIa light-curve fitting within FLRW cosmologies using the SALT2 method. They adopted maximum likelihood estimators (MLE) that take into account the underlying Gaussianity of particular light-curve parameters (March et al. 2011). NGS16 concluded that the significance for cosmic acceleration, as compared to an empty Milne model (or any cosmology with constant expansion), is ‘marginal’ (at ≲ 3σ significance). This conclusion was challenged by Rubin & Hayden (2016, RH16), who introduced 12 additional light-curve parameters to allow for possible unaccounted systematics, concluding that the 2.8σ significance found by NGS16 increased to 3.7σ for a general FLRW model, or to 4.2σ for the spatially flat case. However, RH16 did not consider whether the increased model complexity was justified from a Bayesian standpoint.

In the SALT2 method, each observed supernova redshift is used to determine a theoretical distance modulus,

\begin{equation} \mu \equiv 25 + 5 \log _{10}\left( d_{\scriptstyle \rm L}\over \hbox{Mpc}\right), \end{equation}

(1)

where |$d_{\scriptstyle \rm L}$| is the luminosity distance for each cosmological model. This is then compared to the observed distance modulus, which is related to the supernova light curve by

\begin{equation} \mu _{\scriptstyle \rm SN}=m_{\scriptstyle \rm B}^*-M_{\scriptstyle \rm B}+\alpha x_{\scriptstyle 1}-\beta c, \end{equation}

(2)

where |$m_{\scriptstyle \rm B}^*$| is the apparent magnitude at maximum in the rest-frame B band, |$M_{\scriptstyle \rm B}$| is the corresponding absolute magnitude of the source, x₁ and c are empirical parameters that describe the light-curve stretch and colour corrections for each supernova, while α and β are parameters that are assumed to be constant for all SneIa.

The theoretical distance modulus (1) is determined for a bolometric flux, which is not directly measured. The SALT2 (Guy et al. 2005, 2007) relation (2) can thus be viewed as a model for a band correction, |$\Delta \mu _{\scriptstyle \rm B}$|⁠, that is linear in the variables |$x_{\scriptstyle 1}$| and c,

\begin{equation} \Delta \mu _{\scriptstyle \rm B}\equiv (m-M)-(m_{\scriptstyle \rm B}^*-M_{\scriptstyle \rm B})=\alpha x_{\scriptstyle 1}-\beta c, \end{equation}

(3)

where m and M are the bolometric apparent and absolute magnitudes in the observer and emitter rest frames, respectively.

In the SALT2 method, the light-curve parameters are simultaneously fit together with the free cosmological parameters on the entire data set.

NGS16 assumed that all SneIa in the JLA catalogue (Betoule et al. 2014) are characterized by parameters, |$M_{\scriptstyle \rm B}$|⁠, |$x_{\scriptstyle 1}$| and c, drawn from the same independent global Gaussian distributions, with means |$M_{\scriptstyle \rm B,0}$|⁠, |$x_{\scriptstyle 1,0}$| and |$c_{\scriptstyle 0}$|⁠, and standard deviations |$\sigma _{\scriptstyle M_{\rm B,0}}$|⁠, |$\sigma _{\scriptstyle x_{\rm 1,0}}$| and |$\sigma _{\scriptstyle c_{\rm 0}}$|⁠, respectively. These six free parameters were then fitted along with the light-curve parameters α, β and the cosmological parameters.

RH16 claimed that the mean light-curve stretch and colour parameters, |$x_{\scriptstyle 1,0}$| and |$c_{\scriptstyle 0}$|⁠, of the Gaussian distributions analysed by NGS16 show some redshift dependence. This may be partially due to astrophysical effects in the host population, or – particularly for the colour parameter – may arise from the colour–luminosity relation combined with redshift-dependent detection limits. In other words, Malmquist type biases may not be completely corrected for in the JLA catalogue (Betoule et al. 2014). In the absence of a known astrophysical model for such corrections, RH16 introduced 12 additional empirical parameters by replacing the global Gaussian means according to

\begin{equation} x_{\scriptstyle 1,0}\rightarrow x_{\scriptstyle 1,0,J} + x_{\scriptstyle z,J} z,\qquad \hbox{and}\qquad c_{\scriptstyle 0}\rightarrow c_{\scriptstyle 0,J} + c_{\scriptstyle z,J} z, \end{equation}

(4)

where the index J runs over the four independent subsamples in the JLA catalogue: (1) SNLS (Supernova Legacy Survey); (2) SDSS (Sloan Digital Sky Survey); (3) nearby supernovae; (4) HST (Hubble Space Telescope), with |$x_{\scriptstyle z,4}=0$|⁠, |$c_{\scriptstyle z,4}=0$| on account of limited HST data. The widths |$\sigma _{\scriptstyle x_{\rm 1,0}}$|⁠, |$\sigma _{\scriptstyle c_{\rm 0}}$| were still treated as global parameters.

We perform a Bayesian comparison of fits of the JLA catalogue (Betoule et al. 2014) to the luminosity distance–redshift relation for the spatially flat ΛCDM model, and for the timescape model (Wiltshire 2007a,b, 2009; Duley et al. 2013). We first use the MLE methodology of NGS16 directly, and then investigate the effect of changes to light-curve fitting suggested by RH16. An empty universe with constant expansion rate is also analysed, as a convenient demarcation of accelerating from non-accelerating expansion in the FLRW case.

Details of the theoretical luminosity distances used in (1) are given in Appendix A. The model differences that we are testing are best appreciated by comparing the distance moduli of the timescape and spatially flat ΛCDM models relative to an empty universe, as shown in Fig. 1. The timescape distance modulus, |$\mu _{\scriptstyle \rm TS}(z)$|⁠, is closer to ΛCDM than the empty case. Nonetheless, |$\mu _{\scriptstyle \rm TS}(z)$| is always closer to |$\mu _{\scriptstyle \rm empty}(z)$| than |$\mu _{\scriptstyle \Lambda {\rm CDM}}(z)$| is, a consequence of cosmic acceleration being an apparent effect in the timescape model.

Figure 1.

The residual distance moduli |$\mu _{\scriptstyle \Lambda {\rm CDM}}(z)-\mu _{\scriptstyle \rm empty}(z)$| and |$\mu _{\scriptstyle \rm TS}(z)-\mu _{\scriptstyle \rm empty}(z)$| with the same |$H_{\scriptstyle 0}$|⁠. The best-fitting parameters of Table 2 are assumed: |$\Omega _{\scriptstyle M0}=0.365$| for spatially flat ΛCDM and f_v0 = 0.778 for timescape. For redshifts z < 0.927 over which |$\mu _{\scriptstyle \rm TS}(z)<\mu _{\scriptstyle \Lambda {\rm CDM}}(z)$|⁠, the maximum difference between the ΛCDM and timescape distance moduli is |$\mu _{\scriptstyle \Lambda {\rm CDM}}(z)-\mu _{\scriptstyle \rm TS}(z)=0.031$| at z = 0.372.

Further technicalities about systematic issues in implementing the SALT2 method are discussed in Appendix B.

3 STATISTICAL METHODS

3.1 The likelihood construction

We adopt the likelihood construction (March et al. 2011) used by NGS16. The likelihood, |$\mathcal {L}$|⁠, is the probability density of the observed data – here |$(\hat{z},{\hat{m}}_{\scriptstyle \rm B}^*,\hat{x}_{\scriptstyle 1},\hat{c})_i$|⁠, i = 1, 2, …, N on N supernovae – given a model, M. The likelihood can be written as (March et al. 2011)

\begin{eqnarray} \mathcal {L}&\equiv& \mathcal {P}\left[(\hat{z}, {\hat{m}}_{\scriptstyle \rm B}^*, \hat{x}_{\scriptstyle 1}, \hat{c})_i \,|\, \mathrm{M}\right] \nonumber \\ &=& \int {\rm d}M_{\scriptstyle \rm B}^{N}\,{\rm d}x_{\scriptstyle 1} ^{N}\,{\rm d}c ^{N}\,\mathcal {P} \left[(\hat{z}, {\hat{m}}_{\scriptstyle \rm B}^*, \hat{x}_{\scriptstyle 1}, \hat{c})_i \,|\, (M_{\scriptstyle \rm B}, x_{\scriptstyle 1}, c)_i,\mathrm{M}\right] \nonumber \\ &&\times\, \mathcal {P}\left[ (M_{\scriptstyle \rm B}, x_{\scriptstyle 1}, c)_i \,|\, \mathrm{M}\right], \end{eqnarray}

(5)

where hatted quantities denote measured data values including all experimental noise, and unhatted quantities are intrinsic parameters that characterize the statistical distributions from which the supernovae are drawn. Only the intrinsic parameters satisfy the SALT2 relation (2). The empirical light-curve model (2) and the theoretical distance modulus (1) together constitute the model, M.

The expansion in (5) allows one to explicitly model the intrinsic scatter of the supernovae. For the NGS16 model (I), we assume that the intrinsic parameters of each supernova are drawn from identical independent Gaussian distributions

\begin{eqnarray} \mathcal {P}[ (M_{\scriptstyle \rm B}, x_{\scriptstyle 1}, c)_i \,|\, \mathrm{M}] &=&\prod _i^N \mathcal {P}[ (M_{\scriptstyle {\rm B},i}, x_{\scriptstyle 1,i}, c_{\scriptstyle i}) \,|\,\mathrm{M}] \nonumber \\ &=& \prod _i^N \mathcal {N}(M_{\scriptstyle {\rm B},i}\,,\, M_{\scriptstyle \rm B,0},\sigma _{\scriptstyle M_{{\rm B},0}}) \nonumber \\ &&\times \,\mathcal {N}(x_{\scriptstyle 1,i}\,,\, x_{\scriptstyle 1,0}, \sigma _{x_{\scriptscriptstyle 1,0}})\, \mathcal {N}(c_{\scriptstyle i}\,,\, c_{\scriptstyle 0}, \sigma _{c_{\scriptscriptstyle 0}}), \end{eqnarray}

(6)

where |$\mathcal {N}(y\,;\,y_{\scriptstyle 0},\sigma )\equiv (2\pi \sigma ^2)^{-1/2}\exp \left[-(y-y_{\scriptstyle 0})^2 /(2\sigma ^2)\right]$| for each triple |$\lbrace y,y_{\scriptstyle 0},\sigma \rbrace$|⁠, with 3N × 3N diagonal covariance matrix |$\Sigma _l= \mathrm{diag}\left(\sigma ^2_{\scriptstyle M_{\scriptscriptstyle {\rm B},0}},\sigma ^2_{x_{1,0}}, \sigma ^2_{\scriptstyle c_{\scriptscriptstyle 0}},\sigma ^2_{\scriptstyle M_{\scriptscriptstyle {\rm B},0}}, \dots \right)$|⁠. The experimental part of the likelihood (5), |$\mathcal {P}\left[(\hat{z},{\hat{m}}_{\scriptstyle \rm B}^*,\hat{x}_{\scriptstyle 1},\hat{c})_i \,|\, (M_{\scriptstyle \rm B},x_{\scriptstyle 1},c)_i, \mathrm{M}\right]$| is assumed to be a Gaussian in the intrinsic supernova parameters, and the covariance matrix of experimental statistical and systematic uncertainties is denoted Σ_d. Performing the integral in (5) one obtains the final expression of the likelihood (Nielsen et al. 2016)

\begin{eqnarray} \mathcal {L}&=&\mathcal {P}\left[(\hat{z}, {\hat{m}}_{\scriptstyle \rm B}^*, \hat{x}_{\scriptstyle 1}, \hat{c})_i \,|\, \mathrm{M}\right] = \mathcal {P}\left[({\hat{m}}_{\scriptstyle \rm B}^*- \mu (\hat{z}), \hat{x}_{\scriptstyle 1}, \hat{c})_i \,|\, \mathrm{M}\right] \nonumber \\ &=& \left[(2\pi )^{3N}\det \left(\Sigma _d+A^\intercal \Sigma _l A\right) \right]^{-1/2} \nonumber \\ &&\times\, \exp \left[-{ {1\over 2}}\left(\hat{Z}-Y_{\scriptstyle 0} A\right) \left(\Sigma _d+A^\intercal \Sigma _l A\right)^{-1}\left(\hat{Z}-Y_{\scriptstyle 0} A\right) ^\intercal \,\right] \end{eqnarray}

(7)

where |$\hat{Z}=({\hat{m}}_{\scriptstyle {\rm B},1}^*-\mu (\hat{z}_{\scriptstyle 1}),\hat{x}_{\scriptstyle 1,1}, \hat{c}_{\scriptstyle 1},{\hat{m}}_{\scriptstyle {\rm B},2}^*-\mu (\hat{z}_{\scriptstyle 2}),\dots )$| is a 3N-dimensional row vector containing the distance modulus residual and light-curve data, |$Y_{\scriptstyle 0}=(M_{\scriptstyle \rm B,0},\,x_{\scriptstyle 1,0},\,c_{\scriptstyle 0},\,M_{\scriptstyle \rm B,0},\dots )$| is a 3N-dimensional row vector of the intrinsic Gaussian means and A is the block diagonal matrix that propagates |$Y_{\scriptstyle 0}$| to

\begin{equation*} Y_{\scriptstyle 0} A=(M_{\scriptstyle \rm B,0}-\alpha x_{\scriptstyle 1,0} + \beta c_{\scriptstyle 0},\,x_{\scriptstyle 1,0},\,c_{\scriptstyle 0},\,\dots ). \end{equation*}

Note that the cosmological model enters only explicitly through the conversion |$\mu (\hat{z})$| of the observed redshift to a distance modulus. There can, however, be implicit model dependence in the experimental covariance matrix¹ Σ_d or in corrections made to data prior to the analysis.

To implement the RH16 parametrization (4) we replace (6) by

\begin{eqnarray} \mathcal {P}[ (M_{\scriptstyle \rm B}, x_{\scriptstyle 1}, c)_i \,|\, \mathrm{M}] &=& \prod _{J=1}^4\,\prod _{i=1}^{N_{\scriptscriptstyle J}} \mathcal {N}(M_{\scriptstyle {\rm B},i},\,M_{\scriptstyle \rm B,0},\, \sigma _{\scriptstyle M_{{\rm B},0}}) \nonumber \\ &&\times\, \mathcal {N}(x_{\scriptstyle 1,i},\,x_{\scriptstyle 1,0,J} + x_{\scriptstyle z,J}\hat{z},\, \sigma _{x_{\scriptscriptstyle 1,0}})\,\mathcal\nonumber\\ &&\times\,{N}(c_{\scriptstyle i},\,c_{\scriptstyle 0,J}+c_{\scriptstyle z,J}\hat{z},\,\sigma _{c_{\scriptscriptstyle 0}}). \end{eqnarray}

(8)

We recover (7) with the one difference: in place of three repeated entries, the vector |$Y_{\scriptstyle 0}$| is now partitioned into different pieces for each subsample, |$Y_{\scriptstyle 0}=(M_{\scriptstyle \rm B,0},\,x_{\scriptstyle 1,0,1}+ x_{\scriptstyle z,1}{\hat{z}}_{\scriptstyle 1},\,c_{\scriptstyle 0,1} + c_{\scriptstyle z,1}{\hat{z}}_{\scriptstyle 1},\,\dots , M_{\scriptstyle \rm B,0},\,x_{\scriptstyle 1,0,2}+x_{\scriptstyle z,2}{\hat{z}}_{\scriptstyle i},\,c_{\scriptstyle 0,2} + c_{\scriptstyle z,2}{\hat{z}}_{\scriptstyle i},\dots , \dots ,M_{\scriptstyle \rm B,0},x_{\scriptstyle 1,0,4},c_{\scriptstyle 0,4})$|⁠.

From the likelihood (7), we can define frequentist confidence regions and goodness of fit measures or alternative Bayesian versions of these, following conventional statistical procedures summarized in Appendix C.

In practice, estimating the Bayesian evidence is a computationally intensive task, much more so than what is required to obtain parameter estimates. We use standard Markov Chain Monte Carlo methods to sample parameter space. We estimate the evidence using the publicly available multinest (Feroz, Hobson & Bridges 2009) code,² with python interface pymultinest (Buchner et al. 2014), for the efficient evaluation of the evidence integral (C5) with likelihood (7). The accuracy of the Bayesian evidence estimate is controlled by the number of ‘live’ points, |$n_{\scriptstyle \rm live}$|⁠, with an error |$\sigma \sim O(n^{-1/2}_{\scriptstyle \rm live})$|⁠. In our analysis, we choose 1000 points for the 8 or 9 parameter base model and add 100 more points for each additional parameter.

3.2 Choice of priors

Given the sensitivity of the Bayes factor to priors it is important these are chosen as objectively as possible. The choice of priors is summarized in Table 1.

Table 1.

All nuisance parameters in each model have identical priors. In the timescape model |$\Omega _{\scriptstyle M0}$| is defined in terms of f_v0, hence we take the latter to be the more ‘fundamental’ parameter and assign the prior to it.

Cosmological parameters	Prior distribution	Range
Timescape: fv0	Uniform	[0.500, 0.799]
Flat ΛCDM: |$\Omega _{\scriptstyle M0}=1-\Omega _{\scriptstyle \Lambda 0}$|	Uniform	[0.143, 0.487]
Nuisance parameters	Prior distribution	Range
α	Uniform	[0, 1]
|$x_{\scriptstyle 1,0}$|	Uniform	[−20, 20]
|$\sigma _{\scriptstyle x_{1,0}}$|	Uniform on |$\log _{\scriptstyle 10}\sigma _{\scriptstyle x_{1,0}}$|	[−5, 2]
β	Uniform	[0, 4]
|$c_{\scriptstyle 0}$|	Uniform	[−20, 20]
|$\sigma _{\scriptstyle c_0}$|	Uniform on |$\log _{\scriptstyle 10}\sigma _{\scriptstyle c_0}$|	[−5, 2]
|$M_{\scriptstyle \rm B,0}$|	Uniform	[−20.3, −18.3]
|$\sigma _{\scriptstyle M_{\rm B,0}}$|	Uniform on |$\log _{\scriptstyle 10}\sigma _{\scriptstyle M_{\rm B,0}}$|	[−5, 2]
Additional stretch and colour
parameters for models II–VIII	Uniform	[−20, 20]

3.2.1 Nuisance parameters

All nuisance parameters are common to both timescape and ΛCDM models and we therefore assign the same priors to both models. Where possible,³ we adopt priors that have been used in previous Bayesian studies of the SALT2 method (March et al. 2011; Shariff et al. 2016). The standard deviations |$\lbrace \sigma _{\scriptstyle x_{1,0}}, \sigma _{\scriptstyle c_0}, \sigma _{\scriptstyle M_{\rm B,0}}\rbrace$| are ‘scale’ parameters (of the residuals) and so it is more appropriate to assign a log-uniform prior to these parameters. The priors for the nuisance parameters are wide to ensure the most likely regions of parameter space are supported, and provided they are wide enough, this will have no overall effect on the Bayes factor (as the evidence of each model will be similarly scaled).

3.2.2 Cosmological parameters

Only one free cosmological parameter can be constrained by supernovae: |$\Omega _{\scriptstyle M0}$| for spatially flat ΛCDM or f_v0 for the timescape model. Conventionally, the combination of |$\Omega _{\scriptstyle M0}$| and |$H_{\scriptstyle 0}$| for the standard cosmology is strongly constrained by the CMB acoustic peaks (Ade et al. 2016). Measurements of the Baryon Acoustic Oscillation (BAO) scale in galaxy clustering statistics (Aubourg et al. 2015; Alam et al. 2017) at low redshifts and the Lyman α forest (Aubourg et al. 2015; Delubac et al. 2015) provide independent constraints. In the case of the timescape model, however, our ability to model the CMB is still limited by systematic uncertainties of 8–13 per cent (Nazer & Wiltshire 2015).

We therefore determine priors for f_v0 in the timescape model using best present knowledge. For the CMB, we use results of a model-independent analysis of the acoustic peaks (Aghamousa et al. 2015) with Planck satellite data, and choose a prior from a 95 per cent confidence fit of the angular scale of the sound horizon. To date BAO studies all implicitly assume the FLRW model, and do not yet provide an equivalent model-independent constraint. We therefore adopt a prior using FLRW-model estimates of the angular diameter of the BAO scale, including the full range of values which are currently in tension (Aubourg et al. 2015; Delubac et al. 2015; Alam et al. 2017). We take generous 95 per cent confidence limits determined by assuming that either the low-redshift galaxy clustering results (Aubourg et al. 2015; Alam et al. 2017) or the z = 2.34 Lyman α results (Delubac et al. 2015) are correct. Priors for the spatially flat ΛCDM model are determined by an identical methodology. Further details are given in Appendix D.

4 RESULTS

4.1 Analysis with supernova parameters drawn from global Gaussian distributions

Since there is a degeneracy between the Hubble constant, |$H_{\scriptstyle 0}$|⁠, and the magnitude, |$M_{\scriptstyle \rm B,0}$|⁠, we fix |$H_{\scriptstyle 0}$| for each model. The value of |$M_{\scriptstyle \rm B,0}$| then depends on this choice. We are then left with one free cosmological parameter, the matter density parameter |$\Omega _{\scriptstyle M0}$| in the spatially flat ΛCDM model, and the present epoch void fraction f_v0 in the timescape model. We can alternatively define an effective ‘dressed matter density parameter’ |$\Omega _{\scriptstyle M0}={\textstyle {1\over 2}}(1-f_{\rm v0})(2+f_{\rm v0})$| (Wiltshire 2007a, 2009), which takes similar numerical values to the concordance ΛCDM model, allowing likelihood functions to be plotted on the same scale. (This parameter does not obey the Friedmann equation sum rule, however.) The nine parameters |$\lbrace \Omega _{\scriptstyle M0},\alpha ,x_{\scriptstyle 1,0},\sigma _{\scriptstyle x_{\rm 1,0}},\beta ,c_{\scriptstyle 0},\sigma _{\scriptstyle c_{\rm 0}},M_{\scriptstyle \rm B,0},\sigma _{\scriptstyle M_{\rm B,0}}\rbrace$| are then fit for each model by determining the likelihood function with all parameters other than |$\Omega _{\scriptstyle M0}$| (or f_v0) treated as nuisance parameters. The empty universe has eight parameters since |$\Omega _{\scriptstyle M0}=0$|⁠.

4.1.1 Statistical homogeneity scale cuts

An important systematic issue in the timescape cosmology is the fact that an average expansion law only holds on scales greater than the statistical homogeneity scale (SHS) ≳ 70–120 h⁻¹Mpc (Hogg et al. 2005; Scrimgeour et al. 2012). This corresponds to a CMB rest-frame redshift of order |$z\mathop {\sim }\limits 0.023$|–0.04. In fact, SneIa analyses using the MLCS method have typically excluded SneIa below a cutoff at z = 0.024 (Riess et al. 2007). However, the JLA catalogue (Betoule et al. 2014) includes 53 SneIa, with z < 0.024.

Following Smale & Wiltshire (2011), we determine cosmological model distances in the CMB frame, but make a redshift cut at the SHS, taken at |$\mathop {\sim }\limits 100\,h^{-1}\hbox{Mpc}$|⁠. Furthermore, to examine the effect of the SHS cut on the fit of light-curve parameters, we perform the entire analysis while progressively varying the minimum redshift in the range 0.01 ≤ z_min < 0.1; i.e. up to a redshift three times larger than the SHS. Systematic effects associated with the SHS can then be revealed. Our key results will be quoted for a cut at |$z_{\scriptstyle \rm SHS}=0.033$| in the CMB rest frame. The best-fitting MLE parameters with z_min = 0.033 are presented in Table 2.

Table 2.

Best-fitting MLE parameters corresponding to the likelihood |$\mathcal {L} ({\rm Data} | \mathrm{M})$| with the model M representing the cosmological model, the SALT2 procedure and the intrinsic distributions of SneIa parameters. SneIa at redshifts z < 0.033 (statistical homogeneity scale) are excluded. Confidence limits are given for the one free cosmological parameter. In the timescape case this corresponds to |$f_{\rm v0}=0.778^{+0.063\,(1\sigma )\,0.155\,(2\sigma )}_{-0.056\,(1\sigma )\,0.104\,(2\sigma )}$|⁠. The value of |$M_{\scriptstyle \rm B,0}$| is obtained for h = 0.668 for the timescape, and h = 0.7 for the two FLRW models. The difference of parameters from NGS16 is principally due to the SHS cut at z_min = 0.033, the effect of which is seen in Fig. 2.

Model	|$\Omega _{\scriptstyle M0}$|	α	|$x_{\scriptstyle 1,0}$|	|$\sigma _{\scriptstyle x_{1,0}}$|	β	c0	|$\sigma _{c_0}$|	|$M_{\scriptstyle \rm B,0}$|	|$\sigma _{\scriptstyle M_{\rm B,0}}$|
Timescape	|$0.309^{+0.070\,(1\sigma )\,0.127\,(2\sigma )}_{-0.088\,(1\sigma )\,0.210\,(2\sigma )}$|	0.134	0.1050	0.899	3.13	− 0.0211	0.0689	− 19.1	0.104
Spatially flat ΛCDM	|$0.365^{+0.033\,(1\sigma )\,0.066\,(2\sigma )}_{-0.031\,(1\sigma )\,0.060\,(2\sigma )}$|	0.134	0.1061	0.899	3.14	− 0.0215	0.0688	− 19.0	0.104
Empty universe	−	0.133	0.1013	0.900	3.13	− 0.0204	0.0690	− 19.0	0.106

For the priors given in Table 1, the Bayesian evidence in favour of the timescape model relative to the spatially flat ΛCDM model is ln B = 0.085 ± 0.01 with z_min = 0.033, or lnB = 0.600 ± 0.007 with z_min = 0.024. Since |ln B| < 1 the two models⁴ are statistically indistinguishable. This conclusion is insensitive to O(1σ) changes to the width of the uniform priors on f_v0 and |$\Omega _{\scriptstyle M0}$|⁠, or to variations of the minimum redshift as shown in Fig. 2(e).

Figure 2.

MLE best-fitting parameters, and Bayes factor, for the NGS16 model as z_min is varied: (a) α; (b) β; (c) |$x_{\scriptstyle 1,0}$|⁠; (d) |$c_{\scriptstyle 0}$|⁠; (e) |$\Omega _{\scriptstyle M0}$| (or f_v0) with 1σ bounds; (f) |$\ln B=\ln (E_{\scriptstyle \rm TS}/E_{\rm model})$|⁠. Vertical dotted lines at z_min = 0.024 and z_min = 0.033 indicate the expected rough redshift range of an emerging statistical homogeneity scale.

While the Bayes factors do not show significant variation with z_min, the values of particular best-fitting light-curve parameters show a marked change at the SHS. As shown in Fig. 2, there is a marked 30 per cent jump in |$c_{\scriptstyle 0}$| as z_min is varied from 0.01 up to z ≃ 0.033, when compared to the subsequent fluctuations if z_min is increased up to 0.1. For |$x_{\scriptstyle 1,0}$| there is a similar jump, although a linear trend remains in the range 0.033 < z_min < 0.1. The β parameter shows a small (3 per cent) jump up to the SHS followed by ±1 per cent fluctuations, while α remains relatively constant, fluctuating by ±2 per cent over the whole range.

Since the light-curve parameters are remarkably close for all three cosmologies while showing a jump as the SHS emerges (distinct from the residual |$c_{\scriptstyle 0}$| trend for the empty model with z_min ≳ 0.05), there is clear evidence for some systematic effect at precisely the scale where we expect it.

4.2 Analysis with linear redshift variation for mean stretch and colour parameters

Although RH16 considered four distinct subsamples, the mean stretch parameter actually shows a global increasing trend in the ΛCDM case evident in Rubin & Hayden (2016, fig. 1, left-hand panels). Our procedure of varying the minimum redshift cut on the whole sample also isolates any global trend. Such a trend is indeed evident in Fig. 2(c) beyond the SHS, with x_{1, 0} increasing 40 per cent as the minimum sample redshift increases from z_min = 0.033 to z_min = 0.1. Beyond z_min = 0.034, an equivalent global trend in the mean colour parameter, |$c_{\scriptstyle 0}$|⁠, is not evident in Fig. 2(d), however, except in the case of the empty universe, which shows a 13 per cent decrease in c₀ between z_min = 0.034 and z_min = 0.1. A global shift in x_{1, 0} would seem more consistent with an astrophysical systematic in the source population, rather than sampling biases with different thresholds for different samples.

To fully understand the differences that arise on making the RH16 changes (4), we have also investigated the effect of adding a smaller number of free parameters, by considering linear z relations in just one of the parameters |$x_{\scriptstyle 1,0}$| or |$c_{\scriptstyle 0}$|⁠, and the difference between global linear relations and a split by subsamples. The advantage of our fully Bayesian approach is that not only can we compare the relative Bayesian evidence for different cosmological models with the same light-curve parameters, but we can also compare the merits of different empirical light-curve models. The values of the Bayesian evidence are shown in Table 3, along with a selection of parameters. The changes to the parameters α and β are negligible between models, and are not tabulated. We do not tabulate all additional (up to 12) parameters for the case of the split subsamples, but an average.

Table 3.

Selected parameters fit for z_min = 0.033, with the following empirical models for light-curve parameters: (I) constant |$x_{\scriptstyle 1,0}$|⁠, constant |$c_{\scriptstyle 0}$|⁠; (II) global linear |$x_{\scriptstyle 1,0}$|⁠, constant |$c_{\scriptstyle 0}$|⁠; (III) split linear |$x_{\scriptstyle I,1,0}$|⁠, constant |$c_{\scriptstyle 0}$|⁠; (IV) constant |$x_{\scriptstyle 1,0}$|⁠, global linear |$c_{\scriptstyle 0}$|⁠; (V) constant |$x_{\scriptstyle 1,0}$|⁠, split linear |$c_{\scriptstyle 0,I}$|⁠; (VI) global linear |$x_{\scriptstyle 1,0}$|⁠, global linear |$c_{\scriptstyle 0}$|⁠; (VII) split linear |$x_{\scriptstyle 1,0}$|⁠, split linear |$c_{\scriptstyle 0}$|⁠; (VIII) global linear |$x_{\scriptstyle 1,0}$|⁠, split linear |$c_{\scriptstyle 0}$|⁠.-

Timescape	k	|$\Omega _{\scriptstyle M0}$|	|$\langle x_{\scriptstyle 1,0}\rangle$|	|$\langle c_{\scriptstyle 0}\rangle$|	|$\langle x_{\scriptstyle z}\rangle$|	|$\langle c_{\scriptstyle z}\rangle$|	BIC	|$\Delta _{\scriptstyle \rm BIC}$|	ln E	Cb	|$\ln B_{\scriptstyle 1}$|	fv0
I (NGS16)	9	0.309	0.105	− 0.021			− 185.0		80.38	8.53		0.778
II	10	0.278	− 0.073	− 0.021	0.511		− 199.6	− 14.6	81.67	9.28	− 1.28	0.802
III	15	0.278	− 0.183	− 0.021	0.806		− 169.9	15.1	66.46	14.11	13.93	0.801
IV	10	0.000	0.104	0.002		− 0.065	− 249.9	− 64.9	78.52	9.19	1.86	1.000
V	15	0.010	0.092	0.054		− 0.351	− 157.7	27.3	81.29	14.15	− 0.91	0.993
VI	11	0.000	− 0.071	0.001	0.499	− 0.062	− 189.7	− 4.7	79.35	10.21	1.94	1.000
VII (RH16)	21	0.000	− 0.123	0.054	0.490	− 0.348	− 200.0	− 15.0	65.85	19.47	14.53	1.000
VIII	16	0.000	− 0.085	0.061	0.501	− 0.348	− 229.3	− 44.3	82.03	15.30	− 1.65	1.000
Flat ΛCDM	k	|$\Omega _{\scriptstyle M0}$|	|$\langle x_{\scriptstyle 1,0}\rangle$|	|$\langle c_{\scriptstyle 0}\rangle$|	|$\langle x_{\scriptstyle z}\rangle$|	|$\langle c_{\scriptstyle z}\rangle$|	BIC	|$\Delta _{\scriptstyle \rm BIC}$|	ln E	Cb	|$\ln B_{\scriptstyle 1}$|	|$\ln B_{\scriptstyle 2}$|
I (NGS16)	9	0.365	0.106	− 0.021			− 192.5		80.30	8.93		0.08
II	10	0.353	− 0.069	− 0.021	0.503		− 241.2	− 48.7	81.64	10.01	− 1.34	0.03
III	15	0.353	− 0.186	− 0.021	0.847		− 159.8	32.7	66.62	14.60	13.68	− 0.16
IV	10	0.303	0.106	− 0.002		− 0.057	− 192.9	− 0.4	79.60	9.98	0.70	− 1.08
V	15	0.296	0.093	0.052		− 0.354	− 228.5	− 36.0	83.77	14.87	− 3.47	− 2.47
VI	11	0.292	− 0.069	− 0.002	0.501	− 0.057	− 179.1	13.4	80.87	10.97	2.89	− 1.52
VII (RH16)	21	0.286	− 0.127	0.051	0.534	− 0.352	− 155.3	37.2	68.97	20.58	11.33	− 3.12
VIII	16	0.286	− 0.080	0.059	0.499	− 0.354	− 232.8	− 40.3	84.95	15.89	− 4.65	− 2.92
Empty	k	|$\Omega _{\scriptstyle M0}$|	|$\langle x_{\scriptstyle 1,0}\rangle$|	|$\langle c_{\scriptstyle 0}\rangle$|	|$\langle x_{\scriptstyle z}\rangle$|	|$\langle c_{\scriptstyle z}\rangle$|	BIC	|$\Delta _{\scriptstyle \rm BIC}$|	ln E	Cb	|$\ln B_{\scriptstyle 1}$|	|$\ln B_{\scriptstyle 2}$|
I (NGS16)	8	–	0.101	− 0.020			− 181.5		78.18	8.11		2.20
II	9	–	− 0.078	− 0.019	0.517		− 190.1	− 8.6	79.92	9.02	− 1.74	1.75
III	14	–	− 0.095	− 0.020	0.749		− 218.9	− 37.4	64.43	13.75	13.76	2.03
IV	9	–	0.098	0.002		− 0.054	− 185.7	− 4.2	78.56	9.05	− 0.37	− 0.03
V	14	–	0.087	0.054		− 0.336	− 180.4	1.1	79.85	14.19	− 1.66	1.45
VI	10	–	− 0.072	0.002	0.489	− 0.051	− 186.3	− 4.8	79.62	10.17	0.23	− 0.27
VII (RH16)	20	–	− 0.122	0.054	0.460	− 0.332	− 198.7	− 17.2	64.31	18.68	13.88	1.55
VIII	15	–	− 0.081	0.061	0.482	− 0.334	− 221.4	− 39.9	80.74	14.89	− 2.55	1.30

Notes. k ≡ number of free parameters; quantities |$\langle \Phi \rangle \equiv (\sum N_{\scriptstyle I}\,\Phi _{\scriptstyle I})/(\sum _{I=1}^4 N_{\scriptstyle I})$| denote an average over subsamples with I = 1…4 for |$x_{\scriptstyle 1,0,I}$|⁠, |$c_{\scriptstyle 0,I}$| and I = 1…3 for |$x_{\scriptstyle z,I}$|⁠, |$c_{\scriptstyle z,I}$| for split models or 〈Φ〉 ≡ Φ otherwise; BIC = Bayesian Information Criterion; E = Bayesian evidence; C_b = Bayesian complexity; |$\Delta _{\scriptstyle \rm BIC}=\hbox{BIC}_{\rm model}- \hbox{BIC}_{\scriptstyle \rm I}$| and |$\ln B_{\scriptstyle 1}=\ln (E_{\scriptstyle \rm I}/E_{\rm model})$| are evaluated with cosmological model fixed; |$\ln B_{\scriptstyle 2}=\ln (E_{\scriptstyle \rm TS}/E_{\rm model})$| is evaluated with light-curve model fixed.

4.2.1 Stretch parameter |$x_{\scriptstyle 1,0}$|

Consistent with remarks above, relative to the baseline NGS16 model I, light-curve model II provides positive (but not strong) Bayesian evidence for a global linear trend in |$x_{\scriptstyle 1,0}$| independent of cosmological model, with |$\ln B_{\scriptstyle 1}=1.28$|⁠, 1.34, 1.74 for the timescape, ΛCDM and empty models, respectively. The BIC evidence for the same conclusion is very strong (timescape, ΛCDM models) or strong (empty model). By contrast, model III gives strong evidence |$|\ln B_{\scriptstyle 1}|>13$| against a split linear law in |$x_{\scriptstyle 1,0}$| independent of cosmological model. The Bayesian penalty for introducing new empirical parameters depends on the choice of the priors, but our conclusion is robust to reasonable changes. Furthermore, the frequentist BIC statistic |$\Delta _{\scriptstyle \rm BIC}$| also strongly disfavours model III relative to models I, II in the ΛCDM and timescape cases.

4.2.2 Colour parameter |$c_{\scriptstyle 0}$|

In contrast to the stretch parameter, results involving the colour parameter vary greatly with cosmological model. Despite model IV having the global minimum BIC statistic for timescape, |$\ln B_{\scriptstyle 1}$| shows no significant evidence⁵ for any global linear redshift law. Relative Bayesian evidence for a split linear law in |$c_{\scriptstyle 0}$| is marginal for timescape, positive for the empty universe and strong for ΛCDM, with |$\ln B_{\scriptstyle 1}=-0.91,-1.66$| and −3.47, respectively.

The original RH16 model VII suffers similar problems to model III in terms of Bayesian evidence, evidently on account of the split linear law in |$x_{\scriptstyle 1,0}$|⁠. However, model VIII has the strongest Bayesian evidence of all models. It adds a global linear redshift law in |$x_{\scriptstyle 1,0}$| to model V. The improvement in Bayesian evidence for model VIII relative to model V is marginal for timescape and the empty universe, and positive for ΛCDM, with |$\Delta \ln B_{\scriptstyle 1}=-0.74$|⁠, −0.89 and −1.18, respectively.

Although |$\ln B_{\scriptstyle 2}$| for model VIII gives positive (but not strong) relative evidence for ΛCDM over timescape, any conclusion drawn from this depends on additional empirical light-curve parameters which now depend on the cosmological model⁶. Furthermore, since the maximum likelihoods are comparable, the difference in Bayesian evidence is primarily due to the timescape maximum likelihood being driven to the unphysical limit f_v0 → 1 for any light-curve model with linear variations in |$c_{\scriptstyle 0}$|⁠, as is shown in Fig. 3, which compares likelihoods in |$\Omega _{\scriptstyle M0}$| (or f_v0) for the NGS16 and RH16 models for two choices of z_min.

Figure 3.

Profile likelihoods in |$\Omega _{\scriptstyle M0}$| and f_v0 for model I (NGS16) and model VII (RH16): (a) NGS16, z_min = 0.024; (b) NGS16, z_min = 0.033; (c) RH16, z_min = 0.024; (d) RH16, z_min = 0.033. Model IV, V, VI and VIII results are very similar to model VII.

Very similar results were found by Smale & Wiltshire (2011) in applying SALT2 without the methodology of NGS16, leading to a large discrepancy in the predictions of the SALT2 and MLCS2k2 fitters for timescape. Since direct application of the NGS16 methodology to the JLA catalogue agrees with some previous MLCS2k2 fits to smaller data sets (Leith et al. 2008; Smale & Wiltshire 2011), we conclude that systematics similar to linear redshift variations in |$c_{\scriptstyle 0}$| may be the key to earlier discrepancies.

4.2.3 Cosmological model dependence of linear redshift changes to SALT2 methodology

To understand the origin of such differences consider the Taylor series expansions (A13)–(A15) for the timescape, spatially flat ΛCDM and empty universe models, as given in Appendix A. Leading coefficients for (A13) and (A14) are shown graphically in Fig. 4 as a function of the free cosmological parameter.

Figure 4.

Coefficients in the Taylor series (A13), (A14), |$\mu =\mu _{\scriptstyle 0}(z)+\sum _{n=1}\mu _n z^n$|⁠, of the spatially flat ΛCDM and timescape models, as a function of the free cosmological parameter, |$\Omega _{\scriptstyle M0}$| or f_v0. For timescape, the coefficients |$\mu _{\scriptstyle {\rm TS},n}$| are very slowly varying monotonic functions of f_v0 on the range 0.6 < f_v0 ≤ 1, whereas the coefficients |$\mu _{\scriptstyle \Lambda {\rm CDM},n}$| are polynomials of order n. For each n, |$|\mu _{\scriptstyle {\rm TS},n}|<|\mu _{\scriptstyle \Lambda {\rm CDM},n}|$|⁠, reflecting the ‘flatter’ distance modulus (cf. Fig. 1). Linear changes of |$\Omega _{\scriptstyle M0}$| can become degenerate with empirical light-curve parameters linear in z for parameters close to the minimum of |$\mu _{\scriptstyle \Lambda {\rm CDM},2}$| at |$\Omega _{\scriptstyle M0}=0.296$|⁠. The change in the coefficients between NGS16 model I and models VII/VIII is indicated.

All cosmologies show improvement to a global increase in |$x_{\scriptstyle 1,0}$| with redshift and piecewise decreases in |$c_{\scriptstyle 0}$|⁠, including the empty model which has no free parameter to adjust. However, if we incorporate linear corrections x_zz to |$x_{\scriptstyle 1,0}$|⁠, or c_zz to |$c_{\scriptstyle 0}$|⁠, in the SALT2 relation, then the difference of (2) and (A13)–(A15) gives a potential degeneracy between empirical parameters x_z or c_z and changes in the free cosmological parameter if the linear term in (A13)–(A15) can be changed without greatly altering the next most important O(z²) term. Such a possibility is admitted by ΛCDM but not timescape.

For ΛCDM, the O(z) term in (A14) is linear in |$\Omega _{\scriptstyle M0}$|⁠, and the O(z²) term is quadratic in |$\Omega _{\scriptstyle M0}$| with a minimum at |$\Omega _{\scriptstyle M0}= {\textstyle {8\over 27}}=0.296$|⁠. For model V with split linear redshift laws in |$c_{\scriptstyle 0}$| only, the best-fitting |$\Omega _{\scriptstyle M0}$| coincides precisely with this minimum. The decrease in |$\Omega _{\scriptstyle M0}$| by adding a global linear z dependence to |$x_{\scriptstyle 1,0}$| is approximately the same, |$\Delta \Omega _{\scriptstyle M0}=-0.01$|⁠, in going from model V to VIII, or from model I to II. The difference in (A14) between models I and VII/VIII,

\begin{eqnarray} &&\!\!\!\!&&{ \mu _{\scriptstyle \Lambda {\rm CDM}}(0.286)-\mu _{\scriptstyle \Lambda {\rm CDM}}(0.365)=0.1287\,z- 0.0085\,{z}^{2}} \nonumber \\ &&\,\, - 0.0481\,{z}^{3}+ 0.0249\,{z}^{4}+ 0.0161\,{z}^{5}- 0.0232\,{z}^{6}+\cdots \end{eqnarray}

(9)

is dominated by the linear redshift changes, with negligible changes in the O(z²) term.

By contrast the terms in the Taylor series (A13) for timescape are very slowly varying monotonic functions of f_v0 on the range 0.6 < f_v0 ≤ 1.0 (as shown in Fig. 4), so changes in |$\mu _{\scriptstyle \rm TS}$| are much more constrained. The difference in (A13) between models I and VII/VIII, is

\begin{eqnarray} &&\!\!\!\!&&{ \mu _{\scriptstyle \rm TS}(1.0)-\mu _{\scriptstyle \rm TS}(0.778)=0.0674\,z+ 0.0444\,{z}^{2}} \nonumber \\ &&\,\,\,-0.0242\,{z}^{3} + 0.0190\,{z}^{4}- 0.0193\,{z}^{5}+ 0.0173\,{z}^{6}+\cdots \end{eqnarray}

(10)

A large change in f_v0 is required make changes in |$\mu _{\scriptstyle \rm TS}$| comparable to the ΛCDM case, and the effect of increasing f_v0 increases both the O(z) and O(z²) terms. As seen in Fig. 5, the likelihood is consequently peaked along narrow ridges in the (⁠|$x_{\scriptstyle z},f_{\rm v0}$|⁠) and (⁠|$c_{\scriptstyle z,J},f_{\rm v0}$|⁠) planes, with almost constant values of |$x_{\scriptstyle z}$| and |$c_{\scriptstyle z,J}$| and no upper bound on f_v0.

Figure 5.

Likelihood function contours for model VIII with z_min = 0.033 projected in the planes: (a) (x_z, Ω_M0); (b) (c_{z, 1}, Ω_M0) (SNLS sample, mean redshift 〈z〉 = 0.636); (c) (c_{z, 2}, Ω_M0) (SDSS sample, mean redshift 〈z〉 = 0.199); (d) (c_{z, 3}, Ω_M0) (low z sample with z > 0.033, mean redshift 〈z〉 = 0.0495). 67.3 per cent, 95.5 per cent and 99.7 per cent confidence contours are shown. In panel (a) x_z contours for model II are also shown to demonstrate the effect of adding the c_{z, J} parameters. For spatially flat ΛCDM, the maximum likelihood is driven to the vicinity of the minimum |$\Omega _{M0}= \frac{8}{27}$| of the O(z²) Taylor series term (A14). The timescape Taylor series (A13) consists of slowly varying monotonic functions of f_v0, and the maximum likelihood is driven to the edge of parameter space, f_v0 → 1.

5 DISCUSSION

Our study has a number of important consequences. First, the timescape and spatially flat ΛCDM model luminosity distance–redshift fits to the JLA catalogue are statistically indistinguishable using either the approach of NGS16, or with modifications to only the mean stretch parameter. As shown in Table 3, the Bayesian complexity, C_b, is lower (better) for timescape than for ΛCDM, for every choice of light-curve model.

This completely reframes a debate (Nielsen et al. 2016; Rubin & Hayden 2016; Shariff et al. 2016; Haridasu et al. 2017; Tutusaus et al. 2017) about whether cosmic acceleration is marginal or not, within the confines of a FLRW expansion history. Current supernova data does not distinguish between the standard ΛCDM model and the non-FLRW expansion history of the timescape model, which has non-zero apparent cosmic acceleration but with a marginal amplitude. The apparent deceleration parameter (A12) for the best-fit value of Table 2 is |$q_{\scriptstyle 0}\equiv q(f_{\rm v0})=-0.043^{+0.004}_{-0.000}$|⁠.

Within the class of FLRW models, the significance of cosmic acceleration is often assessed by comparison to the empty universe model. That model is unphysical, since standard nucleosynthesis and recombination can never occur in a universe with a(t) ∝ t regardless of its matter content.⁷ The timescape model has positive |$\ln B_{\scriptstyle 2}$| compared to the empty universe.⁸ Nonetheless, |$|\mu _{\scriptstyle \rm TS}(z)-\mu _{\scriptstyle \rm empty}(z)|<|\mu _{\scriptstyle \Lambda {\rm CDM}}(z)-\mu _{\scriptstyle \rm empty}(z)|$| (cf. Fig. 1) at late epochs, for a simple physical reason. The timescape model is void dominated at z < 1, and the expansion of individual voids is close to an empty universe. While the timescape model has apparent acceleration at late epochs, its expansion law is closer to that of an empty universe than that of the ΛCDM model.

The second important consequence of our study is that allowing linear changes with redshift to the mean colour parameter, |$c_{\scriptstyle 0}$|⁠, produces cosmological model dependence. Since the redshift–distance relation of the timescape model effectively interpolates (Wiltshire 2009, 2014) between those of spatially flat ΛCDM models with different values⁹ of |$\Omega _{\scriptstyle M0}$|⁠, particular care must be taken with piecewise linear relations in redshift.

The improved 16 parameter model VIII (this being a better fit than the original 21 parameter RH16 model) has positive Bayesian evidence for ΛCDM relative to the timescape model. However, this is contingent on degeneracies in the likelihood function between the free cosmological parameter and additional empirical parameters. The RH16 parametrization allows the ΛCDM deceleration parameter |$q_{\scriptstyle 0}=-1+{\textstyle {3\over 2}}\Omega _{\scriptstyle M0}$| contained in the O(z) term of the Taylor series (A14) to be adjusted¹⁰ near the global minimum |$\Omega _{\scriptstyle M0}=0.296$| of the O(z²) term in (A14). However, the same procedure drives the timescape free parameter to an unphysical limit, f_v0 → 1. No fundamental model underlies the empirical parametrization (4). Variations in |$c_{\scriptstyle 0}$| are most plausibly related to selection effects, given we cannot fit them by a global law. However, selection effects would be more correctly modelled by removing the tail of a Gaussian distribution rather than shifting its mean linearly in redshift.

Our results show that NGS16 did not account for every possible selection bias that remains in the JLA catalogue, consistent with some comments of RH16. Nonetheless, NGS16 are correct to point out the possible pitfalls in fitting SneIa data when empirical light-curve parameters are mixed with the cosmological parameters of a single class of cosmological models. If SneIa are to be used to distinguish cosmological models, then systematic uncertainties and selection biases should be corrected in as model-independent manner as possible before the data is reduced.

A related issue which remains to be explored is the extent to which the corrections for selection biases that have already been made in the JLA catalogue depend on the FLRW model. Betoule et al. (2014) follow a procedure of Mosher et al. (2014, section 6.2), who used the snana package to estimate selection biases. Simulated data (using the FLRW model) is used in such estimates. While efforts have been made to consider different dark energy equations of state (Mosher et al. 2014), models which do not satisfy the Friedmann equation fall outside the scope of such analyses.

Whether or not the timescape model is ultimately a better fit than the standard FLRW model, it may provide a useful diagnostic tool in comparing methods for SneIa light-curve reduction purely at the empirical level. In particular, it has an analytic non-FLRW redshift–distance relation which is very close to that of the ΛCDM model, but which is considerably more constrained in the free parameter f_v0 than the ΛCDM model is in |$\Omega _{\scriptstyle M0}$|⁠.

Finally, Fig. 2(b), (c) and (d) shows evidence for a ≃ 100 h⁻¹Mpc statistical homogeneity scale which has an effect on global fits of light-curve parameters – most notably a 30 per cent shift of |$c_{\scriptstyle 0}$| – independent of the cosmological model. These systematics, which occur at a scale relevant from independent observations (Hogg et al. 2005; Scrimgeour et al. 2012), must be explained irrespective of the cosmological model.

Acknowledgements

We thank Chris Blake, Tamara Davis, Ahsan Nazer, David Rubin, Subir Sarkar and Bonnie Zhang for helpful discussions and correspondence, and Thomas Buchert for hospitality at the ENS, Lyon, France. This work was supported by Catalyst grant CSG-UOC1603 administered by the Royal Society of New Zealand.

Code availability: The code and data used in this analysis are available at https://doi.org/10.5281/zenodo.831360.

1

The propagation of the error σ_z to σ_μ depends on the model. However, by (A13)–(A15), to leading order for small z, σ_μ ≃ 5σ_z/(z ln 10) for all cases.

2

This package is based on the Nested sampling algorithm (Skilling 2004).

3

Given the complications introduced by empirical changes (4) to |$x_{\scriptstyle 1,0}$|⁠, |$c_{\scriptstyle 0}$|⁠, we adopt uniform priors for these parameters.

4

Both models have positive relative Bayesian evidence compared to the empty model. Although the evidence is not particularly strong, |$|\ln B_{\scriptstyle 2}|\mathop {\sim }\limits 2.2$| incorporates priors which demand standard recombination for both ΛCDM and timescape. By that criterion the empty model is simply ruled out.

5

The empty universe has marginal evidence, consistent with Fig. 2(d) for z_min > 0.05.

6

Some of the largest differences occur in the SNLS subsample: |$c_{\scriptstyle \rm 0,\scriptstyle \Lambda CDM}=0.0483$| and |$c_{\scriptstyle \rm 0,TS}=0.0565$|⁠, a 17 per cent difference. For the NGS16 model, by contrast, differences are 2 per cent.

7

In particular, the R_h = ct model is unphysical for this reason (Lewis, Barnes & Kaushik 2016).

8

This is true for the NGS16 model I and all light-curve models for which |$\ln B_{\scriptstyle 1}$| shows an improvement independent of cosmology, viz. models II, V, VIII.

9

Note that the Planck best-fitting value |$\Omega _{\scriptstyle M0}=0.3175$| (Ade et al. 2016) is lower than the best-fitting value for the spatially flat ΛCDM model value |$\Omega _{\scriptstyle M0}=0.365$| from Table 2, consistent with the timescape expectation.

10

For the NGS16 model I and models VII/VIII one has best fits |$q_{\scriptstyle 0}=-0.453$| and |$q_{\scriptstyle 0}=-0.571,$| respectively. The respective spatially flat ΛCDM values quoted by Rubin & Hayden (2016), namely |$q_{\scriptstyle 0}=-0.412$|⁠, −0.552, (or |$\Omega _{\scriptstyle M0}=0.392$|⁠, 0.299), differ mostly on account of our SHS cut at z_min = 0.033.

11

Since we do not constrain |$H_{\scriptstyle 0}$|⁠, we do not specifically investigate the relationship between light-curve parameters and determinations of the local Hubble constant, which have been discussed in the past (Conley et al. 2007; Jha et al. 2007; Hicken et al. 2009; Smale & Wiltshire 2011; Wiltshire et al. 2013; McKay & Wiltshire 2016). In the timescape model, higher average values of |$H_{\scriptstyle 0}$| are expected below the SHS.

12

We use the NASA/IPAC Extragalactic Database standard, 371 km  s⁻¹ in the direction (ℓ, b) = (264.14°, 48.26°) (Fixsen et al. 1996).

13

We note that only models II, IV, VI are extensions of the 9 parameter base model, i.e. model I is nested in II, IV and VI, while II and IV are nested models of VI. Model V is nested in model VIII.

14

The AIC and BIC statistics also include a penalty using simple approximations to Bayesian methods which derive from different assumptions about the priors. The factor of two difference in the IC evidence scale (Kass & Raftery 1995) reflects the factor of 2 multiplying |$\ln \left(\mathcal {L}_{\rm max}\right)$| in the definitions (C2), (C3).

15

If we were to use conventional narrower priors for ΛCDM then the timescape model is either unfairly advantaged or disadvantaged, depending on whether the maximum likelihood lies within the range of the narrower prior or not. For the NGS16 model, for example, this is not the case for the spatially flat ΛCDM model, and the timescape model is unfairly advantaged. For model VIII the situation is reversed.

APPENDIX A: LUMINOSITY DISTANCES IN THE FLRW AND TIMESCAPE COSMOLOGIES

We compare SneIa observations to distance moduli (1) for theoretical luminosity distances determined from the FLRW and timescape models. Regardless of the matter content of the universe, the distance modulus for any general FLRW model can be expanded as a Taylor series of derivatives of the scale factor a(t) for small redshifts, z. This leads to a distance modulus (Visser 2004)

\begin{eqnarray} \mu _{ \rm FLRW}&=& 25+5\log _{ \rm 10}\left({\tt c}\,z\over H_{ 0}\,\hbox{Mpc}\right) +{5\over \ln 10}\Bigl \lbrace { {1\over 2}}(1-q_{ 0})z +{ {1\over 24}}\Bigl [9{q_{ 0}}^{2}-10q_{ 0}-7-4j_{ 0}-\Omega _{ k0}\Bigr ]z^2\nonumber\\ && +\,{ {1\over 24}}\Bigl [s_{ 0}+5-10{q_{ 0}}^3-16{q_{ 0}}^2-9q_{ 0}+j_{ 0}\left(8q_{ 0}+7\right) -\,4\Omega _{ k0}\left(q_{ 0}+1\right)\Bigr ]z^3 +\cdots \Bigr \rbrace , \end{eqnarray}

(A1)

where |${\tt c}$| is the speed of light, and |$H_{\scriptstyle 0}$|⁠, |$q_{\scriptstyle 0}$|⁠, |$j_{\scriptstyle 0}$|⁠, |$s_{\scriptstyle 0}$| and |$\Omega _{\scriptstyle k0}$| are the present values of the Hubble, deceleration, jerk, snap and spatial curvature parameters: H(t) ≡ a^{− 1}∂_ta; |$q(t)\equiv -a^{-1}H^{-2}\mathrm{\partial} _t^2 a$|⁠; |$j(t)\equiv a^{-1}H^{-3}\mathrm{\partial} _t^3 a$|⁠; |$s(t)\equiv a^{-1}H^{-4}\mathrm{\partial} _t^4 a$|⁠; |$\Omega _{\scriptstyle k}(t)\equiv -k{\tt c}^2(Ha)^{-2}$|⁠.

The luminosity distance–redshift relation in the ΛCDM model is given exactly by

\begin{eqnarray} d_{\scriptstyle \rm L}&=&\frac{(1+z){\tt c}}{H_{\scriptstyle 0}\sqrt{|\Omega _{\scriptstyle k0}|}}\,\hbox{sinn}\left(\sqrt{|\Omega _{\scriptstyle k0}|}\int \limits ^1 _{1/(1+z)}\frac{{\rm d}y}{\mathcal {H}(y)}\right), \nonumber \\ &&\mathcal {H}(y)\equiv \sqrt{\Omega _{\scriptstyle R0}+\Omega _{\scriptstyle M0}y+\Omega _{\scriptstyle k0}y^2+\Omega _{\scriptstyle \Lambda 0}y^4}, \nonumber \\ && \,\hbox{sinn}(x)\equiv \left\lbrace \begin{array}{@{}l@{\quad }l@{}}\sinh (x), &\Omega _{\scriptstyle k0}>0\\ x,&\Omega _{\scriptstyle k0}=0\\ \sin (x),&\Omega _{\scriptstyle k0}<0\\ \end{array}\right., \end{eqnarray}

(A2)

where |$\Omega _{\scriptstyle R0}$|⁠, |$\Omega _{\scriptstyle M0}$| and |$\Omega _{\scriptstyle \Lambda 0}$| are respectively the present epoch values of the radiation, non-relativistic matter and cosmological constant density parameters, which at all epochs obey the Friedmann equation sum rule |$\Omega _{\scriptstyle R}+\Omega _{\scriptstyle M}+\Omega _{\scriptstyle k}+\Omega _{\scriptstyle \Lambda }=1$|⁠. Since |$\Omega _{\scriptstyle R0}=4.15\times 10^{-5}h^{-2}$|⁠, the radiation term can be neglected on the distance scales of supernovae. Furthermore, for FLRW models |$\Omega _{\scriptstyle k0}$| is constrained to be close to zero by the angular diameter distance of the sound horizon in the CMB. Thus we will restrict attention to the spatially flat ΛCDM model, with two effective free parameters, |$H_{\scriptstyle 0}$| and |$\Omega _{\scriptstyle M0}\simeq 1-\Omega _{\scriptstyle \Lambda 0}$|⁠. We use equation (A2) with |$\Omega _{\scriptstyle k0}=0$|⁠, |$\Omega _{\scriptstyle R0}=0$|⁠, for computations but note that in the Taylor series (A1), |$q_{\scriptstyle 0}=-1+{\textstyle {3\over 2}}\Omega _{\scriptstyle M0}$|⁠, |$j_{\scriptstyle 0}=1$|⁠, |$s_{\scriptstyle 0}=1-{\textstyle {9\over 2}}\Omega _{\scriptstyle M0}$|⁠.

We also consider the FLRW model with linear expansion a(t) ∝ t. This solution is obtained by setting |$\Omega _{\scriptstyle k0}=1$|⁠, |$\Omega _{\scriptstyle R0}=\Omega _{\scriptstyle M0}=\Omega _{\scriptstyle \Lambda 0}=0$| in (A2) or |$\Omega _{\scriptstyle k0}=1$|⁠, |$q_{\scriptstyle 0}=j_{\scriptstyle 0}=s_{\scriptstyle 0}=\dots =0$| in (A1). Following NGS16 we denote this the empty universe, but note any matter content is admissible as long as the luminosity distance is exactly |$d_{\scriptstyle \rm L}={\tt c}\,z(1+{\textstyle {1\over 2}}z)/H_{\scriptstyle 0}$|⁠.

The timescape model (Wiltshire 2007a,b, 2009; Duley et al. 2013), does not evolve by the Friedmann equation, and its distance modulus does not yield a Taylor series that coincides with (A1) beyond the leading Hubble term. Instead observables are determined by conformally matching radial null geodesics of the regional ‘finite infinity’ geometry of observers in gravitationally bound systems to a statistical geometry determined by fitting a spherically symmetric metric to a solution (Wiltshire 2007b, 2009; Duley et al. 2013) of the Buchert equations (Buchert 2000, 2001).

For the purpose of supernova distance analysis, the radiation density parameter (though somewhat differently calibrated to the CMB; Duley et al. 2013) is negligible at the present epoch. To an accuracy of 0.3 per cent, the expansion history at late epochs is then given analytically (Wiltshire 2007b, 2009). The ‘dressed’ luminosity distance, |$d_{\scriptstyle \rm L}$|⁠, and angular diameter distance, |$d_{\scriptstyle \rm A}$|⁠, are given by

\begin{eqnarray} d_{\scriptstyle \rm L}=(1+z)^2d_{\scriptstyle \rm A}, \end{eqnarray}

(A3)

\begin{eqnarray} d_{\scriptstyle \rm A}={\tt c}\,t^{2/3}\int _t^{t_{\scriptscriptstyle 0}} {2\,{\rm d}t^{\prime }\over (2+{f_{\rm v}}(t^{\prime }))(t^{\prime })^{2/3}} ={\tt c}\,{t^{2/3}\left({\cal F}(t_{\scriptstyle 0})-{\cal F}(t)\right)}, \end{eqnarray}

(A4)

\begin{eqnarray} {\cal F}(t)\equiv 2t^{1/3}+{b^{1/3}\over 6}\ln \left((t^{1/3}+b^{1/3})^2\over t^{2/3}-b^{1/3}t^{1/3}+b^{2/3}\right) +{b^{1/3}\over \sqrt{3}}\tan ^{-1}\left(2t^{1/3}-b^{1/3}\over \sqrt{3}\,b^{1/3}\right), \end{eqnarray}

(A5)

where the volume-average time parameter, t, is defined implicitly in terms of the redshift by

\begin{equation} z+1={(2+{f_{\rm v}}){f_{\rm v}}^{1/3}\over 3f_{\rm v0}^{1/3}{\vphantom{\tilde{H}}\bar{H}}_{\scriptstyle \!0}t} ={2^{4/3}t^{1/3}(t+b)\over f_{\rm v0}^{1/3}{\vphantom{\tilde{H}}\bar{H}}_{\scriptstyle \!0}t(2t+3b)^{4/3}}, \end{equation}

(A6)

|$b\equiv 2(1-f_{\rm v0})(2+f_{\rm v0})/(9f_{\rm v0}{\vphantom{\tilde{H}}\bar{H}}_{\scriptstyle \!0})$|⁠, f_v0 is the present epoch value of the void volume fraction,

\begin{equation} {f_{\rm v}}(t)={3f_{\rm v0}{\vphantom{\tilde{H}}\bar{H}}_{\scriptstyle \!0}t\over 3f_{\rm v0}{\vphantom{\tilde{H}}\bar{H}}_{\scriptstyle \!0}t+(1-f_{\rm v0})(2+f_{\rm v0})}, \end{equation}

(A7)

and |${\vphantom{\tilde{H}}\bar{H}}_{\scriptstyle \!0}$| is the ‘bare Hubble constant’ related to the observed Hubble constant by |${\vphantom{\tilde{H}}\bar{H}}_{\scriptstyle \!0}=2(2+f_{\rm v0})H_{\scriptstyle 0}/(4{f_{\rm v0}}^2+f_{\rm v0}+4)$|⁠. The parameter t is related to the time parameter, τ, measured by typical observers in bound structures by

\begin{equation} \tau ={ {2\over 3}} t+{2(1-f_{\rm v0})(2+f_{\rm v0})\over 27f_{\rm v0}{\vphantom{\tilde{H}}\bar{H}}_{\scriptstyle \!0}}\ln \left(1+{9f_{\rm v0}{\vphantom{\tilde{H}}\bar{H}}_{\scriptstyle \!0}t \over 2(1-f_{\rm v0})(2+f_{\rm v0})}\right). \end{equation}

(A8)

The effective dressed scale factor a(τ) is given by

\begin{equation} a\equiv {\vphantom{\tilde{\gamma }}\bar{\gamma }}^{-1}{\vphantom{\tilde{a}}\bar{a}}, \end{equation}

(A9)

where |${\vphantom{\tilde{a}}\bar{a}}$| is the bare or volume-average scale factor and |${\vphantom{\tilde{\gamma }}\bar{\gamma }}$| is the phenomenological lapse function. These have simple analytic forms in terms of the volume-average time parameter, t, namely

\begin{equation} {\vphantom{\tilde{a}}\bar{a}}={{\vphantom{\tilde{a}}\bar{a}}_{\scriptstyle 0}\bigl (3{\vphantom{\tilde{H}}\bar{H}}_{\scriptstyle \!0}t\bigr )^{2/3}\over 2+f_{\rm v0}}\left[3f_{\rm v0}{\vphantom{\tilde{H}}\bar{H}}_{\scriptstyle \!0}t+ (1-f_{\rm v0})(2+f_{\rm v0})\right]^{1/3} \end{equation}

(A10)

and

\begin{equation} {\vphantom{\tilde{\gamma }}\bar{\gamma }}={ {1\over 2}}(2+{f_{\rm v}})={3(t+b)\over (2t+3b)}, \end{equation}

(A11)

respectively (Wiltshire 2007b, 2009). The bare Hubble parameter, |${\vphantom{\tilde{H}}\bar{H}}\equiv \mathrm{\partial} _t{\vphantom{\tilde{a}}\bar{a}}/{\vphantom{\tilde{a}}\bar{a}}$|⁠, and dressed Hubble parameter, H ≡ ∂_τa/a, are given respectively by |${\vphantom{\tilde{H}}\bar{H}}=(2+{f_{\rm v}})/(3t)$| and |$H=\left(4{f_{\rm v}}^2+{f_{\rm v}}+4\right){\vphantom{\tilde{H}}\bar{H}}/[2(2+{f_{\rm v}})]$|⁠. The bare deceleration parameter, |${\vphantom{\tilde{q}}\bar{q}}\equiv -{\vphantom{\tilde{a}}\bar{a}}^{-1}{\vphantom{\tilde{H}}\bar{H}}{}^{-2}\mathrm{\partial} _t^2{\vphantom{\tilde{a}}\bar{a}}$|⁠, is always positive. However, on account of the different time parameters (A8) the dressed deceleration parameter inferred by observers in bound systems, |$q\equiv -a^{-1}H^{-2}\mathrm{\partial} _\tau ^2 a$|⁠, may change sign from positive to negative, indicating apparent acceleration. The dressed deceleration parameter is given by

\begin{equation} q={-\left(1-{f_{\rm v}}\right)(8{f_{\rm v}}^3+39{f_{\rm v}}^2-12{f_{\rm v}}-8)\over \left(4+{f_{\rm v}}+4{f_{\rm v}}^2\right)^2}. \end{equation}

(A12)

The onset of apparent acceleration is determined by a root of the cubic in f_v in the numerator of (A12), and begins when f_v ≃ 0.587.

One may substitute (A3)–(A5) in (1) and then invert (A6) as a series in z using an algebraic computing package to arrive at a low-redshift Taylor series for the distance modulus, |$\mu _{\scriptstyle \rm TS}$|⁠, equivalent to (A1) for the FLRW model. The first terms are given below, along with equivalent expressions for the spatially flat ΛCDM and empty universe models as determined from (A1):

\begin{eqnarray} &&{ \mu _{ \rm TS}=\mu _{ 0}(z)+{5\over \ln 10}\left\lbrace \left[\frac{24\,f_{\rm v0}^{4}-23\,f_{\rm v0}^{3}+99\,f_{\rm v0}^{ 2}+8}{2\left( 4\,f_{\rm v0}^{2}+f_{\rm v0}+4 \right) ^{2}}\right]z\right.},\nonumber \\ &&-\left.\left[ {\frac{1984\,f_{\rm v0}^{8}-4352\,f_{\rm v0}^{7}+16515f_{\rm v0}^{6}+14770\,f_{\rm v0}^{5}+7819\,f_{\rm v0}^{4}-11328f_{\rm v0}^{3}+32080\,f_{\rm v0}^{2}-128\,f_{\rm v0}+960}{24\left( 4f_{\rm v0}^{2}+f_{\rm v0}+4 \right) ^{4}}}\right]z^2 +\dots \right\rbrace , \end{eqnarray}

(A13)

\begin{eqnarray} \mu _{ \Lambda {\rm CDM}} =\mu _{ 0}(z)+{ {5\over \ln 10}}\left\lbrace \left(1-{ {3\over 4}}\Omega _{ M0}\right)z -\left[{ {1\over 2}}+{ {1\over 2}}\Omega _{ M0}-{ {27\over 32}}\Omega _{ M0}^2\right]z^2 +\left[{ {1\over 3}}-{ {1\over 8}}\Omega _{ M0}+{ {21\over 16}}\Omega _{ M0}^2-{ {45\over 32}}\Omega _{ M0}^3\right]z^3 +\cdots \right\rbrace , \end{eqnarray}

(A14)

\begin{eqnarray} \mu _{ \rm empty}=\mu _{ 0}(z)+{ {5\over \ln 10}}\displaystyle\left\lbrace { {1\over 2}}z-{ {1\over 8}}z^2+{ {1\over 24}}z^3 +\cdots \right\rbrace . \end{eqnarray}

(A15)

Here, the term |$\mu _{\scriptstyle 0}(z)\equiv 25+5\log _{\scriptstyle \rm 10}[{\tt c}z/(H_{\scriptstyle 0}\,\hbox{Mpc})] =25+5\log _{\scriptstyle \rm 10}(2997.9\,h^{-1})+5\log _{\scriptstyle \rm 10}z$| is common to all models, the Hubble constant being |$H_{\scriptstyle 0}=100\,h\,\hbox{km}\,\,\hbox{s}^{-1}\,\hbox{Mpc}^{-1}$|⁠.

APPENDIX B: IMPLEMENTATION OF THE SALT2 METHOD

The SALT relation (2) refers to the actual emitter (em) and observer (obs), but the luminosity distance relations (A2) and (A4) refer to ideal observers who determine an isotropic distance–redshift relation. Consequently, the theoretical relations (A2) and (A4) must be transformed to the frame involving the actually measured redshift |$\hat{z}=(\lambda _{\rm obs}-\lambda _{\rm em})/\lambda _{\rm em}$| before the SALT relation is applied. The luminosity distance entering (1) is then

\begin{equation} \hat{d_{\scriptstyle \rm L}}(\hat{z})=\frac{1 + \hat{z}}{1 + z} d_{\scriptstyle \rm L}(z)=(1+\hat{z})D(z), \end{equation}

(B1)

where for each cosmological model, |$D(z)=d_{\scriptstyle \rm L}/(1+z)=(1+z)d_{\scriptstyle \rm A}$| is the (effective) comoving distance, and

\begin{equation} 1+\hat{z}=(1 + z)(1+z^{\rm pec}_{\rm obs})(1+z^\phi _{\rm obs}) (1+z^{\rm pec}_{\rm em})(1+z^\phi _{\rm em}) \end{equation}

(B2)

gives the measured redshift, |$\hat{z}$|⁠, in terms of the cosmological redshift, z, the local Doppler redshifts of observer, |$z^{\rm pec}_{\rm obs}$| and emitter, |$z^{\rm pec}_{\rm em}$|⁠, and gravitational redshifts at the two locations, |$z^\phi _{\rm obs}$| and |$z^\phi _{\rm em}$|⁠.

For our observations, |$\hat{z}$|⁠, is the heliocentric redshift as the Earth's annual motion is averaged to the rest frame of the Sun. In the standard cosmology gravitational potential effects are assumed to be small, and the only relevant terms on the right-hand side of (B2) are assumed to be local boosts of order |$v/c\mathop {\sim }\limits 10^{-3}$|⁠. This leads to 0.1 per cent corrections to the luminosity distance which are often neglected. However, as noted by Calcino & Davis (2017) differences of 0.1 per cent in |$d_{\scriptstyle \rm L}$| lead to order 1 per cent corrections to cosmological parameters, which we have confirmed in our analysis.

In the timescape model, as in any inhomogeneous cosmology, expansion below the |$\mathop {\sim }\limits 100\,h^{-1}\hbox{Mpc}$| SHS will generally differ from that of a global average geometry plus local boosts. Equivalently, very slowly varying time-dependent gravitational potentials also make a contribution to (B2). Such terms encode non-kinematic differential expansion (Bolejko, Nazer & Wiltshire 2016) from inhomogeneities below the SHS. Spatial variations in the term |$z^\phi _{\rm em}$| may have significant consequences for interpreting the local ‘peculiar velocity field’ for sources within ≲ 100 h⁻¹Mpc of our location (Wiltshire et al. 2013; McKay & Wiltshire 2016) but any net anisotropy on SneIa redshifts on larger scales should only make a small correction to the standard boost between the heliocentric and CMB frames. Indeed, it could be a source for a small systematic redshift uncertainty of the type considered by Calcino & Davis (2017). However, we do not investigate that possibility in this paper as the RH16 empirical light-curve models we study are already very complex. The peculiar velocity and gravitational potential terms in (B2) that we are unable to determine will be assumed to contribute to statistical uncertainties in measured redshifts only.

We therefore compute cosmological luminosity distances in the CMB rest frame, exclude data below the SHS and study the effect of different choices for this cutoff.¹¹ Furthermore, we apply the SALT2 relation in the heliocentric frame using the values tabulated in the JLA catalogue (Betoule et al. 2014), and calculate the corresponding CMB rest-frame redshifts ourselves.¹²

We do not use the JLA tabulated CMB frame redshifts (Betoule et al. 2014) since in addition to our |$z^{\rm pec}_{\rm obs}$| correction, these values also already include a correction, |$z^{\rm pec}_{\rm em}$|⁠, for the peculiar velocity field (Hudson et al. 2004; Conley et al. 2011) of galaxies up to z = 0.071, implicitly assuming the FLRW model.

APPENDIX C: MODEL COMPARISON

C1 Frequentist approach

We are interested in the dependence of the likelihood (7) on the model parameters, Θ. We write |$\mathcal {L}(\Theta ) \equiv \mathcal {L}( \mathrm{Data} | \Theta , \mathrm{M})$|⁠. We are interested in a subset of parameter-space |$\Theta _{\scriptstyle 1}\subset \Theta$|⁠, for which we construct a profile likelihood |$\mathcal {L}_p (\Theta _{\scriptstyle 1}) \equiv \max _{\Theta _2}[\mathcal {L} (\Theta )]$|⁠, where maximization is over the nuisance parameters |$\Theta _2 = \Theta \setminus \Theta _{\scriptstyle 1}$|⁠. In our case, |$\Theta _{\scriptstyle 1}$| is usually the free parameter(s) of the cosmological model, and Θ₂ the intrinsic supernovae parameters and the empirical parameters, α, β, of the light-curve fitter.

Confidence regions for the parameters of interest are estimated from the coverage probability p_cov of a region in the k-dimensional slice of parameter space, |$k\equiv \mathop {\rm dim}\Theta _{\scriptstyle 1}$|⁠, given asymptotically by the integral

\begin{equation} p_{\rm cov} (\mathrm{region}) = \int _{0}^{-2 \ln \left(\mathcal {L}_p(\Theta _{\scriptscriptstyle 1}) /\mathcal {L}_{\rm max}\right)} f_{\chi ^2} (\Theta _{\scriptstyle 1}^{\prime }, k){\rm d}\Theta _{\scriptstyle 1}^{\prime }, \end{equation}

(C1)

where |$\mathcal {L}_{\rm max}\equiv \max _{\Theta _{\scriptscriptstyle 1}}[ \mathcal {L}_p (\Theta _{\scriptstyle 1})]$|⁠, and |$f_{\chi ^2}(x,k)$| is the probability density function of a χ² distributed variable with k degrees of freedom. Having constructed confidence intervals from (C1), one can compare nested models.

Since we wish to compare independent non-nested models,¹³ we need to go beyond the procedure of (C1). The Akaike Information Criterion (AIC; Akaike 1974) and Bayesian Information Criterion (BIC; Schwarz 1978) are two widely used measures of the relative information loss for non-nested models, given respectively by

\begin{eqnarray} {\rm AIC}=2k-2\ln \left(\mathcal {L}_{\rm max}\right), \end{eqnarray}

(C2)

\begin{eqnarray} {\rm BIC}=k\ln N-2\ln \left(\mathcal {L}_{\rm max}\right), \end{eqnarray}

(C3)

where k is the number of independent model parameters, and N the number of data points fit. The AIC estimate of relative probability of minimal information loss for two models is |$P_{AIC}\equiv \exp \left[-{\textstyle {1\over 2}} \left({\rm AIC}_{\scriptstyle 1}-{\rm AIC}_{\scriptstyle 2}\right)\right]$|⁠, and similarly for BIC. The BIC gives a greater penalty for introducing additional parameters than AIC if N ≥ 8. Differences of at least 2, 6 and 10 are considered to be respectively ‘positive’, ‘strong’ and ‘very strong’ evidence (Kass & Raftery 1995) for the model with the lower IC value. Both tests reduce to a maximum likelihood comparison when k₁ = k₂, as is the case for the timescape and spatially flat ΛCDM models.

C2 Bayesian approach

The frequentist methods place emphasis on the maximum likelihood, which is of limited use. We therefore perform a fully Bayesian analysis of the JLA data set to determine the relative statistical support for each cosmological model, as well as for the introduction of additional redshift-dependent light-curve parameters.

Given data, D, and a model, M, determined by a set of n parameters Θ = (θ₁, θ₂, …, θ_n), by Bayes theorem the posterior probability distribution, p(Θ|D, M), is given by

\begin{equation} p(\Theta |D,\mathrm{M}) = \frac{\mathcal {L}(\Theta ) \pi (\Theta |\mathrm{M})}{p(D|\mathrm{M})}, \end{equation}

(C4)

where |$\mathcal {L}(\Theta )\equiv p(D|\Theta ,\mathrm{M})$| is the likelihood, π(Θ|M) is the prior distribution and p(D|M) is the Bayesian evidence. The prior represents a subjective initial state of belief in the parameters based on external information or previous experiments, while the evidence is a normalization constant,

\begin{equation} E\equiv p(D|\mathrm{M})=\int \mathcal {L}(\Theta ) \pi (\Theta |\mathrm{M})\,{\rm d}\Theta , \end{equation}

(C5)

to ensure the posterior is a probability distribution. It is independent of parameters and as such does not play a role in parameter estimation, but becomes important for model comparison.

Given two models, M₁ and M₂, for the same data D, the Bayes factor (Kass & Raftery 1995)

\begin{equation} B\equiv \frac{E_1}{E_2}=\frac{p(D|\mathrm{M}_1)}{p(D|\mathrm{M}_2)}, \end{equation}

(C6)

gives a measure for which model is more probable in view of the data. The Bayes factor implicitly applies the principle of Occam's razor¹⁴ with a penalty for adding extra parameters. This makes model selection natural in the Bayesian framework. Values of B > 1 indicate preference for model 1, B < 1 for model 2. On a standard scale, evidence with |ln B| < 1 is ‘not worth more than a bare mention’ (Kass & Raftery 1995) or ‘inconclusive’ (Trotta 2007), while 1 ≤ |ln B| < 3, 3 ≤ |ln B| < 5 and |ln B| ≥ 5 indicate ‘positive’, ‘strong’ and ‘very strong’ pieces of evidence, respectively (Kass & Raftery 1995).

In the Bayesian approach, the nuisance parameters are marginalized over, i.e. integrated out from the posterior p(Θ|D, M). For example, the marginal posterior of θ₁ is obtained from the n-dimensional posterior by

\begin{equation} p(\theta _1|D,\mathrm{M}) = \int p(\theta |D,\mathrm{M})\,{\rm d}\theta _2\,{\rm d}\theta _3\ldots {\rm d}\theta _n, \end{equation}

(C7)

and from this 1-dimensional distribution parameter inferences can be made. The posterior mean value is given by

\begin{equation} \bar{\theta }_1 = \int \theta _1\, p(\theta _1|D,\mathrm{M})\,{\rm d}\theta _1, \end{equation}

(C8)

and more generally

\begin{equation} \bar{f} = \int f(\theta _1)\, p(\theta _1|D,\mathrm{M})\,{\rm d}\theta _1, \end{equation}

(C9)

for some parameter-dependent quantity f. Credible intervals, or uncertainties in parameters, can also be obtained from the marginal posterior. For example, a 68 per cent equal-tailed credible interval is defined in such a way that (1 − 0.68)/2 = 0.16 of the probability lies on either side of the interval.

In cases where the Bayes factor is close to unity giving no clear preference for either model, the Bayesian complexity (Spiegelhalter et al. 2002) can provide a secondary measure to the model selection process. It is defined as

\begin{equation} C_b \equiv -2\left(D_{\scriptstyle \rm KL}(p,\pi )-\widehat{D_{\scriptstyle \rm KL}}\right), \end{equation}

(C10)

where

\begin{equation} D_{\scriptstyle \rm KL}(p,\pi )\equiv \int p(\Theta |D,\mathrm{M})\ln \left[\frac{p(\Theta |D,\mathrm{M})}{\pi (\Theta |\mathrm{M})}\right]{\rm d}\Theta , \end{equation}

(C11)

is the Kullback–Leibler divergence measuring the information gain of the inference, and |$\widehat{D_{\scriptstyle \rm KL}}$| is a point estimator evaluated at the posterior mean |$\bar{\Theta }$| measuring the expected information gain:

\begin{equation} \widehat{D_{\scriptstyle \rm KL}}\equiv p(\bar{\Theta }|D,\mathrm{M})\ln \left[\frac{p(\bar{\Theta }|D,\mathrm{M})}{\pi (\bar{\Theta }|\mathrm{M})}\right] =\ln \mathcal {L}(\bar{\Theta })-\ln p(D|\mathrm{M}), \end{equation}

(C12)

where we have used Bayes theorem in the second equality. As the data may not be able to constrain all parameters, the Bayesian complexity determines the effective number of parameters supported by the data. Thus for models with |ln B| < 1, the model with the lower C_b indicates the simpler model and is therefore preferred. By defining the effective chi-squared|$\chi ^2(\Theta )\equiv -2\ln \mathcal {L}$| and invoking Bayes theorem (C4), we can rewrite (C10) as

\begin{equation} C_{\rm b}=\overline{\chi ^2(\Theta )}-\chi ^2(\bar{\Theta }), \end{equation}

(C13)

with |$\overline{\chi ^2}$| being the posterior mean of χ².

APPENDIX D: COSMOLOGICAL MODEL PRIORS

We construct priors for the timescape model (Wiltshire 2007a,b, 2009; Duley et al. 2013) based on CMB and BAO observations, to the best of our knowledge. We will also construct equivalent priors for the ΛCDM model based on the same assumptions. The resulting priors are wider than in conventionally assumed, but do not unfairly weight a Bayesian comparison by integrating the ΛCDM model likelihood function over a narrow parameter range as compared to the timescape case.¹⁵

D1 CMB acoustic scale constraint

In the case of the CMB, a cosmology-independent analysis of the angular scale and heights of the acoustic peaks has been undertaken by Aghamousa et al. (2015) from the Planck data. We use the information resulting from the angular scale of the acoustic peaks alone. The angular scale depends on the angular diameter distance of the sound horizon alone, which is constrained in the timescape model. By contrast, the relative peak heights are related to the baryon-to-photon ratio, |$\eta _{\scriptstyle B\gamma }$|⁠, and the spectral index, n_s, which are parameters with the largest systematic uncertainties in the timescape case.

A non-parametric fit of the acoustic scale alone gives |$286\le \ell _{\scriptstyle \rm A} \le 305$| at 95 per cent confidence (Aghamousa et al. 2015). Our CMB prior is then determined by demanding that the angular diameter distance of the sound horizon at decoupling matches the corresponding angular scale |$\theta _{\scriptstyle \rm A}=\pi / \ell _{\scriptstyle \rm A}$|⁠; i.e. |$0.01030\le \theta _{\scriptstyle \rm A}\le 0.01098$|⁠. In earlier work (Leith et al. 2008; Smale & Wiltshire 2011; Duley et al. 2013), given that non-parametric fits had not been performed, we had demanded a match to the FLRW parametric estimate of the acoustic scale |$\theta _{\scriptstyle \rm A}=0.01041$| to within 2 per cent, 4 per cent or 6 per cent. The non-parametric fit represents a considerable improvement, particularly since the FLRW model value is not in the mid-range of the non-parametric 95 per cent confidence interval.

To constrain the angular diameter distance of the sound horizon |$d_{\rm A\,dec} ={\bar{D}}_{\rm s}(z_{\rm dec})/\theta _{\scriptstyle \rm A}$| in the timescape model, we determine the redshift of decoupling, z_dec, and the comoving distance of the sound horizon |$\bar{D}_{\rm s}$| at that epoch (Duley et al. 2013; Nazer & Wiltshire 2015), which require the baryon-to-photon ratio to be specified. In the FLRW model, this ratio is very tightly constrained by the ratio of CMB peak heights, as first measured by WMAP (Bennett et al. 2003). However, to achieve a similarly precise constraint in the timescape model we would need to include backreaction in the primordial plasma (Nazer & Wiltshire 2015), which is beyond the scope of current investigations. In previous work (Leith et al. 2008; Smale & Wiltshire 2011; Duley et al. 2013), we used a range of pre-WMAP baryon-to-photon ratios (Tytler et al. 2000), |$4.6<10^{10}\eta _{\scriptstyle B\gamma }<5.6$|⁠, for which all light element abundance measurements are within 2σ, i.e. with no primordial lithium abundance anomaly. In the present analysis, we wish to use the same priors on |$\eta _{\scriptstyle B\gamma }$| for both the timescape and ΛCDM models, and thus need to include the standard model value |$\eta _{\scriptstyle B\gamma }=2.736\times 10^{-8}\Omega _{\scriptstyle M0}h^2=(6.08\pm 0.07)\times 10^{-10}$| for which the primordial lithium abundance is problematic. We therefore adopt the more conservative pre-WMAP range given by Olive, Steigman & Walker (2000), namely |$4.2<10^{10}\eta _{\scriptstyle B\gamma }<6.3$|⁠.

D2 Baryon Acoustic Oscillation constraints

Determinations of the BAO scale from galaxy clustering at low redshifts and Lyman α forest statistics at z = 2.34 provide complementary constraints on the expansion history. In previous work (Leith et al. 2008; Smale & Wiltshire 2011; Duley et al. 2013), we simply demanded that the timescape effective comoving BAO scale match a single estimate determined from the FLRW cosmology to within ±2 per cent, ±4 per cent or ±6 per cent, which was a crude method but the best available given the earlier precision of measurements. The number and precision of measurements has now improved.

For the present investigation, we have considered estimates of the BAO scale at different redshifts (Aubourg et al. 2015; Delubac et al. 2015; Alam et al. 2017) using the best available data from the BOSS survey. Unfortunately the standard FLRW cosmology plays an implicit role in the data reduction, and limits the extent to which bounds can be placed on non-FLRW models. The systematic issues can be most directly understood by noting that the BAO scale is determined separately in the angular and radial directions, by converting angular separations and redshift separations in the galaxy–galaxy correlation function into the displacements

\begin{equation} \alpha _\perp =\frac{\left[d_{\scriptstyle \rm A}(z)/r_{\rm d}\right]}{\left[d_{\scriptstyle \rm A}(z)/r_{\rm d}\right]_{\rm fid}} \qquad \hbox{and}\qquad \alpha _\parallel =\frac{\left[d_{\scriptstyle \rm H}(z)/r_{\rm d}\right]}{\left[d_{\scriptstyle \rm H}(z)/r_{\rm d}\right]_{\rm fid}} \end{equation}

(D1)

where r_d is the present comoving scale of the sound horizon at the baryon drag epoch, d_H(z) ≡ c/H(z), and the subscript ‘fid’ refers to quantities computed in a fiducial FLRW model that is used to convert the raw angular and redshift displacements into 3-dimensional comoving space. (Here, we neglect the effect of redshift–space distortions which are also often modelled with N-body Newtonian simulations based on the ΛCDM model.)

The conversion to 3-dimensional comoving space can be problematic for a non-FLRW model. While use of purely angular results should pose no problems for the timescape cosmology, the conversion of redshift increments to a radial comoving distance involves different assumptions about spatial curvature in the FLRW and timescape models. One could in principle use the values determined by a fiducial ΛCDM model (Delubac et al. 2015; Aubourg et al. 2015; Alam et al. 2017) to recompute the radial comoving distance except for an additional problem: in particular redshift ranges the relative Alcock–Paczyński factor |$\left[\alpha _\perp /\alpha _\parallel \right]_{\scriptstyle \rm \Lambda CDM}/ \left[\alpha _\perp /\alpha _\parallel \right]_{\scriptstyle \rm TS}=\left[H(z)d_{\scriptstyle \rm A}(z)\right]_{\scriptstyle \rm \Lambda CDM}/ \left[H(z)d_{\scriptstyle \rm A}(z)\right]_{\scriptstyle \rm TS}$| between a fiducial ΛCDM model and the timescape model varies over the redshift slices |$\Delta z\mathop {\sim }\limits 0.2$| used in the BOSS survey (Alam et al. 2017) by an amount similar in magnitude to the uncertainty. Consequently, to have any confidence in radial measurements, one really needs to recompute the radial BAO scale from the raw data assuming a fiducial timescape model. That is beyond the scope of this paper.

For the present analysis, we will consequently restrict constraints on the BAO scale to 2σ bounds obtained from the angular estimates of BOSS data (Alam et al. 2017) at low redshifts 0.38 ≲ z ≲ 0.61 and in the Lyman α forest (Delubac et al. 2015) at z = 2.34. In the former case, the radial and angular measurements are actually also somewhat correlated. Consequently, and also in view of the fact that the measurements at low and high redshifts are in tension in the ΛCDM model, we will take bounds that result from the union of the constraints at low and high redshifts, rather than their intersection. In practice, the bounds are mostly set by the Lyman α measurement since it has a much larger uncertainty.

D3 Joint constraints

The joint 2σ confidence regions from applying the CMB constraint and either the low-z or z = 2.34 BAO constraint is shown in Fig. D1. Since a range of possible baryon–to–photon ratios are admitted, with no information from the relative heights of the acoustic peaks used in either model, the width of the allowed regions is larger than in conventional analyses for ΛCDM.

Figure D1.

Cosmological parameter constraint priors from on the timescape model (left-hand panel) and the spatially flat ΛCDM model (right-hand panel). Independent 2σ constraints are determined for: (i) the angular scale of sound horizon in the CMB (contours from top left to bottom right); (ii) the fit of the angular BAO scale from BOSS galaxies at low redshift (contours from bottom left to top right); (iii) the angular BAO scale from one measurement of the Lyman α forest at z = 2.34 (wide contours). A range of possible baryon-to-photon ratios are allowed, with the extremes indicated. The joint confidence region is determined by applying the CMB constraint and allowing one or other BAO constraint.

For timescape the confidence regions are f_v0 ∈ (0.588, 0.765) at 1σ, f_v0 ∈ (0.500, 0.799) at 2σ, f_v0 ∈ (0.378, 0.826) at 3σ. The corresponding effective dressed |$\Omega _{\scriptstyle M0}={\textstyle {1\over 2}}(1-f_{\rm v0})(2+f_{\rm v0})$| is |$\Omega _{\scriptstyle M0}\in (0.325,0.534)$| at 1σ, |$\Omega _{\scriptstyle M0}\in (0.281,0.625)$| at 2σ and |$\Omega _{\scriptstyle M0}\in (0.245,0.740)$| at 3σ. For spatially flat ΛCDM, the corresponding confidence regions are |$\Omega _{\scriptstyle M0}\in (0.162,0.392)$| at 1σ, |$\Omega _{\scriptstyle M0}\in (0.143,0.487)$| at 2σ and |$\Omega _{\scriptstyle M0}\in (0.124,0.665)$| at 3σ. We adopt the 2σ bounds as priors.