Percent-Level Test of Isotropic Expansion Using Type Ia Supernovae

TL;DR

This paper introduces a new non-parametric method to test the isotropy of the universe's expansion using Type Ia supernovae, setting a stringent upper limit on spatial variations in the Hubble parameter.

Contribution

It presents a novel, robust statistical test for isotropy based on supernova data, improving constraints on spatial variations in cosmic expansion.

Findings

Less than 1% rms spatial variation in Hubble parameter at 99.7% confidence

The method effectively accounts for systematic effects

Provides a new upper limit on isotropy deviations

Abstract

We propose and implement a novel, robust, and non-parametric test of statistical isotropy of the expansion of the universe, and apply it to around one thousand type Ia supernovae from the Pantheon sample. We calculate the angular clustering of supernova magnitude residuals and compare it to the noise expected under the isotropic assumption. We also test for systematic effects and demonstrate that their effects are negligible or are already accounted for in our procedure. We express our constraints as an upper limit on the rms spatial variation in the Hubble parameter at late times. For the sky smoothed with a Gaussian with fwhm=60 deg, less than 1% rms spatial variation in the Hubble parameter is allowed at 99.7% confidence.

Equations

$$r_i\equiv \frac{m_i-m_i^{\mathrm{th}}}{{\sigma }_i}$$

$$C_{\ell }=\frac{1}{2\ell +1}\sum _{m=-\ell }^{\ell }|a_{\ell m}{|}^2,$$

$$a_{\ell m}\equiv \frac{\int r(\hat{\mathbf{n}})W(\hat{\mathbf{n}})Y_{\ell m}^{\ast }(\hat{\mathbf{n}})d^2\hat{\mathbf{n}}}{\int W(\hat{\mathbf{n}})d^2\hat{\mathbf{n}}/(4\pi )}$$

$${\chi }^2=\left(C_{\ell }-{\bar{C}}_{\ell }^{\mathrm{boot}}\right)(M^{-1}{)}_{\ell {\ell }^{\prime }}\left(C_{{\ell }^{\prime }}-{\bar{C}}_{{\ell }^{\prime }}^{\mathrm{boot}}\right),$$

$${\delta }_H(\hat{\mathbf{n}})\equiv \frac{\delta H}{H}(\hat{\mathbf{n}})\ll 1.$$

$$⟨(r^{\mathrm{data}}{)}^2⟩={\left(\frac{5}{\ln 10}\right)}^2⟨{\left(\frac{{\delta }_H}{\sigma }\right)}^2⟩+N(r^{\mathrm{iso}})$$

$$\mathrm{Var}(r^{\mathrm{data}})\equiv ⟨(r^{\mathrm{data}}{)}^2⟩=\sum _{\ell =1}^{{\ell }_{\max }}\frac{2\ell +1}{4\pi }{C}_{\ell }$$

$$({\delta }_H{)}_{\mathrm{rms}}\simeq \left(\frac{\ln 10}{5}\right)⟨{\sigma }^{-2}{⟩}^{-1/2}\sqrt{\mathrm{Var}(r^{\mathrm{data}})}$$

$${\left[({\delta }_H{)}_{\mathrm{rms}}\right]}^{\mathrm{FWHM}=60^{\circ }}<0.009\text{mbox}(99.7\%CL),$$

Full text

Introduction. The simplest inflationary-cosmology scenarios Guth (1981); Linde (1982); Mukhanov and Chibisov (1981); Albrecht and Steinhardt (1982) generically predict that the expansion of the universe is isotropic. However, violations of statistical isotropy can certainly be accommodated in models with additional complexity (e.g. Ackerman et al. (2007); Pullen and Kamionkowski (2007); Watanabe et al. (2009)) and even perturbative effects on the expansion rate in the standard cosmological model (which are, however, expected to be small, e.g. Bonvin et al. (2006a, b)). While tests of statistical isotropy of the early universe have typically been carried out by analyzing the cosmic microwave background anisotropy maps (e.g. Hanson and Lewis (2009); Bennett et al. (2011); Copi et al. (2010); Kim and Komatsu (2013); Ade et al. (2016)), it is well worthwhile to investigate the isotropy of the late-time universe. The latter is particularly interesting given the lack of a fundamental understanding of the physical nature of dark energy that powers the accelerated expansion of the universe.

