Introduction

Young, star-forming galaxies are long thought to be the major drivers of the cosmic hydrogen reionization because they should have strong UV emission and are sufficiently abundant at z > 6 (e.g., H. Yan & R. A. Windhorst 2004; R. J. Bouwens et al. 2006;S. L. Finkelstein et al. 2012;B. E. Robertson 2022;H. Atek et al. 2024). Their exact contribution to the ionizing photon background, however, still depends on how effectively their Lyman-continuum (LyC; λ < 912 Å) photons can escape and reach the surrounding intergalactic medium (IGM). Direct measurement of this escape fraction ( f esc ) for star-forming galaxies at z > 6 is impossible, because the line-of-sight IGM H I absorption at such redshifts wipes out any LyC photons. As an alternative, f esc has been measured for analogs at lower redshifts (z  3-4) where the IGM hydrogen is fully ionized (e.g., C. C. Steidel et al. 2001;E. Vanzella et al. 2012;B. Siana et al. 2015;A. Grazian et al. 2016; Y. I. Izotov et al. 2016Izotov et al. , 2018Izotov et al. , 2021;;C. C. Steidel et al. 2018; S. R. Flury et al. 2022;A. Griffiths et al. 2022;A. Saldana-Lopez et al. 2022;A. Citro et al. 2024), and the correlations between f esc and various observables are sought after in order to provide some viable routes to indirectly measure f esc at z > 6. Among all possible correlations, the one with the rest-frame UV slope (commonly denoted as β; f λ ∝ λ β ) is the most promising (e.g., E. Zackrisson et al. 2013;J. Chisholm et al. 2022).

There have been ample studies of the UV slopes of galaxies at z > 6 using the Hubble Space Telescope (HST) deep survey data (e.g., N. P. Hathi et al. 2008; R. J. Bouwens et al. 2010; J. S. Dunlop et al. 2012Dunlop et al. , 2013;;S. L. Finkelstein et al. 2012), and the investigation has been extended to higher redshifts and fainter limits by the James Webb Space Telescope (JWST; e.g., M. W. Topping et al. 2022;F. Cullen et al. 2023;T. Nanayakkara et al. 2023;D. Austin et al. 2024; A. M. Morales et al. 2024). However, the vast majority of those are based on photometric samples of candidates that have not yet been spectroscopically confirmed. Thanks to the growing number of JWST spectroscopic programs, the situation is now quickly changing (e.g., S. Fujimoto et al. 2023;M. Tang et al. 2023;G. Roberts-Borsani et al. 2024;A. Saxena et al. 2024). Nevertheless, there are still few studies using confirmed z > 6 galaxies. In particular, there seems to be a lack of extremely young galaxies (age  30 Myr) in the samples. Such galaxies presumably should have the bluest UV slopes because their UV emission must be dominated by O and B stars. This work presents a case study using JWST data of four z ≈ 8.3-8.5 galaxies behind the MACS J0416.1-2403 cluster (hereafter MACS0416). MACS0416 is one of the six Hubble Frontier Fields (HFF; J. M. Lotz et al. 2017) and one of the targets of the Reionization Lensing Cluster Survey (D. Coe et al. 2019). A few investigations have used HST to search for lensed high-redshift galaxies behind this cluster (e.g., D. Coe et al. 2015;L. Infante et al. 2015;N. Laporte et al. 2015), and some candidates with photometric redshifts z ph  8-9 have been found. Among them, two have reported spectroscopic redshifts at z > 8, both based on Atacama Large Millimeter/submillimeter Array (ALMA) spectroscopy. One is MACSJ0416.1_Y1 (hereafter "Y1" for simplicity), which was first selected as a z ≈ 8 candidate (N. Laporte et al. 2015; z ph = 8.1-8.9 with best-fit value 8.57) and was later confirmed at z = 8.31 through the detections of the [O III] 88 μm line (Y. Tamura et al. 2019) and the [C II] 157.7 μm line (T. J. L. C. Bakx et al. 2020). The other is MACS0416.1-JD (hereafter "JD"), which was also first discovered as a z ≈ 8 candidate by D. Coe et al. (2015;object "FFC2-1151-4540") with z ph = 7.3-8.6 and best-fit z ph = 8.1. JD was independently rediscovered by N. Laporte et al. (2015; their object "MACSJ0416.1_Y2") with z ph = 8.3-8.6 and best-fit z ph = 8.47. (See also L. Infante et al. 2015;N. Laporte et al. 2016.) N. Laporte et al. (2021) reported the detection of the [O III] 88 μm line at z = 9.28 and renamed this object JD, which we have adopted for brevity. As we will show, the new JWST data confirm the redshift of Y1 at z = 8.31 but change that for JD to z = 8.34. In addition, the close neighbor to JD, object "MACSJ0416.1_Y3" as reported by (N. Laporte et al. 2015; z ph = 8.8-9.7 with best-fit z ph = 9.29), is also at z = 8.34. Furthermore, we have identified a new JWST NIRCam F090W dropout, which we name "f090d_018." It turns out to have a similar redshift, z = 8.49. Interestingly, all four objects have very young ages (∼4-18 Myr) inferred from their spectral energy distributions (SEDs), and yet only two of them show blue UV slopes of β < -2.0.

