A Strong He ii λ1640 Emitter with an Extremely Blue UV Spectral Slope at z = 8.16: Presence of Population III Stars?

Cosmic hydrogen reionization and cosmic production of the first metals are major phase transitions of the Universe occurring during the first billion years after the Big Bang; however, these are still underexplored observationally. Using the JWST/NIRSpec prism spectroscopy, we report the discovery of a sub-L * galaxy at z spec = 8.1623 ± 0.0007, dubbed RX J2129–z8He II, via the detection of a series of strong rest-frame UV/optical nebular emission lines and the clear Lyman break. RX J2129–z8He II shows a pronounced UV continuum with an extremely steep (i.e., blue) spectral slope of β=−2.53−0.07+0.06 , the steepest among all spectroscopically confirmed galaxies at z spec ≳ 7, in support of its very hard ionizing spectrum that could lead to a significant leakage of its ionizing flux. Therefore, RX J2129–z8He II is representative of the key galaxy population driving the cosmic reionization. More importantly, we detect a strong He II λ1640 emission line in its spectrum, one of the highest redshifts at which such a line is robustly detected. Its high rest-frame equivalent width (EW = 21 ± 4 Å) and extreme flux ratios with respect to UV metal and Balmer lines raise the possibility that part of RX J2129–z8He II’s stellar population could be Pop III (Pop III)-like. Through careful photoionization modeling, we show that the physically calibrated phenomenological models of the ionizing spectra of Pop III stars with strong mass loss can successfully reproduce the emission line flux ratios observed in RX J2129–z8He II. Assuming the Eddington limit, the total mass of the Pop III stars within this system is estimated to be 7.8 ± 1.4 × 105 M ⊙. To date, this galaxy presents the most compelling case in the early Universe where trace Pop III stars might coexist with metal-enriched populations.


Introduction
The identification and characterization of the first generation of stars (Pop III) are of paramount importance, as this quest provides key insights into the stellar evolution physics in the early Universe (Bromm 2013).The observational signatures expected from early galaxies dominated by these young, metal-free, massive Pop III stars are highly associated with very hard ionizing flux, given the high effective temperature of Pop III stars and a top-heavy initial mass function (IMF; Nakajima & Maiolino 2022;Zackrisson et al. 2023).This flux can ionize the interstellar medium (ISM) and produce strong hydrogen and helium emission lines, such as strong Lyα, Balmer lines (e.g., Hα, Hβ), He II λ1640 (:=He II), and He II λ4686, without the corresponding metal emission lines (Schaerer 2002(Schaerer , 2003;;Zackrisson et al. 2011).It is believed that Pop III galaxies with stellar mass (M * ) ∼ 10 4-5 M e consisting of a pure Pop III stellar population and no metal-enriched ones are extremely hard to detect, even with JWST, unless they are substantially magnified by gravitational lensing (Windhorst et al. 2018;Vikaeus et al. 2022).More massive Pop III galaxies with M *  10 6 M e can be first photometrically selected due to their blue UV continuum and spectroscopically confirmed by their high-ionization nebular emission lines and nebular continuum (Trussler et al. 2023).Yet they should be really rare, with their characteristic nebular emission features fading within a short timescale of ∼10 Myr at z ∼ 8-9 (Katz et al. 2023).
Recent cosmological hydrodynamic simulations show that Pop III stars can also continue to form at later epochs, provided that pockets of primordial gas can be preserved during cosmic evolution.This is possible either in halos that gain their gas from regions not yet polluted by outflows from nearby starforming galaxies or in halos whose progenitors had a suppression of star formation (Ciardi & Ferrara 2005;Venditti et al. 2023).This naturally leads to a type of "hybrid" Pop III galaxies, which have better chances to be identified using current facilities.One important clue to such "hybrid" Pop III stellar populations is the existence of the strong He II line, with a rest-frame equivalent width (EW) as high as 100-150 Å, powered by the hard ionizing background radiation from the Pop III stars.
At low redshifts, such broad nebular He II emission is very rare, and its spectral properties are usually consistent with the ionizing sources of Wolf-Rayet (W-R) stars, stripped stars, X-ray binaries, or active galactic nuclei (AGNs; see, e.g., Berg et al. 2016;Nanayakkara et al. 2019;Senchyna et al. 2019;Saxena et al. 2020).In addition, this He II emission is usually accompanied by other prominent high-ionization lines, such as N V λ1240, C IV λλ1548, 1551 (:=C IV), and C III] λλ1907, 1909 (:=C III]), of similar excitation energy (Sobral et al. 2015;Stark et al. 2015).Therefore, a clear path moving forward to identify Pop III stars in and beyond the epoch of reionization (EoR) is to search for pronounced He II emission without the accompanying rest-frame UV metal lines that are often observed in the low-redshift Universe.
In this work, we present our detailed analysis of a galaxy (dubbed RX J2129-z8He II) spectroscopically confirmed at z spec = 8.1623 ± 0.0007, showing strong He II emission and an absence of any UV metal lines-an excellent candidate likely having a mixture of Population II and Pop III stars.The overall physical properties of RX J2129-z8He II are showcased in Table 1.This identification benefits from the synergy of the powerful JWST spectroscopy and the lensing magnification boost from a foreground massive cluster of galaxies.
This Letter is orchestrated as follows.First, in Section 2, we describe the novel JWST observations analyzed in this work.Then, in Section 3, we briefly explain the lens model used to correct for the cluster magnification of our source.The details of our analysis methods and primary scientific results are given in Section 4. Finally, we summarize our main findings in Section 5. Throughout this Letter, we adopt the standard concordance ΛCDM model with Ω m = 0.3, Ω Λ = 0.7, H 0 = 70 km s −1 Mpc −1 , and the AB magnitude system (Oke & Gunn 1983).

Observations and Data Reduction
The JWST observations analyzed in this work were acquired by a Director's Discretionary program (DD-2767; PI: P. Kelly) targeting the field of the galaxy cluster RXJ 2129.7+0005(henceforth RXJ 2129) at z = 0.234.The primary goal of DD-2767 is to measure the light curves and spectra of a strongly lensed supernova, SN 2022riv, at z = 1.52 that Kelly et al. (2022aKelly et al. ( , 2022b) ) discovered in the RXJ 2129 field by a Hubble Space Telescope (HST) SNAP program (GO-16729; PI: P. Kelly).These observations consist of an imaging component using the Near-Infrared Camera (NIRCam) and a spectroscopic component using the Near-Infrared Spectrograph (NIRSpec), which are described in Sections 2.1 and 2.2, respectively.

JWST/NIRCam Data Reduction and Mosaicking
The imaging exposures were taken on 2022 October 6 (UT dates quoted throughout) using the F115W, F150W, F200W, F277W, F356W, and F444W filters, covering the wavelength range of λ obs ä [1, 5] μm.The exposure times are 4982 s for the F150W/F356W filters and 2061 s for the other four NIRCam filters (F115W, F200W, F277W, and F444W), equivalent to a 5σ limiting depth of ∼29 mag.We reduce these NIRCam images using the standard JWST pipeline version 1.9.4 and making use of the JWST_1040.PMAP context.We follow the standard three stages of the image data reduction pipeline CALWEBB_DETECTOR1, CALWEBB_IMAGE2, and CALWEBB_IMAGE3 to calibrate each image and build the mosaic image of each band.Furthermore, we remove the "snowball"18 and the 1/f noise using the scripts from https://github.com/chriswillott/jwst.The "wisps" features 19 are removed using the default wisps templates. 20We produce the imaging mosaics at the 40 and 20 mas pixel scales, astrometrically aligned to the Gaia DR2 astrometry frame.
In Figure 1, we show a color-composite image produced from these data, where our source of interest (RX J2129-z8He II) is highlighted in the inset zoom-in panel.Our source in fact has two Emission Line Analyzes (Section 4.5) EW He II [Å] 21 Note.The parameter values are given in median posterior constraints followed by a 1σ CI.The lensing magnification is predicted from the up-to-date GLAFIC lens model (see Section 3), which is accounted for in the estimates of these physical properties when necessary.All reported line flux ratios and limits have been corrected for dust extinction.
components, A and B, with A dominating the total flux in all filters in the NIRCam long-wavelength channel.Hereafter, unless otherwise specified, we use RX J2129-z8He II to designate its component A, which is targeted by NIRSpec as described below.We also take advantage of the archival HST imaging obtained by the Cluster Lensing and Supernova Survey21 (CLASH; PI: M. Postman; Postman et al. 2012).CLASH provides the imaging data in eight filters taken with the Advanced Camera for Surveys (ACS)/WFC (F775W, F814W, and F850LP) and WFC3/IR (F105W, F110W, F125W, F140W, and F160W).We compile a photometric catalog combining both the new NIRCam imaging and the existing HST imaging in a selfconsistent manner (see Section 4.1).The image stamps of our target are displayed in the top row of Figure 2.

