Fisher matrix forecasts on the astrophysics of galaxies during the epoch of reionisation from the 21-cm power spectra

We use the L210_N4320 dark matter-only N-body simulation of the Genesis suite of simulations (Power et al. in prep). Containing $4320^{3}$ particles in a $210^{3}~{}h^{-3}$ Mpc³ volume, the simulation has a mass resolution of $\sim 5\times 10^{8}~{}h^{-1}~{}\rm M_{\odot}$. L210_N4320 was run with the SWIFT (Schaller et al., 2018) cosmological code, using the VELOCIraptor (Elahi et al., 2019a) halo identifier and merger trees generated using TreeFrog (Elahi et al., 2019b). The simulation consists of 120 snapshots between redshifts 30 and 5.

In order to resolve haloes down to the atomic cooling limit at $z=20$, the L210_N4320 simulation was augmented using the Monte-Carlo algorithm-based code DarkForest  (Qiu et al., 2020) achieving an effective halo mass resolution of $\sim 2\times 10^{7}~{}h^{-1}~{}\rm M_{\odot}$ (L210_AUG of Balu et al. 2023). DarkForest achieves this by sampling from a conditional mass function based on the extended Press-Schechter theory (Lacey & Cole, 1993; Bond et al., 1991; Bower, 1991) that has been modified to match the halo mass functions from N-body simulations. These new low-mass haloes are introduced into the simulation volume and appended to the merger trees. Importantly, DarkForest also assigns positions to these newly added haloes based on the local halo density field (Ahn et al., 2015) and the linear halo bias (Tinker et al., 2010). For the rest of this work, we use L210_AUG  as our fiducial simulation.

Meraxes (Mutch et al., 2016) contains detailed and physically motivated prescriptions for galaxy formation and evolution. These include, for example, gas infall into dark matter haloes from the IGM followed by its radiative cooling and subsequent star formation, as well as eventual supernova feedback and metal enrichment of the ISM. Active galactic Nuclei (AGN) feedback from central black holes of the galaxies was implemented in Qin et al. (2017) and calculations of the galaxies’ UV luminosity in Qiu et al. (2019). In addition, Meraxes also has a coupled treatment of reionisation based on the 21cmFAST semi-numerical code (Mesinger et al., 2011; Murray et al., 2020, see section 2.3 for more details). In the following sections, we give detailed but non-exhaustive descriptions of the various physics implementations pertinent to our astrophysical model. Table 1 lists the free parameters of our model along with their adopted fiducial values.

The thermal evolution and ionisation state of the IGM follows 21cmFAST (Mesinger et al., 2011; Murray et al., 2020), modified to take advantage of our realistic galaxy population. Meraxes sources the dark matter density and velocity fields from the underlying N-body simulation. We begin by gridding our simulation volume and assigning the galaxies to voxels based on their positions. In the present work, motivated by the typical H ii  bubble sizes during the EoR (Furlanetto et al. 2004, Wyithe & Loeb 2004), we divide our simulation into $1024^{3}$ voxels corresponding to a cell dimension of $\sim 0.21~{}h^{-1}$ Mpc. We give a brief summary of our method in the following subsections.

A radio telescope measures the brightness temperature ($\delta T_{b}$) of a H I cloud, which is the offset of the spin temperature ($T_{\rm S}$) against a background radiation source (assumed to be the CMB) of temperature $T_{\gamma}$ as a function of frequency $\nu$, given by

where $\tau_{\nu_{0}}$ is the optical depth at the 21-cm transition frequency $\nu_{0}$, $x_{\mathrm{H\,{\scriptscriptstyle I}}}$ is the neutral hydrogen fraction, $\delta_{nl}$ is the fluctuations in the underlying dark matter density field, $H$ is the Hubble parameter, and $dv_{r}/dr$ is the line-of-sight component of the radial derivative of the peculiar velocity. All of the terms in equation 18 are evaluated at $z=\nu_{0}/\nu-1$.

The L210_AUG  simulation from Balu et al. (2023) was calibrated against existing observables including UV LFs from Bouwens et al. (2015) and Bouwens et al. (2021), and stellar mass functions from Song et al. (2016) and Stefanon et al. (2021). The reionisation history was calibrated using the Thomson scattering optical depth of free electrons to CMB photons from the Planck Collaboration et al. (2020) (see Fig. 3 & 4 of Balu et al., 2023).

