Exploring reionization and high-z galaxy observables with recent multiredshift MWA upper limits on the 21-cm signal

ABSTRACT

We use the latest multiredshift (z = 6.5−8.7) upper limits on the 21-cm signal from the Murchison Widefield Array (MWA) to explore astrophysical models which are inconsistent with the data. We explore these limits in the context of reionization astrophysics by using 21CMMC to connect the disfavoured regions of parameter space to existing observational constraints on reionization such as high-z galaxy ultraviolet (UV) luminosity functions, the background UV photoionization rate, the intergalactic medium (IGM) neutral fraction, the electron scattering optical depth and the soft-band X-ray emissivity. We find the vast majority of disfavoured models to already be inconsistent with existing observational constraints. These can be broadly classified into two types of models: (i) ‘cold’ reionization and (ii) pure matter density fluctuations in a cold, neutral IGM (i.e. no reionization). Interestingly, a small subsample of models inconsistent with the MWA is consistent with the aforementioned constraints (excluding the X-ray emissivity). This implies that the current MWA limits are already providing unique information to disfavour models of reionization, albeit extremely weakly. We also provide the first limits on the soft-band X-ray emissivity from galaxies at high redshifts, finding 1σ lower limits of ϵ_{X, 0.5−2 keV} ≳ 10^34.5 erg s⁻¹ Mpc⁻³. Finally, we recover 95 per cent disfavoured limits on the IGM spin temperature of |$\bar{T}_{\rm S}\lesssim$| 1.3, 1.4, 1.5, 1.8, 2.1, and 2.4 K at z = 6.5, 6.8, 7.1, 7.8, 8.2, and 8.7. With this, we infer the IGM must have undergone, at the very least, a small amount of X-ray heating. Note, the limits on ϵ_{X, 0.5–2 keV} and |$\bar{T}_{\rm S}$| are conditional on the IGM neutral fraction.

galaxies: high-redshift, intergalactic medium, dark ages, reionization, first stars, diffuse radiation, early Universe, cosmology: theory

1 INTRODUCTION

Prior to the formation of the first stars and galaxies, the Universe is opaque to visible radiation due to the neutral hydrogen fog that pervades. This fog can only be lifted once the cumulative ionizing radiation escaping from the stars and galaxies exceed the recombination rate of the neutral hydrogen. This transition is referred to as the Epoch of Reionisation (EoR), and corresponds to the final major baryonic phase change of the Universe. Observing this phase transition is vitally important as it can reveal insights into the properties of the first astrophysical sources.

Our most promising avenue for detecting the EoR is through observing the 21-cm hyperfine transition of the neutral hydrogen in the intergalactic medium (IGM) (see e.g. Furlanetto, Oh & Briggs 2006; Morales & Wyithe 2010; Pritchard & Loeb 2012; Zaroubi 2013; Barkana 2016). The radiation from these first stars and galaxies leaves an imprint on the thermal and ionization state of the IGM, detectable via this 21-cm signal. Importantly, since it is a line transition, the spatial and frequency (redshift) dependence of the 21-cm signal can yield a three-dimensional (3D) time-lapse of the history of the IGM. From this, we will be able to infer the ultraviolet (UV) and X-ray properties of the astrophysical sources responsible for reionization.

Unfortunately, observing the cosmic 21-cm signal is extremely difficult, owing to the fact that the signal is five orders of magnitude fainter than the astrophysical foregrounds. Nevertheless, this has not deterred numerous experiments from embarking on detecting this elusive signal. These experiments can be broadly classified into two types: (i) global signal experiments which spatially average the signal over the entire sky and (ii) large-scale interferometers sensitive to the spatial fluctuations in the 21-cm signal.

Global signal experiments are conceptually simpler as they typically consist of only a single dipole. Completed or ongoing experiments include, the Experiment to Detect the Global EoR Signature (EDGES; Bowman & Rogers 2010), the Sonda Cosmológica de las Islas para la Detección de Hidrógeno Neutro (SCI-HI; Voytek et al. 2014), the Shaped Antenna measurement of the background RAdio Spectrum (SARAS; Patra et al. 2015), Broad-band Instrument for Global HydrOgen ReioNisation Signal (BIGHORNS; Sokolowski et al. 2015), the Large Aperture Experiment to detect the Dark Ages (LEDA; Greenhill & Bernardi 2012; Bernardi et al. 2016), Probing Radio Intensity at high-Z from Marion (PRIZM; Philip et al. 2019) and the Netherlands–China Low-Frequency Explorer (NCLE¹).

The first-generation of large-scale radio interferometers, the Low-Frequency Array (LOFAR; van Haarlem et al. 2013), the Murchison Wide Field Array (MWA; Tingay et al. 2013; Wayth et al. 2018) and the Precision Array for Probing the Epoch of Reionisation (PAPER; Parsons et al. 2010) have limited sensitivities, requiring long integration times to potentially yield a low signal-to-noise ratio (S/N) detection of the cosmic 21-cm signal. These experiments have additionally informed the development of the next- generation of significantly larger interferometers; the Hydrogen Epoch of Reionization Array (DeBoer et al. 2017) and the Square Kilometre Array (Mellema et al. 2013; Koopmans et al. 2015). With these next-generation instruments, not only should we be able to provide high S/N statistical detections across multiple redshifts, we should also be able to provide the first three-dimensional tomographic images of the EoR.

Other than a reported detection of an absorption feature near z ≈ 17 by EDGES (Bowman et al. 2018a), whose cosmological origins are still heavily questioned (see e.g. Draine & Miralda-Escudé 2018; e.g. Hills et al. 2018; Bowman et al. 2018b; Bradley et al. 2019; Singh & Subrahmanyan 2019), we are left only with upper limits on the 21-cm signal. For example, limits on the sky-averaged signal have been obtained with LEDA (Bernardi et al. 2016), EDGES high-band (Monsalve et al. 2017) and SARAS2 (Singh et al. 2017). Upper limits on the 21-cm spatial fluctuations, measured through the power spectrum (PS), have been measured at many frequencies (redshifts) and spatial scales throughout the literature, making direct comparisons complicated. The first upper limits were achieved with the Giant Metrewave Radio Telescope (GMRT; Paciga et al. 2013) at z ≈ 8.6. Since then, upper limits have been published by LOFAR at z = 9.6−10.6 (Patil et al. 2017) and z = 19.8−25.2 (Gehlot et al. 2019), PAPER at z ≈ 7.5−11 (Cheng et al. 2018; Kolopanis et al. 2019), MWA at z ∼ 7 (Dillon et al. 2015; Beardsley et al. 2016; Barry et al. 2019; Li et al. 2019) and by the Owens Valley Radio Observatory Long Wavelength Array (OVRO-LWA; Eastwood et al. 2019) at z ≈ 18.4.

Recently, both LOFAR (Mertens et al. 2020) and the MWA (Trott et al. 2020) published their current best upper limits on the 21-cm PS. For LOFAR this culminated in a best upper limit at z ≈ 9.1 from 141 h of observations. In the case of the MWA, following careful quality control of the data, deep multiredshift limits were achieved from 298 h of observations at six separate redshifts spanning z = 6.5−8.7. While these new limits still remain a few orders of magnitude above fiducial theoretical models, they are aggressive enough to explore extreme models of reionization. For example, ‘cold’ reionization, where the 21-cm PS amplitude can be in excess of |$\Delta ^{2}_{21}\gtrsim 10^4~{\rm mK}^{2}$| due to large temperature contrasts between the neutral and ionized IGM (e.g. Mesinger, Ewall-Wice & Hewitt 2014; Parsons et al. 2014). This occurs when the neutral IGM undergoes little to no heating and adiabatically cools faster than the cosmic microwave background (CMB) temperature with the expansion of the Universe.

In this work, we shall focus explicitly on exploring the astrophysical models inconsistent with the MWA multiredshift limits from Trott et al. (2020). Similar analyses have already been performed for LOFAR in the context of general IGM properties such as the neutral fraction and spin temperature (Ghara et al. 2020), the excess radio background (Mondal et al. 2020) and high-z galaxy and reionization observables Greig et al. (2020). Here, we follow the same approach as Greig et al. (2020).

Our exploration takes advantage of 21CMMC  ² (Greig & Mesinger 2015, 2017; Greig & Mesinger 2018; Park et al. 2019), a Markov chain Monte Carlo (MCMC) sampler of the seminumerical reionization code 21cmfast  ³(Mesinger & Furlanetto 2007; Mesinger, Furlanetto & Cen 2011). 21CMMC forward models the full 3D cosmic 21-cm signal provided by 21cmfast in a fully Bayesian framework allowing us to compare against observations of the first billion years of the Universe. In particular, we adopt the Park et al. (2019) galaxy model parametrization allowing direct comparison against high-z galaxy UV luminosity functions (LFs).