JWST/NIRSpec Microshutter Array Spectroscopy and Data Reduction
The follow-up JWST spectroscopy was carried out on 2022 October 22 using the NIRSpec instrument in the Multi-Object Spectroscopy (MOS) mode.RX J2129-z8He II was among the NIRSpec targets that were preselected based on their photometric redshifts (Williams et al. 2023).The JWST/NIRSpec observation was carried out using the prism disperser, which offered a very wide, continuous wavelength coverage of λ obs ä [0.6, 5.3] μm with a resolution of R = λ/Δλ ∼ 50-400.In the inset of Figure 1, we show a zoom-in view of RX J2129-z8He II with the position of the MOS slit superposed.The slit is composed of three shutters in the microshutter array (MSA), resulting in a total size of ∼0 2 × 1 4. The observation adopts the standard three-point nodding pattern to facilitate background subtraction and a total exposure time of 4464.2 s.
Figure 1.A color-composite NIRCam image of the RX J2129 galaxy cluster.We use the F115W+F150W, F200W+F277W, and F356W+F444W imaging as the blue, green, and red colors, respectively.The regions of formally infinite magnification at z = 8.16 predicted by the up-to-date GLAFIC lens model (Oguri 2010(Oguri , 2021) ) are represented by the magenta curves.The location of RX J2129-z8He II is highlighted by the red circle.It is estimated to have a magnification factor of μ = 2.26 ± 0.14 (see Section 3).The inset shows a zoom-in view (3″ × 3″) of RX J2129-z8He II.Within this inset, we also show the orientation of the NIRSpec slit (0 2 × 1 4), constructed by three MSA slitlets, through which the prism spectroscopy of RX J2129-z8He II is taken by the JWST Director's Discretionary program (DD-2767; PI: P. Kelly).
We reduce the NIRSpec MSA spectroscopic data following our customized procedures.Our reduction uses the context file JWST_1040.PMAP, which contains up-to-date reference files for NIRSpec.This version of the reference files is considered to be a major update to the NIRSpec S flats, F flats, and readnoise files, which have been calibrated using in-flight data.First of all, the level 1 CALWEBB_DETECTOR1 calibration pipeline is employed to reduce the raw exposures (UNCAL.FITS) into count rate maps (CRMs; i.e., RATE.FITS), during which the detector artifacts and cosmic rays are flagged and removed.The resulting CRMs are visually inspected, and any remaining artifacts and bad pixels are manually masked.We then use the custom reduction software MSAEXP22 to perform the remaining steps of the reduction.It first preprocesses the rate images to equalize the pedestal of each exposure and identify and correct the "snowball" and 1/f noise features.Then it calls some specific modules from the standard level 2 CALWEBB_SPEC2 calibration pipeline to carry out the bulk of the data reduction.These modules include ASSIGNWCS, EXTRACT2DSTEP, FLATFIELDSTEP, PATHLOSSSTEP, and PHOTOMSTEP, which perform world coordinate system initialization, 2D extraction of spectra, slit-level flat-fielding, path-loss correction, and wavelength and flux calibration of each science exposure.
MSAEXP subsequently subtracts background from each exposure using the two associated exposures at the other dithered positions from the nodding sequence and drizzles the background-subtracted science exposures onto a common wavelength grid via inverse-variance weighting, with outlier rejection and bad pixel masking.We adopt an oversampling rate of 2 for the output wavelength grid to Nyquist-sample the line-spread function in order to improve the sampling of the emission line profiles.Finally, MSAEXP obtains the 1D spectrum from the 2D combined spectral traces following an optimal extraction methodology (Horne 1986), which uses the actual light profile of our target as the optimal aperture for spectral extraction.The object light profile along the crossdispersion direction is modeled as a Gaussian function fit to the collapsed 2D spectra (summed along the dispersion direction).
As a result, we obtain the 2D/1D prism spectra, as shown in the middle/bottom panels of Figure 2, covering an uninterrupted wide wavelength range of λ obs ä [0.6, 5.3] μm.For z ∼ 8 galaxies, this corresponds to contiguous coverage of the rest-frame UV and optical spectral energy distribution (SED) at λ rest ä [700, 6000] Å.We detect pronounced emission features of the [O III] doublets and Hβ line complex at λ obs ∼ 4.5 μm and a clear continuum break at λ obs ∼ 1.1 μm with prominent continuum flux redward of the break.The identification of this continuum break as the Lyman break occurring at the Lyα wavelength (λ rest = 1216 Å) is in excellent agreement with the unambiguous [O III]+Hβ emission feature in the prism spectrum, leading to an accurate spectroscopic confirmation of RX J2129-z8He II at z spec = 8.1623 ± 0.0007, when the Universe is ∼613 Myr old.This puts RX J2129-z8He II deep in the EoR, when the intergalactic medium (IGM) was mostly neutral hydrogen (Fan et al. 2006;Stark 2016).
We also perform extensive tests to verify the calibration of NIRSpec prism spectra.On one hand, we cross-check the flux levels of our extracted spectra and that of our broadband photometry from NIRCam imaging in the similar wavelength regime.They show good agreement with each other (within 10%; see Figure 3).On the other hand, we conduct detailed investigation of the wavelength calibration by performing line identifications in low-z galaxies showing multiple emission lines across their entire prism spectra, observed in the same MSA mask and extracted in the same fashion.We find that the wavelength offset between the best-fit centroid and that expected for the mean redshift shows a scatter of <0.004 μm in the observer frame, much smaller than the instrument resolution.So we conclude that our flux and wavelength calibrations are sufficiently accurate.

Strong Lensing Models
We adopt the cluster lens model of RX J2129 constructed using the GLAFIC software (Oguri 2010(Oguri , 2021)), updated from the initial work of Okabe et al. (2020).In total, 22 individual images of seven multiply lensed background galaxies spectroscopically confirmed by Caminha et al. (2019) are used as the strong lensing constraints.The macroscopic mass model consists of one cluster-scale dark matter halo in the elliptical Navarro-Frenk-White profile (Navarro et al. 1997) combined with galaxy-scale halos modeled using the pseudo-Jaffe profile according to the scaling relation of the velocity dispersion and luminosity of cluster member galaxies (Kawamata et al. 2018).The best-fit model is derived using the Markov Chain Monte Carlo (MCMC) sampling process with a χ 2 minimization assuming a positional uncertainty of 0 4. A hundred additional realizations are also created to bootstrap the 1σ statistical uncertainties for magnifications.As a consequence, we obtain the best-fit and 1σ confidence interval (CI) of the magnification estimates of RX J2129-z8He II to be 2.26 and [2.12, 2.40].We also double-check with an independent lens model built by the Lenstool software (Caminha et al. 2019) and derive consistent results.

JWST/NIRCam and HST/ACS-WFC3 Photometry
For comprehensive photometry of RX J2129-z8He II, we not only use the six-band imaging mosaics produced from JWST/ NIRCam data but also include the publicly released eight-band HST imaging mosaics from CLASH (Postman et al. 2012).We first transform the CLASH 0 065 imaging mosaics to mosaics on a 0 04 plate scale and then point-spread function (PSF) match the eight HST ACS/WFC3 filters and six JWST/NIRCam filters to the F444W resolution (FWHM = 0 14).We utilize a window function to remove the high-frequency noise in the Fourier domain.We modify the rms images in order to match the F444W resolution for every other filter using , where RMS F444W,matched and RMS i are the rms after PSF matching and the original rms in filter i, respectively, and W i,F444W is the kernel used to match all PSFs to that of the F444W filter through PSF F444W,matched = PSF i * W i,F444W .
In the JWST/NIRCam F200W image, RX J2129-z8He II shows two components, "A" and "B," which are marked by the arrows in the F200W stamp shown in Figure 2. The JWST/ NIRSpec MOS slit falls on component A, which is also predominant in the total flux in the long-wavelength NIRCam filters.We use the F200W image as the detection image and obtain photometry for the two components separately.To maximize the detection S/Ns, fluxes in each filter are measured within fixed apertures of 2 × FWHM of the F444W image (0 28) in diameter.Furthermore, the flux is dust-corrected for galactic extinction through ´, where A i is the extinction in the ith band (Cardelli et al. 1989).The resultant photometric measurements of the two components of RX J2129-z8He II are presented in Table 2.For filters with no detections, we report the 2σ upper limits.The NIRCam photometry normalized at the individual F200W flux is shown in Figure 4.