The fiducial values of our 8 astrophysical model parameters are given in the fourth column of Table 1. These parameters can be divided into four groups of two in the order they are shown in the table, having direct control on the X-ray properties ($L_{\rm X}/\rm{SFR}$ & $E_{0}$), the ionising UV photon escape fraction ($f_{\rm esc,0}$ & $\alpha_{\rm esc}$), star formation ($\Sigma_{\rm SF}$ & $\alpha_{\rm SF}$), and the supernova feedback ($\varepsilon_{0}$ & $\eta_{0}$) of the galaxies.

The frequency dependence of the 21-cm signal imparts a line-of-sight evolution to the 21-signal. To account for this, Meraxes produces a 21-cm light-cone by stitching together (linearly) interpolated $\delta T_{b}$ co-evolution grids. Following Greig & Mesinger (2018), we subdivide our 21-cm light-cone into cubical grids of a size equivalent to our simulation volume ($1024^{3}$) and calculate the spherically averaged dimensionless 21-cm PS in each (shown as the solid blue curve in Fig. 1). We obtain 11 redshift chunks spanning $z\in[5,24]$, and assign redshift values to them corresponding to the mid-point of each chunk.

The sensitivity of a radio interferometer to the 21-cm PS can be divided into two components: (1) thermal noise, which is important on the small scales and (2) sample (cosmic) variance, prominent on the large scales. The total noise power $\sigma[\Delta_{21}^{2}(k)]$ is given by adding these two in quadrature:

where $\Delta_{\rm N}^{2}(k)$ is the thermal noise given by (Morales, 2005; McQuinn et al., 2006; Parsons et al., 2012):

where $X$ and $Y$ relate the bandwidths and solid angles to the comoving distance to the source, $\Omega^{{}^{\prime}}$ is a factor dependant on the beam of the telescope (see Parsons et al., 2014), $t$ is the integration time for the mode $k$, and $T_{\rm sys}$ is the system temperature given by $T_{\rm sys}=T_{\rm sky}+T_{\rm rec}$ with $T_{\rm sky}$ and $T_{\rm rec}$ being the sky and receiver temperature respectively. The quadrature addition and the form of the cosmic variance noise in equation 19 assume that the errors are Gaussian distributed which is reasonable for the relevant scales in this work (see Qin et al., 2021a, b; Prelogović & Mesinger, 2023).

In this work, we focus on a future observation of the 21-cm PS by the SKA. We limit our attention to the upcoming first phase of SKA – the so-called SKA1-low¹⁰¹⁰10See the official SKA1 System Baseline Design document for further details.. We only include the stations in the ‘Central Area’ of the SKA1-low, resulting in 296 stations of diameter 35 m distributed across a circular area with 1.7 km diameter, we calculate the interferometer’s sensitivity to the 21-cm PS using the 21cmsense ¹¹¹¹11www.github.com/steven-murray/21cmSense python package (Pober et al., 2013, 2014). We assume an observational campaign corresponding to 6 hours per night for 180 days, i.e. a total of $1080$ hours and sky temperature $T_{\rm sky}$ to be dominated by galactic synchrotron emission (Thompson et al., 2017) scaling with frequency $\nu$ as $T_{\rm sky}=60(\nu/300\rm{MHz})^{-2.55}\rm{K}$. Additionally, assuming that the telescope reflects $10$ per cent of the response to the sky (Pober et al., 2014), we set $T_{\rm rec}=40\rm{mK}+0.1T_{\rm sky}$. We combine partially coherent baselines to improve power spectrum sensitivity and use the ‘moderate’ foreground removal scenario of 21cmsense wherein we avoid the modes that are contained within the so-called foreground wedge (Datta et al., 2010) extending $k_{\parallel}=0.1h$ Mpc^-1 beyond the horizon limit. As we ignore the modes that fall within the foreground wedge, we are limited to $k_{\rm min}=0.16$ Mpc^-1 for our analysis. We also have an upper limit of $k_{\rm max}=1.4$ Mpc^-1 which arises from a combination of both the spatial scales that are probed by the SKA1-low (set by the longest baseline we consider in our modelling of the array), as well as Poisson, shot noise from our simulation. Figure 1 shows the PS noise (grey shaded region) that is generated for our fiducial 21-cm PS. The red-shaded region shows $k$-modes that are ignored in this work.

