# Searching for Extremely Blue UV Continuum Slopes at z = 7-11 in JWST/NIRCam Imaging: Implications for Stellar Metallicity and Ionizing Photon Escape in Early Galaxies

Michael W. Topping, Daniel P. Stark, Ryan Endsley, Adele Plat, Lily Whitler, Zuyi Chen, Stéphane Charlot

<figure><img src="image_1.png" /></figure>

## To cite this version:

Michael W. Topping, Daniel P. Stark, Ryan Endsley, Adele Plat, Lily Whitler, et al., Searching for Extremely Blue UV Continuum Slopes at z = 7-11 in JWST/NIRCam Imaging: Implications for Stellar Metallicity and Ionizing Photon Escape in Early Galaxies. The Astrophysical Journal, 2022, 941, 10.3847/1538-4357/aca522. [hal-03974051](https://hal.science/hal-03974051)

HAL Id: hal-03974051

https://hal.science/hal-03974051

Submitted on 7 Feb 2023

HAL is a multi-disciplinary open access archive for the deposit and dissemination of scientific research documents, whether they are published or not. The documents may come from teaching and research institutions in France or abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est destinée au dépôt et à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

<figure><img src="image_2.png" /></figure>
<figure><img src="image_3.png" /></figure>

# Searching for Extremely Blue UV Continuum Slopes at $z=7–11$ in JWST/NIRCam Imaging: Implications for Stellar Metallicity and Ionizing Photon Escape in Early Galaxies

Michael W. Topping¹, Daniel P. Stark¹, Ryan Endsley¹, Adele Plat¹, Lily Whitler¹, Zuyi Chen¹ ID, and Stéphane Charlot² ID

¹Steward Observatory, University of Arizona, 933 N Cherry Avenue, Tucson, AZ 85721, USA; michaeltopping@arizona.edu

²Sorbonne Université, UPMC-CNRS, UMR7095, Institut d'Astrophysique de Paris, F-75014, Paris, France

Received 2022 August 3; revised 2022 October 27; accepted 2022 November 21; published 2022 December 20

## Abstract

The ultraviolet (UV) continuum slope ($\beta$, where $f_{\lambda} \propto \lambda^{\beta}$) of galaxies is sensitive to a variety of properties, from the metallicity and age of the stellar population to dust attenuation throughout the galaxy. Considerable attention has focused on identifying reionization-era galaxies with very blue UV slopes ($\beta < -3$). Not only do such systems provide a signpost of low-metallicity stars, but they also identify galaxies likely to leak ionizing photons from their H II regions as such blue UV slopes require the reddening effect of nebular continuum to be diminished. In this paper we present a search for reionization-era galaxies with very blue UV colors in recent JWST/NIRCam imaging of the Extended Groth Strip field. We characterize UV slopes for a large sample of $z \simeq 7–11$ galaxies, finding a median of $\beta = -2.0$. Two lower luminosity ($M_{UV} \simeq -19.5$) and lower stellar mass ($6–10 \times 10^7 M_{\odot}$) systems exhibit extremely blue UV slopes ($\beta = -2.9$ to $-3.1$) and rest-optical photometry indicating weak nebular line emission. Each system is very compact ($r_e \lesssim 260$ pc) with very high star formation-rate surface densities. We model the spectral energy distributions (SEDs) with a suite of BEAGLE models with varying levels of ionizing photon escape. The SEDs cannot be reproduced with our fiducial ($f_{esc,H II} = 0$) or alpha-enhanced ($Z_{*} < Z_{ISM}$) models. The combined blue UV slopes and weak nebular emission are best-fit by models with significant ionizing photon escape from H II regions ($f_{esc,H II} = 0.5–0.8$) and extremely low-metallicity massive stars ($Z_{*} = 0.01–0.06 Z_{\odot}$). The discovery of these galaxies highlights the potential for JWST to identify large numbers of candidate Lyman continuum leaking galaxies in the reionization era and suggests low-metallicity stellar populations may be common in dwarf galaxies at $z > 7$.

Unified Astronomy Thesaurus concepts: High-redshift galaxies (734); Galaxies (573); Reionization (1383); Metallicity (1031); Chemical abundances (224)

## 1. Introduction

Deep multiwavelength imaging surveys have been at the forefront of studies of galaxy formation at high redshift. The sensitivity and large sample sizes comprising such surveys have utilized models of the integrated spectral energy distributions (SEDs) produced by galaxies to infer many key properties of the galaxy population, including stellar masses, stellar population ages, star formation histories, and strengths of strong nebular emission lines (e.g., Ono et al. 2012; Labbé et al. 2013; Bouwens et al. 2014; Smit et al. 2014; Bouwens et al. 2015; Finkelstein et al. 2015; Roberts-Borsani et al. 2016; Stark 2016; De Barros et al. 2019; Endsley et al. 2021; Stefanon et al. 2022). New data sets from the James Webb Space Telescope (JWST) will soon dramatically advance the available sample sizes and measurement precision of probing galaxies in the earliest epochs.

The power-law slope of the ultraviolet (UV) continuum ($\beta$, where $f_{\lambda} \propto \lambda^{\beta}$) has been a prominent metric for establishing galaxy properties in the high-redshift universe. Great strides were made in the construction of statistical samples at high redshift brought on by near-infrared imaging using the Hubble Space Telescope (HST)/Wide Field Camera 3 (WFC3) camera. The derived UV continuum slopes are found to be blue, with values of

roughly $\beta = -2$ being common at $z \simeq 7$ (e.g., McLure et al. 2011; Finkelstein et al. 2012a; Dunlop et al. 2012; Rogers et al. 2013; Bouwens et al. 2014; Bhatawdekar & Conselice 2021). The UV slopes have been shown to become bluer both at lower UV luminosities and at higher redshifts (e.g., Bouwens et al. 2014), suggesting the $z \simeq 7$ population likely faces less attenuation from dust than is common in galaxies at later times.

Early efforts with HST were driven by the search for the bluest galaxies, with the motivation that these may provide a signpost of stellar populations dominated by very massive stars and extremely low metallicities. A variety of programs have presented measurements claiming a small population of faint galaxies at $z \gtrsim 7$ exhibiting very blue UV slopes ($\beta \simeq -3$; e.g., Bouwens et al. 2010; Labbe et al. 2010; Ono et al. 2010; Jiang et al. 2020). While such measurements face considerable uncertainty owing to the underlying photometric errors, the discovery of these sources have nonetheless pushed the limits of stellar population models (see Bouwens et al. 2010; Wilkins et al. 2011). Many models do suggest that UV continuum slopes approaching $\beta \simeq -3$ can in principle be produced in young and low-metallicity stellar populations, but in practice they are not expected to be observed. At the young stellar population ages where such blue slopes are found, nebular emission from two-photon and free-bound continuum processes (hereafter nebular continuum) significantly redden the UV continuum slopes (e.g., Nussbaumer & Schmutz 1984; Bottorff et al. 2006).

The only way to achieve UV slopes approaching values of $\beta \simeq -3$ is if the nebular continuum contribution is significantly

<figure><img src="image_4.png" /></figure>

Original content from this work may be used under the terms of the Creative Commons Attribution 4.0 licence. Any further

distribution of this work must maintain attribution to the author(s) and the title of the work, journal citation and DOI.
reduced, as would occur if ionizing photons escape from H II regions without being reprocessed into nebular emission (Bouwens et al. 2010; Ono et al. 2010; Raiter et al. 2010; Robertson et al. 2010). While escape of ionizing photons from H II regions is only the first step toward leakage (the radiation also must escape the interstellar matter, interstellar medium (ISM), and circumgalactic medium), it is a necessary precondition, and sources with a large fraction of their ionizing radiation escaping from H II regions (hereafter $f_{esc,H II}$) should be excellent Lyman continuum (LyC) leaker candidates. As such, the detection of galaxies with $\beta \simeq -3$ may identify extremely efficient ionizing agents, with hard spectra from low-metallicity stars and large escape fractions.

The detections described above have motivated attempts to develop the modeling infrastructure necessary to link very blue UV slopes to information on ionizing photon escape fractions. Early papers explored implications for some of the first detected sources with $\beta < -2.8$ (Bouwens et al. 2010; Ono et al. 2010; Wilkins et al. 2011), but constraints on $f_{esc,H II}$ were quite limited given the available data. In more recent studies, efforts have focused on modeling the reduction in strength of nebular emission lines that arises when there is significant ionizing photon escape (e.g., Robertson et al. 2010; Zackrisson et al. 2013, 2017; Plat et al. 2019). Zackrisson et al. (2013) motivate that both UV slopes and H$\beta$ EWs can be used simultaneously to identify galaxies with significant LyC leakage. Meanwhile, observations have begun to develop large samples of LyC leakers at lower redshifts (e.g., Fletcher et al. 2019; Flury et al. 2022a; Pahl et al. 2022), allowing these models to be more directly tested (e.g., Yamanaka et al. 2020; Flury et al. 2022b; Chisholm et al. 2022). And with the commencement of JWST cycle 1, it is rapidly becoming possible to better constrain the extremely blue tail of the distribution of UV slopes at $z > 7$. And unlike in the past, it is now possible to provide robust constraints on the rest-optical emission line strengths of $z > 7$ objects with blue colors, allowing much-improved modeling of the influence of ionizing photon escape on the SED.

In this paper, we present results stemming from a search for galaxies with very blue UV continuum slopes among the $z > 7$ population in the recent NIRCam imaging of the Extended Groth Strip (EGS) field by JWST. We identify two galaxies with both extremely blue UV slopes and rest-optical photometry that is suggestive of relatively weak nebular lines, both indicative of LyC leakage and low-metallicity stellar populations. We fit the SEDs with a suite of models allowing variation of $f_{esc,H II}$ (Plat et al. 2019) and discuss the properties of the sources in the context of known LyC leakers at lower redshifts (e.g., Flury et al. 2022a; Pahl et al. 2022). Section 2 provides a summary of the models that we use in this paper. Section 3 briefly introduces our data products and reduction techniques (Section 3.1) before discussing our search for very blue galaxies in $z \simeq 7-11$ galaxy samples in NIRCam imaging. In Section 4, we investigate these extremely blue objects in the context of models that include the effects of ionizing radiation escape and compare the properties of our candidates to known LyC leakers at lower redshifts. Finally, Section 5 presents a brief summary. Throughout this paper we use AB magnitudes (Oke & Gunn 1983) and assume a cosmology with $\Omega_m = 0.3$, $\Omega_\Lambda = 0.7$, and $H_0 = 70 \text{ km s}^{-1} \text{ Mpc}^{-1}$. We additionally adopt solar abundances from Asplund et al. (2009; i.e., $Z_\odot = 0.014$).

## 2. Stellar Population Models

In this paper, we investigate the presence of extremely blue UV slopes in star-forming galaxies at $z > 7$. Our goals in doing so are both to characterize the stellar populations in the reionization era and to search for galaxies that may be LyC leaking candidates. To interpret these systems, we will compare the observed SEDs to a suite of models that allow for a range of escape fractions using density-bounded H II regions (Plat et al. 2019). The models make use of the Bayesian Analysis of Galaxy SEDs tool (BEAGLE; Chevallard & Charlot 2016; described in more detail below) and utilize the latest version of the Bruzual & Charlot (2003) stellar population synthesis models. BEAGLE self-consistently includes the contribution of nebular emission calculated using the photoionization models described in Gutkin et al. (2016). The density-bounded models that form the basis of our varying escape fraction models have been developed and described in detail in Plat et al. (2019). In this section, we discuss the UV slopes expected in our fiducial ionization-bounded models ($f_{esc,H II} = 0$) and then describe how both the UV slopes and rest-optical emission lines vary as escape fraction is increased.

The sensitivity of $z > 7$ SEDs to the presence of ionizing photon escape is a direct result of the prominence of nebular continuum emission in the SEDs. The shift of the galaxy population toward larger specific star formation rates (sSFRs) at fixed mass between $z \simeq 2$ and $z \simeq 7$ (Topping et al. 2022) implies a corresponding shift toward younger light-weighted ages (Whitler et al. 2022a). With the considerable effort put into characterization of SEDs with HST and Spitzer in the last two decades, we now have a reasonably good idea of the distribution of light-weighted ages up to $z \simeq 7$, with significant evidence that a substantial fraction of the population is dominated by the light from very young stellar populations (Smit et al. 2014; Endsley et al. 2021; Stefanon et al. 2022). At these young ages, nebular continuum emission makes a significant contribution to the total continuum flux from the far-UV (FUV) into the optical.

In Figure 1 we show the fraction of the observed continuum that is of nebular origin as a function of wavelength for our fiducial ionization-bounded BEAGLE models. We show the nebular continuum fraction for three different ages, representative of the range of ages seen at $z \simeq 7$, and we highlight the dependence on metallicity. These have been run with our standard parameterization that assumes a constant star formation history (SFH), a Chabrier (2003) initial mass function (IMF) with an upper-mass cutoff of 300 $M_\odot$, $\log(U) = -2.5$, a solar abundance pattern, and no dust. At young ages, the nebular continuum contributes a sizable fraction of the total flux at all wavelengths. At low metallicities, galaxies with an age of 1 Myr have SEDs that are composed of ~40% nebular continuum at 1500 Å, rising to nearly ~80% at the blue side of the Balmer jump. While the normalization of this continuum fraction decreases with age, at 10 Myr the contribution is still over one-third of the total flux and comprises 10%-20% at 100 Myr. Clearly, this represents a significant effect in galaxies that are common at $z > 7$. This figure additionally illustrates the significant effect that the metallicity has on the nebular contribution. The lower metallicity models (with their associated higher electron temperatures) yield SEDs with a more significant contribution of the nebular continuum at all wavelengths.
<figure><img src="image_5.png" /><figcaption>**Figure 1.** Fraction of the total flux that is the nebular continuum as a function of rest wavelength. We show the results for models of 1, 10, and 100 Myr in panels (a), (b), and (c), respectively. In each panel, we display continuum fractions calculated for five different assumed metallicities. The nebular continuum is strongest at the youngest ages and lowest metallicities and decreases in strength toward older ages and higher metallicities.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**Figure 2.** (a) Rest-UV continuum slope $\beta$ as a function of age and metallicity measured from models constructed using BEAGLE. These models assume a constant SFH, solar abundance patterns, and $\log(U) = -2.5$, and they contain no dust attenuation. Shown here are results from models of stellar light only (dashed lines) and those with self-consistently computed nebular continuum emission (solid lines). As described in Section 2, the addition of nebular emission reddens the UV continuum by $\delta\beta \sim 0.6$ at the youngest ages and lowest metallicities and has a decreasing effect toward older galaxies and higher metallicities. (b) Same as panel (a) except the UV slopes are calculated using models produced using Prospector. We note that while we display UV slopes of models with ages up to 1 Gyr for reference, this extends beyond the age of the universe at the redshifts considered in this study.</figcaption></figure>

It is clear from Figure 1 that the contribution of the nebular continuum is more prominent in the near-UV (NUV) and optical. We thus may expect to see very blue UV slopes in the FUV (where stellar continuum contributes over 50% of the light) with slightly redder slopes in the NUV where the contribution of nebular continuum is greater. As a result, it will be important to consider the UV slopes in both the FUV and NUV if interested in whether the nebular continuum is making a significant contribution. Objects with significant leakage (and hence diminished nebular continuum) would be expected to show blue UV slopes in both the FUV and NUV, whereas

those with very young SEDs and minimal leakage may show a break in their continuum slopes between the FUV and NUV.

To guide our investigations in the following section, we also quantify how the UV slope varies with stellar population parameters for the fiducial BEAGLE models we are using. This has been done previously for other stellar population models (e.g., Bouwens et al. 2010; Dunlop et al. 2013; Zackrisson et al. 2013), but we seek to test whether there are significant differences in the intrinsic UV slopes for the updated models we are considering in this paper. Figure 2 shows how the power-law slope $\beta$ varies as a function of age and metallicity,
calculated as the power-law slope over the range 1280–2580 Å
as in Calzetti et al. (1994). This figure includes UV slopes
calculated using BEAGLE models of the stellar emission only,
as well as models that include the contribution from the nebular
continuum. The physics driving variations in the UV slope with
age and metallicity has been described elsewhere (e.g., Reddy
et al. 2012; Topping & Shull 2015). While we impose a
Chabrier (2003) IMF in this study, we note that different
assumptions, such as a top-heavy IMF, can impact the relative
abundances of the most massive stars, affecting the stellar
continuum shape. Results from Wilkins et al. (2012) suggest
that imposing a top-heavy IMF with a high-mass slope of
$\alpha = -1.5$ could result in a steepening of the UV slope of
$\Delta\beta \sim 0.2$ compared to when a high-mass slope of $\alpha = -2.35$ is
assumed. Jeřábková et al. (2017) further explored how UV
slopes vary under different IMF assumptions as a function of
age. They demonstrated that at all but the youngest ages (i.e.,
$\lesssim 5$ Myr), the choice of IMF played a minor role in the UV
slope but observed a similar steepening to that found by
Wilkins et al. (2012) amounting to a $\Delta\beta \sim 0.2$ at the youngest
ages. These results suggest that a varying IMF plays an
important role in setting the UV slope; however, the leakage of
the ionizing photon is still needed to achieve the bluest slopes
in our sample.

The BEAGLE models asymptote to $\beta = -3.2$ at the youngest ages ($<3$ Myr) and lowest metallicities ($\log(Z/Z_\odot) = -2$) considered. We note similarly blue UV slopes are seen in the Prospector models (Johnson et al. 2021) with the default FSPS stellar templates (Conroy et al. 2009; Conroy & Gunn 2010), which implement the MIST stellar isochrones (Choi et al. 2016; see right panel of Figure 2). Here we have adopted the same model setup as in our fiducial BEAGLE models. These results thus indicate that BEAGLE and Prospector have no problem reproducing UV slopes as blue as $\beta \approx -3$, provided the metallicities are low, ages are young, and stellar continuum dominates the light. In addition to age and metallicity, stellar multiplicity can be an additional factor regulating the shape of the UV continuum. We explore the possibility that the introduction of binary stars could yield significantly bluer UV slopes using stellar templates from BPASS (Eldridge et al. 2017; Stanway & Eldridge 2018). However, the BPASS models reach a minimum $\beta$ of $-3.15$, similar to that found using BEAGLE and Prospector, indicating our fiducial models capture the range of possible UV slopes under a variety of assumptions. The solid lines in Figure 2 demonstrate the reddening effect of the nebular continuum in the BEAGLE and Prospector models. The impact of the nebular continuum is strongest at young ages and low metallicities (see Figure 1), causing the intrinsic UV slope (i.e., before dust attenuation) to approach $\beta \approx -2.6$ at the youngest ages. We note that once the nebular continuum has been added, the intrinsic UV slopes are not strongly dependent on age (for a constant SFH), with older models ($\approx 100$ Myr) fully capable of reproducing very blue UV slopes, provided minimal attenuation is present.

It is clear that according to fiducial models within BEAGLE
and Prospector, selecting galaxies with $\beta < -3$ will pick out a
population with moderately low stellar metallicity and
significant leakage of ionizing photons from the H II regions,
and this threshold motivates our search in Section 3. To self-
consistently model sources we find, we include nebular
emission in conditions with nonzero $f_{\text{esc,H II}}$, adopting models
with varying amounts of ionizing radiation escape calculated in
the density-bounded regime from Plat et al. (2019). These

models span a large range of $f_{\text{esc,H II}}$, up to extreme values of
50%–80%, consistent with constraints of the most extreme
leakers observed at lower redshift (e.g., Vanzella et al. 2018).
To test how imposing varying amounts of $f_{\text{esc,H II}}$ on the model
SEDs impacts the resulting UV slopes and nebular emission
line strengths, we construct grids of models at fixed $f_{\text{esc,H II}}$
spanning the full range of ages (1 Myr–1 Gyr) and metallicities
($\log(Z/Z_\odot) = -2-0.3$) available and assume no dust.

The resulting ranges of $\beta$ and [O III]+H$\beta$ EW found throughout
these variable $f_{\text{esc,H II}}$ grids are shown in Figure 3. As expected,
increasing $f_{\text{esc,H II}}$ allows models to be constructed that yield bluer
UV slopes. Additionally, larger values of $f_{\text{esc,H II}}$ restrict the
maximum EWs that can be observed, in line with the results from
previous studies (e.g., Zackrisson et al. 2013; Plat et al. 2019).
However, we note that strong emission lines are still present in
some cases with significant leakage of ionizing photons. For
example, we find [O III]+H$\beta$ EWs of over 1000 Å for young ages
(1–3 Myr) and moderate metallicities (0.2–0.5 $Z_\odot$) with
significant leakage ($f_{\text{esc,H II}} = 0.8$). This follows naturally from
the fact that emission lines for these young sources are extremely
strong in BEAGLE models, and hence the lines remain
moderately strong even if the flux has been reduced by the
escape of ionizing photons. It is only at lower gas-phase
metallicities and slightly older ages where the lines become
considerably weaker in the high $f_{\text{esc,H II}}$ models. We thus expect
very blue UV slopes to flag sources with young stellar ages and
large values of $f_{\text{esc,H II}}$, and the strength of the nebular lines then
dictates whether the metallicity is very low (weak lines) or
moderate (stronger lines).

The model predictions outlined in Figure 3(b) demonstrate
the importance of assuming no dust attenuation for this
selection. For example, even the presence of dust corresp-
onding to a $\tau_{V,eff} = 0.1$ reddens the UV continuum yielding a
UV slope that is indistinguishable from an unattenuated SED
with $f_{\text{esc,H II}} = 0$ for all but the most extreme values of $f_{\text{esc,H II}}$.
As such, focusing our selection toward the bluest UV slopes
provides the best chance of identifying objects with large
$f_{\text{esc,H II}}$. Additionally, while the increased presence of dust can
lead to weaker emission lines, the amount of attenuation
required will require a significant reddening of the UV
continuum to values outside our selection. Variation to the
strength of these emission lines can also occur due to changes
to the IMF (e.g., Zackrisson et al. 2013). However, assuming
an IMF resulting in bluer UV slopes would likely strengthen
the emission lines, making such objects less attractive as high
$f_{\text{esc,H II}}$ candidates based on our selection.

Based off what we have reviewed in this section, our selection in the NIRCam imaging of the EGS field will focus on galaxies with $\beta \approx -3$, taking care to ensure that the slope is blue in the FUV and NUV. We will then look to the flux excesses in the rest optical to probe the rest-optical line strengths, but we will not rule sources out if they have moderately strong line emission.

**3. Observations and Analysis**

3.1. NIRCam Imaging and Photometry of the EGS Field

We utilized imaging obtained as part of the first data release
from the Cosmic Evolution Early Release Science (CEERS;
Finkelstein et al. 2017) survey.<fnref n="3" /> These data include deep JWST/

<footnote n="3">³ https://ceers.github.io/</footnote>
<figure><img src="image_7.png" /><figcaption>**Figure 3.** (a) Regions of the β-EW parameter space probed by SED models with varying amounts of $f_{esc,H}$ II. These regions were derived by constructing grids of SED models spanning a broad range of ages and metallicities (i.e., log(Age/yr) = 6 to 9); log(Z/Z⊙) = -2 to 0.3). Models with an $f_{esc,H}$ II = 1 will trend toward bluer β, down to a minimum of -3.2 (see Figure 2) but do not produce any line emission (i.e., EW). As such, they are not displayed in this figure. (b) Same as left panel, except calculated from BEAGLE models with a modest attenuation from dust, corresponding to $\tau_{V,eff}$ = 0.1. The addition of dust in these models results in reddened UV slopes and decreased EWs.</figcaption></figure>