Photometric Redshift Estimates
We use the EAzY software (Brammer et al. 2008) to estimate the photometric redshift from the broadband photometry presented in Table 2.All together, we use the HST filters ACS/F775W, ACS/F814W, ACS/F850LP, WFC3/F105W, WFC3/F110W, WFC3/F125W, WFC3/F140W, and WFC3/ F160W, as well as the JWST/NIRCam filters F115W, F150W, F200W, F277W, F356W, and F444W.We adopt the standard set of galaxy SED templates EAZY_V1.1_LINES.SPECTRA.PARAM, which includes star-forming galaxies with strong emission lines.The resulting photometric redshifts of components A and B are z phot = 8.35 ± 0.29 and z phot = 8.93 ± 0.91, respectively (see Figure 5), in good agreement with the spectroscopic redshift that we derive from the NIRSpec prism spectroscopy.As shown in Figure 5, the model and observed photometry for component A of RX J2129-z8He II shows high consistency, except in filters NIRCam/F150W and WFC3/ F160W, indicating strong emission features at λ obs ∼ 1.5 μm.Component B shows a similar excess in observed photometry in the bluest NIRCam bands.

Source Morphology
As shown in the multiple image stamps displayed in the upper panels of Figure 2, the entire RX J2129-z8He II galaxy consists of two components, A and B, with the former dominating the flux in the rest-frame optical probed by the NIRCam long-wavelength channels, and the latter clearly manifesting in the rest-frame UV covered by the shortwavelength channels, especially F150W.Therefore, our morphology analyzes are performed in two NIRCam bands, namely, F150W and F444W, representing the rest-frame UV and optical light profiles of our galaxy, respectively.We use GALFIT (Peng et al. 2002) to model both components in F150W.First, two models are chosen to fit component A simultaneously.One is a Sérsic profile describing the nucleated emission in the center with a Sérsic index (n) of 4, an effective radius (R e ) of 0 048, an axis ratio (b/a) of 0.85, and a position angle (θ) of −48°.8.An exponential disk model is adopted to fit the underlying extended structure with a disk scale length (R s ) of 0 08, b/a of 0.26, and θ of −50°.The inset of Figure 1 shows the position of the NIRSpec MOS slit with respect to RX J2129-z8He II, which primarily covers its component A. As implied by the observed prism spectra shown in the upper panel of Figure 2, the broadband flux of the NIRCam F444W filter is likely dominated by the high-EW [O III]+Hβ nebular emission lines (also see Table 3).Indeed, we observe a more extended light profile in the restframe optical than in the UV.We thus model the image of component A in F444W with an exponential disk and a Sérsic profile.The diffuse and extended structure is well reproduced by an exponential disk model with R s of 0 12, b/a of 0.3, and θ of −45°.8.The nucleated structure is modeled by a Sérsic profile (n = 4) with R e = 0 068, b/a = 0.88, and θ = 5°.9.Component B in F444W is very faint and can be well fitted by an exponential disk model with R s of 0 04, b/a of 0.35, and θ of −40°.9.The resulting model residual is shown in the lower right panel of Figure 6.

Spectrophotometric Analyzes
We employ multiple independent methods to conduct detailed spectrophotometric analyzes of both the broadband photometry and the full NIRSpec spectrum that we obtained for RX J2129-z8He II.Performing full spectrum fitting is highly critical since it is the only appropriate approach to extract the detailed information of stellar population from the lowresolution prism spectroscopy.We first employ the BAGPIPES software (Carnall et al. 2018) to perform SED fitting of both our spectroscopic and photometric data simultaneously.We adopt a widely used, double power-law star formation history (SFH) model to describe the evolution of the cosmic star formation rate (SFR) density (Madau & Dickinson 2014), which has the flexibility to account for both rising and declining star formation activities:   BAGPIPES utilizes the nested sampling algorithm to perform efficient Bayesian inference and obtain posterior distribution of parameters.We take the [16th, 50th, 84th] percentiles of the parameter posteriors as the median and 1σ CI and show the 1D posterior distributions of some key parameters in Figure 7.The overview of the physical properties of RX J2129-z8He II is presented in Table 1.After correcting for magnification (μ = 2.26 ± 0.14), we obtain the following physical picture of RX J2129-z8He II: it is a very low-mass (M * ∼ 10 7.8 M e ), young (t age ∼ 216 Myr) galaxy, actively forming stars (SFR ∼ 9.6 M e yr -1 ) with subsolar metallicity ( * ( ) Z Z log 0.9 ~- ) and little dust (A V ∼ 0.1).We recover a fast-rising SFH of RX J2129-z8He II, with a predominant young stellar population of 10-50 Myr old.The SFH recovered by BAGPIPES is shown in Figure 8, with the red curve and shaded region marking the best-fit and 1σ confidence range, respectively.
We follow the prescription in Alavi et al. (2014) to estimate the absolute UV magnitude, where μ is the lensing magnification factor expressed in magnitude units.At z = 8.16, the WFC3/F140W filter covers the rest-frame λ rest ä [1300, 1800] Å and is therefore appropriate to use.The UV absolute magnitude is thus calculated to be M UV ∼ −19.6 mag, which places RX J2129-z8He II at the level of 40% L * at z ≈ 8 (Bouwens et al. 2015).Importantly, the bulk of the ionizing UV photons that caused the cosmic reionization are thought to come from sub-L * systems (Yan & Windhorst 2004;Finkelstein et al. 2019).RX J2129-z8He II shows pronounced detection of the far-UV continuum.We follow the standard formalism to derive its UV spectral slope (β).We fit f λ ∝ λ β in the wavelength range of λ rest ä [1300, 2600] to the model galaxy SEDs produced by BAGPIPES with all possible emission features properly masked.The resulting UV slope is measured to be 2.53 0.07 . We obtain consistent results if measuring β using the smoothed prism spectrum directly observed.This classifies RX J2129-z8He II among the spectroscopically confirmed galaxies in the EoR having the steepest UV continuum slope-a strong implication of significant leakage of its ionizing radiation to the IGM (Bouwens et al. 2010;Zackrisson et al. 2013).The completion of reionization by z ∼ 6 requires that the absolute escape fraction of the Lyman continuum (LyC) photons from galaxies should be f 10% esc LyC  , on average (Finkelstein et al. 2012;Robertson et al. 2015).In the EoR, a negative correlation between M UV and β has been seen from photometric analyzes, in support of We therefore obtain an estimate of f 0.16 0.03 esc LyC =  for RX J2129-z8He II.
Figure 9 summarizes the currently available measurements of β and M UV of galaxies at z  7, among which the spectroscopically confirmed sources are highlighted in color (Watson et al. 2015;Zitrin et al. 2015;Hashimoto et al. 2018;Curtis-Lake et al. 2023;Morishita et al. 2023;Roberts-Borsani et al. 2023;Schaerer et al. 2022;Williams et al. 2023;Bunker et al. 2023), and the photometric redshift (photo-z) selected candidates are represented in gray (Tacchella et al. 2022;Castellano et al. 2022;Naidu et al. 2022;Finkelstein et al. 2022;Topping et al. 2022;Cullen et al. 2023).While steeper slopes (β ∼ − 3) have been reported for photo-z and Lymanbreak selected galaxies (see also Bouwens et al. 2010), RX J2129-z8He II has by far the steepest UV spectral slope and thus the largest LyC escape fraction among all the spectroscopically confirmed galaxies at z spec 7. Note that throughout this work, β stands for β obs in the nomenclature used by Chisholm et al. (2022) without dust correction.After correcting for dust, Chisholm et al. (2022) find that the intrinsic UV spectral slopes (β int ) of the LzLCS galaxies observed with HST/COS are clustered in the range of [−2.8, −2.6].
We also double-check the derived stellar population properties of RX J2129-z8He II by modeling its broadband photometry and spectra using the PROSPECTOR software (Johnson et al. 2021).
PROSPECTOR is built upon a fully Bayesian framework that makes it possible to fit the SED of high-redshift galaxies with complex SFH models.Regarding the basic fitting setup, we adopt the Flexible Stellar Population Synthesis code (Conroy et al. 2009;Conroy & Gunn 2010) using the MIST stellar isochrone libraries (Choi et al. 2016;Dotter 2016) and the MILES stellar spectral libraries (Falcón-Barroso et al. 2011).During the modeling, we use the MCMC sampling code dynesty (Speagle 2020) based on the nested sampling technique.We assume the Kroupa IMF (Kroupa 2001) for consistency with the IMF used in the nebular continuum from Byler et al. (2017) and the line emission model that we adopt in our SED modeling.We adopt the IGM absorption model from Madau (1995).Following Charlot & Fall (2000), we assume a two-component dust attenuation model, where the dust attenuation of nebular emission and young stellar populations and of old stellar populations are treated differently.
We test our results with two different parameterizations of the galaxy's SFH.First, we assume a piecewise, nonparametric form composed of six lookback time bins (i.e., the time prior to the time of observation), where SFR is fitted as a constant in each bin.Among the six lookback time bins, the first two bins are fixed to be 0-5 and 5-10 Myr in order to capture recent star formation activities; the last bin is assumed to be t(z = 20) − t H , where t H is the Hubble time at z = 8.1623; and the remaining three bins are evenly spaced in the logarithmic space of 10 Myr t(z = 20).A similar setup of nonparametric SFH has been used in recent literature about high-redshift galaxy SED fittings (Robertson et al. 2023;Tacchella et al. 2023).During the nonparametric SFH fitting, we adopt a continuity prior that has been demonstrated to work well across various galaxy types using mock observations from cosmological simulations (Leja et al. 2019).The recovered SFH is represented by the black histogram shown in Figure 8 with the 1σ confidence region marked by the dashed histograms.We note that the continuity prior is equivalent to a symmetric prior on stellar age and a constant SFH prior with SFR(t) = M * /t H .The fact that the fitting converges to a monotonically rising SFH suggests that the data strongly favor a very young stellar age of RX J2129-z8He II, hence fully in line with its strong He II detection.
Second, we assume a parametric delayed-τ SFH model, i.e., ( ) ( ) t µ ´twith τ denoting the peak time of star formation sampled in the range of τ ä (0, t H ). The resulting SFH is shown as a blue curve and shaded regions in Figure 8 as well.Using two independent techniques (i.e., BAGPIPES and PROSPECTOR) under three sets of different SFH assumptions (i.e., double power-law, nonparametric, and delayed-τ models), we achieve a selfconsistent conclusion as follows.The SFH of RX J2129-z8He II is sharply rising, and there exists a predominant stellar population 10-50 Myr old.This young stellar population is highly likely to be responsible for the strong He II line identified in the rest-frame UV spectra of RX J2129-z8He II.