In addition to exploring the possible constraints on our fiducial model from the 21-cm PS we also consider the improvements that are achievable with a joint analysis of both the 21-cm PS and the UV LFs.

In addition to updates on the SNe feedback (see section 2.2.2), Qiu et al. (2019) explored different implementations of dust attenuation in Meraxes . Here, we adopt the parameterisation for the UV optical depth that depends on the dust-to-gas ratio (Charlot & Fall, 2000) of the galaxies. The model differentiates between a short-lived birth cloud of the stars as well as the interstellar medium (ISM) of the galaxy. Photons are absorbed by both the birth cloud (albeit for a short while until it gets depleted by star formation) and the ISM. The free parameters of Meraxes – particularly the ones impacting the galaxy physics – have been calibrated (see Qiu et al., 2019; Balu et al., 2023) to their fiducial values against infrared excess (IRX)-$\beta$, UV LFs, and stellar mass functions at $z>5$.

The blue curves in Fig. 2 are the UV LFs from our calibrated simulation. The error bars are determined by multiplying our simulated UV LF data by the corresponding fractional uncertainty on each data point from the Bouwens et al. (2021) observational data. In addition to the 21-cm PS, we thus have 6 UV LFs from $z\in[5,10]$.

In Fig. 3 we show the forecasts on our astrophysical parameter set from the 21-cm PS corresponding to a $\sim 1000$ h observation with SKA1-low. The dark (light) contours correspond to the 1 (2)-$\sigma$ 2-D marginalised confidence intervals for each astrophysical parameter and the diagonal subplots show the 1-D normalised marginal probability distribution function. The fifth column of table 1 lists the 1-$\sigma$ uncertainty on our parameters. The fractional uncertainties, $|\sigma(\theta_{i})/\theta^{fid}_{i}|$, are given within brackets and are shown as blue squares in Fig. 5.

We obtain tight constraints ($\lesssim 10$ per cent) for the X-ray parameters (luminosity $L_{\rm X}/\rm{SFR}$ & minimum X-ray photon energy escaping galaxies $E_{0}$), escape fraction normalisation ($f_{\rm esc,0}$ ) and star formation efficiency ($\alpha_{\rm SF}$) while the SNe ejection efficiency ($\epsilon_{0}$) is constrained to $\sim 25$ per cent. On the other hand, the critical mass normalisation ($\Sigma_{\rm SF}$) and the efficiency of SNe reheating ($\eta_{0}$) remain relatively unconstrained with $\sim 234$ per cent and $\sim 117$ per cent fractional $1\sigma$ uncertainties respectively. The relatively poor constraints on $\Sigma_{\rm SF}$ and $\eta_{0}$ are primarily because these parameters have negligible impact on the 21-cm PS.

As previously mentioned, our parameter set forms four groups of two with each group controlling a different aspect of galaxy properties.

1. 1.

   X-ray parameters: The 21-cm signal is very sensitive to the amount and energy of the X-ray photons in the early Universe. The 21-cm PS, therefore, provides tight constraints on the X-ray parameters, $E_{0}$ ($\sim 11$ per cent) and $L_{\rm X}/\rm{SFR}$ ($\sim 10^{-2}$ per cent), in our model. These forecasts are consistent with similar works in the literature (e.g. Greig & Mesinger, 2017; Park et al., 2019; Mason et al., 2022), which is unsurprising given they have the same X-ray implementation despite the differences in the source modelling (i.e. SAM).

2. 2.

   UV escape fraction: $f_{\rm esc}$ is more sensitive to the normalisation $f_{\rm esc,0}$ as opposed to the redshift power-law exponent ($\alpha_{\rm esc}$). This is also reflected in the derivative of the 21-cm PS with respect to these parameters (green and pink curves in Fig. 6). We thus forecast tighter constraints on $f_{\rm esc,0}$ (6.55 per cent) compared to $\alpha_{\rm esc}$ (66.93 per cent).