The outline of this work is as follows. In Section 2, we summarize the 21cmfast astrophysical model before outlining the 21CMMC setup in Section 3. In Section 4, we discuss our main results before providing our conclusions in Section 5. Unless otherwise stated, all quantities are in co-moving units with the following adopted cosmological parameters: (⁠|$\Omega _\Lambda$|⁠, Ω_M, Ω_b, n, σ₈, H₀)  = (0.69, 0.31, 0.048, 0.97, 0.81, 68 km s⁻¹ Mpc⁻¹), consistent with recent results from the Planck mission (Planck Collaboration XIII 2016).

2 MODELLING THE 21-CM SIGNAL

The cosmic 21-cm signal is modelled using the seminumerical simulation code, 21cmfast (Mesinger & Furlanetto 2007; Mesinger et al. 2011). In particular, we use the Park et al. (2019) astrophysical parametrization, which allows the star formation rate (SFR) and ionizing escape fraction to depend on the mass of the host dark matter halo. This, following some simple conversions, enables 21cmfast to produce UV LFs that are able to be compared against observed high-z galaxy LFs. Additionally, we also include a recipe for an on-the-fly ionizing photon non-conservation correction (Park et al., in preparation) to account for the fact that 3D excursion set approaches that track ionizations are not photon conserving (e.g. McQuinn et al. 2005; Zahn et al. 2007; Paranjape & Choudhury 2014; Paranjape, Choudhury & Padmanabhan 2016; Hassan et al. 2017; Choudhury & Paranjape 2018; Hutter 2018; Molaro et al. 2019). Below, we shall outline the important aspects of 21cmfast which lead to the modelling of the 21-cm signal, and defer the reader to the aforementioned publications for further details.

2.1 Galaxy UV properties

The typical stellar mass of a galaxy, M_*, is assumed to be directly related to its host halo mass, M_h (e.g. Kuhlen & Faucher-Giguère 2012; Dayal et al. 2014; Behroozi & Silk 2015; Mitra, Choudhury & Ferrara 2015; Mutch et al. 2016; Ocvirk et al. 2016; Sun & Furlanetto 2016; Yue, Ferrara & Xu 2016; Hutter et al. 2020):

$$\begin{eqnarray*} M_{\ast }(M_{\rm h}) = f_{\ast }\left(\frac{\Omega _{\rm b}}{\Omega _{\rm m}}\right)M_{\rm h}, \end{eqnarray*}$$

(1)

where f_* is the fraction of galactic gas in stars given by

$$\begin{eqnarray*} f_{\ast } = f_{\ast , 10}\left(\frac{M_{\rm h}}{10^{10}\, {\rm M}_{\odot }}\right)^{\alpha _{\ast }}, \end{eqnarray*}$$

(2)

with f_{*, 10} being its normalization and a power-law⁴ index, α_*.

The SFR is then estimated by dividing the stellar mass by a characteristic time-scale,

$$\begin{eqnarray*} \dot{M}_{\ast }(M_{\rm h},z) = \frac{M_{\ast }}{t_{\ast }H^{-1}(z)}, \end{eqnarray*}$$

(3)

where t_* is a free parameter allowed to vary between zero and unity and H⁻¹(z) is the Hubble time.

A galaxy’s UV ionizing escape fraction, f_esc, is equivalently allowed to vary with halo mass,

$$\begin{eqnarray*} f_{\rm esc} = f_{\rm esc, 10}\left(\frac{M_{\rm h}}{10^{10}\, {\rm M}_{\odot }}\right)^{\alpha _{\rm esc}}, \end{eqnarray*}$$

(4)

with normalization set by f_{esc, 10} and a power-law index, α_esc.

Finally, to account for inefficient cooling and/or feedback processes which can prevent small-mass haloes from hosting active, star-forming galaxies, a duty-cycle is included to suppress their contribution:

$$\begin{eqnarray*} f_{\rm duty} = {\rm exp}\left(-\frac{M_{\rm turn}}{M_{\rm h}}\right). \end{eqnarray*}$$

(5)

This results in a fraction, (1 − f_duty), of host haloes which do not host star-forming galaxies, whose suppression scale is controlled by M_turn (e.g. Giroux, Sutherland & Shull 1994; Shapiro, Giroux & Babul 1994; Hui & Gnedin 1997; Barkana & Loeb 2001; Springel & Hernquist 2003; Mesinger & Dijkstra 2008; Okamoto, Gao & Theuns 2008; Sobacchi & Mesinger 2013a, b).

In summary, we have six free parameters, f_{*, 10}, f_{esc, 10}, α_*, α_esc, M_turn and t_*, which describe the galaxy UV properties.

2.2 Galaxy X-ray properties

In the early Universe, X-rays escaping from the first galaxies, likely from stellar remnants, are thought to be responsible for the heating of the IGM. X-ray heating in 21cmfast is modelled by calculating the cell-by-cell angle-averaged specific X-ray intensity, |$J(\boldsymbol{x}, E, z)$|⁠, (in erg s⁻¹ keV⁻¹ cm⁻² sr⁻¹). This is obtained by integrating the co-moving X-ray specific emissivity, |$\epsilon _{\rm X}(\boldsymbol{x}, E_e, z^{\prime })$| back along the light-cone:

$$\begin{eqnarray*} J(\boldsymbol{x}, E, z) = \frac{(1+z)^3}{4\pi } \int _{z}^{\infty } {\rm d}z^{\prime } \frac{c {\rm d}t}{{\rm d}z^{\prime }} \epsilon _{\rm X} {\rm e}^{-\tau }, \end{eqnarray*}$$

(6)

where e^−τ accounts for attenuation by the IGM. In the emitted frame, E_e = E(1 + z^′)/(1 + z), the co-moving specific emissivity is

$$\begin{eqnarray*} \epsilon _{\rm X}(\boldsymbol{x}, E_{\rm e}, z^{\prime }) = \frac{L_{\rm X}}{\rm SFR} \left[ (1+\bar{\delta }_{\rm nl}) \int ^{\infty }_{0}{\rm d}M_{\rm h} \frac{{\rm d}n}{{\rm d}M_{\rm h}}f_{\rm duty} \dot{M}_{\ast }\right], \end{eqnarray*}$$

(7)

where |$\bar{\delta }_{\rm nl}$| is the mean, non-linear density in a shell around |$(\boldsymbol{x}, z)$|⁠, |$\frac{{\rm d}n}{{\rm d}M_{\rm h}}$| corresponds to the halo mass function (HMF)⁵ and the quantity in square brackets is the SFR density along the light cone.

This expression is normalized by the specific X-ray luminosity per unit star formation escaping the host galaxies, L_X/SFR (erg s⁻¹ keV|$^{-1}\, {\rm M}^{-1}_{\odot}$| yr). The X-ray luminosity is assumed to be a power law with respect to photon energy, |$L_{\rm X} \propto E^{- \alpha _{\rm X}}$|⁠, which is attenuated below a threshold energy, E₀, to account for the absorption of low-energy X-rays by a neutral interstellar medium within the host galaxy. Finally, this specific luminosity is then normalized to an integrated soft-band (<2 keV) luminosity per SFR (in erg s|$^{-1}\, {\rm M}^{-1}_{\odot}$| yr), which is a free parameter in the model:

$$\begin{eqnarray*} L_{{\rm X}\lt 2\, {\rm keV}}/{\rm SFR} = \int ^{2\, {\rm keV}}_{E_{0}} {\rm d}E_{\rm e} ~ L_{\rm X}/{\rm SFR} ~. \end{eqnarray*}$$

(8)

A limit of |$2\, {\rm keV}$| corresponds to a mean-free path of the order of the Hubble length at high redshifts, implying harder X-ray photons do not contribute to IGM heating (e.g. McQuinn 2012).

In summary, there are three free model parameters describing the X-ray properties of the first galaxies, |$L_{{\rm X}\lt 2\, {\rm keV}}/{\rm SFR}$| , E₀ and α_X.

2.3 Ionization and thermal state of the IGM

At any redshift within 21cmfast the evolved IGM density and velocity fields are obtained following second-order Lagrange perturbation theory (e.g Scoccimarro 1998) from an initial high-resolution linear density field. The ionization field is then determined from the evolved density field by the excursion-set approach (Furlanetto, Zaldarriaga & Hernquist 2004), whereby the balance between the cumulative number of ionizing photons and the number of neutral hydrogen atoms plus cumulative recombinations are tracked within spheres of decreasing radii. A cell is considered ionized when,

$$\begin{eqnarray*} n_{\rm ion}(\boldsymbol{x}, z | R, \delta _{R}) \ge (1 + \bar{n}_{\rm rec})(1-\bar{x}_{\rm e}), \end{eqnarray*}$$

(9)

where |$\bar{n}_{\rm rec}$| is the cumulative number of recombinations (e.g. Sobacchi & Mesinger 2014) and |$(1-\bar{x}_{\rm e})$|⁠, corresponds to ionizations by X-rays, which are expected to contribute at the ∼10 per cent level (e.g. Ricotti & Ostriker 2004; Mesinger, Ferrara & Spiegel 2013; Madau & Fragos 2017; Ross et al. 2017; Eide et al. 2018). This first term, n_ion, is the cumulative number of IGM ionizing photons per baryon inside a spherical region of size, R and corresponding overdensity, δ_R,