Emission Line Fitting and Diagnostics
We utilize the PPXF software (Cappellari 2023) to perform accurate emission line analyzes, fixing the source redshift to z = 8.1623.To be self-consistent, we rely on the BC03 stellar  b The intrinsic line fluxes and upper limits after the corrections of dust extinction and lensing magnification (μ = 2.26 ± 0.14).We adopt A V = 0.12 ± 0.04 measured from our spectrophotometric analyzes to deredden the line fluxes since we only achieve upper limits on the nebular dust extinction.
population library (Bruzual & Charlot 2003), generated using the Chabrier IMF (Chabrier 2003) and dust-corrected with the Calzetti reddening curve (Calzetti et al. 2000).The best-fit continuum model is shown as the red curve in the upper panel of Figure 3.We then fit Gaussian profiles to nebular emission features at known wavelength positions, given a fixed redshift.During the emission line profile fitting, the range of their intrinsic width is set to [0, 300] km s −1 .The flux ratio of the [O III] λλ4959, 5007 doublet is fixed to the theoretical value of 0.33.For the weak emission lines in the UV bands, we bind their shift of line centers to those stronger emission lines (e.g., He II) to improve the fitting robustness.We avoid fitting for any line emission blueward of λ rest ∼ 1400 Å due to strong attenuation by the highly neutral IGM (Dijkstra 2014).
The best-fit emission line profiles are demonstrated with orange curves in the upper panel of Figure 3.The resulting emission line fluxes (observed and intrinsic) and EWs23 are listed in Table 3.Here the statistical uncertainties of the line fluxes and EWs are obtained from bootstrapping the fitting process to 50 mock spectra generated using the best-fit model spectrum and flux errors.
At an S/N threshold of 3, we detect six emission lines: ), and He II, among which [O III], Hβ, and He II are detected with an S/N of 37, 7, and 5, respectively.In addition, we give 2σ upper limits for other lines in Table 3, including Figure 8. SFH of RX J2129-z8He II, obtained from three independent methods.The SFH given by our default BAGPIPES full spectrum fitting under the assumption of a double power-law model is represented by the red curve (best-fit value) with shaded regions (1σ confidence ranges).The SFHs given by the PROSPECTOR SED fitting assuming nonparametric and delayed-τ models are shown in black histograms and blue curves, respectively.All three SFHs returned by different assumptions and techniques reach excellent consensus that the SFH of RX J2129-z8He II is rapidly rising, in support of the existence of a very young stellar component likely associated with the Pop III signatures seen in the JWST/NIRSpec spectra.
He II λ4686, [O III] λ4363, C III], and C IV.We do not see any signature of Lyα emission and put a 2σ upper limit of EW Lyα < 10 Å.The absence of a Lyα line is likely due to the damping wing opacity caused by the diffuse neutral IGM at z ∼ 8 and/or dense self-shielding H I gas clouds inside the large ionized H II bubbles (Dijkstra 2014).
We apply a well-established Bayesian forward-modeling inference framework to the measured line fluxes summarized in Table 3.We constrain jointly three key properties of the ISM: gasphase metallicity ( ( ) 12 log O H + ), nebular dust extinction (A V N ), and the instantaneous SFR converted from the dereddened Hβ flux ( f Hβ ) corrected for the magnification factor of μ = 2.26 ± 0.14.Note that ( ) 12 log O H + and A V N are unaffected by lensing magnification.The likelihood function is defined following Wang et al. (2017Wang et al. ( , 2019Wang et al. ( , 2020Wang et al. ( , 2022aWang et al. ( , 2022b) ) as Here f EL i and EL i s correspond to the intrinsic line ([O III], Hβ, Hγ, Hδ, and [O II]) flux and uncertainty with extinction corrected using the Cardelli dust attenuation law (Cardelli et al. 1989) with R V = 3.1 following Valentino et al. (2017).
and R i s represent the expected flux ratio and its intrinsic scatter of each line with respect to Hβ, given by the widely used Maiolino strong line calibrations (Maiolino et al. 2008) and the default Balmer decrements assuming case B recombination conditions.Thus, the instantaneous SFR can be estimated from the Balmer line luminosity with the Kennicutt calibration (Kennicutt 1998) and the Chabrier IMF (Chabrier 2003), assuming case B recombination:


We adopt the EMCEE software to perform the MCMC Bayesian parameter sampling, with 100 walkers each sampling 1000 iterations.After a burn-in of 200 for each walker, we have sampled the parameter space 80,000 times.We show the resultant 1D and 2D parameter constraints in Figure 10.The median constraints and 1σ CIs are shown in Table 1.To understand the energy source that powers the photoionization of the nebular emission observed in RX J2129-z8He II, we rely on the massexcitation (MEx) diagram (Juneau et al. 2014;Coil et al. 2015), which is an empirical relation between M * and the [O III]/Hβ flux ratio with demarcation schemes separating the H II region and AGN.We confirm nondetection of our galaxy from a total of ∼36k second exposures of Chandra X-ray observations (ObsID 9370, PI: Allen; and ObsID 553, PI: Garmire) covering the galaxy cluster center field of RX J2129.Assuming a fiducial incident spectrum with an absorbed power law with a photon index of 1.7 and an absorption column density of 10 21 cm −2 , we reach a 3σ upper limit of 1.5 × 10 45 erg s -1 in the energy range of [0.5, 8] keV, lower than the measured X-ray luminosity of X-raybright AGNs at z > 5 (see, e.g., Li et al. 2021).It is likely that RX J2129-z8He II is a star-forming galaxy with negligible contamination from AGN ionization, although we caution that AGNs need not be X-ray bright and the MEx diagram may differ substantially for very low metallicities.The definitive exclusion of AGN contamination should come from much deeper JWST/ NIRSpec spectroscopy with sufficient wavelength resolution to cast stringent limits on the high-ionization emission lines in the rest-frame UV and optical.The color-coded circles represent the spectroscopically confirmed galaxies at z spec > 7 with β measurements (Watson et al. 2015;Zitrin et al. 2015;Hashimoto et al. 2018;Williams et al. 2023;Morishita et al. 2023;Schaerer et al. 2022;Curtis-Lake et al. 2023;Bunker et al. 2023).Among this cohort, RX J2129-z8He II (indicated as the star) has the steepest slopes.The gray diamonds show the galaxy candidates with z phot  8 (Tacchella et al. 2022;Castellano et al. 2022;Naidu et al. 2022;Finkelstein et al. 2022;Topping et al. 2022;Cullen et al. 2023).The dotted and dashed lines denote the population-averaged correlation of M UV and β at z ∼ 7 (Bouwens et al. 2010(Bouwens et al. , 2014) ) and z ∼ 8 (Bhatawdekar & Conselice 2021), respectively.The shaded regions correspond to the parameter space where the escape fraction of the ionizing radiation ( f esc LyC ) is above 5%, 10%, and 20%, respectively, derived by Chisholm et al. (2022) from the LzLCS data set.This strongly suggests that RX J2129-z8He II is hitherto the most promising spectroscopically confirmed galaxy that enables significant LyC leakage into the IGM to contribute to the cosmic reionization.
Since we only achieve an upper limit on nebular dust extinction from our emission line diagnostics (see Figure 10), we utilize A V = 0.12 ± 0.04 given by our spectrophotometric analyzes in Section 4.4 to derive dereddened line fluxes and uncertainties.After correcting for dust attenuation, we obtain a high intrinsic flux ratio of [O III]/[O II] = 11.7 ± 4.2, further supporting considerable LyC escape (Izotov et al. 2018) x - (Tang et al. 2019), in agreement with the expectation that at z  8, sub-L * galaxies are the major sources of reionization (Yan & Windhorst 2004;Finkelstein et al. 2019).
Notably, the He II emission is the only rest-frame UV line clearly detected in RX J2129-z8He II.This line has intrinsic (corrected for both lensing magnification and dust extinction) flux f He II = 120 ± 22 × 10 −20 erg s −1 cm −2 and EW He II = 21 ± 4 Å.This is hitherto the highest-redshift He II line potentially powered by star formation. 24The inset in the upper panel of Figure 3 shows a zoom-in view of our fitting result to the source spectrum at λ rest ä [1600, 1700], where He II and the oxygen auroral lines O III] λλ1661, 1666 (:= O III]) are marked in magenta vertical dashed lines.We verify that our wavelength calibration of the NIRSpec prism spectroscopy is sufficiently accurate with a 1σ uncertainty being ∼0.004 μm in the observer frame, corresponding to ∼4 Å in the rest frame, distinguishing between He II and O III] at a ∼5σ CI.On the other hand, if the line were indeed O III], the only possible scenario that could power the doublets with such high EW ∼ 20 Å is AGN ionization, which would result in much brighter C IV and C III] (Hirschmann et al. 2019).This reinforces our conclusion that the detection of the He II line is robust against negligible contamination from the oxygen auroral lines of O III].We also perform detailed checks of the relative width and profile of the He II line, in comparison to those of [O III]-Hβ.We do not see any signs of broadening of He II and find that the normalized profile and the relative width of He II match reasonably well those of Hβ and [O III] λ5007, strongly indicating that the He II line observed in RX J2129-z8He II also originates from the nebular H II regions.
The nebular He II line requires a hard ionizing background radiation, which is usually attributed to W-R stars, stripped stars, X-ray binaries, or AGNs (Nanayakkara et al. 2019;Saxena et al. 2020).The location of RX J2129-z8He II in the MEx diagram (Juneau et al. 2014;Coil et al. 2015), the lack of variability from archival HST imaging, and the nondetection in deep Chandra exposures altogether disfavor a supermassive black hole as the cause.As opposed to He II, RX J2129-z8He II shows an absence of the UV carbon lines, setting it apart from other galaxies at various redshifts with He II detection (Berg et al. 2016;Sobral et al. 2018;Senchyna et al. 2019;Nanayakkara et al. 2019), as shown in the right panel of Figure 11.Interestingly, some sources in the MUSE HUDF sample at 2 < z < 4 (Nanayakkara et al. 2019) approach our measured limits of f He II /f C IV and f He II /f C III] .However, their He II EWs are reported to be less than 10 Å, not as high as that measured in our galaxy, and their He II emission is likely caused by stellar binarity and density-bounded H II regions with subsolar ISM metallicities (Plat et al. 2019).We do not detect any UV metal lines (e.g., C III], C IV, O III], N IV] λλ1483, 1487, N V λ1240) either, inconsistent with the signatures of W-R and stripped stars that are often seen in local galaxies (Morris et al. 2008;Leitherer 2020), although we caution that the strength of these wind lines for W-R and stripped stars is metallicity-dependent and is likely to be very weak at extremely low metallicities (e.g., much lower than onetenth solar).
Alternatively, the He II line could be due to high-mass, metal-free Pop III stars that have exceedingly hard ionizing fields capable of sufficient He + ionization.While the ISM of RX J2129-z8He II is already metal-enriched to roughly onetwelfth solar ( ( ) 12 log O H + ∼ 7.63), its dereddened flux ratio of f He II /f Hβ = 1.7 ± 0.4 is several orders of magnitude larger than that powered by "normal" metal-poor O/B stellar populations (see the left panel of Figure 11).We check that the UV and optical photometry of this potential "Pop III" component of our galaxy is consistent with the expected colors of Pop III galaxies reported in Trussler et al. (2023).In the next section, we investigate in detail the feasibility of reproducing the nebular emission observed in this system using Pop III stellar evolution and photoionization models.