NIRCam imaging over the EGS observed during 2022 June in three short-wavelength (SW) bands (F115W, F150W, and F200W), in addition to four long-wavelength (LW) bands comprising three wide and one medium band (F277W, F356W, F410M, and F444W). A detailed description of the reduction procedure is provided in Whitler et al. (2022b) and Endsley et al. (2022); however, we give a brief summary here. Calibrated coadded mosaics for each band were produced by first processing individual detector exposures through the JWST pipeline<fnref n="4" /> after implementing 1/f noise subtraction and background subtraction. We adopt photometric zero-points derived by Gabe Brammer,<fnref n="5" /> calculated using observations of a flux calibration standard in NIRCam for a subset of filters in module B in addition to NIRISS and NIRCam observations of the LMC astrometric calibration field for the remaining filters and module A. These zero-points have been shown to produce reasonably consistent color-magnitude diagrams from the JWST Resolved Stellar Populations Early Release Science program (Boyer et al. 2022). The coadded images of individual filter, detector, and pointing combinations were then matched to the Gaia astrometric frame by running TWEAKREG using the Gaia-aligned HST/WFC3 F160W Complete Hubble Archive for Galaxy Evolution (CHArGE) mosaic (see below) as a reference image (see Chen et al. 2022 for more details). The resulting alignment achieves an rms offset of ~6–15 mas compared to the Gaia-aligned HST imaging. The calibrated and aligned NIRCam images for each detector and pointing were then coadded into final mosaics for each band and resampled onto a common pixel grid with 0.03/pixel spacing.