In this paper we present a novel test of the isotropy of cosmic expansion and apply it to current type Ia supernovae (SNIa) data. While an investigation of the isotropy of the universe using SNIa data has been carried out by numerous previous works Antoniou and Perivolaropoulos (2010); Colin et al. (2011); Campanelli et al. (2011); Cai and Tuo (2012); Mariano and Perivolaropoulos (2012); Appleby and Linder (2013); Kalus et al. (2013); Cai et al. (2013); Zhao et al. (2013); Beltran Jimenez et al. (2015); Appleby et al. (2015); Bengaly et al. (2015, 2017); Lin et al. (2016); Javanmardi et al. (2015); Huterer et al. (2017); Wang and Wang (2018); Andrade et al. (2018a); Sun and Wang (2018); Andrade et al. (2018b); Deng and Wei (2018), our methodology (described below) extends these efforts. Our approach is parameter-free, robust, and explicitly independent of assumptions about the distribution of the data. We now describe the data we use, our methodology, and present the results along with estimates of the effects of systematic errors.

Data. For our analysis, we use the “Pantheon” compilation of SNIa Scolnic et al. (2018). The Pantheon sample combines 279 SNIa $(0.03<z<0.68)$ from the Pan-STARRS1 Medium Deep Survey with SNIa from Sloan Digital Sky Survey (SDSS), SuperNova Legacy Survey (SNLS), and various low-z and Hubble Space Telescope samples to produce a SNIa sample of 1048 objects in the redshift range $0.01<z<2.26$. The Pantheon sample was produced using the PS1 Supercal process Scolnic et al. (2015), which determined a global calibration solution to combine 13 different SNIa samples. The latter analysis also corrects for expected biases in light-curve fit parameters and their errors using the method outlined in Ref. Kessler and Scolnic (2017).

Methodology. The SNIa data consists of individual magnitude measurements $m_{i}\equiv m(z_{i},\hat{\mathbf{n}})$, where $z_{i}$ is the redshift of a supernova in the cosmic microwave background frame and corrected for peculiar velocities Scolnic et al. (2018), and $\hat{\mathbf{n}}$ is its location on the sky. The individual SNIa magnitude errors $\sigma_{i}$ are generalized, in modern SNIa analyses, to the full covariance matrix, $\mathbf{C}=\mathbf{S}+\mathbf{N}$, where $\mathbf{S}$ and $\mathbf{N}$ are the signal and noise matrices, respectively. The noise matrix encodes statistical magnitude measurement errors and covariances due to unknown fit parameters that correlate the measurements, such as the color and stretch. The signal matrix is nonzero at low redshift because, roughly speaking, nearby SNIa are pulled by the same structures, resulting in correlated peculiar velocities. At high redshift the signal matrix is nonzero mainly because of the effects of lensing on SNIa magnitudes.

We work with SNIa magnitude residuals divided by individual statistical errors

[TABLE]

where $m^{\rm th}_{i}=5\log_{10}[H_{0}d_{L}(z_{i},\Omega_{M})]+\mathcal{M}$ is the theoretically expected magnitude for an object at redshift $z_{i}$ and a given value of the matter density relative to critical $\Omega_{M}$, in the flat $\Lambda$CDM universe which we assume. Here $\mathcal{M}$ is the nuisance parameter that combines the absolute magnitude of SNIa with the Hubble constant $H_{0}$. The distribution of the residuals on the sky is shown in Fig. 1; the inset in the Figure shows that the residuals are approximately Gaussian-distributed.

The goal of this study is to put an upper limit on the value of the “signal” in the distribution of the SNIa magnitude residuals caused by violations of statistical isotropy. Since we are searching for the excess signal with respect to statistical noise, we choose to normalize the magnitude residuals in Eq. (1) by the diagonal statistical error $\sigma_{i}$, and not elements of the full (signal plus noise) covariance $\mathbf{C}$. Note that the statistical measurement error constitutes the majority of the contribution to the diagonal of the noise covariance, $\sigma_{i}^{2}\simeq 0.99N_{ii}$, thus we are effectively dividing by the square roots of the latter. The total signal in the clustering of residuals also has guaranteed contributions from the peculiar velocities of the SNIa (which are correlated because the velocity’s origin is the gravitational pull of the nearby large-scale structures), and from the systematic uncertainty in the precise values of the cosmological parameters $\Omega_{M}$ and $\mathcal{M}$ which are required to calculate $m^{\rm th}(z)$. We demonstrate below that these additional contributions to the signal are very small compared to the noise level in the SNIa data and can be ignored.