This paper is organized as follows. Section 2 describes the JWST data, and Section 3 presents the photometric and spectroscopic results. The analysis is given in Section 4, and Section 5 provides a brief summary. We adopt a flat Λ cold dark matter cosmology with h = 0.7, Ω Λ = 0.7, and Ω M = 0.3. All magnitudes are in the AB system, and all coordinates are in the ICRS frame (equinox 2000).

JWST NIRCam Observations and Data Reduction

NIRCam Imaging

MACS0416 is one of the targets observed by the JWST GTO program Prime Extragalactic Areas for Reionization and Lensing Science (PEARLS; PI: R. Windhorst; PID 1176; R. A. Windhorst et al. 2023). It was imaged by NIRCam in eight filters, F090W, F115W, F150W, and F200W in the short wavelength (SW) channel and F277W, F356W, F410M, and F444W in the long wavelength (LW) channel. The observations were carried out in three epochs: 2022 October 7 (Ep1), 2022 December 29 (Ep2), and 2023 February 10 (Ep3). The observations were described by Yan et al. (2023, their Table 1), who gave details of the data reduction. In H. Yan et al. (2023), the images were stacked on a per-epoch basis for the transient studies. For this work, we combined the data in all three epochs, which reached total integration times of 8761 s in F150W, F200W, F277W, and F356W, 10,909 s in F115W and F410M, and 11,338 s in F090W and F444W. In Ep3, Y1 is contaminated by a strong diffraction spike of a bright star, and therefore we only combined the Ep1 and Ep2 data for the study of this object. The total integration times are 5841 s in F150W, F200W, F277W, and F356W, and 7559 s in F090W, F115W, F410M, and F444W. As in H. Yan et al. (2023), we created stacks at scales of both 0 06 (60 mas) and 0 03 (30 mas), which have magnitude zero-points of 26.581 and 28.087, respectively. These images are all aligned to the HFF HST images. Figure 1 shows NIRCam image cutouts around Y1, JD, and the new object f090d_018.

NIRCam Wide-field Slitless Spectroscopy

MACS0416 was also observed by the NIRCam instrument in its wide-field slitless spectroscopy (WFSS) mode. This mode has two settings, Grism R (GR) and Grism C (GC), which disperse light in either the direction of detector rows (R) or columns (C) in the LW channel, both having spectral resolution R ≈ 1600 at ∼4 μm. The field was observed by the programs PID 2883 ("MAGNIF: Medium-band Astrophysics with the grism of NIRCam in Frontier Fields;" PI: F. Sun) and PID 3538 ("Unveiling the properties of high-redshift low/intermediatemass galaxies in Lensing fields with NIRCam Wide Field Slitless Spectroscopy;" PI: E. Iani).

The MAGNIF observations of MACS0416 were carried out on 2023 August 20 with GC using the F480M filter with a total integration time of ∼3100 s. These observations covered only JD and f090d_018 with Y1 ∼3″ outside the field. The observations of PID 3538 in this field were carried out on 2023 December 22-24 and 2024 January 17 with both GR and GC. Four filters were used, F300M, F335M, F410M, and F460M, with total integration times in each band of 1761 and 1546 s, respectively, for the two grisms.