Pop III Star Photoionization Model, Formation Rate, and Total Mass
One crucial aspect of this work is to explain the flux ratios measured between the strong emission lines observed in the rest-frame UV and optical spectra of RX J2129-z8He II using ), and instantaneous SFR already corrected for lensing magnification (μ = 2.26 ± 0.14) of RX J2129-z8He II, obtained from our Bayesian forwardmodeling emission line diagnostic method.On top of each 1D posterior panel is shown the resulting median with 1σ uncertainties for each parameter.Yet notice that we can only draw a 1σ upper limit on A V N as presented in Table 1.The panel in the upper right shows the −χ 2 /2 values for all 80,000 parameter sampling iterations.They clearly reach a global minimum of χ 2 , implying a good convergence of our Bayesian inference.state-of-the-art photoionization models and reasonable assumptions on the input ionizing spectra from Pop III stellar populations.We carry out such analysis with the state-of-theart MAPPINGS V photoionization code (Sutherland & Dopita 1993).MAPPINGS V has comprehensive consideration of microphysics in the ISM, providing accurate predictions of fluxes of over 80,000 cooling and recombination lines.By propagating a simple toy model of Pop III star ionizing spectra calibrated by Schaerer (2002Schaerer ( , 2003) ) through the MAPPINGS V photoionization code, we verify that the observed line flux ratios can be reproduced under reasonable assumptions of the ionizing spectral hardness and spectral shapes that are believed to be typical in Pop III stellar populations.We assume the following form of broken power laws for this toy model: where β = −2.53 is determined from our full spectrum fitting analysis and α represents the power-law index in the He 0 and He + ionization energy range.We experiment with two indices of α = −1.8 and α = −2.5 in constructing this phenomenological model of Pop III star extreme-UV spectra and arrive at the conclusion that α does not affect the resulting line ratios as long as the spectral hardness is maintained (see Table 4).The C 0 coefficient is computed to match the measured rest-frame EW of EW He II 21 Å, assuming that observationally, the far-UV continuum is mostly coming from Pop III stars. 25The other coefficients (i.e., C 1,2,3 ) are determined to match the hardness ratios of the ionizing spectra reported in Schaerer (2002), given by the ratios of the ionizing photon flux (Q(H), Q(He 0 ), and Q(He + )).
We take the time-averaged Q values calculated from the Pop III stellar evolution tracks with M ini = 500 M e as specified in Tables 4 and 5 of Schaerer (2002) to account for the effects of stellar evolution, which typically decreases the spectral hardness by a factor of 2. Conventionally, Pop III stars with very high mass are preferred because of the lack of efficient cooling channels (due to the absence of the CNO elements; see reviews by Bromm & Larson 2004;Bromm 2013).Both scenarios with and without mass loss are considered to investigate the effects of mass loss due to stellar winds and core-collapse supernovae on ratios of the emergent rest-frame UV/optical emission lines (Schaerer 2003).
The resulting phenomenological models of Pop III star ionizing spectra calibrated against realistic Pop III stellar evolution tracks and atmosphere models are represented by the black and green lines in Figure 12.For comparison, we also calculate the empirical Crucially, we find that high [O III]/[O II] alone is not necessary a good indication of the hard ionizing spectra, potentially from Pop III stars.Right: UV line ratios measured from galaxies that have reported He II detections (S/N  2.5).These measurements are color-coded in EW He II reported in the respective work.The arrows denote 2σ lower limits.RX J2129-z8He II shows greatly elevated He II/C III] and He II/C IV flux ratios, as compared to local galaxies (Berg et al. 2016;Senchyna et al. 2019), galaxies at the cosmic noon epoch (Nanayakkara et al. 2019), and CR7 at z = 6.6 (Sobral et al. 2018).We also overlay three sets of evolutionary tracks at different stellar metallicities given by the BPASS v2.1 stellar population synthesis models (Xiao et al. 2018).Among these three tracks, the measured flux ratio lower limits of RX J2129-z8He II clearly prefer that with extremely low metallicity.The gray diamonds represent the predictions of AGN narrow-line regions calculated using the CLOUDY code (Ferland et al. 1998), with larger symbol size corresponding to higher metallicity (i.e., 0.03, 0.1, 0.5 and 1× solar).
ionizing spectra of single metal-poor stars using the TLUSTY non-LTE plane-parallel code modeling massive star atmospheres (Hubeny & Lanz 2011;Hubeny et al. 2021), with three characteristic temperatures of 55 kK, 47.5 kK, and 45 kK shown as the blue, red, and cyan curves, respectively, in Figure 12.Their stellar phase metallicity is taken to be [Fe/H] = −3 (onethousandth solar).We assume three conditions of the gas-phase metallicity for all input stellar spectra: Z ISM /Z e = 0.5, 0.1, and 10 −3 in our photoionization calculations.
The MAPPINGS V photoionization code takes these ionizing spectra as input and produces the emergent rest-frame UV and optical spectra via self-consistent Monte Carlo radiative transfer.We adopt the plane-parallel geometric configuration in radiative transfer calculations, since we believe that the main ionizing sources (Pop III stars) are located ∼1 kpc away (i.e., in component B) from the ISM, where our NIRSpec slit spectroscopy is taken (i.e., in component A).The Pop III star scenarios, on the other hand, offer sufficient extreme-UV ionization fields giving rise to strong He II emission and can reproduce the observed ratio of He II/ Hβ = 1.7 ± 0.4 in RX J2129-z8He II.We find that the Pop III ionizing spectra both with and without mass loss lead to comparable ionizing radiation, with the former producing slightly elevated He II/Hβ ratios, due to increased spectral hardness (Schaerer 2002(Schaerer , 2003)).The observed ratio of [O III]/ [O II] = 11.7 ± 4.2 better favors the Pop III spectral models with strong mass loss.Finally, we also see good agreement between the observed [O III]/Hβ = 5.5 ± 0.8 and our calculations setting ISM metallicity to one-tenth solar (see Table 4 and circles in the left panel of Figure 11), which lends support to the constraint on as shown in Table 1), we follow the prescriptions of Raiter et al. (2010) and Cai et al. (2015).Under the fiducial conditions of case B recombination and an ionization-bounded nebula with constant density and temperature, the He II line luminosity generated by Pop III stars can be expressed as a function of He + ionizing photon flux (Q(He + )) and Pop III SFR, i.e.,  where f 1640 is the theoretical He II line luminosity for constant star formation models normalized to SFR = 1 M e yr −1 , given in Table 7 of Schaerer (2002; also see Table 5).From our emission line analyzes, we measured the He II line flux corrected for both lensing magnification and dust extinction to be f He II = 120 ± 22 × 10 −20 erg s −1 cm −2 , corresponding to a He II luminosity of L He II = 9.7 ± 1.8 × 10 41 erg s −1 .Therefore, we can estimate SFR Pop III under various Pop III stellar evolution models as shown in Table 5 from the observed L He II and the escape fraction of f esc = 16% inferred from the UV spectral slope of β = −2.53.We caution that there is a large uncertainty attached to the conversion from β to f esc , and the escape fractions of the ionizing flux in different energy regimes (e.g., f esc (>13.6 eV) versus f esc (>54.4 eV)) can be quite different.In fact, the galaxies for which the β-f esc relation has been derived empirically (Equation (3)) are not systems that host Pop III stars (Chisholm et al. 2022).Applying this relation to potential Pop III star host galaxies might be expected to underestimate f esc , since the Pop III stellar models with extremely top-heavy IMFs predict relatively red UV slopes due to increased nebular continuum (Cameron et al. 2023).Nevertheless, with discretion, we estimate the ratio between SFR Pop III and SFR O/B calculated using Equation (6) to be in the range of 3%-72%.Assuming the Eddington limit and L He II being ∼1% of the total bolometric luminosity, we can also derive the total mass of the Pop III stars to be 7. ) and the total mass of Pop III stars within this system reside in reasonable ranges.Nevertheless, we caution readers about the caveat that some assumptions are still highly uncertain and demand more work in this fast-rising field.First and foremost, the theoretical framework of Schaerer (2002Schaerer ( , 2003) ) and Raiter et al. (2010) is based upon the top-heavy Salpeter (Salpeter 1955) IMF. 26 A factor of 8 difference in the inferred SFR is caused by different low-mass cutoffs in the IMF mass range.Second, we see that models incorporating strong mass loss  decrease the SFR by a factor of 3.However, the detailed physical mechanisms (e.g., rotation, pair-instability supernovae, multiplicity, etc.) driving such high mass loss are still under investigation (Ekström et al. 2008).Lastly, we caution that the He II line emission coefficient c 1640 = 5.67 × 10 −12 erg is computed for case B with nebular electron density n e = 100 cm −3 and T e = 30 kK, under the assumption that the radiation field is optically thick.Realistically, there can exist steep density/temperature gradients and significant deviations from the fiducial case B conditions around extremely hot Pop III stars.In particular, the relatively high LyC escape fraction ( f 16% esc LyC = ) estimated in RX J2129-z8He II, indicated by its large O32 flux ratio ( f [O III] /f [O II] = 11.7) and blue UV continuum slope (β = −2.53), is in favor of the ionizationbounded H II regions, rather than the density-bounded ones assumed in case B. This density-bounded condition can suppress low-ionization emission lines such as Hβ and potentially artificially boost ratios between high-ionization and low-ionization lines such as He II/Hβ (Plat et al. 2019).Nonetheless, our work presents the critical first step forward to apply physically motivated Pop III stellar evolution models to JWST high-quality NIRSpec spectroscopy of galaxy candidates revealing Pop III star signatures.