3. 3.

   Star formation: Among the star formation parameters, the efficiency of star formation ($\alpha_{\rm SF}$) is highly constrained relative to the critical mass normalisation ($\Sigma_{\rm SF}$). The primary impact of $\Sigma_{\rm SF}$ is to establish a critical mass threshold needed for the stars to start forming. This parameter thus has only a secondary role in the stellar mass content of galaxies as once a galaxy has accreted enough mass, the amount of stellar material is controlled by $\alpha_{\rm SF}$ (see section 2.2.1). Equation 2 also shows that these two parameters enter our model as a multiplicative term of the form $(1-\Sigma_{\rm SF})\alpha_{\rm SF}$ adding to their degeneracy.

4. 4.

   Supernovae feedback: Compared to the star formation parameters, we find that the SNe parameters are less constrained, with $\sim 25$ per cent for the SNe ejection efficiency $\epsilon_{0}$ while the efficiency of reheating by the SNe $\eta_{0}$ is only constrained to $\sim 117$ per cent. Similar to $\Sigma_{\rm SF}$, the low constraints on these parameters are due to the fiducial model values being such that the 21-cm signal is relatively insensitive to these processes.

We next add in the Fisher matrix corresponding to the UV LFs from 6 redshifts $\in[5,10]$. As the UV LFs functions are sensitive to parameters that do not have a pronounced impact on the ionisation morphology, the addition of the UV LFs into the analysis helps improve the overall constraints on the astrophysical model (Park et al., 2019). We show the impact of varying the two least constrained parameters when using just the 21-cm PS, $\Sigma_{\rm SF}$ and $\eta_{0}$, on the UV LFs in Fig. 2. The green and pink shaded regions are from varying $\eta_{0}$ and $\Sigma_{\rm SF}$ respectively, by an order of magnitude about their fiducial values (we point out that these two parameters are positive by definition). The variation in the UV LFs is larger than the 1-$\sigma$ errors emphasising that including the UV LFs in our Fisher analysis will provide additional constraining power. As previously mentioned, we do not vary the X-ray parameters, $L_{\rm X}/\rm{SFR}$ & $E_{0}$, while computing the UV LF Fisher matrix.

Figure 4 shows the constraints from the joint analysis of the 21-cm PS and UV LFs, and the last column of Table 1 lists the forecasted 1-$\sigma$ uncertainty along with the fractional uncertainties for our model parameters. The green stars in Fig. 5 are the fractional uncertainties for the joint analysis of both the 21-cm PS and the UV LFs, and the orange crosses represent the same from just the UV LFs¹⁶¹⁶16In Fig. 8 we show the constraints from the UV LFs alone.. Figures 3 & 5 demonstrate that combining the UV LFs into the analyses helps to improve the constraints. This is most notable for parameters having a direct impact on the stellar mass of the galaxies i.e. the ones related to star formation and supernova feedback as UV LFs are more sensitive to these parameters compared to the 21-cm PS. This results in significant improvements in $\Sigma_{\rm SF}$ (from $\sim 234$ to $\sim 40$ per cent) and $\eta_{0}$ (from $\sim 117$ to $\sim 8$ per cent) for the joint case of both the 21-cm PS and the UV LFs as the degeneracies between stellar properties and escape fraction can be weakened. This then translates into improvements in other parameters.

The relatively larger improvement in $\Sigma_{\rm SF}$ compared to $\alpha_{\rm SF}$ following the inclusion of information from the UV LFs stems from the breaking of parameter degeneracies. In Fig. 3 (21-cm PS only), $\Sigma_{\rm SF}$ is poorly constrained (highlighted by the large uncertainties), and $\Sigma_{\rm SF}$ and $\alpha_{\rm SF}$ are anti-correlated. In Fig. 8 (UV LFs only), $\Sigma_{\rm SF}$ is now more tightly constrained and positively correlated with $\alpha_{\rm SF}$, indicating the majority of $\Sigma_{\rm SF}$’s constraining power comes from the UV LF. Thus, the inclusion of UV LFs provide considerable improvements to the constraining power of $\Sigma_{\rm SF}$. For $\alpha_{\rm SF}$ on the other hand, it is almost equally well constrained both by the 21-cm PS (Fig. 3) and the UV LFs (Fig 8). Thus combining the two observables provides modest improvements (dominated by the change in correlation).