To reduce the grism data, we retrieved the Level 1b "uncal" files from the Mikulski Archive for Space Telescopes and ran them through the calwebb_detector1 routine of the JWST data reduction pipeline (version 1.13.4 in the context of jwst_1223.pmap) to obtain the "rate.fits" files. We then followed the procedures described by F. Sun et al. (2023) to further process the data. All single exposures were registered to Gaia Data Release 3 astrometry. There is a systematic offset between the HFF astrometry and that of Gaia in this field, which can be corrected by R.A.(Gaia) = R.A.(HFF) +0 21, and decl.(Gaia) = decl.(HFF) -0 08. This was taken into account when extracting the spectra.

Data Analysis

Overview

From their high-resolution (beam size 81.1 mas × 112.2 mas) ALMA [O III] 88 μm image of Y1, Y. Tamura et al. (2023) resolved the system into three knots ("O1," "O2," and "O3"). Using this as a guide, these authors separated their HST WFC3 image into three components ("E," "C," and "W" from east to west) that coincide with the three [O III] 88 μm knots. These three components are clearly distinguished in the NIRCam SW images (Figures 1 and 2), and the entire system extends ∼0 73 along the long axis.

JD has a close companion object 0 91 away to its northeast, which was also selected by N. Laporte et al. (2015) as a Y-band dropout (their "MACSJ0416.1_Y3"). As we will show below, its colors meet the requirements of a J-band dropout as well, and its redshift is almost the same as that of JD. We refer to it as "JD-N" in this work to indicate its association with JD.

The object f090d_018 is among a sample of F090W dropouts that we selected in this field (B. Sun et al. 2024, in preparation). While this object is very close to a bright member galaxy of the M0416 cluster (at z = 0.401; see V. Kokorev et al. 2022), we can still obtain reliable photometry for this object.

Photometry and Spectral Energy Distributions

To construct the SEDs of our targets, we carried out photometry on the 60 mas images. In addition to the NIRCam data, we also incorporated the HST WFC3 data in F105W, F125W, F140W, and F160W obtained by the HFF program. 25 Although these WFC3 images are in the wavelength range covered by NIRCam, their inclusion improves the wavelength sampling. All the NIRCam and WFC3 images were convolved to match the point-spread function (PSF) of the F444W image, which has the coarsest spatial resolution (PSF FWHM = 0 145). WebbPSFs (M. D. Perrin et al. 2014Perrin et al. , 2015) ) were used in the convolution of the NIRCam images. For the WFC3 images, we used the empirical PSFs available at the WFC3 instrumentation site. 26 Matched-aperture photometry was done by running SEXTRACTOR (E. Bertin & S. Arnouts 1996) in the dual-image mode with F444W as the detection band, and we adopted the MAG_ISO magnitudes. The magnitude uncertainties were calculated using rms maps derived with ASTRORMS. 27 Two details affected the photometry: (1) JD-N has a faint neighbor that is not visible in F090W, but the neighbor does not pass our selection for F090W dropouts. The photometry for JD-N was obtained after masking this neighbor. (2) f090d_018 could be affected by the outskirts of the z = 0.401 foreground galaxy to the southeast, and therefore we modeled this foreground galaxy using GALFIT (C. Y. Peng et al. 2002Peng et al. , 2010) ) and subtracted it from all images before doing photometry of f090d_018.

We also obtained the SEDs for Y1ʼs three components by decomposing them using GALFIT on the 30 mas images (not convolved for the PSF sizes). Figure 3 demonstrates the decomposition, which was done only for the NIRCam images but not the WFC3 ones because the latter do not have sufficient spatial resolution. Briefly, we used a single Sérsic profile (J. L. Sérsic 1963) to model the light distribution of each component and fit the three components simultaneously. The centers of the three components were determined on the F200W image, where they are most clearly separated, and the fit was done with the centers fixed. Table 1 summarizes the photometry for all sources. For Y1 as a whole and the other three objects, the 2σ upper limits in F090W were calculated on the rms map at the source positions within a circular aperture whose size matches the MAG_ISO aperture or r = 0 2, whichever is larger. For the three components of Y1, these upper limits were calculated within an r = 0 2 circular aperture.

Spectroscopic Identifications