Conclusions and Discussion
In this Letter, through the comprehensive analysis of the novel JWST/NIRSpec prism spectroscopy and NIRCam imaging, we present RX J2129-z8He II, an intriguing galaxy in the EoR.Fitting simultaneously the entire spectroscopic (continuum + lines) and photometric data sets, we confirm RX J2129-z8He II as an extreme emission line galaxy with strong continuum emission at z spec = 8.1623 ± 0. ) among all spectroscopically confirmed galaxies at z  7 reported to date, implying that it has a large LyC escape fraction, f 16% esc LyC ~, the largest among its cohort.Therefore, RX J2129-z8He II is a promising representative of the predominant galaxy populations that are capable of producing and leaking their extreme-UV photons to the IGM, causing the bulk of the neutral hydrogen in the IGM to be ionized.Furthermore, its prominent He II line and the large flux ratios between the He II and UV metal and Balmer lines make it the best candidate to date where this very galaxy could have a significant Pop III stellar constituent within its stellar populations.
Using the MAPPINGS V photoionization code and a series of phenomenological models of Pop III star spectra calibrated against theoretical Pop III stellar evolution models, we can successfully reproduce the measured ratios of He II/Hβ and [O III]/[O II], assuming strong mass loss and one-tenth ISM metallicity.Crucially, we also find that high [O III]/[O II] alone is not a good indication of the presence of Pop III stars, providing important guidelines for future JWST surveys with the aim of identifying Pop III signatures in metal-poor systems in the EoR.All this, together with the large EW He II ∼ 21 Å, suggests that RX J2129-z8He II likely has a mixture of metalpoor O/B and Pop III stars, the latter of which are the energy source of the He II line emission.With the caveat that our f esc estimate based on observed UV spectral slope can be an underestimation due to nebular continuum, we estimate the ratio between the formation rates of Pop III and O/B stars in the range of [3, 72]% and derive the total mass of Pop III stars within RX J2129-z8He II to be 7.8 ± 1.4 × 10 5 M e , consistent with the predictions based on numerical simulations.
The putative Pop III stars could either be spatially mixed with the enriched populations or segregated.Indeed, recent cosmological hydrodynamic simulations suggest that Pop III stars can continue to form until the end of the EoR, triggered by sustained accretion of pristine gas from the cosmic web onto galaxy disks, resulting in "hybrid" Pop III galaxies (Venditti et al. 2023).We note that the prism spectroscopy is taken on component A of RX J2129-z8He II only; its component B resides ∼0 2 (i.e.,∼1 kpc proper) to the southeast (see Section 4.3).We carry out photometry for A and B separately and find that they have comparable photo-z estimates, z is a factor of ∼2 lower than A's, indicating that B's ISM is more chemically pristine.B's SED shape might allow for a He II line with higher EW, but we caution that B's photo-z carries substantial uncertainties.Simulations predict that Pop III stars can continue to form if there still are clumps of primordial gas not yet polluted by metal-enriched outflows from nearby star-forming galaxies or metal-loaded stellar winds from massive stars close by (Ciardi & Ferrara 2005;Venditti et al. 2023).A deep NIRSpec/integral field unit observation with high spectral resolution covering both components can provide a more definitive answer to the origin of this prominent He II emission observed in RX J2129-z8He II.