$$\begin{eqnarray*} n_{\rm ion} = \bar{\rho }^{-1}_{\rm b}\int ^{\infty }_{0}{\rm d}M_{\rm h} \frac{{\rm d}n(M_{h}, z | R, \delta _{R})}{{\rm d}M_{\rm h}}f_{\rm duty} \dot{M}_{\ast }f_{\rm esc}N_{\gamma /b}, \end{eqnarray*}$$

(10)

where |$\bar{\rho }_b$| is the mean baryon density and N_γ/b is the number of ionizing photons per stellar baryon.⁶ The final term in equation (9), |$(1-\bar{x}_{\rm e})$|⁠, corresponds to the contribution to the number of ionizations by X-rays, which are expected to contribute at around the ∼10 per cent level (e.g. Ricotti & Ostriker 2004; Mesinger et al. 2013; Madau & Fragos 2017; Ross et al. 2017; Eide et al. 2018).

To calculate the thermal state of the neutral (and partially ionized) IGM, 21cmfast self-consistently calculates the heating and cooling rates from structure formation and Compton scattering with the CMB photons, heating from partial ionizations as well as the heating and ionizations caused by the X-rays (as discussed earlier). With these, the IGM spin temperature, T_S, is determined as the weighted mean of the kinetic gas, T_K, Lyman-α (Lyα) radiation, T_α, and CMB, T_CMB, temperatures:

$$\begin{eqnarray*} T^{-1}_{\rm S} = \frac{T^{-1}_{\rm CMB} + x_{\alpha }T^{-1}_{\alpha } + x_{\rm c}T^{-1}_{\rm K}}{1 + x_{\alpha } + x_{\rm c}}, \end{eqnarray*}$$

(11)

where x_α describes the Wouthuysen–Field coupling coefficient (Wouthuysen 1952; Field 1958) and x_c describes the collisional coupling coefficient between the free electrons and the CMB photons.

This background of Lyα radiation is obtained from two processes (see e.g. Mesinger et al. 2011, for further details). First, excitation of the neutral hydrogen atoms by the X-ray photons (J_{α, X}), which is proportional to the X-ray heating rate provided energy injection is balanced by photons redshifting out of Lyα resonance (Pritchard & Furlanetto 2007). The second contribution arises from the direct stellar emission⁷ of Lyman band photons by the first sources, where Lyman-n resonance photons absorbed by the IGM will cascade through Lyα producing a background, J_{α, *}. This is produced by the cumulative sum over all Lyman resonances (e.g. Barkana & Loeb 2005). Note, there are no additional free parameters used to calculate the Lyα background other than those used to calculate the X-ray and UV ionizing photons.

2.4 Ionizing photon non-conservation correction

Excursion-set approaches for tracking ionizations in three dimensions, like that employed by 21cmfast, do not conserve ionizing photons. This is driven by the surplus of ionizing photons remaining within a cell after exceeding the ionization criteria (equation 9) that are not propagated further into the IGM. In effect, this results in an effective bias on the ionizing escape fraction, f_esc. This behaviour has been studied extensively in the literature and has been shown to result in a loss of ∼10−20 per cent of the ionizing photons (e.g. McQuinn et al. 2005; Zahn et al. 2007; Paranjape & Choudhury 2014; Paranjape et al. 2016; Hassan et al. 2017; Choudhury & Paranjape 2018; Hutter 2018; Molaro et al. 2019).

Recently, explicit photon conserving algorithms for seminumerical simulations have been introduced by Choudhury & Paranjape (2018) and Molaro et al. (2019). However, despite being orders of magnitude faster than full radiative-transfer simulations they remain too slow when forward modelling the 21-cm signal in the high-dimensional parameter spaces required to characterize the ionizing, soft UV, and X-ray properties of the first galaxies.

Alternatively, Park et al. (in preparation), introduce an approximate correction to the effective bias on f_esc by analytically solving for the evolution of the ionization fraction given a source model, assuming no correlations between the sources and sinks.⁸ By comparing this analytic expression against a calibration curve (generated from 21cmfast including only ionizations) the delay in the resultant reionization history (owing to the loss of photons) can be corrected for by modifying the redshift at which ionizations are determined. For the duration of the EoR, this correction results in a shift in redshift of Δz ∼ 0.3 ± 0.1.

2.5 21-cm brightness temperature

Finally, the observed 21-cm signal is measured as a brightness temperature contrast relative to the CMB temperature, T_CMB (e.g. Furlanetto et al. 2006):

$$\begin{eqnarray*} \delta T_{\rm b}(\nu) &=& \frac{T_{\rm S} - T_{\rm CMB}(z)}{1+z}\left(1 - {\rm e}^{-\tau _{\nu _{0}}}\right)~{\rm mK}, \end{eqnarray*}$$

(12)

where |$\tau _{\nu _{0}}$| is the optical depth of the 21-cm line,

$$\begin{eqnarray*} \tau _{\nu _{0}} &\propto & (1+\delta _{\rm nl})(1+z)^{3/2}\frac{x_{\mathrm{H\, {\scriptscriptstyle I}}{}}}{T_{\rm S}}\left(\frac{H}{{\rm d}v_{\rm r}/{\rm d}r+H}\right). \end{eqnarray*}$$

(13)

Here, |$x_{\mathrm{H\, {\scriptscriptstyle I}}{}}$| corresponds to the neutral hydrogen fraction, |$\delta _{\rm nl} \equiv \rho /\bar{\rho } - 1$| is the gas overdensity with ρ and |$\bar{\rho}$| corresponding to the local and mean gas densities, respectively, H(z) is the Hubble parameter, dv_r/dr is the line-of-sight component of the velocity gradient and T_S is the gas spin temperature. All quantities, obtained following the method described in the previous sections are evaluated at redshift z = ν₀/ν − 1, where ν₀ is the 21-cm frequency and we drop the spatial dependence for brevity. Redshift space distortions along the line of sight are additionally included as outlined in Mao et al. (2012), Jensen et al. (2013), and Greig & Mesinger (2018).

3 21CMMC SETUP

21CMMC (Greig & Mesinger 2015, 2017; Greig & Mesinger 2018; Park et al. 2019) is the publicly available, massively parallel MCMC sampler of the 3D seminumerical reionization simulation code 21cmfast (Mesinger & Furlanetto 2007; Mesinger et al. 2011). It is based on the publicly available Python module CosmoHammer (Akeret et al. 2013). MCMC sampling is performed using the Emcee Python module (Foreman-Mackey et al. 2013), which is an affine invariant ensemble sampler from Goodman & Weare (2010). For each proposed parameter set in the MCMC, 21CMMC performs an independent 3D realization of the 21-cm signal, allowing any quantity to be sampled (e.g. the 21-cm PS). Below, we outline the astrophysical parameter set to be explored within this work (Section 3.1), the interpretation of the MWA upper limits (Section 3.2) and the 21CMMC setup (Section 3.3).

3.1 Astrophysical parameter set

For the adopted model described in Section 2, we have nine astrophysical parameters. Below, we summarize each parameter as well as the corresponding assumed parameter ranges (which we have adopted based on previous works, e.g. Greig et al. 2017; Park et al. 2019). Additionally, we summarize these parameter ranges in the top row of Table 1. Throughout, we assume flat priors on each astrophysical parameters.

• f_{*, 10}: The fraction of galactic gas in stars evaluated at a halo mass of 10¹⁰M_⊙. The log quantity is allowed to vary as, log₁₀(f_{*, 10}) ∈ [ − 3, 0].

• α_*: the power-law index for the halo mass-dependent star formation. This is varied between, α_* ∈ [ − 0.5, 1].

• f_{esc, 10}: The ionizing UV escape fraction evaluated at a halo mass of 10¹⁰ M_⊙. The log of this quantity varies between, log₁₀(f_{esc, 10}) ∈ [ − 3, 0].

• α_esc: The power-law index for halo mass-dependent ionizing UV escape fraction. This is varied between, α_esc ∈ [ − 1, 0.5].

• t_*: The star formation time-scale as a fraction of the Hubble time, which is varied in the range, t_* ∈ (0, 1].

• M_turn: The characteristic halo mass below which the abundance of active star-forming galaxies are exponentially suppressed by a duty cycle (see equation 5). We allow this to vary between, log₁₀(M_turn) ∈ [8, 10]. This range corresponds to star formation by atomically cooled hydrogen, where the lower limit corresponds to the threshold for atomic cooling and the upper limit is consistent with host halo mass estimates from observations of Lyman break galaxies at z ∼ 6 − 8 (e.g. Kuhlen & Faucher-Giguère 2012; Barone-Nugent et al. 2014).