Figure 2 shows color composites of these two systems with the WFSS dispersion directions indicated. The NIRCam WFSS mode is very suitable for the detection of emission lines, which is the main focus here. To optimize the line detection, we subtracted a continuum estimated (following D. Kashino et al. 2023) by running a median filter 1 × 51 pixels in size with a 9 pixel "hole" at the center along each row or column.

As Y1 consists of three components, it would be ideal to extract the spectra for each component. In practice, however, the individual extraction could only split the signal into two parts, which largely coincide with Y1-E and Y1-W, respectively. This is due to both the faintness of Y1-C at 3-5 μm and the coarse resolution at these wavelengths. For simplicity, we refer to the separate extractions as Y1-E and Y1-W even though both could include some contribution from Y1-C. As shown in Figure 2, the PID 3538 GR dispersion direction is very close to the spatial extension of the Y1 system, causing severe spectral contamination between the two components. For this reason, we only used the GC data. The extracted 2D and 1D spectra are shown in Figure 4.

For both Y1-E and Y1-W, the [O III] λλ4959, 5007 lines are clearly seen in the F460M data. The detection of the [O II] λ3727 line in the F335M data is marginal in the 1D spectrum but convincing in the 2D spectrum. To determine the redshifts, we fitted a Gaussian profile to the [O III] λ5007 line, which is the strongest, to obtain its observed central wavelength. The redshifts thus derived are z = 8.309 ± 0.002 for Y1-E and z = 8.312 ± 0.002 for Y1-W. We also fitted Gaussians to the weaker lines, using the initial guesses of the line centers at where they should be under these redshifts. The results based on the best-fit central wavelengths thus obtained are all in excellent agreement with the redshifts based on the [O III] λ5007 line.

For JD and JD-N, we extracted the spectra using the data from both MAGNIF and PID 3538. However, we discarded the GC data from the latter because the dispersion direction is almost parallel to the position angle of the two objects (Figure 2). The spectra are shown in Table 2 summarizes the line measurements. The line widths were converted to velocities by Δv = cΔλ/λ c , where λ c and Δλ are the mean and FWHM of the Gaussian fit, respectively, and c is the speed of light. The total line intensities were obtained by integrating the fitted Gaussian profile within a 4 × FWHM wavelength range centered at λ c . The associated errors were estimated using the nonsmoothed 1D spectra. When calculating the observed equivalent widths (EWs), the F410M magnitudes in Table 1 were used as the continuum flux density. To account for continuum contributions from Y1-C in the spectra of Y1-E and Y1-W, half the F410M flux density of Y1-C was added to the

Note.

The measurements were all obtained by fitting a Gaussian profile to the emission lines (not corrected for the magnifications). The spectroscopic redshifts are based on the fitted central wavelengths of the [O III] λ5007 line, which is the strongest. The integrated line fluxes I are in units of 10 -17 erg cm -2 s -1 ; the line widths Δv were calculated using the FWHM values of the lines and are in units of km s -1 . The observed EWs (in units of Å) were calculated assuming a flat f ν continuum at the level determined by the F410M magnitude. The subscripts "5007," "4959," and "3727" represent the [O III] λ5007, [O III] λ4959, and [O II] λ3727 lines, respectively. For the nondetections, the 2σ f ν and f λ upper limits are in units of μJy and 10 -20 erg s -1 cm -2 Å -1 , respectively.

Y1-E and Y1-W continua. Not surprisingly, these objects all have very large [O III] EW values, which can largely account for the brightening of their SEDs in F444W. For example, the measured [O III] λ5007 EW of JD implies that this line alone increases the object's F444W brightness (as compared to that in F410M) by -0.83 ± 0.06 mag, which is almost the same as the observed m 410m 444 = -0.83 ± 0.09 mag.

Discussion

Redshifts, Magnification, and Environment

The redshifts of Y1-E and Y1-W agree with each other within the uncertainties, and their rest-frame relative velocity is only 97 ± 64 km s -1 . The average is z = 8.311 ± 0.003, which is in very good agreement with the previously reported z = 8.3118 ± 0. 2019) is 0.66 Jy km s -1 , which corresponds to 7.92 × 10 -18 erg cm -2 s -1 at the observed frequency (364.377 GHz). Therefore, the line intensity ratio of [O III] λ5007 and 88 μm is 2.79.