Figure 2 .
Figure 2. JWST+HST imaging and JWST/NIRSpec spectroscopy of RX J2129-z8He II.Top: RX J2129-z8He II's postage stamp images in multiple filters, cut from the HST CLASH (PI: M. Postman) and JWST/NIRCam (DD-2767; PI: P. Kelly) imaging mosaics.These stamps are oriented such that north is up and east is to the left, which is different from the orientation in Figure 1.The strong flux deficit in the ACS/F814W and NIRCam/F115W bands as compared to the other NIRCam filters indicates z  8.The rightmost panel in the top row zooms in on the [O III]-Hβ line complex taken from the full NIRSpec prism 1D spectrum shown in the bottom row.Middle and bottom: NIRSpec prism spectroscopy of RX J2129-z8He II.In the middle trace, we show the unsmoothed 2D prism spectra combined from the three individual slit-dithering sequences using the up-to-date reduction software MSAEXP (see Section 2.2).We detect prominent emission features of [O III] and Hβ at λ obs ∼ 4.5 μm and a strong continuum break blueward of Lyα at λ obs ∼ 1.1 μm.In the bottom row, we show the 1D spectrum optimally extracted from the 2D spectral trace in gray histograms with the 1σ uncertainty in cyan.We mark the locations of multiple rest-frame UV and optical emission lines using the magenta dotted lines, and the Lyman break is the blue dashed line.The strong detections of [O III]+Hβ and the Lyman break pinpoint a secure spectroscopic redshift of RX J2129-z8He II to be z spec = 8.1623 ± 0.0007.Bottom inset: a zoom-in view on the 2D spectrum around the rest-frame UV wavelength range, where the location of the He II line is marked by the red circle.The 2D spectrum is shown in units of S/N (i.e., f obs obs s n

Figure 3 .
Figure3.Spectrophotometric and emission line analyzes of the full JWST-HST observations of RX J2129-z8He II.Top: the optimally extracted 1D NIRSpec prism spectrum is shown in blue histograms with a 1σ error spectrum in cyan bands.The broadband photometric results from NIRCam and WFC3+ACS are represented by the red circles and green diamonds, respectively.We perform full spectrum analyzes and fit Gaussian profiles to nebular emission features to estimate line fluxes and upper limits.The resulting best-fit continuum model and emission lines are represented by the red and orange curves.The top inset panel zooms in on the He II line, with significant detection.Bottom: the ratio between the continuum-subtracted f λ flux and the measured 1σ error spectrum.For the wavelength ranges without significant emission line flux contribution, this ratio follows a ( ) 0, 1 distribution, implying that the noise properties are well defined.As clearly shown here, the most pronounced emission features in the rest-frame UV and optical wavelengths are the He II and [O III]-Hβ lines, respectively.The gray shaded region marks the wavelength range affected by strong attenuation by the highly neutral IGM at z ∼ 8.
It has a Jeffery's prior on the exponents, a, b ä (0.01, 1000), and a flat prior on the peak time of star formation, τ ä (0, t H ), where t H is the Hubble time at the observed redshift.BAGPIPES relies on the BC03(Bruzual & Charlot 2003) stellar population synthesis model and the nebular emission model created by the CLOUDY photoionization code(Ferland et al. 2017).For other key assumptions, we choose the Chabrier IMF(Chabrier 2003), the Calzetti dust attenuation law(Calzetti et al. 2000) with( ) A 0, 3 V S Î, a stellar metallicity range of Z/Z e ä (0, 2), and a conservative Gaussian redshift prior of z ∼ N(8.16, 0.01).

Figure 4 .
Figure 4. We compare the photometry for both components of RX J2129-z8He II, normalized at their individual F200W flux.In each filter, the data points are slightly offset in wavelength for clearer visualization.The photometry suggests that component B has a comparably blue UV continuum yet possibly with substantially fainter oxygen lines in the rest-frame optical (probed by F444W), indicating that component B is more metal-poor than A. Also, component B's SED shape might accommodate He II with similar EW, albeit at lower confidence.

Figure 5 .
Figure 5. EAZY fit to the broadband photometry of RX J2129-z8He II, performed on components A and B separately.The red points with error bars correspond to observed magnitudes with 1σ measurement uncertainties, whereas the triangles represent 2σ upper limits, all given in Table 2.The blue curve in each panel represents the best-fit SED model, with the open circles showing the model photometry in respective filters.The inset panel shows the probability distribution function of the photometric redshift, i.e., P(z).The photometric redshift constraints of the two components are both in good agreement with z spec = 8.1623 determined from the NIRSpec prism spectroscopy of component A. This strongly suggests that component B is also part of the galaxy RX J2129-z8He II.

Figure 6 .
Figure6.The rest-frame UV and optical morphology of the entire RX J2129-z8He II galaxy, comprising two components: A and B. We use GALFIT to perform morphological modeling of both components in JWST/NIRcam imaging taken at filters F150W (upper panels) and F444W (lower panels).At z spec = 8.16, F150W and F444W correspond to rest-frame UV and optical wavelength ranges, respectively.In each row, we show the original observed image, GALFIT models of both components, and the residual removing their GALFIT models.On the observed images, we overlay the NIRSpec MSA slit, represented by the blue box.The size of all stamps is 3″ × 3″, with north up and east to the left.As indicated by the residual panels, we achieve reasonable morphological models of the entire galaxy.

Figure 9 .
Figure9.UV spectral slope (β) as a function of absolute UV magnitude (M UV ).The color-coded circles represent the spectroscopically confirmed galaxies at z spec > 7 with β measurements(Watson et al. 2015;Zitrin et al. 2015;Hashimoto et al. 2018;Williams et al. 2023;Morishita et al. 2023;Schaerer et al. 2022;Curtis-Lake et al. 2023;Bunker et al. 2023).Among this cohort, RX J2129-z8He II (indicated as the star) has the steepest slopes.The gray diamonds show the galaxy candidates with z phot  8(Tacchella et al. 2022;Castellano et al. 2022;Naidu et al. 2022;Finkelstein et al. 2022;Topping et al. 2022;Cullen et al. 2023).The dotted and dashed lines denote the population-averaged correlation of M UV and β at z ∼ 7(Bouwens et al. 2010(Bouwens et al. , 2014) ) and z ∼ 8 (Bhatawdekar & Conselice 2021), respectively.The shaded regions correspond to the parameter space where the escape fraction of the ionizing radiation ( f esc LyC ) is above 5%, 10%, and 20%, respectively, derived byChisholm et al. (2022) from the LzLCS data set.This strongly suggests that RX J2129-z8He II is hitherto the most promising spectroscopically confirmed galaxy that enables significant LyC leakage into the IGM to contribute to the cosmic reionization.