We additionally utilized deep observations from HST collated as part of the CHArGE. Specifically, this additional data comprises HST/Advanced Camera for Surveys (ACS)

F435W, F606W, and F814W, obtained as part of the All-Wavelength Extended Groth Strip International Survey (AEGIS; Davis et al. 2007), Cosmic Assembly Near-Infrared Deep Extragalactic Legacy Survey (CANDELS; Grogin et al. 2011; Koekemoer et al. 2011), and UVCANDELS (PI: Teplitz<fnref n="6" />) surveys. Imaging data from these programs were combined using the software package GRIZLI (Brammer & Matharu 2021) and then subsequently registered to the Gaia coordinate system and sampled at a consistent pixel scale of 0.04 pixel. Finally, to provide consistent photometric measurements across all filters, the JWST/NIRCam LW bands were point-spread function (PSF)-matched to those of WFC3/F160W mosaics (also from CHArGE; 0.08 pixel), while the JWST/NIRCam SW and HST/ACS images were adjusted so that their PSF matched those of the ACS/F814W image. Details regarding this PSF-matching procedure largely follow that of Endsley et al. (2021) and are discussed in more detail in Endsley et al. (2022).

Photometry on this combined HST/ACS and JWST/NIRCam data set were calculated using elliptical apertures constructed to capture the full galaxy emission in each of the PSF-homogenized SW and LW bands. Errors on the photometry were estimated by placing many iterations of the same-sized aperture among blank regions of the sky surrounding each object. Corrections for contamination from nearby bright objects were made by fitting models to the flux distribution of the contaminant and subtracting off the excess emission. Further details regarding the photometric measurements are provided in Whitler et al. (2022b) and Endsley et al. (2022).