For JD, our redshift z = 8.341 ± 0.002 disagrees with the previously reported z = 9.28 (N. Laporte et al. 2021), which was based on the claimed signal-to-noise ratio ≈ 6 detection of the [O III] 88 μm line. That redshift could not explain the emission line observed at 4.68 μm. Therefore, we believe that their identification of the [O III] 88 μm line is incorrect and that the claimed detection is likely due to a false positive. There have been a few cases in the literature reporting seemingly significant millimeter line detections indicative of high redshifts but being subsequently refuted (e.g., G. Popping 2023; Y. Harikane et al. 2024, gave recent examples), and JD reinforces the need for caution.

According to the lens model of J. M. Diego et al. (2024), the magnification factors for Y1, the JD-JD-N pair, and f090d_018 are μ = 1.21, 2.26, and 1.32, respectively. Based on their magnitudes in F150W (sampling the rest-frame UV range of ∼1430-1790 Å), their intrinsic absolute UV magnitudes after correcting for lensing magnifications are M UV = -20.83, -19.60, -19.77, and -19.57 mag for Y1, JD, JD-N, and f090_018, respectively. The observed size of Y1 as a whole is only marginally affected by magnification, and its corresponding physical size along the long axis is ∼3.4 kpc. The separation of JD and JD-N, on the other hand, is affected significantly by lensing. Their separation in the source plane (see Table 1) is ∼0 51, which corresponds to 2.4 kpc. The velocity offset between JD and JD-N is only 160 ± 64 km s -1 . Therefore, JD and JD-N should form an interacting pair.

After correcting for lensing magnification, the projected distance between Y1 and the JD/JD-N pair is only ∼165 kpc, and their comoving radial distances differ by only ∼8.2 Mpc. Therefore, these three objects are likely part of an overdensity at z ≈ 8.3. This is similar to the z ≈ 8.2 overdensity recently identified by J. M. Helton et al. (2024) in the GOODS-S field. Furthermore, Y1 + JD/JD-N might signpost an even large structure: f090d_018 at z = 8.49 is separated from Y1 by ∼77 kpc (projected) and ∼48.7 Mpc in the radial direction, suggesting that it could be a member of a filament-like structure stretching from the Y1 + JD/JD-N complex. A future paper (N. Foo et al. 2024, in preparation) will present a detailed analysis of this possible structure.

Stellar Populations

To understand the stellar populations of these galaxies, we fitted their SEDs using the BAGPIPES software (A. C. Carnall et al. 2018), which utilizes the stellar population synthesis models of G. Bruzual & S. Charlot (2003) with the P. Kroupa (2001) initial mass function. The fitting was done at fixed redshifts of 8.31, 8.34, and 8.49 for Y1, JD/JD-N, and f090d_018, respectively. An exponentially declining star formation history (SFH) in the form of SFR ∝e -t/ τ was assumed. The option to include nebular emission lines was enabled, and we used the Calzetti dust-extinction law (D. Calzetti et al. 1994;D. Calzetti 2001) with A V ranging from 0 to 2.0 mag. The metallicity was allowed in the range of 0 Z * /Z e 2.5, and the ionization parameter could vary in ( ) U 2.5 log 0.5 -< <-. The results of the SED fitting are summarized in Table 3, and Figure 5 shows the 16th-84th percentile range of the posterior spectra superposed on the SEDs. Corner plots showing the posteriors for the fitted parameters are given in Figure 6.

Overall, these four objects share many similarities, and the three components of Y1 are also quite similar to each other. In particular, Y1, JD and JD-N all have ages  10 Myr. To our knowledge, these are the youngest galaxy ages ever reported in the literature. 28 F090d_018 is only about twice as old, albeit Note. The quoted value for each parameter is the median of the posterior distribution, while the lower and upper error bars correspond to the 16th and 84th percentile values. The SFRs are the instantaneous values at the latest age bins, which have width 5773 yr. These values, as well as the stellar masses, are not corrected for the lensing magnifications, which (based on the J. M. Diego et al. 2024 lens model) are μ = 1.21 for Y1, μ = 2.26 for JD and JD-N, and μ = 1.32 for f090d_018.

with a large uncertainty. The most distinct difference among them is in their SFRs, which span 2 orders of magnitude: after correcting the magnifications, the SFRs range from 3.6 (f090d_018) to 200 M e yr -1 (Y1).