• E₀: the minimum energy threshold above which X-ray photons can escape their host galaxy. We allow this to vary between, E₀ ∈ [0.2, 1.5] keV, which corresponds of an integrated column density of, |${\rm log_{10}}(N_{\mathrm{H\, {\scriptscriptstyle I}}{}}/{\rm cm}^{2}) \in [19.3,23.0]$|⁠.

• |$L_{{\rm X}\lt 2\, {\rm keV}}/{\rm SFR}$|⁠: The soft-band X-ray luminosity per unit star formation from the E₀ − 2 keV energy band. We vary this between, |${\rm log_{10}}(L_{{\rm X}\lt 2\, {\rm keV}}/{\rm SFR}) \in [30, 42]$|⁠. This lower bound is considerably lower than what is typically adopted (i.e. Park et al. 2019) in-order to explore extreme ‘cold’ reionization scenarios.

• α_X: The power-law index of the X-ray source spectral energy distribution (SED), which we allow to vary between α_X ∈ [ − 1, 3].

Table 1.

Summary of the 68th and 95th percentile limits on the disfavoured regions on the nine astrophysical parameters included in 21CMMC using the z = 6.5 − 8.7 upper limits from the MWA (Trott et al. 2020).

	log10(f*, 10)	α*	log10(fesc, 10)	αesc	log10(Mturn)	t*	|${\rm log_{10}}\left(\frac{L_{{\rm X}\lt 2{\rm keV}}}{\rm SFR}\right)$|	E0	αX
					(M⊙)		|$({\rm erg\, s^{-1}\, M_{\odot }^{-1}\, yr})$|	(keV)
Prior ranges	[−3.0, 0.0]	[−0.5,1.0]	[−3.0, 0.0]	[−1.0,0.5]	[8.0,10.0]	(0.0,1.0]	[30.0,42.0]	[0.2,1.5]	[−1.0,3.0]
68th percentile limits	[−2.40, −0.75]	[−0.5, 0.75]	[−3.0, −1.75]	[−1.0, 0.21]	[8.40, 10.0]	(0.0, 0.55]	[30.0, 35.9]	[0.35,1.5]	[−1.0,3.0]
95th percentile limits	[−2.88, −0.11]	[−0.5, 0.96]	[−3.0, −0.80]	[−1.0, 0.45]	[8.06, 10.0]	(0.0, 0.91]	[30.0, 37.8]	[0.14,1.5]	[−1.0,3.0]

3.2 The latest MWA 21-cm upper limits

Recently, Trott et al. (2020) published deep, multiredshift upper limits on the 21-cm signal from the EoR using the MWA. These limits are obtained at k = 0.07−3.0 h Mpc⁻¹ across six redshift bins from z = 6.5−8.7 using 298 h of carefully selected clean observations over four observing seasons. At z = 6.5, these correspond to the lowest yet available upper limits on the reionization epoch.

Unfortunately, these upper limits are still too large to begin to provide statistical constraints on astrophysical models of reionization. However, we can instead use these upper limits to explore astrophysical models that exceed (i.e. are disfavoured by) the existing observational limits placed by the MWA. This follows the same approach to that applied to the recent upper limits on the 21-cm signal at z ≈ 9.1 achieved by LOFAR (Ghara et al. 2020; Greig et al. 2020; Mondal et al. 2020).

In order to explore astrophysical models disfavoured by the MWA, we construct a likelihood function of the following form:

$$\begin{eqnarray*} \mathcal {L}(\theta) \propto \prod ^{n}_{i,j} \frac{\int ^{(1 + \sigma)\Delta ^{2}_{21,{\rm mod.}}}_{(1 - \sigma)\Delta ^{2}_{21,{\rm mod.}}} p_{\rm ex.}(\Delta ^{2}_{21,{\rm mod.}}(k_{i},z_{j},\theta)) {\rm d}\Delta ^{2}_{21,{\rm mod.}}}{\int ^{(1 + \sigma)\Delta ^{2}_{21,{\rm mod.}}}_{(1 - \sigma)\Delta ^{2}_{21,{\rm mod.}}} {\rm d}\Delta ^{2}_{21,{\rm mod.}}}, \end{eqnarray*}$$

(14)

where θ corresponds to the astrophysical parameter set (i.e. model), |$\Delta ^{2}_{21, {\rm mod.}}(k_{i},z_{j},\theta)$| is the model 21-cm PS (from 21cmfast) as a function of k and z and |$p_{\rm ex.}(\Delta ^{2}_{21,{\rm mod.}}(k_{i},z_{j},\theta))$| is the probability of the model 21-cm PS exceeding the MWA upper-limit (see Appendix  A for further details). The integral over the range ±σ accounts for uncertainties on the amplitude of the modelled 21-cm PS.

The probability of the 21-cm PS to be in excess of the upper limits is obtained from the true, measured probability density functions (PDFs) of the 21-cm PS amplitude used to construct the observed upper limits from Trott et al. (2020). In effect, |$p_{\rm ex.}(\Delta ^{2}_{21,{\rm mod.}}(k_{i},z_{j},\theta))$| goes to zero for decreasing 21-cm PS amplitudes below the MWA upper limits and approaches unity for model 21-cm PS amplitudes far in excess of the upper limit.

In this work, we combine the data from all six redshift bins (z = 6.5, 6.8, 7.1, 7.8, 8.2, and 8.7) with the first four Fourier bins, k = 0.14, 0.21, 0.28, and 0.35 h Mpc⁻¹. Note that for smaller scales (i.e. larger k), the upper limits (and associated errors) become sufficiently large that they provide very little additional information (see e.g. Fig. 1). It is important to note that we only consider the disfavoured regions obtained from using the full data set. As such, the disfavoured regions are going to be preferentially skewed towards regions of the parameter space where the ionizing source populations can produce a 21-cm PS that can more readily exceed the lowest amplitude 21-cm upper limits (i.e. at the lowest redshifts). Importantly, at this point with still relatively extreme upper limits, we caution against taking too strong a physical meaning from the disfavoured regions of parameter space.

Figure 1.

A random sample of 500 21-cm PS as a function of redshift disfavoured by the latest MWA upper limits (Trott et al. 2020). The four panels correspond to the first four Fourier modes (k = 0.14, 0.21, 0.28, and 0.35 h Mpc⁻¹). The upper limits (yellow) correspond to the 2σ MWA upper limits (Trott et al. 2020), whereas the blue (red) curves correspond to the 21-cm PS colour coded by the resultant neutral fraction at z = 5.9.

3.3 Simulations

Now, with our likelihood described as above, we outline the simulation setup adopted for 21CMMC. In order to sufficiently sample the upper limit at the largest spatial scale (k = 0.14 h Mpc⁻¹), we perform 3D realizations of the 21-cm signal in a volume with length 250 Mpc and 128 voxels per sidelength. For our total uncertainty, σ, on the modelled 21-cm PS, we sum in quadrature a conservative 20 per cent multiplicative modelling uncertainty⁹ on the sampled 21-cm PS with the estimated sample variance from our simulation setup. Note that for this simulation setup, the estimated sample variance on the 21-cm PS at k = 0.14 h Mpc⁻¹ is ≲ 10 per cent.

For this exploration of the disfavoured regions of astrophysical parameter space with 21CMMC, we use 400 walkers each performed for 300 iterations. In total, this results in 1.2 × 10⁵ realizations of 21cmfast and have confirmed that under this setup it is sufficient to reach a converged result.

4 RESULTS

4.1 Disfavoured 21-cm PS

Before exploring the implications of the latest MWA upper limits (Trott et al. 2020) on the individual astrophysical parameters describing reionization, it is first illustrative to explore the 21-cm PS which are disfavoured by the observational limits. In Fig. 1, we show the evolution in the 21-cm PS as a function of redshift for the first four Fourier modes (k = 0.14, 0.21, 0.28, 0.35 h Mpc⁻¹). The thin blue and red curves correspond to 500 randomly sampled 21-cm PS which are separated by the resultant neutral fraction at z = 5.9 in order to distinguish between two types of models.

Typically, when presenting the 21-cm PS evolution as a function of redshift, one expects to observe three distinctive peaks at large-scales in the PS amplitude (e.g. Pritchard & Furlanetto 2007; Baek et al. 2010; Mesinger et al. 2011). These are, in order of increasing redshift, driven by (i) reionization, (ii) X-ray heating and (iii) Wouthuysen–Field (WF) coupling. In the case of these extreme ‘cold’ reionization models (i.e. absence of X-ray heating), we expect no secondary peak, with a single large-amplitude peak during reionization along with a WF-coupling peak at higher redshifts.¹⁰

This reionization peak is clearly exhibited by the blue curves in Fig. 1, which corresponds to IGM neutral fractions below 90 per cent¹¹ at z = 5.9. Here, the 21-cm PS amplitude peaks at z = 6 − 7 before rapidly dropping in amplitude as reionization completes. These models are indicative of these ‘cold’ reionization models, with the fluctuations being driven by the patchy nature of reionization and the amplitude as a result of the large temperature contrasts in the cold, neutral IGM. These astrophysical models clearly exceed the recent 2σ upper limits presented by the MWA (Trott et al. 2020), highlighting that the MWA can already begin to rule out such extreme models.