## 3.2. Selection of $z > 7$ Galaxies with Very Blue UV Slopes

The deep observations provided by JWST/NIRCam allow robust selections within several redshift windows at $z > 7$. Our

<footnote n="4">⁴ https://jwst-pipeline.readthedocs.io/en/latest/index.html</footnote>

<footnote n="5">⁵ https://github.com/gbrammer/grizli/pull/107</footnote>

<footnote n="6">⁶ https://archive.stsci.edu/hlsp/uvcandels</footnote>
<figure><img src="image_8.png" /><figcaption>**Figure 4.** Rest-frame coverage of available NIRCam filters throughout the rest-UV continuum. Bands are shown at three different redshifts spanning the different photometric selections described in Section 3.1. The positions of Lyα and the Balmer break are displayed as the vertical lines. Filters that extend past these boundaries, and thus are less sensitive to the UV continuum, are displayed at a reduced opacity.</figcaption></figure>

group has developed selections of galaxies at redshifts spanning $z \sim 6.5-11$, comprising two photometric dropout selections at $z \sim 6.5-8$ and $z \sim 8.5-11$. A detailed discussion of these selections are provided in Whitler et al. (2022b) and Endsley et al. (2022) for the lower and higher redshift windows, respectively, and we provide a brief description here. Initial selection of $z \sim 6.5-8$ galaxies (Endsley et al. 2022) was obtained by objects that satisfy the following conditions:

$$ 
\begin{align*}
(F200W < 28) \land (F606W - F115W > 1.5) \land \\
(F814W - F115W > 1.5) \land (F115W - F200W < 0.5) \land \\
(F814W - F115W > F115W - F200W + 1.5).
\end{align*}
$$

For the purposes of these color cuts, the flux in the dropout bands (i.e., F606W and F814W) are set to their $1\sigma$ upper limit in cases of nondetections (signal-to-noise ratio, S/N < 1). We additionally require nondetections (S/N < 2) in ACS/F435W, F606W, and F814W and optical $\chi^2$ (Bouwens et al. 2015) of $\chi^2_{opt} < 5$ calculated using the three ACS bands. A S/N > 3 detection is required in at least two JWST/LW bands, with a constraint of an S/N > 5 detection in F200W.

A selection of $z \sim 9-11$ galaxies (Whitler et al. 2022b) is obtained through the following photometric conditions ensuring placement of the Lyα break:

$$ 
\begin{align*}
(F115W - F150W > 1.0) \land (F150W - F277W < 0.4) \land \\
(F115W - F150W > 0.8 \times F150W - F277W + 1.5),
\end{align*}
$$

with an added requirement of nondetection (S/N < 2) in the ACS bands and a detection at the $>7\sigma$ level in F200W and $>3\sigma$ for all LW bands. Similar to the approach above, the flux in the F115W dropout band is set to its $1\sigma$ upper limit in the case of nondetection.

This high-redshift sample is supplemented in Whitler et al. (2022b) by a selection at slightly lower redshift ($z \sim 8.5-9$), where the observed-frame Lyα break yields a partial decrease in the emission within F115W and may not suppress the entirety of the flux in this band. Specifically, galaxies in this redshift window are identified by the following criteria:

$$ 
\begin{align*}
(F115W - F150W > 0.6) \land (F150W - F277W < 0.4) \land \\
(F115W - F150W > 1.5 \times F150W - F277W + 0.6).
\end{align*}
$$

As before, we additionally require a $>7\sigma$ detection in F200W and a $>3\sigma$ detection in all LW bands, with the further requirement of a $<28$ magnitude in F200W and nondetection ($<2\sigma$) in the ACS bands. We characterize photometric redshift probability distributions using BEAGLE, removing those sources with large probabilities ($>30\%$) of being at $z < 8.2$.

We derive properties for the 123 galaxies in our $z \sim 7-11$ sample using BEAGLE (Chevallard & Charlot 2016). We broadly follow the approach described in Section 2. For galaxies selected as part of the $z \sim 7$ selection, we assume uniform priors in log space on both metallicity and ionization parameter, sampled over the full range allowed by the models (i.e., $\log(Z/Z_\odot) \in [-2.2, 0.3]$; $\log(U) \in [-4, -1]$). Our models assume an SMC dust law (Pei 1992), which has been shown to reproduce the observed properties of high-redshift galaxies (e.g., Bouwens et al. 2016; Reddy et al. 2018). All models are constructed using a constant SFH, and we impose a log-uniform prior on the assumed age ranging from 1 Myr to the age of the universe at the given redshift. We define these as our fiducial models and will consider the influence of varying the escape fraction using density-bounded models in Section 4. A detailed discussion comprising the full set of properties derived for the $z \sim 7$ and $z \gtrsim 8.5$ samples is provided in Endsley et al. (2022) and Whitler et al. (2022b), respectively.

To measure $\beta$ from the observed photometry, we utilize available NIRCam filters that are sensitive to the rest-frame UV continuum. Crucially, we avoid filters that may include emission from Lyα or may be biased by including wavelengths in the Lyα forest. We additionally exclude filters that extend redward of the Balmer break in the rest frame, as they may be biased due to an older stellar population or other effects. Figure 4 illustrates the rest-frame wavelengths probed by the available filters, in addition to examples of the selection of broadband filters used at several redshifts. Example images in the observed broadband filters for two objects in our sample are provided in Figure 5, demonstrating a clear detection at rest-UV wavelengths. UV slopes are then calculated for each object by a power-law fit to the observed photometry that is confined to wavelengths bounded by Lyα and the Balmer break. Uncertainties on $\beta$ are calculated using Monte Carlo simulations, where both the observed photometry is perturbed by their corresponding errors and the photometric redshifts are sampled from the derived posterior probability distribution. With each iteration of these simulations, the set of filters considered to probe the UV continuum are reestablished due to the changing redshift.

Figure 6 provides derived UV slopes for our sample as a function of their absolute UV magnitude. On average, we find that the galaxies in our sample have blue UV slopes, represented by a median value of $-2.0$ for the $z \sim 7$ sample and $-1.9$ for the $z \sim 8.5-11$ sample, with the reddest galaxies reaching values of $\beta = -0.7$. These results are consistent with other measurements of $z \gtrsim 8$ galaxies derived from NIRCam imaging, which cover a similar range of UV luminosity (e.g., Cullen et al. 2022). We note that these are averages of our
<figure><img src="image_9.png" /><figcaption>**Figure 5.** Image cutouts in all NIRCam filters and combined color image of the two extremely blue objects in the sample. Postage stamps are all 1" on a side centered at the location of the object. All of the objects present clear detections in multiple filters, indicating they are not artifacts.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>**Figure 6.** UV slopes as a function of $M_{UV}$ for the full sample of $z \sim 7$ (blue circles; Endsley et al. 2022) and $z \sim 8.5-11$ (red diamonds; Whitler et al. 2022b) galaxies. The two objects identified with extremely blue UV slopes and described in Tables 1 and 3 are shown as the yellow stars.</figcaption></figure>