To estimate the signal power spectrum, we pixelate the sky using the HEALPix Gorski et al. (2005) resolution parameter Nside. Our fiducial analysis is done at $\textsc{Nside}=16$, though we also compare with results at $\textsc{Nside}=$ 8, 32, and 64, finding good agreement. Each pixel has a side of roughly $60^{\circ}/\textsc{Nside}$. Because more than one SNIa may lie in a given pixel, we choose to take the average of the residuals in the pixel. Thus, the value of the $j^{\rm th}$ pixel is given by $p_{j}=(\sum_{r_{i}\in p_{j}}r_{i})/n_{j}$, where $n_{j}$ is the number of SNIa in that pixel and the sum goes over the residuals located in the pixel. We now outline how the data power spectrum is computed and then discuss our noise estimation approach.

To calculate the angular power spectrum of the map of residuals, $C_{\ell}$, we employ a pseudo-$C_{\ell}$ estimator which, given the small sky coverage of the SNIa (at $\textsc{Nside}=16$, the fractional sky coverage is only $f_{\rm sky}\simeq 0.07$), is much more practical than the maximum-likelihood estimators which try to recover the full-sky signal. To get the pseudo-$C_{\ell}$ we adopt the function anafast in HEALPix . Given the significant variation of pixel occupancy by SNIa, it is crucial to weight each pixel by the number of objects in it; this guarantees that large-angle (low-$\ell$) $C_{\ell}$ will not depend on the pixelation as long as the latter is finer than the scales we wish to probe. Our angular power spectrum is given by the usual pseudo-$C_{\ell}$ formula

[TABLE]

with the harmonic expansion of the residuals that applies weight to the pixels

[TABLE]

where the integral is typically discretized as the sum over the pixels whose centers are in directions ${\bf\hat{n}}$ and which have areas $d^{2}{\bf\hat{n}}$. Here $r({\bf\hat{n}})$ is the mean residual in a pixel in the direction ${\bf\hat{n}}$, while the weight $W({\bf\hat{n}})$ is given by the number of objects in that pixel. The denominator in Eq. (3) evaluates to $f_{\rm sky}\langle W_{\rm pix}\rangle$, where $\langle W_{\rm pix}\rangle$ is the average number of SNIa per pixel. Note that the overall normalization of the $C_{\ell}$ is not important for comparing the angular spectrum of our sky to that of the bootstrapped sample since the two have the same normalization, but it is important when we quote limits on the statistical isotropy of the expansion rate.

Having obtained the angular power spectrum of the SNIa, we then produce the distribution of angular power spectra that would be expected in an isotropic universe, which in turn is given by the clustering noise and no signal. Because the number of SNIa varies significantly from pixel to pixel, the computation of noise in the angular clustering of the SNIa is analytically intractable, and would be so even under the simplified assumption of Gaussian noise. To get a reliable upper bound on the background anisotropy it is crucial to estimate the clustering noise, as well as its uncertainty, directly from the data. To address this, we employ a simple, non-parametric bootstrap approach whose principal advantage is that the noise level can be directly estimated from data. Additionally, this approach does not make assumptions about the statistical distribution of residuals, e.g. whether it is Gaussian. At the location of each SNIa we draw, with replacement, a residual from the distribution of the residuals of all SNIa. Having done this for all SNIa, we have a randomized realization of the residuals — hence one isotropic-universe bootstrap — which we refer to as $r^{\rm iso}$. To estimate the statistical uncertainty due to the finite number of SNIa, we repeat this procedure, which is illustrated in Fig. 2, 1,000,000 times. At every 1,000-th bootstrap, we also draw the cosmological parameters $\Omega_{M}$ and $\mathcal{M}$ from their posterior distribution obtained using our cosmological analysis of the Pantheon SNIa, and then reevaluate the residuals. We do this to account for the imperfect knowledge of these parameters. Finally, we calculate the angular power spectrum for each bootstrapped realization. This procedure produces a range of values of $C_{\ell}$ expected in an isotropic universe for a discrete realization of SNIa.