Note, however, that the 21-cm PS can only exceed the latest MWA upper limits at the largest scale (k = 0.14 h Mpc⁻¹). For larger k (smaller scales) the MWA upper limits increase significantly to render all models consistent with the data. Thus, almost all of the constraining power arises from the redshift evolution of the large-scale power. This highlights that, in order to produce the largest gains in ruling out regions of astrophysical parameter space, in the short-term observational efforts should focus on improving the 21-cm PS upper limits at large-scales, where this ‘bump’ in the redshift evolution occurs as a result of reionization.

For the second type of models, corresponding to an IGM neutral fraction above 90 per cent at z = 5.9 (denoted by the red curves in Fig. 1), we observe no distinctive features. The 21-cm PS amplitude exhibits a smooth increase in amplitude with decreasing redshift. Again, the large amplitudes are driven by the large temperature contrasts in the cold, neutral IGM, however, the high neutral fractions (as indicated by the large IGM neutral fraction at z = 5.9) are indicative of models where reionization has yet to commence, or has only just begun. Thus, these smooth, featureless 21-cm PS are simply driven by fluctuations in the underlying matter density field. While these models are already strongly ruled out given existing limits on the IGM neutral fraction at z = 5.9, this is the first-time astrophysical models driven by pure matter density fluctuations have been found to be disfavoured by limits from the 21-cm signal.

4.2 Disfavoured astrophysical parameters

In Fig. 2, we present the marginalized one and two-dimensional posterior distributions recovered from 21CMMC for the nine astrophysical parameters disfavoured by the multifrequency upper limits on the 21-cm PS from the MWA (Trott et al. 2020). It is important to note that these posteriors are obtained using only the MWA upper limits, with no other existing observation constraints used. Dashed (solid) white contours on the 2D posteriors correspond to the 68th (95th) percentile limits. The one-dimensional marginalized 68th and 95th percentile limits on each parameter are summarized in Table 1. Note that the marginalized posteriors for several of these astrophysical parameters show strong degeneracies and multimodal features (in particular f_{*, 10}, α_*, f_{esc, 10}, α_esc, and M_turn). These occur owing to the two classes of models¹² that the existing MWA upper limits are capable of disfavouring (see Section 4.1), which can inhabit different regions of the parameter space.

Figure 2.

Marginalized one and two-dimensional posterior distributions for the nine astrophysical parameters which are disfavoured by the multifrequency upper limits on the 21-cm PS from the MWA (Trott et al. 2020). White dashed (solid) contours correspond to the 68th (95th) percentiles. In the top right-hand panels, we provide 500 randomly sampled LFs which are drawn from the posterior of astrophysical models that exceed at least one of the MWA upper limits and are compared against existing constraints at z = 6−8 (Bouwens et al. 2015, 2017) and z = 10 (Oesch et al. 2018). In the middle right-hand panel, we compare the bounds on the reionization histories disfavoured by MWA against all current observational constraints on the IGM neutral fraction (see text for further details).

At this point, it is important to remember that these marginalized posteriors only contain information about the astrophysical models which can exceed the recent MWA upper limits on the 21-cm signal. These marginalized posteriors do not imply that these regions of astrophysical parameter space are ruled out at any statistical significance. Rather, it is indicative of the types of extreme astrophysical models that observational experiments such as the MWA are already starting to disfavour.¹³ Further, it is a demonstration of the necessity of forward modelling tools such as 21CMMC for exploring and interpreting observational data. Nevertheless, below we discuss some of the implications on the EoR that can be made from the recent MWA upper limits.

Additionally, it is also worth noting that the results presented in this work are specific to the underlying astrophysical model assumptions. For example, the astrophysical parametrization used in 21cmfast assumes only a single ionizing source population and further, ignores any explicit redshift dependence on the escape fraction or stellar mass. However, there is currently no evidence for a more complex parametrization as the simple Park et al. (2019) model has been shown to be consistent with galaxy UV LFs, seminumerical semi-analytic models and hydrodynamical simulations.

The main observation, we can take from Fig. 2 are the limits on the soft-band X-ray luminosity, |$L_{{\rm X}\lt 2\, {\rm keV}}/{\rm SFR}$|⁠. The primary role of this X-ray luminosity is to control the amount of heating the IGM incurs between the dark ages and the EoR as a result the thermal energy deposited by the X-rays. Low values of this quantity are responsible for producing both (i) the extreme ‘cold’ reionization scenarios represented by the blue curves in Fig. 1 and (ii) the large-amplitude 21-cm PS due only to the matter density fluctuations (red curves in Fig. 1). With these latest MWA upper limits we recover, at 95 per cent confidence, disfavoured limits of |${\rm log_{10}}\left(L_{{\rm X}\lt 2\, {\rm keV}}/{\rm SFR}\right) \lesssim 37.8$|⁠, implying the MWA disfavours models with lower soft-band X-ray luminosities. Presently, these disfavoured limits sit well below our current expectations from observations of analogue low-redshift star-forming galaxies (Mineo, Gilfanov & Sunyaev 2012), stacked Chandra observations (Lehmer et al. 2016) and predictions at high-redshift from population synthesis models (Fragos et al. 2013) (see Section 4.3.5 for more details). In contrast, for the remaining two galaxy X-ray properties, E₀ and α_X, we recover no meaningful limits, owing to the fact we are disfavouring low-amplitude soft-band X-ray luminosities. In the absence of X-ray heating, the shape of the X-ray SED does not matter.

Unfortunately, as we are only exploring models in excess of the 21-cm upper limits, it is not fair to perform direct comparisons between these results and those presented by LOFAR at z ≈ 9.1. The primary reason for this is that this approach is sensitive only to models that maximize the 21-cm PS (in order to exceed the upper limits on the 21-cm signal). In this work, the MWA upper limits are most sensitive to reionization models that produce 21-cm PS that peak at z < 7.5 (see Fig. 1), whereas the current LOFAR limits are sensitive to only a single redshift at z ≈ 9.1. Further, the existing LOFAR limits are more sensitive to the larger spatial scales (i.e. k = 0.075 h Mpc⁻¹). Thus the disfavoured regions of each experiment are sensitive to different reionization models. Nevertheless, one could, in principle, combine the recent limits from both the MWA and LOFAR, to improve our understanding of these disfavoured regions, however, we postpone that to future work.

In terms of the galaxy UV properties, the strongest disfavoured limits are on the fraction of galactic gas in stars, f_{*, 10}, the escape fraction, f_{esc, 10}, and the SFR time-scale, t_*. Again, these limits are driven by seeking to maximize the 21-cm PS amplitude between z = 6.5−8.7 (the reionization ‘peak’ at z = 6−7 in Fig. 1). Qualitatively speaking, this occurs when reionization is roughly close to its midpoint (i.e. |$\bar{x}_{\mathrm{H\, {\scriptscriptstyle I}}{}} \sim 0.5$|⁠; Mellema et al. 2006; Lidz et al. 2008).¹⁴ For example, both f_{*, 10} and f_{esc, 10}, which are highly degenerate in the absence of constraints from UV LFs (e.g. Park et al. 2019), control the timing of reionization. Thus, limits on these quantities are driven to very low values (minimizing the number of ionizing photons produced), pushing reionization to occur at lower redshifts (i.e. z = 6−7). However, rather than producing a limit for f_{*, 10}, we see a clear peak for f_{*, 10} at log₁₀(f_{*, 10}) ∼ −2.0. This is driven by the overlap of the resultant posteriors for the two distinct models outlined in the previous section (i.e. ‘cold’ reionization and the matter density fluctuations).

For t_*, the disfavoured limits are driven by their degeneracy with the soft-band X-ray luminosity.¹⁵ In this model, the number of X-ray photons produced is inversely proportional to the SFR time-scale. For extremely short SFR time-scales, this can result in a large number of X-ray photons, subsequently heating the IGM and decreasing the amplitude of the resultant 21-cm PS. As a result, this places a disfavoured limit on short star formation time-scales.

4.3 Comparison against existing observations

Now that we have explored the disfavoured astrophysical parameters above, we now shift our focus towards how these models compare against existing observational constraints on the reionization epoch.

4.3.1 Reionization history