measured data and should not be used as averages of the population. Such analysis requires comprehensive simulations to account for systematics of the selection function (e.g., Dunlop et al. 2012; Rogers et al. 2013; Bouwens et al. 2014), which are beyond the scope of the present work. Two galaxies within these high-redshift selections stand out as having UV slopes significantly bluer than the majority of the sample. These two systems have derived photometric redshifts of $z_{phot} = 6.7$ and $z_{phot} = 7.2$ and show UV slopes of $\beta = -2.9_{-0.3}^{+0.1}$ and $\beta = -3.1_{-0.1}^{+0.3}$, respectively (see Table 1). Each of these systems matches the selection motivated in Section 2 and is a good candidate for LyC leakers with reasonably low metallicities. For both galaxies, solutions at low redshift ($z < 6$) comprise a minimal fraction (<5%) of the total redshift probability distribution. We note that galaxies with UV slopes as blue as $\beta \sim -3$ have been found in other samples of high-redshift galaxies selected from NIRCam imaging (e.g., Atek et al. 2022; Castellano et al. 2022; Furtak et al. 2022).

The observed SEDs are presented in Figure 7, and the corresponding broadband fluxes are provided in Table 2. In each object, the blue UV slope is clearly visible across at least three NIRCam filters, extending across both the FUV and

<table><caption>**Table 1**

Properties of the Extremely Blue Sample</caption><thead><tr><th>Object ID</th><th>R.A.</th><th>Decl.</th><th>z<sub>phot</sub></th><th>F200W</th><th>β</th></tr></thead><tbody><tr><td>EGS-9466</td><td>214.77676</td><td>+52.84293</td><td>6.7<sup>+0.1</sup><sub>-0.1</sub></td><td>27.9</td><td>-2.9<sup>+0.1</sup><sub>-0.3</sub></td></tr><tr><td>EGS-16168</td><td>214.79312</td><td>+52.87012</td><td>7.2<sup>+0.1</sup><sub>-0.2</sub></td><td>27.4</td><td>-3.1<sup>+0.3</sup><sub>-0.1</sub></td></tr></tbody></table>

**Note.** Photometric redshifts are derived using BEAGLE models. As described in Section 3.2, UV slopes are calculated by fitting a power-law slope to the photometry residing within the rest UV.

NUV. This helps reduce the likelihood that the UV slopes arise due to random scatter, while also satisfying expectations for minimal contribution from the nebular continuum, as would be expected if there is significant ionizing photon leakage (see Section 2). The SEDs of these reveal weak-to-moderate photometric excesses in F356W where strong nebular emission lines (i.e., [O III]+Hβ) contribute to the flux. The inferred strength of the rest-optical emission lines, in combination with the exceedingly blue UV slopes, are preferentially fit with SED models that allow for the escape of ionizing radiation. The models describing these galaxies, as well as implications for $f_{esc,H II}$ and metallicity, are described in detail in the next section.

While these extremely blue sources possess moderately weak [O III]+Hβ EWs, a population of very young (e.g., <10 Myr) galaxies with even weaker line emission has emerged at $z \sim 6.5-8$ (Endsley et al. 2022). While no objects in this population, which comprise ~20% of the $z \sim 6.5-8$ selection, exhibit extremely blue UV slopes, the weak line emission is a characteristic of LyC photon escape. While nonzero $f_{esc,H II}$ models are not strongly preferred by BEAGLE, they do still provide an accurate description of this photometry. Specifically, when $f_{esc,H II}$ is fixed to zero, the weak lines are best reproduced by extremely low metallicities (~1%). In this regime, the SED fits arrive at a typical minimum $\chi^2$ of 4.5. When nonzero $f_{esc,H II}$ is assumed, there is more flexibility in the inferred properties while still being able to fit the lines. At high values of $f_{esc,H II}$, we still achieve a typical minimum $\chi^2$ of 5.0, indicating fits with similar fidelity to the $f_{esc,H II} = 0$ models. The BEAGLE fits indicate the presence of dust, corresponding to a $\tau_{V,eff} = 0.13$, preventing UV slopes from reaching values of ~-3. So while these objects do not possess the combined extremely blue UV slopes and low emission line strengths to confidently suggest high levels of $f_{esc,H II}$, such
<figure><img src="image_11.png" /><figcaption>**Figure 7.** Best-fit SEDs for the two extremely blue objects in our sample. The best-fit SED calculated from BEAGLE using models that allow for the escape of ionizing radiation is displayed in the left panel, while the model with fixed $f_{esc,HII} = 0$ is given on the right. In both panels, the observed photometry is displayed as the cyan diamonds, while the predicted photometry from the models is displayed as colored squares. The minimum $\chi^2$ found for each set of models, as well as the best-fit $f_{esc,HII}$, are given in the lower right of each panel.</figcaption></figure>

<table><caption>**Table 2**
Flux Measurements and Errors for the Two Galaxies Identified with Very Blue UV Slopes in nJy</caption><thead><tr><th>Object ID</th><th>F435W</th><th>F606W</th><th>F814W</th><th>F115W</th><th>F150W</th><th>F200W</th><th>F277W</th><th>F356W</th><th>F410M</th><th>F444W</th></tr></thead><tbody><tr><td>EGS-9466</td><td>-0.1 ± 4.0</td><td>0.6 ± 3.6</td><td>2.7 ± 5.4</td><td>60.1 ± 1.9</td><td>49.7 ± 2.2</td><td>38.4 ± 2.2</td><td>32.5 ± 1.9</td><td>51.0 ± 2.1</td><td>28.2 ± 3.9</td><td>35.5 ± 4.0</td></tr><tr><td>EGS-16168</td><td>-0.5 ± 6.2</td><td>4.7 ± 4.5</td><td>2.4 ± 7.5</td><td>36.4 ± 2.7</td><td>27.9 ± 3.7</td><td>24.2 ± 3.2</td><td>16.0 ± 3.8</td><td>24.9 ± 3.9</td><td>41.4 ± 7.3</td><td>26.3 ± 4.8</td></tr></tbody></table>

**Note.** The derivation of these fluxes is described in Section 3.1.

properties may still be present in these systems. A complete discussion of the potential characteristics of this population is provided in Endsley et al. (2022).

## 4. The Physical Properties of $z > 7$ Galaxies with Extremely Blue UV Slopes

In the previous section, we identified two objects in $z \simeq 7-11$ galaxy selections that have extremely blue ($\beta \simeq -3$) UV continuum slopes. In this section we explore the properties of these systems in more detail. We first consider properties derived from our fiducial (ionization-bounded) BEAGLE

models. We follow the same model setup as described in Section 3 for our initial photometric redshift runs to characterize the properties of our galaxies. Namely, we apply a uniform prior in redshift restricting the allowed values to $z = [5, 12]$, impose a log-uniform prior on metallicity in the range $log(Z/Z_\odot) \in [-2.2, 0.3]$, and allow the ionization parameter to vary within $log(U) \in [-4, -1]$ with a uniform prior. We also fit each object using ionization-bounded models that approximate the effects of alpha enhancement, such that the metallicity of the stars (set by iron and iron-peak elements) is lower than that of the gas (set by oxygen and other alpha
<table><caption>**Table 3**

Properties of the Extremely Blue Objects Derived Using Best-fit BEAGLE Models Allowing for Variations in $f_{esc,H \text{ II}}$, as Described in Section 4</caption><thead><tr><th>Object ID</th><th>log(Z/Z<sub>☉</sub>)</th><th>τ<sub>V</sub></th><th>log(M/M<sub>☉</sub>)</th><th>SFR/M<sub>☉</sub> yr<sup>-1</sup></th><th>Age (Myr)</th><th>f<sub>esc,H II</sub></th><th>sSFR (Gyr<sup>-1</sup>)</th><th>[O III]+Hβ EW (Å)</th><th>log(U)</th></tr></thead><tbody><tr><td>EGS-9466</td><td>-1.6<sup>+0.8</sup><sub>-0.4</sub></td><td>0.01<sup>+0.02</sup><sub>-0.01</sub></td><td>7.8<sup>+0.4</sup><sub>-0.1</sub></td><td>5<sup>+7</sup><sub>-3</sub></td><td>12<sup>+18</sup><sub>-7</sub></td><td>0.74<sup>+0.18</sup><sub>-0.44</sub></td><td>79<sup>+121</sup><sub>-46</sub></td><td>440<sup>+80</sup><sub>-110</sub></td><td>-1.8<sup>+0.3</sup><sub>-0.2</sub></td></tr><tr><td>EGS-16168</td><td>-1.5<sup>+0.7</sup><sub>-0.4</sub></td><td>0.02<sup>+0.02</sup><sub>-0.02</sub></td><td>8.0<sup>+0.4</sup><sub>-0.3</sub></td><td>3<sup>+5</sup><sub>-2</sub></td><td>40<sup>+40</sup><sub>-28</sub></td><td>0.49<sup>+0.40</sup><sub>-0.32</sub></td><td>30<sup>+50</sup><sub>-20</sub></td><td>410<sup>+190</sup><sub>-140</sub></td><td>-2.2<sup>+0.4</sup><sub>-0.5</sub></td></tr></tbody></table>

Note. The given uncertainties on these values correspond to the 16th and 84th percentiles of the probability distributions.

elements). And we consider the density-bounded models with varying $f_{esc,H \text{ II}}$ described in Section 2 (see Plat et al. 2019). An overview of the derived properties of these two galaxies is provided in Table 3 and discussed below.