Results. Figure 3 shows the angular power spectrum of the pixelated SNIa average residuals; the black error bars indicate the effect of the uncertain knowledge of cosmological parameters, corresponding to the uncertainty in theoretical magnitudes $m^{\rm th}_{i}$. We also show the distribution of $C_{\ell}$ of bootstrap rearrangements of residuals on the sky: the dark thin red curve shows the mean value of the $C_{\ell}$ due to noise and calculated from our million bootstraps, while the yellow region around it shows the 68% uncertainty in this distribution. Here and in what follows we show results for $\textsc{Nside}=16$. We have checked that our results are basically unchanged for $\textsc{Nside}=8$, 32, and 64.

Figure 3 indicates that the angular power spectrum of the Pantheon data appears consistent with the isotropic assumption given by the bootstraps, which encode the theoretical expectation of clustering noise and no signal. A simple chi-squared test confirms this; we calculate the quantity

[TABLE]

where the sum over the multipoles is implied and where $C_{\ell}$ and $\bar{C}_{\ell}^{\rm boot}$ correspond to the data and the mean of the bootstraps respectively. The coupling matrix $M_{\ell\ell^{\prime}}=\langle(C_{\ell}^{\rm boot}-\bar{C}_{\ell}^{\rm boot})(C_{\ell^{\prime}}^{\rm boot}-\bar{C}_{\ell^{\prime}}^{\rm boot})\rangle$ is calculated directly from the bootstraps, and is non-diagonal because the SNIa do not cover the full sky. We find that $\chi^{2}/{\rm dof}=1.41$ for a total of 48 degrees of freedom, being a little under the 2-$\sigma$ level (the significance is even smaller for other pixelations we looked at, $\textsc{Nside}=8,32$, and 64). Therefore, the null hypothesis of an isotropic expansion rate cannot be rejected, and this preliminary test of the isotropy is passed. We now turn to a more quantitative interpretation of our results.

We would like to get additional insight on how well our data constrain the isotropy in an as model-independent way as possible. To that effect, we consider a (redshift-independent) fractional variation in the expansion rate at $z\lesssim 1$

[TABLE]

Propagating it through to the magnitude and labeling the perturbed magnitudes with a tilde, it follows that $\widetilde{m}_{i}({\bf\hat{n}})=m_{i}({\bf\hat{n}})+(5/\ln 10)[1-\delta_{H}({\bf\hat{n}})]$. The variance of the residuals subject to such isotropy breaking, $\langle(r^{\rm data})^{2}\rangle$, is then

[TABLE]

where the noise term is given by the variance expected in an isotropic universe due to chance statistical fluctuations (as well as any systematic uncertainties), $N(r^{\rm iso})\equiv\langle(r^{\rm iso})^{2}\rangle$. Here we have assumed no correlation between the random fluctuations in the isotropic residuals $r^{\rm iso}({\bf\hat{n}})$ and the isotropy breaking $\delta({\bf\hat{n}})$, which is justified given their completely different origins.

Conveniently, our angular power spectrum measurement can be converted to the variance of the residuals on the sky via

[TABLE]

where we define the sum in the range that we measured the multipoles, $\ell\in[1,{\ell_{\rm max}}]$.

Because both the data and bootstrap variance increase with maximum multipole, each going roughly as $({\ell_{\rm max}})^{2}$ for a flat power spectrum, the results depend on the resolution of the map. Physically this makes sense, as a finer-resolution map allows for additional, smaller-scale anisotropic spatial modes that can lurk just below the detection level, and hence leads to a weaker overall constraint on the breaking of isotropy. To address this in a way that is both quantitative and physically motivated, we smooth the residual maps. It is sufficient to do this in multipole space; we apply a Gaussian beam $b_{\ell}$ that depends on the desired smoothing FWHM; the effect on the angular power spectrum is then $C_{\ell}^{\rm smoothed}=b_{\ell}^{2}C_{\ell}$.

We are finally ready to recast our results as limits on the variation of the Hubble parameter at the present time. From Eq. (6), we first cast the variance in the observed residuals as the rms Hubble parameter variation. If the noise term could be neglected, the rms Hubble parameter variation would be, according to Eq. (6)

[TABLE]

where $\langle\sigma^{-2}\rangle^{-1/2}\simeq 0.13$ is the inverse-square-weighted intrinsic dispersion of Pantheon SNIa. Given (see below) that we do not observe the evidence for a “signal” — a larger variance of the residuals than that expected in the statistically isotropic universe — the variance evaluated on isotropic bootstraps alone, $N(r^{\rm iso})$, will serve to produce an upper limit on $(\delta_{H})_{\rm rms}$ as per Eq. (8).