Consistent with their young ages, the galaxies have low stellar masses of 0.6-7.7 × 10 8 M e (after correcting for lensing magnifications). These can be achieved by the derived SFRs in the derived ages. (See also Figure 6.) Interestingly, f090d_018, which has the lowest intrinsic stellar mass (5.7 × 10 7 M e ) and SFR (3.6 M e yr -1 ), also has the smallest amount of dust extinction, A V ≈ 0.06 mag.

The existence of the oxygen lines in these objects (as well as the [C II] 157.7 μm line in Y1) indicates that they have already acquired substantial amounts of metals. This is supported by the metallicities derived from the SED analysis. Given the young ages for Y1, JD, and JD-N (Z/Z e ≈ 0.3-0.4), the only route for their metal enrichment was through core-collapse supernovae. On the other hand, the metallicity estimate for f090d_018 (

 Z Z 0.04 0.02 0.03 = -+

) is rather low, which appears consistent with its low SFR.

Ultraviolet Slopes and Lyman-continuum Photon Escape Fractions

Using the photometry reported in Table 1, we calculated the rest-frame UV slope β for these sources. The calculation was done using two different filter pairs, (F150W, F200W) and (F150W, F277W), which approximately correspond to using the rest-frame wavelength pairs of (1600 Å, 2100 Å) and (1600 Å, 3000 Å), respectively. As an alternative, we also used the best-fit template spectrum (corresponding to the 50th percentile of the posterior distribution) of the SED and fitted a power law between 1300 and 1850 Å in the rest frame to obtain this slope. These results are listed in Table 4. The slopes derived using the best-fit template spectra (β PL ) are all significantly flatter than those obtained using the photometry in the blue passband pairs (β 150, 200 and β 150, 277 ), which cautions that care must be taken when interpreting UV slopes obtained by different methods (e.g., J. S. Dunlop et al. 2013;D. Austin et al. 2024).

The UV slopes are related to the LyC photon escape fractions ( f esc ), which are critical in understanding the sources of the cosmic hydrogen reionization. For each β value, we calculated f esc following Equation (11) of J.

Chisholm et al. (2022), ( ) • ( ) ( ) f 1.3 0.6 10 10 . 1 esc 4 1 . 2 2 0 . 1 =  ´b --

The derived values are reported in Table 4. Somewhat surprisingly, only JD and f090d_018 have blue UV slopes β  -2.0. These two objects' inferred f esc values are between ∼6% and 31% depending on the adopted β. We believe that this can be attributed to their much less dust extinction (A V ∼ 0.3 and 0.06 mag) as compared to those of Y1 and JD-N (A V ∼ 0.9 and 1.1 mag). The existence of two objects with high extinction and low escape fraction suggests that a significant fraction of very young galaxies at high redshifts might not contribute to reionization because they formed a considerable amount of dust that dispersed throughout their bodies in only a few megayears and which blocks UV emissions.

Summary

New JWST data have identified four very young galaxies at z = 8.31-8.49 behind the lensing cluster MACS0416. Three of them are strong [O III] λ5007 emitters, and their line EWs can easily explain their apparent flux excess in F444W. Two galaxies, Y1 and JD, had prior spectroscopic redshifts from ALMA detections of the [O III] 88 μm line and/or [C II] 157.7 μm line. We confirm the redshift 8.31 for Y1 but find z = 8.34 for JD as opposed to the previously claimed z = 9.28. Y1 is a merging system that is resolved into three components, extending ∼3.4 kpc along the long axis. Another object, JD-N, is a previously discovered z > 8 candidate that is now confirmed to have the same redshift as JD; the two are only ∼0 51 (∼2.4 kpc) apart in the source plane and therefore very likely merging or at least interacting. The final object, f090d_018 at z = 8.49, is newly identified. These four objects spread over a projected distance of ∼165 kpc in the source plane and a radial distance of ∼49 Mpc, and it is likely that they are part of an overdensity on a filament-like structure.