In Figure 7, we display the best-fit SED using the models that allow for the escape of ionizing radiation from H II regions, as well as the results when we fix $f_{esc,H \text{ II}} = 0$. In both objects, the models with varying $f_{esc,H \text{ II}}$ provide much better fits to the observed photometry. This agreement is apparent both visually and based on $\chi^2$, such that for EGS-9466 and EGS-16168 the varying $f_{esc,H \text{ II}}$ models provide a $\chi^2 = 3.8$ and $\chi^2 = 5.9$, while the $f_{esc,H \text{ II}} = 0$ models only achieve $\chi^2 = 10.3$ and $\chi^2 = 8.8$, respectively. As pointed out in the previous sections, simultaneously reproducing the shape of the UV continuum and relative weakness of the emission requires a reduction in the amount of nebular emission contributing to the flux. The density-bounded models are able to reproduce both the shape of the UV continuum and inferred strength of the emission lines simultaneously, which is not possible under the constraint of $f_{esc,H \text{ II}} = 0$. In the fiducial (ionization-bounded) models, the UV slopes are clearly very poorly matched owing to the contribution of the nebular continuum.

The increased ability to explain the observed photometry with these updated models allows us to explore the stellar populations within these extreme objects at high redshift. In addition to allowing $f_{esc,H \text{ II}}$ to vary within the model galaxy spectra, other exotic populations, such as Population III stars, have been introduced as potential solutions to poorly fit extremely blue SEDs. We find that the observed photometry of these galaxies is well fit through only a high $f_{esc,H \text{ II}}$, and such exotic populations are not necessary. We additionally fit the SEDs of these objects using models that are $\alpha$ enhanced (i.e., $Z_* < Z_{ISM}$), corresponding to conditions observed in high-redshift galaxies (e.g., Steidel et al. 2016; Sanders et al. 2020; Topping et al. 2020; Strom et al. 2022). Within these models, the stellar and gas-phase metallicities are decoupled, thus enabling a hard ionizing spectrum (from low-metallicity stars) to irradiate ISM with moderate metallicities, thereby yielding stronger lines. Specifically, we combined Bruzual & Charlot (2003) stellar templates with nebular emission calculated with Ferland et al. (2017) and tested the case where the stellar metallicity is fixed to 0.2 $Z_{ISM}$ (corresponding to the theoretical limit of Type II supernova enrichment) and also where the stellar metallicity is fixed to $Z_* = 10^{-4}$ and $Z_{ISM}$ is allowed to vary. However, models allowing for alpha enhancement, but without the escape of ionizing radiation, are not able to explain the observed photometry. The high level of inferred $f_{esc,H \text{ II}}$ implies that the UV continuum is increasingly sensitive to continuum emission from the most massive stars. This sensitivity in the absence of nebular contamination provides an improved view into the properties of the massive star population. The BEAGLE SED models overwhelmingly prefer solutions with low stellar metallicities, such that we find a

value of $\log(Z/Z_{☉}) = -1.6$ and $\log(Z/Z_{☉}) = -1.4$ for EGS-9466 and EGS-16168, respectively. As suggested in Section 2, these low metallicities are required to reproduce the observed color excesses in the LW filters, which correspond to [O III]+Hβ EWs of 440 Å and 410 Å for EGS-9466 and EGS-16168, respectively. High-mass stars at such extremely low metallicities would be intense sources of ionizing radiation that, in combination with large values of $f_{esc,H \text{ II}}$, would make these objects highly efficient ionizing agents at this epoch.

We further explore the physical properties of these systems derived from the SED modeling approach with varying $f_{esc,H \text{ II}}$. These properties are summarized in Table 3. Unsurprisingly, we find that both objects have minimal dust content, with V-band optical depths ranging from $\tau_V = 0.01$ to $\tau_V = 0.02$. Furthermore, both objects are relatively low-mass systems, with similar inferred stellar masses of $\log(M/M_{☉}) = 7.8$ and $\log(M/M_{☉}) = 8.0$. By assuming a constant SFH, we find these two objects span a range of young ages of 12 and 40 Myr, consistent with constraints at the young end of the age distribution among galaxies near this epoch (e.g., Endsley et al. 2021; Whitler et al. 2022a). We constrain the rest-UV sizes of these two galaxies using the NIRCam SW bands and following the procedure outlined in Chen et al. (2022). Both galaxies are consistent with being unresolved, indicating that they are highly compact systems, with effective rest-UV radii below 0.26 kpc (see Figure 5). The galaxy properties provided by the SED fitting, in combination with the small size constraints, indicate these systems have high SFR surface densities ($\Sigma_{SFR} \equiv \frac{\text{SFR}/2}{\pi r_{\text{eff}}^2}$). The unresolved nature of these sources means we can only place lower limits on $\Sigma_{SFR}$ of $>11.1 M_{☉} \text{ yr}^{-1}$ kpc$^{-2}$ and $>7.3 M_{☉} \text{ yr}^{-1}$ kpc$^{-2}$, respectively, for EGS-42501 and EGS-47688, implying their true values may be much higher.

These measured and derived properties allow us to place the two $z > 7$ high-$f_{esc,H \text{ II}}$ candidates in context with samples of known LyC leakers found at lower redshifts, where direct detection of the ionizing radiation escape is possible due to relatively low attenuation from the intergalactic medium. In particular, we compare to the Low-Redshift Lyman Continuum Survey (LzLCS; Flury et al. 2022a, 2022b) comprising a sample of objects at $z \sim 0.3$ and the Keck Lyman Continuum Survey (KLCS; Steidel et al. 2018; Pahl et al. 2022) composed of $z \sim 3$ Lyman break galaxies. The comparison of our two objects with these lower-redshift samples is provided in Figure 8. We note that within these comparisons, the $f_{esc}$ values of the low-redshift samples are global escape fractions based on direct detection, while our results comprise candidates of $f_{esc,H \text{ II}}$ defined locally within the H II regions of the high-redshift systems.

In Figure 8(a) we compare the SFRs and sizes of objects within low-redshift samples to our objects and indicate low-redshift galaxies with confirmed $f_{esc}$. It is apparent that galaxies
<figure><img src="image_12.png" /><figcaption>**Figure 8.** Comparison of our two high-*f*<sub>esc,H II</sub> candidates with low-redshift samples. (a) SFR vs. effective radius for the two objects in our sample (yellow stars), the KLCS sample (red circles; Steidel et al. 2018; Pahl et al. 2022), and the LzLCS (blue circles; Flury et al. 2022a). Lines of constant SFR surface density are shown for reference. (b) Distribution of sSFR surface density (Σ<sub>sSFR</sub>) for the LzLCS sample. The values for the two objects in our sample are indicated by the vertical yellow lines. (c) *f*<sub>esc</sub> as a function of *β*. We note that the *f*<sub>esc</sub> of the LzLCS represents a direct measure of the global escape fraction, while the *f*<sub>esc,H II</sub> for our objects is a measure local to the H II regions. Furthermore, *β* for these two samples is measured over a different wavelength range; however, the two measures of the UV slope can be compared in a relative sense, such that objects would be “blue” using both measures.</figcaption></figure>