Even though Fisher matrices implicitly make simplifying assumptions, such as Gaussian likelihoods and Gaussian errors, they nevertheless give useful insights. In this section, we explore the correlations and degeneracies among the model parameters. For brevity, we focus only on the results from the joint analysis of the UV LFs and the 21-cm PS (Fig 4), noting that the trends are similar for the 21-cm PS alone (Fig. 3).

1. 1.

   X-ray parameters: $L_{\rm X}/\rm{SFR}$ anti-correlates with the UV escape fraction normalisation ($f_{\rm esc,0}$ ) and the SNe ejection efficiency ($\epsilon_{0}$). $E_{0}$ shows an anti-correlation with the star formation efficiency ($\alpha_{\rm SF}$) and positively correlates with the SNe reheat efficiency ($\eta_{0}$).

   These correlations can be understood from the impact of the X-rays on the EoR. X-rays can contribute, albeit a small fraction relative to the UV photons, to the ionisation of the local IGM of a galaxy leading to photoionisation regulation of the amount of neutral gas available for accretion onto the galaxy and hence subsequent star formation. Thus increasing $L_{\rm X}/\rm{SFR}$ can be compensated by a decrease in the SNe feedback and/or the ionising UV escape fraction. $E_{0}$ sets the lowest energy for escaping X-rays from the galaxies. An increase in $E_{0}$, therefore, implies less IGM ionisation by the X-rays.

2. 2.

   UV escape fraction: The correlation among $\alpha_{\rm esc}$ and $f_{\rm esc,0}$ is expected given the definition of the UV escape fraction ($f_{\rm esc}$; see equation 10). The UV escape fraction, $f_{\rm esc}$, directly influences the local IGM ionisation state, and hence the amount of baryonic infall into the galaxy for star formation. The amount of cold gas available for star formation is influenced by the SNe ejection efficiency which sets the amount of gas removed from the galaxy. An increase in $f_{\rm esc}$ therefore should be accompanied by a decrease in the SNe feedback.

   In light of this, the relatively strong positive correlation of $f_{\rm esc,0}$ with the supernova ejection efficiency $\epsilon_{0}$ is interesting. This can be understood from the behaviour of the UV escape fraction (i.e. the combination of $f_{\rm esc,0}$ and $\alpha_{\rm esc}$) and the other physical processes. The positive correlation amongst the SNe parameters ($\eta_{0}$ & $\epsilon_{0}$), as can be seen from equations (3 - 7) and Fig. 8, also plays a contributing factor. $\alpha_{\rm esc}$ is also anti-correlated with both $\eta_{0}$ and $\epsilon_{0}$. It is therefore a combination of the very weak correlation of $f_{\rm esc,0}$ with $\eta_{0}$ and the anti-correlation with $\alpha_{\rm esc}$ that results in its positive correlation with $\epsilon_{0}$.

   The comparatively tighter correlation of the escape fraction parameters with the SNe parameters compared to the X-ray parameters reflects the relative importance of the UV photons over the X-ray photons in the IGM ionisation state.

3. 3.

   Star formation: As noted, we report correlations between the star formation efficiency $\alpha_{\rm SF}$ and the minimum X-ray photon energy $E_{0}$.

   The strong correlation between $\alpha_{\rm SF}$ and SNe reheat efficiency $\eta_{0}$ is expected as they are two of the parameters impacting the stellar mass in a galaxy.

   On the other hand, the correlation between $\alpha_{\rm SF}$ and the critical mass normalisation $\Sigma_{\rm SF}$ is surprisingly weak given equation (2). This is because of the relatively small value of $\Sigma_{\rm SF}$. $\Sigma_{\rm SF}$ fixes the critical mass of the cold gas needed for star formation to commence in a galaxy. The small value of the parameter is justified on the physical grounds that a larger value will result in lower available gas for star formation.

4. 4.

   Supernovae feedback: For completeness, we again point out the correlations between the SNe ejection ($\epsilon_{0}$) and reheat ($\eta_{0}$) efficiencies, the strong correlations of $\epsilon_{0}$ with $f_{\rm esc,0}$ & $E_{0}$, and that of $\eta_{0}$ with $\alpha_{\rm SF}$.