In the middle right-hand panel of Fig. 2, we explore the reionization histories of these disfavoured astrophysical models. In the shaded region, we present the full range of reionization histories inconsistent with the latest MWA upper limits on the 21-cm signal. Overlaid on this, we highlight all existing constraints on the IGM neutral fraction. These include limits from the dark pixel statistics of high-z quasars (QSOs; McGreer, Mesinger & D’Odorico 2015), the Lyα fraction (Mesinger et al. 2015), the clustering of Lyα emitters (LAEs; Sobacchi & Mesinger 2015), the Lyα equivalent width distribution of Lyman-break galaxies (LBGs; Mason et al. 2018; Hoag et al. 2019; Mason et al. 2019), the neutral IGM damping wing imprint from high-z QSOs (Greig et al. 2017; Davies et al. 2018; Greig, Mesinger & Bañados 2019) and the midpoint of reionization (z_Re) from Planck (Planck Collaboration VI 2018). Note that this shaded region represents only those models disfavoured by (i.e. exceed) the recent MWA upper limits. The vast majority of astrophysical models consistent with the existing observational constraints produce 21-cm PS which are below the existing upper limits from MWA.

Importantly, since the range of disfavoured reionization histories is completely consistent with existing observational constraints, this implies that the latest MWA upper limits on the 21-cm signal are already sufficient to disfavour models which would otherwise be consistent with existing observations. That is, the MWA is already providing unique constraining power on the astrophysics of reionization. However, owing to the still large-amplitude limits, this constraining power is extremely weak. This differs from the picture presented in the equivalent analysis of the LOFAR upper limits (Greig et al. 2020), where the disfavoured models were already inconsistent with constraints on the IGM neutral fraction. The primary reason for this difference is the fact that the recent MWA limits are at multiple redshift (significantly broadening the range of disfavoured reionization histories). For example, the models currently disfavoured by LOFAR exceed the upper limits at z ≈ 9.1, whereas the models currently disfavoured by the MWA exceed the upper limits anywhere between z = 6.5−8.7.

4.3.2 UV luminosity functions

In the top right-hand panel of Fig. 2, we show 500 randomly drawn UV LFs (thin black curves) from the posterior of models disfavoured by the MWA upper limits from Trott et al. (2020) against observations of unlensed UV LFs at z = 6−8 (orange circles; Bouwens et al. 2015, 2017) and at z = 10 (pink squares; Oesch et al. 2018). This shows that the vast majority of UV LFs disfavoured by the current MWA limits are inconsistent with existing UV galaxy limits by several orders of magnitude. In fact, most of these would be strongly ruled out by existing UV LF observations owing to the relatively small observational uncertainties. Nevertheless, there are still some model UV LFs disfavoured by the current MWA limits consistent with these existing observations, again highlighting that the Trott et al. (2020) limits are providing additional (albeit very weak) constraining power.

It is at this point, we emphasize the importance of using 21cmfast in our Bayesian forward modelling approach. The in-built parametrization of the ionizing sources is able to directly output UV LFs which are capable of being compared to existing observational data. Additionally, this further highlights the synergy between observations of the cosmic 21-cm signal and galaxy UV LFs. Existing galaxy observations already place relatively tight constraints on the UV LF bright end, whereas observing the cosmic 21-cm signal can provide limits on both the very faint end of the underlying UV galaxy LFs and to much higher redshifts. However, we note that interpreting results from the 21-cm signal into galaxy UV LFs is entirely model-dependent.

4.3.3 Electron scattering optical depth, τ_e

In Fig. 3, we compare the electron scattering optical depth, τ_e, from the models disfavoured by the latest MWA upper limits (Trott et al. 2020) to the latest constraints measured by Planck (τ_e = 0.054 ± 0.007; Planck Collaboration VI 2018). The solid vertical line denotes the mean value from Planck, with dashed vertical lines (and shaded regions) denoting ±1σ, ±2σ, and ±3σ from the mean value. The red curve represents a histogram of the τ_e calculated from all models in excess of at least one of the current MWA upper limits.

Figure 3.

Histogram (red curve) of τ_e from all astrophysical models found to be in excess of at least one of the latest MWA upper limits on the 21-cm PS (Trott et al. 2020). Shaded bands correspond to the statistical uncertainty on τ_e as measured by Planck (τ_e = 0.054 ± 0.007; Planck Collaboration VI 2018).

The vast majority of the models disfavoured by the latest MWA upper limits are inconsistent with existing observations at ≳ 2σ. These models prefer very low τ_e values, consistent with the preference for reionization occurring as late as possible to exceed the MWA upper limits at z = 6.5 (as discussed in Section 4.2). It is these models that produce the upper envelope of the shaded region for the disfavoured reionization histories shown in Fig. 2. Interestingly, the τ_e histogram exhibits an extremely long tail toward larger τ_e, crossing the observational constraints from Planck. Again, this implies there are reionization models consistent with existing observational constraints that are already being disfavoured by the existing MWA upper limits (albeit extremely weakly).

4.3.4 Mean UV photoionization rate, |$\bar{\Gamma }_{\rm UVB}$|

Finally, in Fig. 4, we compare against the observational constraints on the mean UV background photoionization rate (red data points) extracted from the proximity zones of z > 6 QSOs (e.g. Calverley et al. 2011; Wyithe & Bolton 2011). Again, we present 500 models (thin black curves) randomly drawn from the posterior of models disfavoured by the latest MWA upper limits (Trott et al. 2020). Once again, the vast majority of these disfavoured models are already ruled out by existing observational constraints, as highlighted by the large scatter (of several orders of magnitude) in the mean UV background photoionization rate obtained from the 21cmfast simulations. Further, the disfavoured models tend to lie on average below the existing observational constraints, implying a reduced output of ionizing photons. This is consistent with what we have established from the previous sections. That being, in order to be able to exceed the current MWA upper limits, reionization is preferred to occur at lower redshifts (or not at all), resulting in a smaller mean photoionization background due to the lower numbers of ionizing photons being produced.

Figure 4.

A comparison of the mean UV background radiation (⁠|$\bar{\Gamma }_{\rm UVB}$|⁠) from a random sample of 500 astrophysical models (black curves) found to be in excess of at least one of the recent MWA upper limits on the 21-cm PS (Trott et al. 2020) against observed constraints from the proximity zones of high-z quasars (Calverley et al. 2011; Wyithe & Bolton 2011).

4.3.5 X-ray emissivity, ϵ_{X, 0.5 − 2 keV}

Previously, we have shown that the models disfavoured by the existing MWA upper limits correspond to low X-ray luminosities (i.e. |${\rm log_{10}}\left(L_{{\rm X}\lt 2\, {\rm keV}}/{\rm SFR}\right) \lesssim 37.8$| at 95 per cent confidence). This implies that the neutral IGM is cold, regardless of whether the fluctuations are driven by the patchy EoR or the underlying matter density field (blue and red curves in Fig. 1, respectively). As a result, we can place limits on the soft-band X-ray emissivity (ϵ_{X, 0.5−2 keV}) of high-redshift galaxies (or indeed any source of energy injection into the neutral IGM). At 1σ, we obtain lower limits on the soft band X-ray emissivity from the latest MWA upper limits of log₁₀(ϵ_{X, 0.5−2 keV}) = 34.72, 34.67, 34.61, 34.46, 34.38, and 34.27 at z = 6.5, 6.8, 7.1, 7.8, 8.2, and 8.7. Note that as these limits are drawn from models that are disfavoured by the current MWA data, thus these limits are conditional on the resultant IGM neutral fraction.

In Fig. 5, we compare the resultant lower limits on the soft band (0.5–2 keV) X-ray emissivity implied by the latest MWA upper limits against existing observational constraints. The observational constraints were obtained from large samples of low redshift (z ≤ 1) galaxies obtained with Chandra; Tzanavaris & Georgantopoulos (2008) (red squares) and Lehmer et al. (2016) (black circles). All curves are obtained from Lehmer et al. (2016). The solid curve corresponds to the theoretical best-estimates for the evolution of the total X-ray emissivity, with contributions from low-mass X-ray binaries (LMXBs), high-mass X-ray binaries (HMXBs), and a hot interstellar medium. The corresponding 1σ shaded region accounts for uncertainties in the measurements of the SFR densities and stellar mass densities (e.g. Madau & Dickinson 2014) as well as uncertainties in the XRB SED due to absorption. The blue-dashed curve corresponds to the HMXB component, which is expected to dominate at high redshifts.

Figure 5.

A comparison of the soft band (0.5–2 keV) X-ray emissivity, ϵ_{X, 0.5−2 keV}, for high redshift galaxies. The lower limits at z = 6.5−8.7 correspond to the 1σ limits on the X-ray emissivity inferred from the recent MWA upper limits on the 21-cm PS (Trott et al. 2020). Observational data is obtained with Chandra; Tzanavaris & Georgantopoulos (2008) (red squares) and Lehmer et al. (2016) (black circles). All curves are obtained from (Lehmer et al. 2016). The solid curves correspond to theoretical best-estimates for the evolution of the total X-ray emissivity, with contributions from LMXBs, HMXBs and a hot interstellar medium. The blue dashed curve corresponds to the HMXB contribution, obtained by scaling Model 269 from (Fragos et al. 2013) to estimates of the stellar mass and SFR density.