in the low-redshift universe that leak ionizing radiation span a wide range of SFRs, sizes, and thus SFR surface densities. While our sources only contain upper limits on Σ<sub>SFR</sub>, it is clear that they are consistent with having Σ<sub>SFR</sub> at the high end of the comparison samples. Significant escape of ionizing radiation may be expected among high-Σ<sub>SFR</sub> systems, as the intense SFR activity can aid in the creation of channels that allow ionizing radiation to escape (e.g., Ma et al. 2016). However, high Σ<sub>SFR</sub> does not appear to be a strict requirement for nonzero f<sub>esc</sub> in these low-redshift systems. The specific SFR surface density (Σ<sub>sSFR</sub> ≡ Σ<sub>SFR</sub>/M<sub>*</sub>) provides an alternative metric that has added sensitivity to the depth of the galaxy potential, which may be an important ingredient in regulating the covering fraction and thus, escape of LyC photons (e.g., Reddy et al. 2022). Strikingly, our two high-f<sub>esc,H II</sub> candidates lie at the highest values probed by the low-redshift sample (see Figure 8(b) supporting this physical picture); however, better statistics will be required to further confirm the prevalence of galaxies with high inferred f<sub>esc,H II</sub> among the high-Σ<sub>sSFR</sub> population. Finally, we examine f<sub>esc</sub> as a function of β for our two objects in the context of results from the LzLCS sample in Figure 8(c). It must be pointed out that due to the available data, both f<sub>esc</sub> and β have different definitions between the two samples, and thus concrete statements regarding their comparison are challenging. However, it is clear from this figure, and discussed in detail using a slightly alternative approach in Chisholm et al. (2022), that the abundance of high-f<sub>esc</sub> objects is qualitatively greater among galaxies with bluer β, suggesting a potential link between the two quantities.

The low number of such blue galaxies makes it challenging to robustly constrain the abundance of these potentially highly efficient LyC leakers in the epoch of reionization. Despite their apparent low abundance, they may become more common at the faint end of the UV luminosity function and indeed may comprise a significant contribution of the total ionizing photon budget at these earlier epochs. Future work using JWST will be needed to confirm the presence of this population and constrain their typical properties. The newly available spectroscopic capabilities from JWST will further cement the properties of

these extreme objects, providing both precise redshifts as well as constraints of internal ISM properties.

**5. Summary and Outlook**

In this paper we have described a search for galaxies with extremely blue UV continuum slopes (β ≃ −3) and moderate-to-weak rest-optical nebular line emission (via color excesses) at z ≃ 7–11 in the CEERS NIRCam imaging of the EGS field. Objects that satisfy these criteria are likely to be at relatively low metallicity with a large fraction of their ionizing photons escaping from their H II regions. We summarize our main conclusions below.

1. We characterize the distribution of UV slopes for a sample of 123 galaxies selected to lie at z = 7 − 11, finding a median value of β = −2.0. We emphasize that this is merely the observed average, and insight into the intrinsic UV slope distribution (e.g., Finkelstein et al. 2012b; Dunlop et al. 2012; Rogers et al. 2013; Bouwens et al. 2014; Bhatawdekar & Conselice 2021) requires detailed simulations that are beyond the scope of this study. Nonetheless, these measurements suggest that our parent sample is comprised of moderately blue UV colors on average, consistent with expectations for these redshifts.

2. Among the sample of high-redshift galaxies, we identify two UV-faint (M<sub>UV</sub> = −19.7 to −19.3) galaxies that exhibit extremely blue UV slopes of β = −2.9 and β = −3.1 with relatively weak-to-moderate color excesses. Collectively, these observations suggest diminished nebular emission (continuum and lines) and low metallicities. We find that the SEDs cannot be fit well by our fiducial ionization-bounded (f<sub>esc,H II</sub> = 0) BEAGLE models. We also consider BEAGLE models that approximate the effects of alpha enhancement (allowing stellar metallicity to be lower than the gas-phase metallicity), finding similarly poor fits. Much better solutions are found when comparing to density-bounded BEAGLE models that include a self-consistent treatment of nebular emission while allowing for the escape of
ionizing radiation (developed in Plat et al. 2019). We find
that simultaneously reproducing the very blue UV slope
and relatively weak color excesses requires density-
bounded models with escape fractions in the range $f_{\text{esc,H II}}$
=0.49–0.74. These objects thus appear to be candidates
for significant LyC leakage at $z > 7$.

3. We investigate the properties of these extremely blue objects and find they are low mass ($\log(M/M_{\odot}) \sim 7.8-8.0$), young (12–40 Myr), and are consistent with being dust free ($\tau_{V} \sim 0.01–0.02$). Each galaxy is unresolved in the UV with NIRCam, with very compact sizes implied ($r_{e} \lesssim 270$ pc). This allows us to put lower limits on the SFR surface densities, constraining them to large values ($> (7.3–11) M_{\odot} \text{ yr}^{-1} \text{ kpc}^{-2}$). These systems thus stand out as dwarf galaxies undergoing recent bursts of star formation, leading to relatively high SFR surface densities.

4. The density-bounded models of these two objects imply they have very low metallicities, with best-fit values of $\log(Z/Z_{\odot}) = -1.6$ and $\log(Z/Z_{\odot}) = -1.5$. These low metallicities are implied by the strength of observed photometric excesses resulting from nebular emission lines, such that we infer $[\text{O III}]+\text{H}\beta$ EWs inferred from the models, where we obtain values of 440 Å and 410 Å, for EGS-9466 and EGS-16168, respectively.

5. We place these two objects (candidates for LyC leakers) in context with low-redshift samples containing galaxies with directly constrained $f_{\text{esc}}$. Their SFR surface densities are consistent with those found for leaking systems within low redshift samples and have comparable $\Sigma_{\text{sSFR}}$ to the highest values found in the comparison samples (Pahl et al. 2022). Finally, we find qualitative agreement with the observed trend between $\beta$ and $f_{\text{esc}}$ found at low redshift (Chisholm et al. 2022), such that the strongest leakers typically have very blue UV slopes. Further work will be required to calibrate these relations at lower redshifts to help interpret samples of LyC leaking candidates at $z > 7$ (such as those in this work) that JWST will soon discover.

The identification of these extremely blue galaxies (coupled
with weak-to-moderate nebular line emission) hints at the
possibility of identifying large numbers of LyC leaking
galaxies in the epoch of reionization, where direct observation
of the ionizing continuum is not possible. Objects with large
$f_{\text{esc,H II}}$, and thus low nebular contribution to their SED, offer a
unique window on the shape of the stellar continuum
(indicating in these objects very low metallicities). Coupling
these constraints on the massive star populations with spectro-
scopic observations that are sensitive to the ISM and gas-phase
properties will yield unique insight into early galaxies.

We thank the anonymous referee for their helpful comments. D.P.S. acknowledges support from the National Science Foundation through the grant AST-2109066. R.E. acknowledges funding from NASA JWimST/NIRCam contract to the University of Arizona, NAS5-02015. L.W. acknowledges support from the National Science Foundation Graduate Research Fellowship under grant No. DGE-2137419. The authors thank Jacopo Chevallard for use of the BEAGLE tool used for much of our SED fitting analysis. We thank Gabe Brammer for providing the optical imaging of the EGS as part of CHArGE program. This material is based in part upon High Performance Computing (HPC) resources supported by the University of Arizona TRIF, UITS, and Research, Innovation, and Impact (RII) and maintained by the UArizona Research Technologies department.

Data Availability

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

Posterior Distributions of the $f_{\text{esc,H II}}$ Candidate Galaxies

In order to investigate possible degeneracies between free
parameters assigned in the SED fitting, we provide corner plots
displaying the posterior probability distributions of each fitted
parameter and the covariance between free parameters. These
results from the BEAGLE SED fits are provided in Figures 9
and 10 for EGS-9466 and EGS-16168, respectively.
<figure><img src="image_13.png" /><figcaption>**Figure 9.** Posterior probability distributions for the parameters inferred for EGS-9466 using the BEAGLE model setup described in Section 4.</figcaption></figure>
<figure><img src="image_14.png" /><figcaption>Figure 10. Same as Figure 9 but for EGS-16168.</figcaption></figure>

ORCID iDs

Zuyi Chen https://orcid.org/0000-0002-2178-5471

Stéphane Charlot https://orcid.org/0000-0003-3458-2275

References

Asplund, M., Grevesse, N., Sauval, A. J., & Scott, P. 2009, *ARA&A*, **47**, 481

Atek, H., Shuntov, M., Furtak, L. J., et al. 2022, *MNRAS*, in press

