1. Introduction
The model of the Universe’s evolution through the Epoch of Reionisation (EoR) is directly affected by our limited understanding of the first galaxies. Debate of their ionising capabilities persists due to their relatively unconstrained ionising photon production and escape fractions (Robertson et al. Reference Robertson, Ellis, Dunlop, McLure and Stark2010; Madau & Haardt Reference Madau and Haardt2015; Naidu et al. Reference Naidu2019; Finkelstein et al. Reference Finkelstein2019). While AGN-quasars are unlikely ionisation source candidates due to their infrequency beyond
$z\gt 3$
(Kulkarni et al. Reference Kulkarni, Worseck and Hennawi2018), bright, highly star-forming galaxies could have produced copious ionising photons for reionisation (Naidu et al. Reference Naidu2019). Alternatively, numerous faint galaxies could provide sufficient ionising flux over a longer history to reionise the intergalactic medium (Bruton et al. Reference Bruton, Lin, Scarlata and Hayes2023).
Operation of the James Webb Space Telescope (JWST) has, in a short time, revolutionised our understanding of galaxies in the EoR and this landscape is constantly changing with new revelations (Bunker et al. Reference Bunker2023; Bagley et al. Reference Bagley2023; Paris et al. Reference Paris2023). Early empirical studies had noticed potentially large [Oiii]5007 + H
$_{\unicode{x03B2}}$
emissions at higher redshift (Schaerer & de Barros Reference Schaerer and de Barros2009; Raiter et al. Reference Raiter, Fosbury and Teimoorinia2010). Now with recent redshift evolution models (Zhai et al. Reference Zhai, Benson, Wang, Yepes and Chuang2019) and JWST/NIRCam number density studies
$z\gt 5.3$
(Matthee et al. Reference Matthee2023) this is evidently ubiquitous in the early universe. Probing galaxies within the EoR reveals almost no dust attenuation as expected, particularly approaching the highest redshift limits (Robertson et al. Reference Robertson2023; Hsiao et al. Reference Hsiao2023; Tacchella et al. Reference Tacchella2023; Haro et al. Reference Haro2023b). However, galaxies with unexpectedly high stellar masses have been found (Boylan-Kolchin Reference Boylan-Kolchin2022; Labbe et al. Reference Labbe2022) bordering the possible limits set by
$\Lambda$
CDM and suggesting extremely efficient star formation. Exceptional findings such as a highly quenched, dusty galaxy at
$z\sim5$
(Donnan et al. Reference Donnan2022; Harikane et al. Reference Harikane2022; Naidu et al. Reference Naidu2022; Haro et al. Reference Haro2023a; Yung et al. Reference Yung, Somerville, Finkelstein, Wilkins and Gardner2023) and spectroscopic confirmation of a
$z\sim13.2$
(Robertson et al. Reference Robertson2022) galaxy have both challenged and affirmed ideas of the timescale over which the Universe evolves.
The LyC is the component of stellar emission with sufficient energy to ionise neutral hydrogen atoms and is a key parameter in theoretical models of EoR. However, deriving physical properties relating to the production of LyC radiation and their escape from the interstellar medium of galaxies comes with both technical and physical challenges. Direct detection of LyC photons is impeded by the opacity of neutral hydrogen in EoR galaxies, instead requiring models based on lower redshift analogues (Izotov et al. Reference Izotov2018). LyC leakers at low redshift exhibit high [Oiii]5007 EWs and high O
$_{32}$
([Oiii]5007/[Oii]3727) ratios albeit with large scatter (Cardamone et al. Reference Cardamone2009; Izotov et al. Reference Izotov2018). For galaxies
${z}\lt 0.5$
, a Lyman-
${\unicode{x03B1}}$
peak separation of less than 200 km/s is the strongest predictor of significant LyC escape at
${z}\sim0$
(
$f_{esc}\gt 0.1$
up to 0.8) (Izotov et al. Reference Izotov2022). This is the basis of EoR-Lyman-
${\unicode{x03B1}}$
correlated
$f_{esc}$
predictions, with
$\gt 10$
Gyr between the two epochs.
It is possible that these low redshift
$f_{esc}$
probes may have a different relationship with LyC escape at high redshifts due to alternate escape methods. Ji et al. (Reference Ji2020) find no Lyman-
${\unicode{x03B1}}$
emission (only absorption) despite strong LyC detections in a
$z\sim3.8$
galaxy; favouring a model of
$f_{esc}$
through LyC-transparent holes rather than through an optically thin interstellar medium (ISM). The assumption that
$z\lt 0.5$
and
$z\gt 5$
galaxies are comparable is in contention with our understanding of these epochs (Madau & Dickinson Reference Madau and Dickinson2014). Studies using one epoch to analyse another must contend with their fundamental differences and what their comparisons really reveal. For example, the Lyman-
${\unicode{x03B1}}$
probe or UV slope may have a different relative dependence on dust extinction and the stellar population due to morphological differences (Hayes et al. Reference Hayes2011; Meng & Gnedin Reference Meng and Gnedin2020).
With the advent of the JWST, recent observations targeting the H
$_{\unicode{x03B1}}$
and H
$_{\unicode{x03B2}}$
Balmer lines have been used to estimate the production efficiency of hydrogen ionising photons (
$\xi_{ion}$
) for many
$z\gt 6$
sources (Tang et al. Reference Tang2023; Simmonds et al. Reference Simmonds2023). However, with the existing uncertainty in dust content of high redshift sources, emission line corrections may not accurately reproduce the original luminosity as intended. Furthermore, limited spectroscopic availability means statistically significant studies come to rely on photoionisation models (Ferland et al. Reference Ferland1998; Levesque et al. Reference Levesque, Kewley and Larson2010) to estimate
$\xi_{ion}$
(Bouwens et al. Reference Bouwens2015) and compare its correlation to other ISM probes such as [Oiii]5007+H
$_{\unicode{x03B2}}$
EW (Kewley et al. Reference Kewley2015; Tang et al. Reference Tang, Stark, Chevallard and Charlot2019, Reference Tang2023)
The lack of photometric measurements/coverage between the UV and NIR diminishes constraints on the modelled physical properties at
${z}\gt 6$
. Quantities such as metal abundances (Vincenzo & Kobayashi Reference Vincenzo and Kobayashi2018; Torrey et al. Reference Torrey2018; Berg et al. Reference Berg, Erb, Henry, Skillman and McQuinn2019), the escape fraction
$f_{esc}$
or the number of ionising photons
$\dot{n}$
(Anderson et al. Reference Anderson, Governato, Karcher, Quinn and Wadsley2017; Iyer & Gawiser Reference Iyer and Gawiser2017) rely heavily on star formation history (SFH) models to infer the absorbed light, which depends strongly on properly sampling the spectrum.
Having a well sampled SED for galaxies with well constrained redshifts allows for the determination of a photometric UV slope (
${\unicode{x03B2}}_p$
). This probes the empirical emission ‘blueness’ which is used in dust corrections (Calzetti et al. Reference Calzetti, Kinney and Storchi-Bergmann1994; Reddy et al. Reference Reddy2018) and as a proxy for the stellar ages and metal/dust content. Studies of local universe starbursts find that galaxies with a bluer, more negative model
${\unicode{x03B2}}$
(determined from
$f_\lambda$
) have lower dust obscuration, similar to
$z\sim2$
studies (Takeuchi et al. Reference Takeuchi, Yuan, Ikeyama, Murata and Inoue2012; Sklias et al. Reference Sklias2013). The relationship of
${\unicode{x03B2}}$
and dust attenuation is degenerate with metallicity and star formation history (Bouwens et al. Reference Bouwens2016) which complicates the relative contributions of each component. However, correlations to dust content are robust at fixed UV luminosity. Reddy et al. (Reference Reddy2018) and Nanayakkara et al. (Reference Nanayakkara2022) find that galaxies at
$4\lt z\lt 7$
are consistently blue and have very low dust attenuation.
The approach this paper takes to these issues in studying the potential influence of the brightest emitters on the EoR is to use an established analogue sample with the stellar and [Oiii]5007 emission properties of higher redshift galaxies at the slightly lower redshift range of
$2.5\lt{ z}\lt 4$
(11.1–12.2 Gyr ago). The selected redshift range is temporally closer to the EoR (12.7 Gyr ago) than many low redshift comparisons, while still having an accessible LyC region for follow up observational
$f_{esc}$
constraints. The sample is taken from the FourStar Galaxy Evolution Survey (ZFOURGE) (Straatman et al. Reference Straatman2016), which combines deep imaging with medium and narrow-band near-IR filters with multi-wavelength observations from several public surveys to accurately determine the redshifts of
$\sim$
$70\,000$
objects to
$\sigma_z\sim 2\% $
accuracy (Nanayakkara et al. Reference Nanayakkara2016; Tran et al. Reference Tran2020).
We will be exploring a subsample of these galaxies presented by Forrest et al. (Reference Forrest2018) as having extreme H
$_{\unicode{x03B2}}$
+[Oiii]5007 EW (
$\gt$
$400$
) and therefore considered analogues of EoR galaxies (Tang et al. Reference Tang2023). Follow up spectroscopy with the KMOS /VLT as part of the Mutli-Object Spectroscopic Emission Line (MOSEL) survey confirms the photometric selection for a subset of these as EELGs (Tran et al. Reference Tran2020; Gupta et al. Reference Gupta2022). Of 19 galaxies targeted, 16 had bright emission lines, where 14 of these had [Oiii]5007 as the brightest line and two had H
$_{\unicode{x03B1}}$
. This sample (40 filters) will be compared to direct EoR galaxy measurements made by the CEERS (14 filters) (Yang et al. Reference Yang2020; Bagley et al. Reference Bagley2023; Yang et al. Reference Yang2023) and JADES (23 filters) (Bunker et al. Reference Bunker2023; Eisenstein et al. Reference Eisenstein2023; Hainline et al. Reference Hainline2023; Rieke et al. Reference Rieke2023) surveys using legacy HST and the recent JWST observations to determine the overlap in their observed and the SED model-derived properties. For consistency, we use the Multi-wavelength Analysis of Galaxy Physical Properties (MAGPHYS) SED fitting code (da Cunha et al. Reference da Cunha, Charlot and Elbaz2008, Reference da Cunha2015) across all samples.
This work is formatted so that Section 2 describes the data and selection philosophies used to identify Extreme Emission galaxy analogues while Section 3 delves into our methods of analysis and the quantities we use to compare the samples. Results of how representative the analogues are and what their high energy emission properties look like are presented in Section 4. We summarize the work in Section 5.
Throughout this paper we assume a flat universe with
$\Omega_m$
= 0.3,
$\Omega_\lambda$
= 0.7 and
$H_0$
= 70 km/s/Mpc for our models.
2. Data
2.1 The analogues
In this section we break down the selection criteria and processing of the
$2.5\lt z\lt 4$
ZFOURGE control sample Section 2.1.1 as well as the EELG analogues Section 2.1.2. See Fig. 1 for the sample processing.
Figure 1. Summary of samples and their constraining quantities as well as some processing information.
2.1.1 ZFOURGE control
ZFOURGE is a survey combining legacy photometric UV to NIR data from 3 well studied regions (CDFS, COSMOS and UDS) with the FOURSTAR instrument upon the Magellan telescope (Straatman et al. Reference Straatman2016); its J, H and K medium band filters spanning the 1–1.8
$\mu {\rm m}$
range. It creates a photometrically well sampled survey and with the inclusion of two narrowband filters to optimally constrain the 4 000 Å break, enables robust 2% accuracy redshifts for over 70 000 galaxies (Nanayakkara et al. Reference Nanayakkara2016). In this work the ZFOURGE survey data was reduced to only retain CDFS sources as this field had the best narrowband filter depth (NB118 and NB209) which is important in isolating the [Oiii]5007 line flux for
$z\sim3$
sources. Galaxies were selected so that their K-band SNR
$\gt$
20 and were between the 2.5–4 redshift range. We used the use = 1 flag which eliminated catastrophic FAST (Kriek et al. Reference Kriek2009) and EAZY (Brammer et al. Reference Brammer, van Dokkum and Coppi2008) SED fits, interloper stars, AGN and sources too close to bright stars. The errors for IRAC bands were set to a floor of 25% (if not already above this) as the large PSF requires that the errors be large at high redshift. (Straatman et al. Reference Straatman2016). These were otherwise unphysically small and were found to significantly reduce the number of galaxies with well constrained physical stellar parameters. This created a total sample of 682 galaxies which includes both the control sample (606 galaxies) and the Extreme Emission Line Galaxy sample (76 galaxies) (see Section 2.1.2).
2.1.2 Filtering ZFOURGE EELGs at
$z \sim 3$
ZFOURGE EELGs were identified in Forrest et al. (Reference Forrest2018) using a stacked superposition of similar galaxy SEDs. The EELGs were determined as having an [Oiii]5007 EW of at least 400 Å and were selected between
$2.5\lt{ z}\lt 4$
. This creates a subsample of 76 EELGs among the 682 total galaxies within the CDFS region, described in the Section 2.1.1. We also use 18 of the 19 reliable spectroscopic redshifts from the MOSEL survey subsample which replace the corresponding ZFOURGE photometric redshifts.
2.2 Epoch of reionisation sample
In this section, we look at the selection criteria and processing of the
$z\gt 5.5$
control samples. This is broken down into the photometric redshift sample taken from the CEERS survey (Section 2.2.1) as well as the subsample of this which have spectroscopic redshifts. This is followed by the JADES survey sample which is only has spectroscopic redshifts (Section 2.2.3).
2.2.1 CEERS
The Cosmic Evolution Early Release Science Survey (Yang et al. Reference Yang2020; Bagley et al. Reference Bagley2023; Yang et al. Reference Yang2023) combines JWST NIRCAM photometry with legacy Extended Groth Strip (EGS) field HST data with the science goal of scouting the emergence of galaxies at cosmic dawn. We use the September 2022 data release processed with GrizliFootnote a where photometric redshifts were derived with EAZY-py. We select galaxies within
$5.5\lt {z}\lt 14$
and eliminate sources with a 95% confidence interval of the
$\chi^2$
value above 0.2 Gyr (
$t_{z[97.5]}-t_{z[2.5]} \lt 0.2$
Gyr), which similarly constrains the low and high ends of the redshift bounds. We further require that the determination of the UV slope be based on at least 3 filters in the UV window (Calzetti et al. Reference Calzetti, Kinney and Storchi-Bergmann1994) and remove sources which do not fit this criteria. This selects only the galaxies with a well constrained photometric redshifts and UV regions, giving us a sample of 461 galaxies within the EoR to compare to the ZFOURGE galaxies, particularly the EELG subsample.
2.2.2 Spectroscopic redshifts
For a small subsample of the CEERS catalogue, we have collected published spectroscopic redshifts from Tang et al. (Reference Tang2023) to minimise fitting errors within MAGPHYS and compare with the results of the publication. These will be referred to as the Tang23 sample throughout the paper.
2.2.3 JADES
The JWST Advanced Deep Extragalactic Survey (JADES, Bunker et al. Reference Bunker2023; Eisenstein et al. Reference Eisenstein2023; Hainline et al. Reference Hainline2023; Rieke et al. Reference Rieke2023) combines JWST NIRCAM photometry with legacy HST Ultra Deep Field (UDF) data. We only select spectroscopically confirmed galaxies at
${z}\gt 5.5$
, deriving a sample of 130 galaxies between
$5.5\lt{z}\lt 13.2$
. The JADES datasets depth and additional filters help to better constrain the physical parameters, and thus we compare mostly to this sample.
3. Methods
In this section we discuss the derivation of parameters used to describe the ‘blueness’ (UV slope,
${\unicode{x03B2}}$
) of a galaxy’s emission profile and its likelihood to contribute LyC to the IGM from both photometric data. We describe the SED models used, their limitations and any modifications we have made.
Figure 2. Observed frame MAGPHYS SED of
$z=3.473$
ZFOURGE EELG indicating the models fit to the photometric filter fluxes (grey circles) and their residuals (black circles, arrows if residual
$\sigma_{res}\gt 1$
). Blue region represents the LyC and is integrated under the unattenuated uncontaminated line (blue) to derive log
$_{10}(\xi_{ion}/(Hz\ {\rm erg}^{-1}))$
photometrically. Orange region is the UV region from which the UV slope is determined from either the attenuated SED
${\unicode{x03B2}}_a$
or from the photometry
${\unicode{x03B2}}_p$
. Gold filters represent the [Oiii]5007 contaminated filters and the difference these make to the attenuated SED is shown in the continuum luminosity difference between the ‘uncontaminated’ (black) and ‘contaminated’ (pink) lines. Residuals of the contaminated bestfit are represented by pink circles where
$\sigma_{res} \gt 1$
are represented by pink arrows.
3.1 MAGPHYS
Multi-wavelength Analysis of Galaxy Physical Properties (MAGPHYS) is an SED fitting package which derives the physical properties of a galaxy from the supplied photometry (da Cunha et al. Reference da Cunha, Charlot and Elbaz2008, Reference da Cunha2015). It does this by assembling a library of dust and stellar models at a predetermined redshift (if not allowed to vary) and then approaches the closest model using a marginalised likelihood distribution of each physical parameter (for more information see the documentation of MAGPHYS). We use the BC03 stellar models (Bruzual & Charlot Reference Bruzual and Charlot2003), (Charlot & Fall Reference Charlot and Fall2000) dust models, (Hildebrand Reference Hildebrand1983) grey body dust emission, (Chabrier Reference Chabrier2003) IMF, (Madau Reference Madau1995) IGM attenuation model and an exponentially declining star formation history model. The SFR timescale is over the past 100 Myr.
In our analysis, we use a modified version of the high_z version of the code with a lower dust prior, an increased range of available redshifts (beyond
$z\gt 10$
) and an approximately 10-fold increased model sampling at higher redshifts (
$z\sim8$
). Dust attenuation is expected to be low at high redshifts due to the early galactic stars being mostly composed of hydrogen and not having yet seeded a significant metallicity content into the ISM, of which dust is made (Shapley et al. Reference Shapley, Sanders, Reddy, Topping and Brammer2023; Cameron et al. Reference Cameron2023). While some dusty galaxies have been discovered in early epochs (Donnan et al. Reference Donnan2022), the original prior is based on local universe observations which are more consistently dusty and so required lowering for improved fits to the higher redshift sample. The introduction of additional models at high redshift allows us to more finely sample the parameter space and derive a better fit to the given photometry. This was found to be particularly necessary to prevent the highest redshift surveys with low sampling from sharing similar ill-fitting models across many different galaxies.
In order to determine the intrinsic production efficiency of ionizing photons we also remove the IGM absorption component of the SED after the other parameters are fit such that their determination is unaffected. The integration of the
$\lambda_{rest}\lt 912$
Å region is then used to derive the flux of the Lyman Continuum region. The physical parameters are still modelled using the Madau (Reference Madau1995) prescription and determined by their likelihood distributions. While MAGPHYS simultaneously determines dust and stellar components, we do not report the dust mass or luminosity as our photometric sampling does not constrain the mid-far IR range.
[Oiii]5007 ‘contaminated’ filters also limit the accuracy of the SED fit, due to these bright lines distorting the flux output in their containing filters. For a discussion of the [Oiii]5007 contamination and the method by which we account for this see Section 3.4. Fig. 2 breaks down a MAGPHYS SED fit to a ZFOURGE EELG photometry which has been sectioned into wavelength regions of interest and identifies the ‘contaminated’ filters.
3.2 Ionizing photon production efficiency
$\log_{10}(\xi_{ion})$
The contribution of a source to reionisation is described by the equation
\begin{equation} \dot{n}_{ion} = f_{esc}\times \xi_{ion} \times {\unicode{x03C1}}_{UV}\end{equation}
where
$\dot{n}_{ion}$
is the production rate of ionizing radiation,
$f_{esc}$
is the escape fraction of LyC light, log
$_{10}(\xi_{ion}/({Hz\ {\rm erg}^{-1}}))$
is the production efficiency of ionizing radiation in
${\rm Hz\ erg}^{-1}$
and
${\unicode{x03C1}}_{UV}$
is the comoving UV luminosity density in
${\rm erg\ s^{-1}\ Hz^{-1}\ Mpc^{-3}} $
.
The production efficiency of ionizing radiation log
$_{10}(\xi_{ion}/(Hz\ {\rm erg}^{-1}))$
is a key determinant of a galaxy’s potential contribution to reionisation; describing how much of the integrated spectrum is in the LyC region relative to the non ionizing UV component which represents the young stellar population. Spectroscopically it is determined by the Balmer line luminosities (
${\rm H}_{\unicode{x03B1}}$
,
${\rm H}_{\unicode{x03B2}}$
or
${\rm H}_{\unicode{x03B3}}$
) line luminosities in ratio with the luminosity of a non-ionizing UV wavelength such as 1 500 Å log
$_{10}(\xi_{ion}/(Hz\ {\rm erg}^{-1})) = \frac{N_{H^0}}{L_{UV}\times c_{rec}}$
where the
$N_{H^0}$
uses the Leitherer & Heckman (Reference Leitherer and Heckman1995) conversion
${\rm N(H^0)}[{\rm s}^{-1}] = \frac{1}{1.36}\times 10^{12}L(H_{\unicode{x03B1}})[{\rm erg\ s}^{-1}]$
from luminosity to a production rate of ionizing photons, assuming that the recombination and ionisation of the nebula are balanced.
$c_{rec} = 2.89$
refers to the case B recombination constant (Osterbrock & Ferland Reference Osterbrock and Ferland2006) for hydrogen that allows use of either the H
$_{\unicode{x03B1}}$
or H
$_{\unicode{x03B2}}$
lines depending on what is available observationally. It should be noted that a dust correction is to be applied to the
${\rm L_{UV}}$
as this is the intrinsic value. Photometrically however, we instead use equation 2 from Wilkins et al. (Reference Wilkins2016):
\begin{equation} (\xi_{ion}/(Hz\ {\rm erg}^{-1})) = \int^{c/912}_{\infty} L_\nu ({\rm h}\nu)^{-1}/L_\nu(1\,500) {\rm d}\nu\end{equation}
where c is the speed of light in Angstroms/s and h is the Planck constant in
${\rm erg\ Hz^{-1}}$
. This essentially divides the Lyman Continuum region luminosity density by a non ionizing luminosity constant, indicating the relative production of ionizing radiation. This is a value strongly dependent on the stellar population synthesis model used requiring a consistent method for meaningful comparison, so we explore only the MAGPHYS SED code for all samples. The unattenuated (dust corrected) luminosity is used for this calculation as it should reflect the intrinsic production efficiency.
3.3 UV slope
${\unicode{x03B2}}$
The UV slope (
${\unicode{x03B2}}$
where
$f_\lambda \propto \lambda^{\unicode{x03B2}}$
) is defined as the slope of the SED between rest frame UV wavelengths selected to exclude the LyC /Lyman-
${\unicode{x03B1}}$
features and to set an upper limit still considered to be in the UV region. It is commonly used as an indicator of dust attenuation and the age of the stellar population (Wilkins et al. Reference Wilkins2013; Williams et al. Reference Williams2018; Reddy et al. Reference Reddy2018; Nanayakkara et al. Reference Nanayakkara2022) which are particularly important for the analogue analysis in this work.
We further define two versions of the UV slope
${\unicode{x03B2}}$
parameter that reflect each the blueness of the direct (not dust corrected) photometry
${\unicode{x03B2}}_P$
and that of the model derived attenuated SED
${\unicode{x03B2}}_A$
between rest wavelengths 1 300–2 600 (the orange highlighted region of Fig. 2) (see Calzetti et al. (Reference Calzetti, Kinney and Storchi-Bergmann1994) for a discussion of the chosen wavelength range). These are determined using the scipy.optimize.curve_fit package and under the constraint that at least 3 filters between this wavelength range have a non-zero flux.
The direct photometry is less commonly considered in these analyses (Rogers et al. Reference Rogers, McLure and Dunlop2013) as the model dependence of the redshift determination creates a degeneracy with the star formation history model parameters that determine the attenuated UV fit. However, our data includes both spectroscopic redshifts and highly accurate photometric redshifts (
$\lt$
2% errors), therefore we can estimate the attenuated UV slope while considering the errors in the observed photometry without relying on the underlying models. A discussion on the value of this parameter in contrast to the model derived attenuated slope can be found in Dunlop et al. (Reference Dunlop2013). Fig. 3 highlights the key differences between these.
Figure 3. Attenuated model (
${\unicode{x03B2}}_a$
vs direct photometry
${\unicode{x03B2}}_p$
UV slopes for each subsample, with the ZFOURGE EELGs as blue stars, ZFOURGE control in cyan pluses, and the EoR samples in pink (JADES), gold (CEERS) and purple (Tang23) circles. The 1 to 1 line is shown in blue. General agreement is found with the exception of the CEERS sample around the
${\unicode{x03B2}}_p \lt -2.6$
region, where the model attempts to constrain the relative number of blue stars that could have formed in the galaxy’s lifetime, limiting the final ‘blueness’ of the slope.
The depth of the CEERS data as well as the limited photometric sampling could significantly distort the value derived from the direct photometric method, however, the model slope is not a solution to this. SED models are typically optimised for typical star forming galaxies at
$z\lt 6$
and hence the modelled UV slope may not be representative of the observed photometry for extreme cases. The model-independent UV slope is robust as it is estimated directly from the photometry and reflects the associated errors. Thus, it provides a consistent, comparative tool between the analogues and EoR galaxies.
3.4 Mass correction in EELGs [Oiii]5007 contamination
The [Oiii]5007 emission line flux poses as a significant contaminant to the photometric filters in which it resides if not taken into consideration (Forrest et al. Reference Forrest2018). The filter containing the emission line adds this to the total measured continuum flux, which is then read by the SED fitting code to determine physical parameters. Significantly bright emission lines raise the measured continuum flux above the real continuum, resulting in a flux excess and causing the physical parameters associated with that portion of the SED to be distorted if not correctly accounted for, as is the case with the current version of MAGPHYS.
Our sample selection is based on the [Oiii]5007 EW, though we entirely remove the [Oiii]5007 contaminated filters while fitting the SED of galaxies with MAGPHYS. Fig. 2 shows the effects of this on an example SED while Fig. 4 shows the effect of contamination on the stellar mass of the ZFOURGE galaxies and the EELG subsample. We note that EELG galaxies around
$z\sim 3.5$
are most significantly affected by this. The mass determination is dependent largely on the optical region of the SED and the heightened [Oiii]5007 flux in this subsample is expected to contribute to the four filters in the K-band (
$\sim$
2 micron) more-so than in most galaxies. The inclusion of emission lines in the SED fitting models themselves is another way of accounting for this effect which is employed in other SED fitting codes (e.g. Chevallard & Charlot Reference Chevallard and Charlot2016). The same removal method was applied to the JADES sample but not the CEERS sample, as its further limited photometric coverage would significantly reduce the model reliability. While this may result in a moderately overestimated mass in the CEERS sample due to the expected high average [Oiii]5007 EW, this alone would not account for the extreme values of some of these galaxies which require better photometric constraints (see Section 4.1 for further discussion).
Figure 4. Stellar Mass with and without inclusion of the [Oiii]5007 contaminated filters for the ZFOURGE control and ZFOURGE EELG samples. We also note a trend that
$2.9\lt z\lt 3.5$
galaxies tend to be more affected by the contamination and our EELG subsample is on the higher end of the discrepancy. This is due to the tendency for the [Oiii]5007 line to fall into multiple filters at this redshift, and this is represented by the colour bar.
4. Results and discussion
The physical parameters relating to the stellar population and dust attenuation are analysed in the ZFOURGE EELG sample and compared to the control samples within the EoR (JADES, CEERS, Tang23) and at similar redshift (ZFOURGE control). The comparisons selected explore the consistency with which the
$2.5\lt z\lt 4$
[Oiii]5007 analogues match the physical parameters and internal processes of EoR galaxies.
4.1 Main sequence
We find that the ZFOURGE EELG subsample is consistently on the lower mass, higher specific star formation rate end of this parameter space when compared to the ZFOURGE control sample (Fig. 5). While ZFOURGE EELG stellar masses are generally an order of magnitude above the EoR counterparts, the sSFR median and 25–75 quartile region is comparable in the ZFOURGE EELG subsample and the JWST survey samples between redshifts 5.5–14 (Table 1). Our ZFOURGE control sample at
$z\sim3.4$
lies close to the stellar mass-redshift-sSFR relation developed by Popesso et al. (Reference Popesso2022) using 27 other studies between
$0\lt z\lt 6$
and
$8.5\lt\log_{10}(M_*)\lt 11.5$
. The sSFR of the ZFOURGE EELG sample, however, is 0.2 dex above this relation (Popesso et al. Reference Popesso2019; Leslie et al. Reference Leslie2020) and more closely related to the semi-empirical samples at
$z\sim6$
(Grazian et al. Reference Grazian2015; Davidzon et al. Reference Davidzon, Ilbert, Faisst, Sparre and Capak2018). This reflects the strength of the [Oiii]5007 EW selection technique at low redshift in collating the most highly star forming galaxies. We note the appearance of striations visible in the CEERS sample. This is due to the relatively lower SNR and fewer photometric filters across the spectrum; preventing the code from fitting models to the more intricate variations in the continuum and thus recycling similar models. This effect continues through the other physical parameters to which the CEERS sample is fit and should be considered accordingly.
Figure 5. Specific star formation rate vs Stellar Mass plot for each of our samples. Colour-scheme is the same as Fig. 3 The median values and interquartile ranges are indicated by the coloured errorbars with orange:CEERS, pink:JADES, dark blue :ZFOURGE EELGs, teal:ZFOURGE control. The upper boundary is caused by the chosen SFR timescale. Popesso et al. (Reference Popesso2022) z-sSFR-
$M_*$
relation selected at
$z\sim3.4$
and
$M_*/M_\odot=10^{10}$
shown as a red cross, and Davidzon et al. (Reference Davidzon, Ilbert, Faisst, Sparre and Capak2018) at
$z=2.22$
and
$z=6$
shown as a blue and green diamonds for comparison.
The stellar mass and star formation rates of both the ZFOURGE control and EELG samples are uniformly higher than that of the EoR samples, with the EELG sample being slightly less massive and less star forming. While the ZFOURGE EELG mass is higher than the EoR samples, it is still
$\sim$
25% less massive than the average population at
$z\sim3$
. The emission of LyC radiation is also not directly dependent on this parameter, so the functionality as an analogue is maintained. The heightened SFR is related to the emission, however, when normalizing by mass (sSFR) we will see that the discrepancy between EELG analogues and EoR galaxies is minimised.
4.2 Dust extinction star formation relation
Dust and star formation are inextricably linked galactic processes. Dust cools and catalyses molecular gas collapse while stellar death synthesizes more dust and transports it through the galaxy (Popping et al. Reference Popping, Somerville and Galametz2017). The correlation of these properties across different epochs could be indicative of their evolving relationship as galaxies age and their dust content increases. Studies such as Zahid et al. (Reference Zahid, Yates, Kewley and Kudritzki2012) link the correlation between
$A_V$
and SFR to the processes that quench star formation. A positive Spearman rank correlation between the dust extinction and star formation rate is observed by Li et al. (Reference Li2019) in nearby galaxies, while Zahid et al. (Reference Zahid, Yates, Kewley and Kudritzki2012) found a mass dependence on the correlation. They found more massive galaxies to be positively correlated while finding an anticorrelation below
$10^{10}\,{\rm M}_\odot$
when the quiescent high mass sample was removed, indicative of differing internal processes governing this relation between the low and high mass samples. Interestingly, the Balmer decrement (
$A_V$
tracer) stellar mass relation shows minimal evolution with redshift for
$z\lt 2$
sources (Battisti et al. Reference Battisti2022).
For galaxies to be analogues of EoR galaxies, they should exhibit similar correlations across various physical parameters, owing to similar internal physical environments. Here we compare our findings across different epochs with our self-consistent, single SED fitting method to determine if the [Oiii]5007 selected ZFOURGE EELG sample has a similar relationship between these properties to EoR galaxies (see Fig. 6). We find that ZFOURGE EELGs have a consistently lower dust extinction value when compared to the ZFOURGE control sample Table 1. The median value is half that of the ZFOURGE control and is well within the 25–75th percentiles of the EoR JADES sample.
We find a positive correlation between SFR and dust attenuation across all our samples. The correlation appears strongest in the ZFOURGE control sample, which tend to have higher
$A_V$
and SFR than the others (See Table 2). This finding is similar to Sakurai et al. (Reference Sakurai, Takeuchi, Yuan, Buat and Burgarella2013) who also find a stronger correlation between
$A_V$
and SFR for galaxies with SFR
$\gt 20\,{\rm M}_\odot/{\rm yr}$
. However, the exact values of our correlations may be biased by the intrinsic scatter and parameter ranges. For example, the weaker correlations observed in ZFOURGE EELGs and EoR (JADES, CEERS) samples could be due to the lack of dusty, low star-forming galaxies in the samples.
Table 1. Parameter medians with 25–75 percentiles for the ZFOURGE control, ZFOURGE EELGs, JADES and CEERS samples.
Table 2. Spearman rank correlations for each subsample (ZFOURGE control, ZFOURGE EELGs, JADES and CEERS) and the combined dataset (Total). Blank spaces reflect correlations above a significance factor of 0.05.
Even after accounting for stellar mass by comparing sSFR and
$A_V$
, we still find significant positive correlations across the four samples. The weaker correlations in the ZFOURGE and CEERS samples could be due to their relatively large parameter ranges compared to the ZFOURGE EELG and JADES samples. The similar correlation in the EELG and JADES samples in contrast to the ZFOURGE control suggests that the physical characteristics of ZFOURGE EELGs are more comparable to the
$z\gt 5.5$
JADES sample than to the control
$z\sim3.4$
sample.
4.3 UV slope (
${\unicode{x03B2}}_p$
) relations
We compare the different samples parameter space relations in Fig. 7 and identify the correlations between each subsample in Table 2. The UV slope interquartile range of the ZFOURGE EELGs matches with the upper quartile of the JADES sample and is well below the ZFOURGE control (Fig. 7 top panel). This identification of the ZFOURGE EELG sample as being similarly blue and dust free as the direct EoR samples while having a higher specific star formation rate than the ZFOURGE control sample (see Table 1) indicates a strong potential for these to be good EoR analogues.
Figure 6. Star formation rate vs Dust attenuation plot for each of our samples. Marker sizes correlate to the Stellar Mass of the galaxy, with larger markers representing more massive galaxies consistent across all samples. Colour scheme represents the same samples as Fig. 3. median values and interquartile ranges are indicated by the coloured errorbars following the same colour-scheme as Fig. 5.
The
${\unicode{x03B2}}_P$
vs
$M_*$
relation (Fig. 8) appears to be discontinuous between the low and high redshift samples. The observed correlation between stellar mass and slope agrees with the conclusions of Pannella et al. (Reference Pannella2009) and with the consensus that increasing stellar mass correlates with increasing dust attenuation and therefore a redder (observed) slope (Bouwens et al. Reference Bouwens2016). This is supported by the lower dust extinction (Table 1) of the high redshift samples, as well as explains the tight parameter space occupancy of the ZFOURGE EELG sample in this figure. The correlation between UV slope and stellar mass is not present in the ZFOURGE EELG sample (Table 2) as it occupies a very constrained portion of the parameter space compared to the other samples (1.8 dex in stellar mass for ZFOURGE EELGs vs 3.6 dex for JADES for example). The UV slope appears to be sensitive to the mass with a strong overall correlation of 0.67 and strong correlations with the EoR and control samples. This suggests that the ZFOURGE EELG sample is consistently blue despite the mass being an order of magnitude greater than the EoR JADES median.
Figure 8. MAGPHYS derived physical parameters vs
${\unicode{x03B2}}_P$
(uncorrected for dust, top panel) vs photometric log
$_{10}(\xi_{ion}/(Hz\ {\rm erg}^{-1}))$
(dust corrected Charlot & Fall Reference Charlot and Fall2000, bottom panel) for each of the samples following the same colour scheme as previous figures. Median values and correlations can be found in Tables 1, 2. Model limitations manifest as ‘walls’ in
$\xi_ion$
and sSFR.
The correlation between the UV slope and dust extinction is readily identifiable, being strongest in the ZFOURGE control and EELG samples, which is similar to the conclusion drawn in Wilkins et al. (Reference Wilkins2013). The parameter space occupied by the ZFOURGE EELGs is the same as that of the JADES and Tang23 samples which could indicate the ‘blueness’ of a galaxy being attributed to the absence of attenuating dust in the model, though better dust constraints are required.
We also note the negative correlation of the
${\unicode{x03B2}}_P$
with
$\xi_{ion}$
plot (Fig. 8). The
$\xi_{ion}$
does not evolve significantly between the different redshift samples, however it does appear to negatively correlate with UV slope in both the EoR samples. This suggests that at least for the EoR samples, bluer galaxies tend to have a higher ionizing to non ionizing photon ratio though further investigation of the model dependence of the
$\xi_{ion}$
parameter are necessary.
4.4 Ionizing photon production
This section explores the correlations between the ionizing photon production efficiency
$\xi_{ion}$
and the physical parameters determined by our model (Fig. 8). We refer to Table 2 for the correlation coefficients and their significance.
The
$\xi_{ion}$
anticorrelates with stellar mass in the low redshift samples (−0.48 ZFOURGE control, −0.36 ZFOURGE EELGs), suggesting that a lower mass corresponds to a higher production efficiency, though this is not clear from its visual appearance (Fig. 8). This weak dependence on stellar mass agrees with the findings of Emami et al. (Reference Emami2020) and Shivaei et al. (Reference Shivaei2018) for this mass range (
$10^{7}$
–
$10^{11} {\rm M}_\odot$
) at a similar redshift range to our study, as well as with Lam et al. (Reference Lam2019) for the EoR samples. This weak anticorrelation with stellar mass potentially implicates smaller systems at early epochs as the powerhouses of ionizing photon production.
The production efficiency correlates most strongly with the sSFR in the ZFOURGE control and its EELG subsample. The assertion that a high specific star formation rate would be reflected in higher LyC emission is known (Castellano et al. Reference Castellano2023) between
$2\lt z\lt 5$
. While a significant correlation could not be determined for the JADES sample due to the small dynamical range of sSFR values, the CEERS sample does reveal a similarly strong relationship in this parameter space.
A notable correlation between
$\xi_{ion}$
and dust extinction is found only in the ZFOURGE EELG sample. The result suggests a relationship between the dust within [Oiii]5007 selected galaxies and their efficiency of LyC photon production. This is likely due to the underlying degeneracy correlating
$A_V$
, age, and SFR. Young galaxies with emerging O star populations haven’t seeded as much dust when compared to their older counterparts, and it is these massive stars that are the driving force for the strong [Oiii]5007 emission lines in the ZFOURGE EELG sample.
5. Summary and conclusion
This paper uses SED modelling to derive the physical parameters of EELGs with analogous [Oiii]5007 emissions to EoR galaxies at
$z\sim3$
and compares them to both a control sample as well as the EoR galaxies observed at
$z\gt 5.5$
with JWST. In this section we summarize the findings of our analysis and the potential avenues further research can take.
-
1. The combined high sSFR (−8.28
${\rm yr}^{-1}$
), comparatively low stellar mass (9.75
${\rm M}_\odot$
), low dust extinction (
$A_V$
= 0.25) and blue UV slope (
${\unicode{x03B2}}_P$
= −1.84) suggest that the [Oiii]5007 selection technique is a good analogue selector. The high SFR and low dust are critical for their analogue status as they must produce copious ionizing radiation that is not attenuated by dust. We need to further investigate the escape of LyC light by studying the
$f_{esc}$
to complete the picture. -
2. We confirm the strength of the [Oiii]5007 selection technique for finding EoR analogues with lower stellar mass, high relative star formation rates for our sample between
$2.5\lt z\lt 4$
. These galaxies are an order of magnitude higher in mass than their EoR counterparts but maintain similar sSFR despite this. -
3. We find that the [Oiii]5007 selected EELGs have a more similar correlation between their dust and star formation parameters to EoR galaxies than to the control sample between
$2.5\lt z\lt 4$
, suggestive of similar internal processes relating to the dust distribution and star formation between the analogues and EoR galaxies. -
4. We find a correlation between the UV slope and dust extinction/attenuation which is strongest in our low redshift samples. We also see a similar correlation with the stellar mass, particularly in our control sample and the CEERS data. Experiencing similar attenuation between these samples indicates low dust obscuration however this is only supported by UV photometry. Further studies using deep ALMA observations are required to determine if the steep UV slopes are due to the stellar population or an absence of dust.
-
5. We find that the median
${\unicode{x03B2}}_P$
of our ZFOURGE EELGs matches more closely with the higher redshift counterparts than with the
$2.5\lt z\lt 4$
control sample, suggesting that our selection technique has recovered galaxies with a similar starbursty nature as the EoR galaxies. -
6. We find the most significant correlation of log
$_{10}(\xi_{ion}/(Hz\ {\rm erg}^{-1}))$
to be with its sSFR which coincides with the findings of Castellano et al. (Reference Castellano2023) for bright
$M_UV\leq 20$
,
$2\lt z\lt 5$
galaxies. This correlation appears the strongest with our low redshift samples, likely due to the wider parameter space.
Acknowledgement
This research was partly supported by the Australian Research Council Centre of Excellence for All Sky Astrophysics in 3 Dimensions (ASTRO 3D), through project number CE170100013. CMT was 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. AH acknowledge support from the ERC Grant FIRST- LIGHT and Slovenian national research agency ARRS through grants N1-0238 and P1-0188.
RJ would like to thank Jessica E Thorne for helpful discussions.
Data Availability
All data and software mentioned is publicly available.
Open access