Although our limits on the soft band X-ray emissivity are still ∼3 orders of magnitude below the expected values, these are the first such limits at these high redshifts. Further, these are obtained from the 21-cm signal, indicating the wealth of information available from the 21-cm signal on the physical properties of the first astrophysical sources. As the lower limits on the 21-cm signal continue to improve, these limits on the X-ray emissivity will approach the theoretical estimates, enabling constraints on the evolution of the X-ray sources to be inferred.

4.3.6 Are any models disfavoured by the MWA consistent with all existing constraints?

In the previous subsections, we have found that it is possible to have models disfavoured by the existing MWA limits to be consistent with individual existing constraints. Here, we explore whether these disfavoured models are consistent with all existing constraints simultaneously. Note, for the purpose of this exploration, we do not consider the soft band X-ray emissivity (Section 4.3.5) as we have already shown the theoretical extrapolations of low-z data out to high-z to be strongly inconsistent with the disfavoured models.

We find that there are a small subset of models disfavoured by the MWA upper limits capable of being simultaneously consistent within ∼2−3σ of (i) the reionization history, (ii) observed UV LFs, (iii) τ_e, and (iv) |$\bar{\Gamma }_{\rm UVB}$|⁠. This implies that already the existing MWA upper limits are capable of providing unique constraining power towards disfavouring models of reionization that would otherwise be consistent with all existing observational constraints on the reionization epoch. Although, currently this additional constraining power is extremely weak. This has important implications in the near future for experiments such as the MWA and LOFAR. As the amplitude of these limits on the 21-cm PS continues to decrease (and for a larger range of observing frequencies) the volume of models disfavoured by these interferometers that are otherwise consistent with existing constraints will continue to increase. Thus, soon we may be able to disfavour more physically reasonable models in addition to the current extreme scenarios.

4.4 Disfavoured IGM properties

After exploring the astrophysical models disfavoured by the latest MWA upper limits (Trott et al. 2020) on the 21-cm signal against existing observations, we now investigate the inferred limits on the globally averaged IGM spin temperature. Note that the IGM spin temperature is not a free parameter of the model, rather it is self-consistently calculated within each simulation voxel within 21cmfast (see e.g. Mesinger et al. 2011).

In Fig. 6, we present the 2D marginalized posteriors for the IGM spin temperature, |$\bar{T}_{\rm S}$| and the IGM neutral fraction, |$\bar{x}_{\mathrm{H\, {\scriptscriptstyle I}}{}}$| for each of the six redshift limits presented in Trott et al. (2020). These 2D posteriors are constructed after marginalizing the output simulation data by the full posterior of astrophysical model parameters. Dashed (solid) contours correspond to the 68th (95th) percentiles. The vertical dot–dashed line corresponds to the value of the adiabatically cooled neutral IGM calculated at mean density, obtained from (recfast; Seager, Sasselov & Scott 1999, 2000).

Figure 6.

2D marginalized posteriors for the IGM neutral fraction, |$\bar{x}_{\mathrm{H\, {\scriptscriptstyle I}}{}}$|⁠, and the IGM spin temperature, |$\bar{T}_{\rm S}$| for the six different MWA redshift bins presented by Trott et al. (2020). Note, these are obtained from our likelihood which combines all six redshift bins simultaneously. Dashed (solid) contours correspond to the 68th (95th) percentile limits.

We recover two distinct islands for the |$\bar{T}_{\rm S}-\bar{x}_{\mathrm{H\, {\scriptscriptstyle I}}{}}$| constraints, driven by the two models outlined in Section 4.1. First, we have an island at |$\bar{x}_{\mathrm{H\, {\scriptscriptstyle I}}{}} \gt 0.9$|⁠, which corresponds to fluctuations solely from the matter density field (i.e. no reionization). Secondly, an island whose neutral fraction decreases for decreasing redshift as reionization progresses (the ‘cold’ reionization scenario). Both islands provide disfavoured limits on the IGM spin temperature, which are driven entirely by astrophysical models that produce little to no heating of the IGM (i.e. |${\rm log_{10}}\left(L_{{\rm X}\lt 2\, {\rm keV}}/{\rm SFR}\right) \lesssim 37.8$| at 95 per cent confidence from Section 4.2). In the absence of any heating source, the neutral gas in the IGM adiabatically cools with the expansion of the Universe (represented by the vertical dot–dashed line). Note that in Fig. 6, the disfavoured regions extend below the limit set by recfast, however, this is driven by the fact that 21cmfast computes the IGM spin temperature on-the-fly in a non-uniform IGM. Due to non-linear structure evolution more of the simulation volume is contained within voids, resulting in the volume-averaged IGM spin temperature dropping below the limit set for the neutral IGM at mean density (the horizontal dot-dashed line).

In order to recover disfavoured limits on just the IGM spin temperature, we marginalize our 2D posteriors over the IGM neutral fraction. In Fig. 7, we present the disfavoured limits on the IGM spin temperature for all six redshifts presented in Trott et al. (2020). Here, the black (red) arrows denote the 68th (95th) percentiles of the disfavoured limits on the IGM spin temperature from the latest MWA upper limits and the black dot–dashed curve corresponds to the value for a neutral IGM at mean density obtained from recfast.

Figure 7.

Marginalized 1D posteriors on the disfavoured limits on the IGM spin temperature, |$\bar{T}_{\rm S}$| for the six different MWA redshift bins presented in Trott et al. (2020). These are obtained from our likelihood, which combines all six redshift bins simultaneously. Black (red) arrows denote the 68th (95th) percentiles on the disfavoured values of the IGM spin temperature, |$\bar{T}_{\rm S}$|⁠. The black dot–dashed line corresponds to the value for the neutral IGM at mean density obtained from recfast.

For the six redshift bins with upper limits on the 21-cm signal provided in Trott et al. (2020), we recover disfavoured limits on the IGM spin temperature |$\bar{T}_{\rm S}\lesssim$| 1.12 (1.32), 1.20 (1.38), 1.28 (1.48), 1.54 (1.82), 1.70 (2.09), 1.90 (2.43) K at 68th (95th) per cent confidence. As the disfavoured limits extend beyond the adiabatically cooled value for a neutral IGM (dot–dashed curve), this implies that the IGM must have undergone a small amount of heating by X-rays. This is broadly consistent with the interpretation of the recent upper limits from LOFAR at z ≈ 9.1, which also prefer a small amount of X-ray heating (disfavoured limits of |$\bar{T}_{\rm S} \lesssim 2.6$| for Greig et al. 2020, |$\bar{T}_{\rm S} \lesssim 2.9$| for Ghara et al. 2020 and |$\bar{T}_{\rm S} \lesssim 10.1$| for Mondal et al. 2020 at 95 per cent confidence). Note, however, that the limits at higher redshifts (i.e. z > 7.5) are completely driven by the 21-cm PS upper limits at lower redshifts (i.e. z = 6.5−7.1), where the 21-cm upper limits are lower in amplitude. This explains the broadening disfavoured regions out to higher redshift, where the variation in the redshift evolution of the IGM spin temperature increases with increasing distance from the better constraining lower redshift upper limits.

5 CONCLUSION

The MWA recently published deep, multiredshift upper limits at z = 6.5−8.7 on the 21-cm PS in the reionization epoch (Trott et al. 2020). These were obtained from 298 h of carefully excised data from four observing seasons, to produce the best upper limits on the 21-cm signal to date at z < 7.5. At present, these upper limits are still too large to begin to rule out regions of astrophysical parameter space. Instead, following a similar approach to that of the recent analyses of the LOFAR upper limits at z ≈ 9.1 (Ghara et al. 2020; Greig et al. 2020; Mondal et al. 2020), we explore regions of astrophysical parameter space that are inconsistent with the observational data. We then extend this further to explore how these disfavoured astrophysical models compare against existing observational constraints on the reionization epoch. We perform this analysis by directly forward modelling the 3D cosmic 21-cm signal using 21CMMC, an MCMC sampler of 3D reionization simulations.

We find two classes of astrophysical models disfavoured by the recent MWA upper limits. These are (i) ‘cold’ reionization models, whereby reionization proceeds in a cold IGM owing to the lack of a heating source (e.g. no X-ray heating) or (ii) pure matter density fluctuations (i.e. no reionization). This is the first work to disfavour a signal driven solely by matter density fluctuations.

With respect to the astrophysical parameters, we find that the latest MWA upper limits primarily restrict the soft-band X-ray luminosity of the first galaxies. At 95 per cent confidence, we recover disfavoured limits of |${\rm log_{10}}\left(L_{{\rm X}\lt 2\, {\rm keV}}/{\rm SFR}\right) \lesssim 37.8$|⁠. These limits sit below our current expectations from observations of analogue low-redshift star-forming galaxies (Mineo et al. 2012), stacked Chandra observations (Lehmer et al. 2016), and predictions at high-redshift from population synthesis models (Fragos et al. 2013).