Bhatawdekar, R., & Conselice, C. J. 2021, ApJ, 909, 144

Bottorff, M. C., Ferland, G. J., & Straley, J. P. 2006, PASP, **118**, 1176

Bouwens, R. J., Aravena, M., Decarli, R., et al. 2016, ApJ, 833, 72

Bouwens, R. J., Illingworth, G. D., Oesch, P. A., et al. 2010, ApJL, 708, L69

Bouwens, R. J., Illingworth, G. D., Oesch, P. A., et al. 2014, ApJ, 793, 115

Bouwens, R. J., Illingworth, G. D., Oesch, P. A., et al. 2015, ApJ, 803, 34

Boyer, M. L., Anderson, J., Gennaro, M., et al. 2022, *RNAAS*, **6**, 191

Brammer, G., & Matharu, J. 2021, ggrammer/grizli: Release v2021, Zenodo, doi:10.5281/zenodo.5012699

Bruzual, G., & Charlot, S. 2003, *MNRAS*, **344**, 1000

Calzetti, D., Kinney, A. L., & Storchi-Bergmann, T. 1994, ApJ, 429, 582

Castellano, M., Fontana, A., Treu, T., et al. 2022, ApJL, 938, L15

Chabrier, G. 2003, PASP, 115, 763

Chen, Z., Stark, D. P., Endsley, R., et al. 2022, *MNRAS*, in press

Chevallard, J., & Charlot, S. 2016, *MNRAS*, **462**, 1415

Chisholm, J., Saldana-Lopez, A., Flury, S., et al. 2022, *MNRAS*, **517**, 5104

Choi, J., Dotter, A., Conroy, C., et al. 2016, ApJ, 823, 102

Conroy, C., & Gunn, J. E. 2010, ApJ, 712, 833

Conroy, C., Gunn, J. E., & White, M. 2009, ApJ, 699, 486

Cullen, F., McLure, R. J., McLeod, D. J., et al. 2022, arXiv:2208.04914

Davis, M., Guhathakurta, P., Konidaris, N. P., et al. 2007, ApJL, 660, L1

De Barros, S., Oesch, P. A., Labbe, I., et al. 2019, *MNRAS*, **489**, 2355

Dunlop, J. S., McLure, R. J., Robertson, B. E., et al. 2012, *MNRAS*, **420**, 901

Dunlop, J. S., Rogers, A. B., McLure, R. J., et al. 2013, *MNRAS*, **432**, 3520

Eldridge, J. J., Stanway, E. R., Xiao, L., et al. 2017, PASA, **34**, e058

Endsley, R., Stark, D. P., Chevallard, J., & Charlot, S. 2021, *MNRAS*, **500**, 5229

Endsley, R., Stark, D. P., Whitler, L., et al. 2022, arXiv:2208.14999

Ferland, G. J., Chatzikos, M., Guzmán, F., et al. 2017, RMXAA, 53, 385

Finkelstein, S. L., Ryan, R. E., Papovich, C., et al. 2015, ApJ, 810, 71

Finkelstein, S. L., Papovich, C., Salmon, B., et al. 2012a, ApJ, 756, 164

Finkelstein, S. L., Papovich, C., Ryan, R. E., et al. 2012b, ApJ, 758, 93

Finkelstein, S. L., Dickinson, M., Ferguson, H. C., et al. 2017, JWST Proposal ID 1345. Cycle 0 Early Release Science
Fletcher, T. J., Tang, M., Robertson, B. E., et al. 2019, ApJ, 878, 87

Flury, S. R., Jaskot, A. E., Ferguson, H. C., et al. 2022a, ApJS, 260, 1

Flury, S. R., Jaskot, A. E., Ferguson, H. C., et al. 2022b, ApJ, 930, 126

Furtak, L. J., Shuntov, M., Atek, H., et al. 2022, arXiv:2208.05473

Grogin, N. A., Kocevski, D. D., Faber, S. M., et al. 2011, ApJS, 197, 35

Gutkin, J., Charlot, S., & Bruzual, G. 2016, MNRAS, 462, 1757

Jeřábková, T., Kroupa, P., Dabringhausen, J., Hilker, M., & Bekki, K. 2017, A&A, 608, A53

Jiang, L., Cohen, S. H., Windhorst, R. A., et al. 2020, ApJ, 889, 90

Johnson, B. D., Leja, J., Conroy, C., & Speagle, J. S. 2021, ApJS, 254, 22

Koekemoer, A. M., Faber, S. M., Ferguson, H. C., et al. 2011, ApJS, 197, 36

Labbé, I., Gonzalez, V., Bouwens, R. J., et al. 2010, ApJL, 716, L103

Labbé, I., Bouwens, R. J., Illingworth, G. D., et al. 2013, ApJL, 777, L19

Ma, X., Hopkins, P. F., Kasen, D., et al. 2016, MNRAS, 459, 3614

McLure, R. J., Dunlop, J. S., de Ravel, L., et al. 2011, MNRAS, 418, 2074

Nussbaumer, H., & Schmutz, W. 1984, A&A, 138, 495

Oke, J. B., & Gunn, J. E. 1983, ApJ, 266, 713

Ono, Y., Ouchi, M., Mobasher, B., et al. 2012, ApJ, 744, 83

Ono, Y., Ouchi, M., Shimasaku, K., et al. 2010, ApJ, 724, 1524

Pahl, A. J., Shapley, A., Steidel, C. C., Reddy, N. A., & Chen, Y. 2022, MNRAS, 516, 2062

Pei, Y. C. 1992, ApJ, 395, 130

Plat, A., Charlot, S., Bruzual, G., et al. 2019, MNRAS, 490, 978

Raiter, A., Schaerer, D., & Fosbury, R. A. E. 2010, A&A, 523, A64

Reddy, N. A., Oesch, P. A., Bouwens, R. J., et al. 2018, ApJ, 853, 56

Reddy, N. A., Pettini, M., Steidel, C. C., et al. 2012, ApJ, 754, 25

Reddy, N. A., Topping, M. W., Shapley, A. E., et al. 2022, ApJ, 926, 31

Roberts-Borsani, G. W., Bouwens, R. J., Oesch, P. A., et al. 2016, ApJ, 823, 143

Robertson, B. E., Ellis, R. S., Dunlop, J. S., McLure, R. J., & Stark, D. P. 2010, Natur, 468, 49

Rogers, A. B., McLure, R. J., & Dunlop, J. S. 2013, MNRAS, 429, 2456

Sanders, R. L., Shapley, A. E., Reddy, N. A., et al. 2020, MNRAS, 491, 1427

Smit, R., Bouwens, R. J., Labbe, I., et al. 2014, ApJ, 784, 58

Stanway, E. R., & Eldridge, J. J. 2018, MNRAS, 479, 75

Stark, D. P. 2016, ARA&A, 54, 761

Stefanon, M., Bouwens, R. J., Illingworth, G. D., et al. 2022, arXiv:2204.02986

Steidel, C. C., Bogosavljević, M., Shapley, A. E., et al. 2018, ApJ, 869, 123

Steidel, C. C., Strom, A. L., Pettini, M., et al. 2016, ApJ, 826, 159

Strom, A. L., Rudie, G. C., Steidel, C. C., & Trainor, R. F. 2022, ApJ, 925, 116

Topping, M. W., Shapley, A. E., Reddy, N. A., et al. 2020, MNRAS, 495, 4430

Topping, M. W., & Shull, J. M 2015, ApJ, 800, 97

Topping, M. W., Stark, D. P., Endsley, R., et al. 2022, MNRAS, 516, 975

Vanzella, E., Nonino, M., Cupani, G., et al. 2018, MNRAS, 476, L15

Whitler, L., Endsley, R., Stark, D. P., et al. 2022b, MNRAS, in press

Whitler, L., Stark, D. P., Endsley, R., et al. 2022a, arXiv:2206.05315

Wilkins, S. M., Bunker, A. J., Stanway, E., Lorenzoni, S., & Caruana, J. 2011, MNRAS, 417, 717

Wilkins, S. M., Gonzalez-Perez, V., Lacey, C. G., & Baugh, C. M. 2012, MNRAS, 424, 1522

Yamanaka, S., Inoue, A. K., Yamada, T., et al. 2020, MNRAS, 498, 3095

Zackrisson, E., Binggeli, C., Finlator, K., et al. 2017, ApJ, 836, 78

Zackrisson, E., Inoue, A. K., & Jensen, H. 2013, ApJ, 777, 39