The principal results are shown in Fig. 4. Here we show the 68%, 95%, and 99.7% upper limits on the expansion rate variance $\langle\delta_{H}^{2}\rangle$ from the bootstraps, along with the variance computed from the Pantheon SNIa data, both as a function of the smoothing scale. The grey band shows the effect on the variance of the data due to the uncertainty in the cosmological-parameter values. As mentioned above, coarser smoothing implies more stringent constraints and vice versa. The dotted line shows the guaranteed signal given by the peculiar velocities of the SNIa in the Pantheon sample, calculated using the formalism in Huterer et al. (2017); the effect is very small and can be neglected. Note also that our limits include the marginalization over the uncertainty in the cosmological parameters since, in generating the bootstraps, we draw $\Omega_{M}$ and $\mathcal{M}$ (and hence values of theoretical magnitudes $m^{\rm th}_{i}$ and the corresponding residuals $r_{i}$) from the posterior distribution of these two parameters obtained from the cosmological analysis of the Pantheon SNIa dataset.

Figure 4 indicates no evidence for breaking of the assumption of isotropy and leads to quantitative limits on its breaking. For example, for the FWHM = 60 degree smoothing, the 99.7 percentile upper limit calculated from the bootstraps is $\mathrm{Var}(r^{\rm iso})=0.0215$, and thus

[TABLE]

or a $\lesssim$ 1% constraint on isotropy of the expansion at large angular scales. To stress-test the dependence of our constraints on data selection we have repeated the analysis with only SNIa out to maximum redshift ${z_{\rm max}}\in\{0.1,0.2,0.5,1.0\}$. While both the variance in the data and the noise limits increase with decreasing ${z_{\rm max}}$, reproduction of Fig. 4 in these cases reveals results qualitatively similar to our fiducial analysis, with no evidence for breaking of the assumption of statistical isotropy.

Conclusions. We have proposed and carried out a non-parametric test of the statistical isotropy of the late-time universe. Our test utilizes the Pantheon set of just over a thousand type Ia supernovae, whose clustering we measure by evaluating the angular power spectrum of the SNIa residuals relative to the best-fit cosmological model. We use a novel — to these tests of isotropy — and simple method of estimating the noise that describes the clustering expected in the isotropic case by bootstrapping the spatial distribution of the SNIa residuals.

To further quantify and summarize our findings, we evaluate the variance of the residuals i.e. calculate their zero-lag correlation function, and express the results in terms of constraints on the rms spatial variation of the expansion rate $(\delta_{H})_{\rm rms}\equiv\langle(\delta H({\bf\hat{n}})/H)^{2}\rangle^{1/2}$, where the latter is constrained at $z\sim 0.3$ where SNIa have the most constraining power. Because this quantity increases as smaller spatial scales are probed we explicitly smooth the angular power spectrum, evaluating the rms variation as a function of the smoothing scale.

Our results show no evidence for breaking of statistical isotropy in the Pantheon sample and, for the first time to our knowledge, constrain it at better than the 1% level at large spatial scales (smoothing FWHM $\gtrsim 60^{\circ}$); see Fig. 4.

We pay particular attention to the control and understanding of systematic errors. Our analysis choices ensure that our results do not depend on the pixelation of the map of SNIa residuals. We explicitly account for the uncertainty in the values of the cosmological parameters used to calculate the residuals and for the fact that there is a guaranteed signal of anisotropy due to the peculiar velocities of nearby objects; both effects are small and we explicitly marginalize over the former. We also find no qualitative change in our results when we restrict the range of redshifts of the SNIa in the Pantheon sample.

Our analysis does not assume the Gaussianity of the SNIa residuals, although the latter does approximately hold. We do assume the $\Lambda$CDM cosmological model; this is justifiable given the lack of evidence for its extensions (e.g. Aghanim et al. (2018)).

Our test therefore constrains the isotropy of the expansion rate at $z\lesssim 1$ at the $\sim 1\%$ level at the largest angular scales and complements the corresponding (though 2-3 orders-of-magnitude stronger) tests in the early universe.

Acknowledgments. JS and DH are supported by NASA under contract 14-ATP14-0005; DH is also supported by DOE under Contract No. DE-FG02-95ER40899. AF is supported by a McWilliams Postdoctoral Fellowship.