Figure 10 .
Figure 10.Corner plot showing the 1D/2D posterior constraints of the fitting parameters, i.e., gas-phase metallicity ( ( ) 12 log O H + ), nebular dust content (A V N), and instantaneous SFR already corrected for lensing magnification (μ = 2.26 ± 0.14) of RX J2129-z8He II, obtained from our Bayesian forwardmodeling emission line diagnostic method.On top of each 1D posterior panel is shown the resulting median with 1σ uncertainties for each parameter.Yet notice that we can only draw a 1σ upper limit on A V N as presented in Table1.The panel in the upper right shows the −χ 2 /2 values for all 80,000 parameter sampling iterations.They clearly reach a global minimum of χ 2 , implying a good convergence of our Bayesian inference.

Figure 11 .
Figure 11.Emission line flux ratios from observations and models.Left: flux ratios of [O III]/[O II] vs.He II/Hβ predicted by the MAPPINGS V photoionization models(Sutherland & Dopita 1993) and measured in RX J2129-z8He II (as in Table1) represented by the magenta star with 1σ error bars.The cyan, red, and blue symbols are calculated using the empirical spectra of metal-poor stars with atmospheric temperatures of 45 kK, 47.5 kK, and 55 kK, respectively.The diamonds, circles, and squares correspond to gas-phase metallicity (Z ISM ) being half, one-tenth, and one-thousandth solar, respectively.It is obvious that "normal" O/B-type stellar populations cannot produce sufficient He II emission, in spite of their high [O III]/[O II] at low Z ISM .We construct physically calibrated phenomenological models of Pop III star spectra, with different assumptions of mass loss and the power-law index (α) in the helium-ionizing energy range (see text).Encouragingly, all these models successfully reproduce the observed value of He II/Hβ.The mass-loss models are better favored considering the observed value of[O III]/[O II].Crucially, we find that high [O III]/[O II] alone is not necessary a good indication of the hard ionizing spectra, potentially from Pop III stars.Right: UV line ratios measured from galaxies that have reported He II detections (S/N  2.5).These measurements are color-coded in EW He II reported in the respective work.The arrows denote 2σ lower limits.RX J2129-z8He II shows greatly elevated He II/C III] and He II/C IV flux ratios, as compared to local galaxies(Berg et al. 2016;Senchyna et al. 2019), galaxies at the cosmic noon epoch(Nanayakkara et al. 2019), and CR7 at z = 6.6(Sobral et al. 2018).We also overlay three sets of evolutionary tracks at different stellar metallicities given by the BPASS v2.1 stellar population synthesis models(Xiao et al. 2018).Among these three tracks, the measured flux ratio lower limits of RX J2129-z8He II clearly prefer that with extremely low metallicity.The gray diamonds represent the predictions of AGN narrow-line regions calculated using the CLOUDY code(Ferland et al. 1998), with larger symbol size corresponding to higher metallicity (i.e., 0.03, 0.1, 0.5 and 1× solar).
The resulting line flux ratios of [O III]/Hβ, [O II]/Hβ, He II/Hβ, and [O III]/ [O II] from these empirical and phenomenological stellar spectral models are given in Table 4.The left panel of Figure 11 shows the model predictions of [O III]/[O II] versus He II/Hβ, compared against our observations.Our detailed photoionization modeling clearly suggests that normal O/B stellar populations, albeit with high atmospheric temperatures and high [O III]/[O II], do not suffice to power prominent highionization He II emission, as they would only result in insignificant He II/Hβ ratios lower than our observation by 4-5 orders of magnitude.This strongly indicates that although [O III]/[O II] is widely used as a key diagnostic of the ionization state of the gas in extragalactic H II regions (Kewley et al. 2013), it is not necessarily a good proxy by itself for identifying potential Pop III stars in metal-poor systems due to serious contamination from normal O/B stellar populations.
line diagnostics shown in Figure10.To estimate the SFR of Pop III stars (SFR Pop III ) and compare with the instantaneous SFR converted from the dereddened Hβ line luminosity

Figure 12 .
Figure 12.Realistic stellar spectra input to our photoionization modeling framework.The cyan, red, and blue curves correspond to the empirical spectra of metal-poor ([Fe/H] = −3) single stars with high atmospheric temperatures of 45 kK, 47.5 kK, and 55 kK, respectively, representing the typical populations of low-metallicity O/B stars.These empirical star model atmospheres are constructed using the non-LTE plane-parallel code TLUSTY (Hubeny & Lanz 2011; Hubeny et al. 2021).The black and green lines show our broken power-law models of Pop III stellar spectra, defined in Equation (7), with different assumptions of mass loss and the power-law index α (see text for more details).The ionizing sections of these Pop III spectra are calibrated by the time-averaged hardness ratios-described by Q(He 0 )/Q(H) and Q(He + )/Q(H) tabulated in Schaerer (2002)-with and without mass loss.The nonionizing sections are determined based upon our actual measurement of the UV spectral slope (β = −2.53)and the emission line EWs (i.e., EW He II = 21 Å).The three vertical dotted gray lines mark the energy levels of 13.6, 24.6, and 54.4 eV, corresponding to the ionization of H, He 0 , and He + ions.
0007.It has a sub-L * M UV , high-EW rest-frame optical lines, and the steepest UV spectral slope (

Table 1
Physical Properties of RX J2129-z8He II

Table 2
Detailed Photometry of RX J2129-z8He II Separating Components A and B from HST/ACS Optical, WFC3 Infrared, and JWST/NIRCam Infrared Imaging

Table 3
Rest-frame UV and Optical Emission Line Properties of RX J2129-z8He II The observed line fluxes and EWs are measured from our emission line analyzes of the NIRSpec/MSA prism spectrum with the PPXF software (see text for details).We list measurements for 2σ detections, or the 2σ upper limits.The error bars correspond to 1σ standard deviations.
aThe observed line fluxes and upper limits before the corrections of dust extinction and lensing magnification.
, consistent with our estimate based on β.From our measurements of EW [O III] , we estimate the ionizing photon production

Table 4
Emission Line Flux Ratios Calculated by Our MAPPINGS V Photoionization Modeling Framework, Assuming Empirical and Phenomenological Stellar Spectra of Normal O/B and Pop III StarsNote.The last row shows our actual measurements with 1σ uncertainties of these key line ratio diagnostics in RX J2129-z8He II using JWST/NIRSpec prism spectroscopy.
(Venditti et al. 2023)onsistent with the predictions based on numerical simulations(Venditti et al. 2023).It is encouraging to see that the SFR ratio (

Table 5
(Kennicutt 1998)X J2129-z8He II at z = 8.1623 Inferred from Equation (8) Note that f 1640 stands for the theoretical He II luminosity normalized to SFR = 1 M e yr −1 for the four stellar evolution models(Schaerer 2002).The actual observed He II luminosity of RX J2129-z8He II is L He II = 9.7 ± 1.8 × 10 41 erg s −1 .The last column represents the ratio between the Pop III SFR (SFR Pop III ) and the instantaneous SFR measured from the intrinsic Hβ luminosity using Equation (6) (SFR O/B ).Note that the Kennicutt calibration(Kennicutt 1998)used here to calculate SFR O/B is generally applicable to galaxies with strong recombination lines powered by the ionizing flux predominantly from O/B-type stars that typically have M * > 10 M e and lifetimes of t age < 20 Myr.Here SFR Pop III has been converted to the value compatible with the Chabrier IMF (Chabrier 2003) to facilitate a direct comparison with SFR O/B .
spectroscopic redshift measured for A (see Section 4.2).Our photometry for component B shows a blue UV continuum comparable to A's but with a much smaller F444W/F200W flux density ratio (see Section 4.1).Assuming that the different color is driven by EW [O III] , we estimate that B's EW [O III]