In terms of galaxy UV properties, the strongest disfavoured limits are for low values of the parameters controlling the normalization of the halo mass-dependent power laws for the fraction of galactic gas in stars, f_{*, 10} and the ionizing escape fraction, f_{esc, 10} and also the star formation time-scale, t_*. In order to exceed the existing MWA upper limits, the model 21-cm PS needs to be near its peak, which qualitatively occurs roughly around the mid-point of reionization. Since the lowest amplitude upper limits are achieved at (i.e. z = 6.5), reionization must be at its mid-point near z = 6.5, which can only be achieved by minimizing the number of ionizing photons available for reionization (i.e. low f_{*, 10} and f_{esc, 10}). For t_*, the limits arise due to the degeneracy with the X-ray luminosity. The number density of X-ray photons is inversely proportional to the star formation time-scale, thus in-order to minimize the amount of X-ray heating to produce large-amplitude 21-cm signals, we strongly disfavour short star formation time-scales.

Next, we extended our exploration to compare the astrophysical models disfavoured by the latest MWA upper limits to existing observation constraints on the reionization epoch. We compared against (i) a census of constraints and limits on the IGM neutral fraction, (ii) observed UV galaxy LFs, (iii) the electron scattering optical depth, (iv) the mean UV background photoionization rate, and (v) the soft-band X-ray emissivity. For all, we found that the vast majority of astrophysical models disfavoured by the existing MWA upper limits were already inconsistent with existing constraints. However, we found a small sample of models which were consistent with existing constraints. This implies that the MWA is already bringing unique constraining information to the astrophysics of reionization, albeit extremely weakly. Using these latest MWA upper limits on the 21-cm signal, we were able to infer the first-ever limits on the X-ray properties of galaxies at high redshifts. For the soft-band X-ray emissivity, conditional on the IGM neutral fraction, we recover 1σ lower limits of log₁₀(ϵ_{X, 0.5 − 2 keV}) = 34.7, 34.7, 34.6, 34.5, 34.4, 34.3 at z = 6.5, 6.8, 7.1, 7.8, 8.2, 8.7.

Finally, we explored the conditional limits on the IGM spin temperature, T_S with the IGM neutral fraction. At 95 per cent confidence, we recover disfavoured limits of |$\bar{T}_{\rm S}\lesssim$| 1.3, 1.4, 1.5, 1.8, 2.1, and 2.4 K at z = 6.5, 6.8, 7.1, 7.8, 8.2, and 8.7, respectively. These limits are found to be above the limits placed by a neutral IGM at mean density undergoing purely adiabatic expansion, implying that the IGM must have undergone some level of X-ray heating. This picture is consistent with that found by Ghara et al. (2020) and Greig et al. (2020), for the recent LOFAR upper limits.

This exploration showcases the value of tools such as 21CMMC, which forward model the cosmic 21-cm signal in a fully Bayesian framework. In doing so, we are able to infer information about the astrophysics of reionization from observations of the 21-cm signal. However, note that the astrophysical interpretations within this work are specific to our astrophysical parametrization and underlying model assumptions. In the near future, as the upper limits on the 21-cm signal continue to improve, we will soon be able to begin to rule out currently viable regions of astrophysical parameter space. Specifically, the MWA in the near term should focus on reducing the amplitude of the 21-cm power at the largest scales, (i.e. k = 0.14 h Mpc) and lowest redshifts (z ∼ 6.5), to maximize the total disfavoured parameter volume.

ACKNOWLEDGEMENTS

We thank the anonymous referee for their feedback which has improved this manuscript. We would like to that thank Andrei Mesinger for comments on an early version of this draft and Bret Lehmer for providing the X-ray emissivity data. Parts of this research were supported by the Australian Research Council Centre of Excellence for All Sky Astrophysics in 3 Dimensions (ASTRO 3D), through project number CE170100013. CMT is supported by an ARC Future Fellowship under grant FT180100321. The International Centre for Radio Astronomy Research (ICRAR) is a Joint Venture of Curtin University and The University of Western Australia, funded by the Western Australian State government. Parts of this work were performed on the OzSTAR national facility at Swinburne University of Technology. OzSTAR is funded by Swinburne University of Technology.

DATA AVAILABILITY

The data underlying this paper will be shared on a reasonable request to the corresponding author.

Footnotes

1

https://www.isispace.nl/projects/ncle-the-netherlands-china-low-frequency-explorer/

2

https://github.com/BradGreig/21CMMC

3

https://github.com/andreimesinger/21cmFAST

4

A power-law dependence between M_* and M_h at z ≳ 5 is consistent with the mean behaviour of both semi-analytic model predictions (e.g Mutch et al. 2016; Yung et al. 2019; Hutter et al. 2020) and semi-empirical fits to observations (e.g. Harikane et al. 2016; Tacchella et al. 2018; Behroozi et al. 2019).

5

Throughout this work, we adopt the Sheth–Tormen HMF (Sheth & Tormen 1999).

6

We assume N_γ/b = 5000, corresponding to a Salpeter initial mass function (Salpeter 1955); however, note that this is highly degenerate with the fraction of galactic gas in stars, f_*.

7

Within 21cmfast the UV spectrum of the stellar emission component is fixed. Also, it does not consider alternative sources of soft UV photons such as quasars (see e.g. Qin et al. 2017; Ricci et al. 2017; Mitra, Choudhury & Ferrara 2018).

8

This assumption only breaks down for the final ∼10 per cent of the EoR (Sobacchi & Mesinger 2014).

9

This modelling uncertainty is motivated by approximations adopted in seminumerical simulations relative to radiative-transfer simulations (e.g. Zahn et al. 2011; Ghara et al. 2018; Hutter 2018).

10

In some instances, depending on the ionizing source parameters (e.g. M_turn) there may only be a single peak. Whereby the WF-coupling peak combined with the strong heating fluctuations from the cold IGM can produce a peak in the 21-cm PS capable of exceeding the existing MWA upper limits. In these models, it is also possible to have reionization commencing.

11

Note this choice of 90 per cent is purely arbitrary, and is used solely to denote whether or not reionization has commenced in these models. This choice is not used for delineating whether or not these models are consistent or inconsistent with existing limits on the IGM neutral fraction.

12

For reference the maximum likelihoods from these two classes of models are |$-{\rm log}(\mathcal {L}) = -22.2$| and −22.6, respectively.

13

These extreme models occupy only a very small sub-volume of the total allowed parameter volume, thus they are not significant enough to result in actual constrains once we marginalize over the full parameter space

14

However, note that this dependence between the peak of the 21-cm PS and |$\bar{x}_{\mathrm{H\, {\scriptscriptstyle I}}{}}$| is strongly model-dependent.

15

Again, this degeneracy only exists in the absence of other constraints on t_*, such as from observed UV galaxy LFs (e.g. Park et al. 2019).

APPENDIX A: PROBABILITY TO BE EXCLUDED BY THE MWA UPPER LIMITS

In this work, the probability for a model 21-cm PS to be disfavoured by the MWA upper limits (e.g. p_ex.) is based directly on the data output from the MWA data-reduction pipeline (e.g. Trott et al. 2020). Specifically, we refer to the output PS from the CHIPS (the Cosmological HI PS estimator) pipeline (Trott et al. 2016). This, takes as input the calibrated observational data cube (in terms of the angular u, v visibilities, and frequency, f) and returns an estimate of the 3D 21-cm power (i.e. in k_x, k_y, k_z).

We then construct a histogram for the 3D 21-cm power for each observing frequency (redshift) and spherically averaged Fourier k-bin. In Fig. A1, we provide an example of this at k = 0.14 h Mpc⁻¹ and z = 6.5.

Figure A1.

A histogram of the 3D 21-cm power at k = 0.14 h Mpc⁻¹ and z = 6.5 recovered following the data reduction pipeline for the MWA (e.g. Trott et al. 2020).

Next, from this histogram, we calculate our p_ex., the probability for the modelled 21-cm PS amplitude (i.e. |$\Delta ^{2}_{21,{\rm mod}}$|⁠) to be disfavoured given the observational data (i.e. |$\Delta ^{2}_{21,{\rm d}}$|⁠). In Fig. A2, we show p_ex. at k = 0.14 h Mpc⁻¹ and z = 6.5. To estimate this probability, for a given |$\Delta ^{2}_{21,{\rm mod}}$|⁠, we calculate the probability for a randomly drawn observation (⁠|$\Delta ^{2}_{21,{\rm d}}$|⁠) to have an amplitude less than the modelled value. This probability is estimated by calculating the cumulative probability distribution. Thus, the probability to be in excess of the existing MWA upper limits approaches unity for increasing model 21-cm PS amplitudes

Figure A2.

The probability, |$p_{\rm ex.}(\Delta ^{2}_{21})$| that a model 21-cm PS amplitude (⁠|$\Delta ^{2}_{21}$|⁠) is disfavoured by the existing MWA upper limits at k = 0.14 h Mpc⁻¹ and z = 6.5. Note, for each redshift and k-bin we construct a new p_ex..