The four objects are magnified by less than a magnitude and have intrinsic M UV ranging from -19.57 to -20.83 mag. Our SED analysis shows that they are all very young systems in the making, with ages of ∼4-18 Myr, arguably the youngest galaxies ever reported at z > 6. These infant galaxies have stellar masses on the order of 10 7-8 M e , and their SFRs range from a few to over a hundred solar masses per year. However, only JD and f090d_018 have blue rest-frame UV slopes (β ranging from -2.2 to -2.8 depending on how it is derived) indicative of high LyC photon escape fractions ( f esc possibly as high as ∼31%). This is largely due to their small dust extinction (A V = 0.06-0.29 mag). Interestingly, these two objects are the least massive (M * = 5.7 and 5.8 × 10 7 M e ) and least active (SFR = 3.6 and 13.6 M e yr -1 ) ones among the four. In contrast, Y1 and JD-N, despite their much higher SFRs, have β  -1.8 to -1.0 (implying f esc  2%) because of their higher dust extinction (A V = 0.92-1.10 mag). This suggests that even very young, very actively star-forming galaxies at high z could have a negligible contribution to the ionizing background if they form dust throughout their bodies too quickly (over a few megayears timescale). Y1 and JD-N are examples that dust formation and pollution processes at z > 8 could indeed be very fast. Note. The UV slope β is defined via f ν ∝ λ β+2 . The numerical superscripts on f esc indicate the passband pairs (F150W, F200W) and (F150W, F277W), respectively, while "PL" indicates the power-law fit to the best-fit template spectrum.

Figure A1. SED fitting results for f090d_018 with redshift as a free parameter to derive z ph . The plots are similar to those in Figures 5 and 6.

Table B1 compares the best-fit stellar mass and age parameters obtained using the τ model with those obtained using these alternative SFHs. For both the τ and the delayed τ models, the age parameter is among the direct outputs. This is not the case when using the lognormal or the double power-law models, however, because there is no clear definition of age in either SFH. To obtain an estimate in these two cases, we derived an age proxy denoted as Age N , which is the time for the galaxy to gain 100% of its stellar mass starting from the time when it had N% of its total stellar mass. For the demonstration purpose here, N was set to 10, 50, and 90. As shown by this comparison, these alternative SFHs resulted in comparable, young ages for all three objects.

Table B1 Comparison of Best-fit Galaxy Properties with Different Star Formation Histories Exponential Delayed τ Lognormal Double PL Y1 ( ) ☉ M M log m 8.97 0.01 0.03 -+ 8.90 0.10 0.07 -+ 8.84 0.07 0.05 -+ 8.85 0.07 0.05 -+ Age/Myr 4.76 0.35 0.28 -+ 3.51 1.82 1.68 -+ ... ... Age 10 /Myr ... ... 4.10 3.33 Age 50 /Myr ... ... 1.39 1.29 Age 90 /Myr ... ... 1.00 1.00 JD ( ) ☉ M M log m 8.12 0.11 0.12 -+ 8.13 0.11 0.18 -+ 8.13 0.10 0.09 -+ 8.15 0.08 0.08 -+ Age/Myr 4.28 1.17 1.06 -+ 5.89 2.60 2.36 -+ ... ... Age 10 /Myr ... ... 6.42 5.13 Age 50 /Myr ... ... 2.05 2.01 Age 90 /Myr ... ... 1.00 1.00 JD-N ( ) ☉ M M log m 8.95 0.04 0.05 -+ 8.99 0.06 0.05 -+ 9.02 0.06 0.06 -+ 8.99 0.05 0.05 -+ Age/Myr 8.69 2.19 2.46 -+ 12.77 4.42 7.24 -+ ... ... Age 10 /Myr ... ... 10.86 7.85 Age 50 /Myr ... ... 3.36 2.82 Age 90 /Myr ... ... 1.00 1.00 f090d_018 ( ) ☉ M M log m 7.88 0.17 0.19 -+ 7.89 0.17 0.21 -+ 7.84 0.13 0.18 -+ 7.91 0.17 0.21 -+ Age/Myr 17.79 7.56 15.50 -+ 33.11 15.59 37.20 -+ ... ... Age 10 /Myr ... ... 21.50 15.84 Age 50 /Myr ... ... 6.62 8.81 Age 90 /Myr ... ... 1.05 1.81

Note. The ages for the lognormal and double power law (PL) SFHs are not available from BAGPIPES direct outputs because the parameter has no clear meaning in either SFH. We define Age N for these two cases, where N = 10, 50, and 90. See text for details.