The JWST Advanced Deep Extragalactic Survey: Discovery of an Extreme Galaxy Overdensity at z = 5.4 with JWST/NIRCam in GOODS-S

We report the discovery of an extreme galaxy overdensity at z = 5.4 in the GOODS-S field using James Webb Space Telescope (JWST)/NIRCam imaging from JADES and JEMS alongside JWST/NIRCam wide-field slitless spectroscopy from FRESCO. We identified potential members of the overdensity using Hubble Space Telescope+JWST photometry spanning λ = 0.4–5.0 μm. These data provide accurate and well-constrained photometric redshifts down to m ≈ 29–30 mag. We subsequently confirmed N = 81 galaxies at 5.2 < z < 5.5 using JWST slitless spectroscopy over λ = 3.9–5.0 μm through a targeted line search for Hα around the best-fit photometric redshift. We verified that N = 42 of these galaxies reside in the field, while N = 39 galaxies reside in a density around ∼10 times that of a random volume. Stellar populations for these galaxies were inferred from the photometry and used to construct the star-forming main sequence, where protocluster members appeared more massive and exhibited earlier star formation (and thus older stellar populations) when compared to their field galaxy counterparts. We estimate the total halo mass of this large-scale structure to be 12.6≲log10Mhalo/M⊙≲12.8 using an empirical stellar mass to halo mass relation, which is likely an underestimate as a result of incompleteness. Our discovery demonstrates the power of JWST at constraining dark matter halo assembly and galaxy formation at very early cosmic times.


INTRODUCTION
In the local Universe, galaxy clusters represent the largest and most massive gravitationally bound structures, consisting of up to thousands of individual galaxies contained within a virialized or virializing dark matter halo, and representing the most extreme matter overdensities allowed by the standard cosmological paradigm of hierarchical structure formation (White & Rees 1978).In the early Universe, the structures that eventually evolved into the galaxy clusters seen today are referred to as "protoclusters", which consist of fewer individual galaxies contained within more complex dark matter halos that are yet to be virialized (for a review of protoclusters, see Overzier 2016).
Observations have suggested that the majority of the stellar mass in the local Universe resides in massive elliptical galaxies, which are preferentially found within galaxy clusters (Dressler 1980).Additionally, the average formation timescales for galaxies in clusters is shorter than that for analogous galaxies in the field (Webb et al. 2020).These results suggest that the physical processes associated with extreme matter overdensities induce earlier star formation, earlier stellar mass assembly, and earlier quenching.However, massive clusters at relatively high redshift (z = 1 − 2) have also been observed to have large amounts of star formation, on par with field populations (Alberts et al. 2014(Alberts et al. , 2016(Alberts et al. , 2021)).Quantifying these effects across the protocluster to cluster boundary remains an important open problem for extragalactic astronomy (Wang et al. 2013).
The impact of environment on galaxy formation and evolution is best understood in the local Universe, where we can observe the effects of transformational processes such as dynamical relaxation, tidal interactions, and mergers (Zabludoff et al. 1996).However, these processes make it difficult to ascertain important evidence (e.g., both the initial relative positions and velocities of the constituent galaxies) related to the early formation and evolution of the most massive gravitationally bound structures.For this reason, searching for protoclusters in the early Universe offers our best chance of understanding the initial formation and subsequent evolution of galaxy clusters today (e.g., Li et al. 2022;Brinch et al. 2023;Morishita et al. 2023).
In this paper, we present the discovery of an extreme galaxy overdensity at z = 5.4 in the GOODS-S field using data from the Near Infrared Camera (NIRCam; Rieke et al. 2005Rieke et al. , 2023a) ) on JWST.The powerful com-bination of deep imaging and wide field slitless spectroscopy (WFSS) provided by JWST/NIRCam allows us to identify this overdensity, characterize the stellar populations of galaxies both inside and outside this largescale structure, and estimate the dark matter halo mass associated with this protocluster.These observations provide important insights into the impact of environment on galaxy formation and evolution immediately after the epoch of reionization (EoR; z > 6) when the Universe was approximately a billion years old.
This paper proceeds as follows.In Section 2, we describe the various data and observations that are used in our analysis, including the photometric redshift determination and emission line detection.In Section 3, we present our analysis and results, including the stellar population modeling and halo mass inference.In Section 4, we summarize our findings and their implications for galaxy evolution in the early Universe.All magnitudes are in the AB system (Oke & Gunn 1983).Uncertainties are quoted as 68% confidence intervals.Throughout this work, we report wavelengths in vacuum and adopt the standard flat ΛCDM cosmology from Planck18 with H 0 = 67.4km/s/Mpc and Ω m = 0.315 (Planck Collaboration et al. 2020).

DATA & OBSERVATIONS
In this work, we use deep optical imaging from the Advanced Camera for Surveys (ACS) on HST alongside deep infrared imaging and WFSS from JWST/NIRCam in the Great Observatories Origins Deep Survey South (GOODS-S; Giavalisco et al. 2004) field.The imaging data and photometry from HST/ACS and JWST/NIRCam are described in Section 2.1.The photometric redshifts and sample selection are described in Section 2.2.The spectral data and line detection from JWST/NIRCam WFSS are described in Section 2.3.
A detailed description of the JWST/NIRCam imaging data reduction and mosaicing will be presented in a forthcoming paper from the JADES Collaboration (Tacchella et al., in preparation).We briefly summarize here the main steps of the reduction and mosaicing process.The data are initially processed with the standard JWST calibration pipeline1 .Customized steps are included to aid in the removal of "1/f" noise, "wisp" artifacts, "snowball" artifacts, and persistence from previous observations (see also Rigby et al. 2022).The JWST Calibration Reference Data System (CRDS) context map jwst 1008.pmap is used, including the flux calibration for JWST/NIRCam from Cycle 1.The background from the sky is modeled and removed using the BackGround2D class from photutils (Bradley et al. 2022).Finally, the image mosaics for each of the fourteen JWST/NIRCam filters are registered to the Gaia DR3 frame (Gaia Collaboration et al. 2023) and resampled onto the same world coordinate system (WCS) with a 30 mas/pixel grid.Assuming a circular aperture with a diameter of 0.3 ′′ , the 5σ point-source detection limit in the F200W filter is m ≈ 30.0 AB mag and m ≈ 29.0 AB mag for the deep and medium regions, respectively.
A detailed description of the JWST/NIRCam source detection was outlined in Robertson et al. (2023) and will be presented in detail in another forthcoming paper from the JADES Collaboration (Robertson et al., in preparation).We briefly summarize here the main steps of the source detection process.Six image mosaics (F200W, F277W, F335M, F356W, F410M, and F444W) are initially stacked using the corresponding error images and inverse-variance weighting to produce a single detection image.These filters were chosen in order to avoid biasing our catalog against SW dropouts (e.g.dropouts in F090W, F115W, or F150W).In this detection image, we construct a source catalog by selecting contiguous regions of greater than five pixels with signal-to-noise ratios S/N > 3 and applying a standard Source Extractor (SExtractor; Bertin & Arnouts 1996) deblending algorithm with parameters nlevels = 32 and contrast = 0.001 using photutils (Bradley et al. 2022).Finally, we perform forced convolved photometry at the source centroids in all HST/ACS and JWST/NIRCam photometric bands, assuming elliptical Kron apertures with parameter = 1.2 (i.e.Kron small) and parameter = 2.5 (i.e.Kron large).To correct for potential missing light, we rescale the Kron small photometry by the flux ratio of Kron large to Kron small in the F444W filter.Using model point spread functions (PSFs) from the TinyTim (Krist et al. 2011) package for HST/ACS and the WebbPSF (Perrin et al. 2014) package for JWST/NIRCam, we apply aperture corrections assuming point source morphologies.Uncertainties are estimated by placing random apertures across regions of the image mosaics to compute a flux variance (e.g., Labbé et al. 2005;Quadri et al. 2007;Whitaker et al. 2011), which are summed in quadrature with the associated Poisson uncertainty for each detected source.

Photometric Redshifts & Sample Selection
Using the previously described photometry, we measure photometric redshifts with the template-fitting code EAZY (Brammer et al. 2008).A more detailed description of this procedure will be discussed in a forthcoming paper from the JADES Collaboration (Hainline et al., in preparation).We briefly summarize here the main steps of the photometric redshift process.EAZY uses a chi-square (χ 2 ) minimization technique to model the broadband spectral energy distributions (SEDs) for galaxies using linear combinations of galaxy templates.It is designed to be both fast and flexible, and has been used extensively in the literature to model the photometric redshifts of galaxies (e.g., Newman et al. 2013;Skelton et al. 2014;Bouwens et al. 2015).
We fit all of the available photometry for each object with EAZY, assuming the rescaled Kron small photometry described in Section 2.1.For objects in the UDF, this includes photometry in nineteen filters.For objects in the JADES deep region but not in the UDF, this includes photometry in fourteen filters.For objects in the JADES medium region, this includes photometry in thirteen filters.We utilize sixteen templates in total to perform the fitting, which includes the nine EAZY "v1.3" templates, two additional templates for simple stellar populations with ages of 5 Myr and 25 Myr, and five more templates with strong nebular continuum emission that were created using the Flexible Stellar Population Synthesis code (FSPS; Conroy et al. 2009;Conroy & Gunn 2010).These templates span a large range of stellar population properties and include contributions from both nebular continuum and line emission, as well as obscuration from dust.
While the photometric calibration of JWST/NIRCam has improved significantly over the course of the last few months, there is still some uncertainty with these calibrations that needs to be taken into account.To this end, we iteratively calculate the photometric offset from the EAZY templates compared to the true JWST/NIRCam photometry, using a sample of galaxies with signal-to-noise ratios between 5 and 20 in F200W.These photometric offsets are relatively small (on the order of a few percent for both HST/ACS and JWST/NIRCam) and are subsequently applied to the entire photometric catalog.We choose not to adopt any apparent magnitude priors, but we do make use of the template error file "template error.v2.0.zfourge".
The primary measurements used here are the EAZY "z a " and "z peak " redshifts.The former corresponds to the fit where the likelihood is maximized (χ2 is minimized), while the latter corresponds to the fit where the likelihood is maximized (χ 2 is minimized) after taking into account the template error file.We allow EAZY to fit across the redshift range of z = 0.2 − 22 with a redshift step size of ∆z = 0.01 (1 + z).To test the accuracy of these photometric redshifts, we compare these predictions with existing spectroscopic redshifts in the GOODS-S field from the Multi Unit Spectroscopic Explorer (MUSE; Inami et al. 2017;Urrutia et al. 2019).While the MUSE spectroscopic redshifts are biased toward the brightest objects detected by JWST at z < 7, we found catastrophic outlier fractions of only 5% (Rieke et al. 2023b) when comparing to the highest quality spectroscopic redshifts available with MUSE.The catastrophic outlier fraction is defined to be the fraction of objects that satisfy Equation 1: To perform an accurate and efficient targeted emission line search within the available spectroscopic data, we require a sample of relatively bright objects, since these are the only objects for which we expect to detect an emission line (see Section 2.3).We also require these objects to have tight photometric redshift constraints, which allow for spectroscopic redshift confirmation us-ing only a single line detection.Most objects with tight photometric redshift constraints have emission lines that fall in one of the medium band filters (e.g., F410M) which allows for tight constraints when paired with the broad band filter coverage of JADES (e.g., F444W).Our selection criteria for the final photometric catalog consist of the following: m < 28.5 AB mag in F444W assuming elliptical Kron apertures with parameter = 2.5, 4.5 < z a < 9.5, 4.5 < z peak < 9.5, ∆z 1 < 1, and ∆z 2 < 2. The first EAZY confidence interval (∆z 1 ) is defined to be the difference between the 16th and 84th percentiles of the photometric redshift posterior distribution and is roughly twice the standard deviation.The second EAZY confidence interval (∆z 2 ) is defined to be the difference between the 5th and 95th percentiles of the photometric redshift posterior distribution and is roughly four times the standard deviation.

Spectroscopic Data & Emission Line Detection
Our spectroscopic data consist of WFSS observations taken with JWST/NIRCam in the F444W filter (λ = 3.9 − 5.0 µm).These data were obtained by the First Reionization Epoch Spectroscopic COmplete Survey (FRESCO; PI: Oesch; PID: 1895) in November of 2022.The FRESCO observations cover a 8.2 ′ × 8.6 ′ area using the row-direction grisms on both modules of JWST/NIRCam (Grism R; R ≈ 1600).The total overlapping area between the JADES and FRESCO footprints is ≈ 41 square arcminutes.The total spectroscopic observing time for FRESCO in GOODS-S is ≈ 16 hours with a typical on-source time of ≈ 2 hours.The 3σ unresolved emission line detection limit around 4.2 µm in the F444W filter is ∼ 1.2 × 10 −18 erg/s/cm 2 , which corresponds to a star formation rate (SFR) detection limit of ∼ 2.1 M ⊙ /yr at z = 5.4 using the conversion factor from Kennicutt & Evans (2012).
A detailed description of the JWST/NIRCam grism data reduction can be found in Sun et al. (2023a).We briefly summarize here the main steps of the reduction process.The data are initially processed with the standard JWST calibration pipeline 2 .We assign WCS to the rate files, perform flat-fielding, and subtract out the sigma-clipped median sky background from each individual exposure after the "ramp-to-slope" fitting in the calibration pipeline.Because we are interested in conducting a targeted emission line search, and we do not expect any of our sources to have a strong continuum due to their general faintness (m = 27 − 28 AB mag), we utilize a median-filtering technique to subtract out any remaining continuum or background on a row-byrow basis following the methodology of Kashino et al. (2023).This produces emission line maps for each grism exposure that are void of any continuum.Although this median-filtering technique is able to properly remove continuum contamination, it sometimes over-subtracts signal in the spectral regions immediately surrounding the brightest emission lines (e.g., emission from [N II] on either side of Hα).The full widths at half maximum (FWHMs) of these emission lines are relatively small and typically on the order of a few pixels, which means that the Hα line flux is preserved by the nine pixel central gap of the median filter and is not affected by the aforementioned over-subtraction.We further remove the "1/f" noise using the tshirt/roeba algorithm3 in both the row and column directions.
We extract two-dimensional (2d) grism spectra using the reduced emission-line maps for all of the objects that are part of the final photometric catalog discussed in Section 2.2.Short wavelength (SW) parallel observations were conducted in two photometric bands (F182M and F210M) and are used for both astrometric and wavelength calibration of the long wavelength (LW) spectroscopic data.We use the spectral tracing and grism dispersion models (Sun et al. 2023a) that were produced using the JWST/NIRCam commissioning data of the Large Magellanic Cloud (LMC; PID: 1076), which are also outlined in Wang et al. (2023).We additionally use the flux calibration models that were produced using JWST/NIRCam Cycle-1 absolute flux calibration observations (PID: 1536/1537/1538).
Using the already extracted 2d spectra, we further extract one-dimensional (1d) grism spectra using a boxcar aperture, assuming a height of five pixels (0.31 ′′ ).We subsequently identify > 3σ peaks automatically in the 1d spectra, assuming various bin sizes (integer units of nm from one to eight) and fit these detected peaks with Gaussian profiles.For each line that is detected with S/N > 3, we tentatively assign a line identification of either Hα or [OIII]λ5008, whichever one minimizes the difference between the best-fit photometric redshift and the tentative spectroscopic redshift.For example, if a line were detected at λ = 4.2 µm and the best-fit photometric redshift is z phot = 5.8, then the initial line identification would be Hα, since the predicted wavelength of this line would be at λ = 4.5 µm, which is closer to the observed wavelength than the predicted wavelength of [OIII]λ5008 (λ = 3.4 µm).Visual inspection is performed on each of these tentative spectroscopic redshift solutions to remove spurious detections caused by either noise or contamination.For sources that pass our visual inspection and have secure line detections, we optimally re-extract the 1d spectra using the F444W surface brightness profile (Horne 1986) and once again fit these detected peaks with Gaussian profiles.According to the grism wavelength calibration uncertainty, the typical absolute uncertainties of our spectroscopic redshifts are ∆z spec = 0.001.
Our final spectroscopic sample includes N = 81 objects at z = 5.2 − 5.5 with > 3σ detections of Hα from the FRESCO spectra.This redshift range was chosen to ensure that Hα would fall in the F410M filter, which is the only medium-band filter for which we have uniform coverage, providing a sanity check for the derived emission line fluxes through a comparison with the F410M excess relative to F444W.Our final spectroscopic sample represents a subset of a larger spectroscopic sample of galaxies from both GOODS fields across a much broader redshift range (Sun et al., in preparation.).For the majority of galaxies in our final spectroscopic sample, neither of the [N II] lines were detected, partially as a result of the aforementioned median-filtering technique (Kashino et al. 2023) but primarily because the line ratio [N II]/Hα is typically low at these redshifts (e.g., Cameron et al. 2023).The NIRCam cutout images alongside the 2d and 1d extracted spectra for these objects are shown in Appendix A. Figure 1 shows the distribution of spectroscopic redshifts for these N = 81 objects while Figure 2 shows the on-sky distribution in angular units.These distributions enabled us to visually identify an overdensity of galaxies around z = 5.4.

ANALYSIS & RESULTS
Using the data and observations from Section 2, we perform various analyses on the N = 81 galaxies in our final spectroscopic sample and present the results.Identification of the extreme galaxy overdensity is described in Section 3.1.Detailed physical modeling of the stellar populations, presentation of the star-forming main sequence, and comparison of inferred SFRs are described in Section 3.2.Determining the dynamic state, estimating the dark matter halo mass, and predicting the future evolution of the overdensity are described in Section 3.3.Placing this overdensity in context with previous works is described in Section 3.4.

Overdensity Identification
Following the technique described in Calvi et al. (2021), we use a Friends-of-Friends (FoF) algorithm to identify the overdensity after looking for threedimensional (3d; two spatial, one spectral) structural groupings (see also Huchra & Geller 1982;Eke et al. 2004;Berlind et al. 2006).This algorithm iteratively selects groups, which consist of one or more galaxies that have projected separations and line-of-sight (LOS) velocity dispersions below the adopted linking parameters (projected separation d link = 500 kpc, chosen to be roughly the virial radius of a typical galaxy cluster; LOS velocity dispersion σ link = 500 km/s, chosen to be roughly the velocity dispersion of such a cluster).These groupings do not depend strongly on the adopted linking parameters, producing similar results when varying either the projected separation or the LOS velocity dispersion by a factor of a few.We identify one large-scale structure consisting of N = 39 galaxies out of the N = 81 galaxies that are part of our final spectroscopic sample.Throughout the rest of this work, we refer to the remaining N = 42 galaxies at z = 5.2 − 5.5 as field galaxies, which consist of (1) isolated galaxies and (2) those in smaller groups as determined by the FoF algorithm4 .The average spectroscopic redshift of the large-scale structure is z = 5.386, spanning a relatively narrow redshift range of 5.374 < z < 5.398.The maximum on-sky separation of the clustered galaxies is roughly 3.6 ′ , corresponding to a physical separation of 8.6 cMpc.
Following the methodology of Chiang et al. (2013), we calculate 3d galaxy overdensities (δ gal = n gal /⟨n gal ⟩−1), where n gal is the number density of galaxies and ⟨n gal ⟩ is the ensemble average number density of galaxies, for each of the N = 81 galaxies that are part of our final spectroscopic sample.To calculate these values, we assume a tophat-weighted spherical window which has a comoving volume equal to (15 cMpc) 3 .The average and standard deviation of the constituent overden- sity values for the large-scale structure identified here is ⟨δ gal ⟩ = 9.2 ± 1.0.This structure is an extreme overdensity at z = 5.4, ∼ 10 times more dense in 3d than the ensemble average at z = 5.2 − 5.5.For comparison, Chiang et al. (2013) found that an average overdensity value of ⟨δ gal ⟩ = 3.04 identifies structures within cosmological simulations as protocluster candidates with 80% confidence.This value is for the z = 5 SFR-limited sample (SFR > 1 M ⊙ /yr), which is the sample that is most similar to our own in terms of selection.Throughout this work, we define a protocluster to be a structure that will eventually collapse into a galaxy cluster at z = 0. Throughout the rest of Section 3, the distinction between the N = 39 "confirmed members of overdensity" and the N = 42 "confirmed members of field" will be used.The confirmed members of the overdensity are shown in Figure 1 (Figure 2) by the turquoise histograms (turquoise pluses) while the confirmed members of the field are shown by the grey histograms (grey points).The median redshift or position of the overdensity is given by the solid magenta line or the magenta point.In both of these figures, the overdensity members ap-pear much more clustered when compared to the field members.The confirmed members of the overdensity are additionally shown in the left panel of Figure 3, color-coded by their spectroscopic redshift.We identify a spatial and kinematic bimodality within the overdensity at z = 5.4, which we return to in Section 3.3.To assign objects between the two components of the bimodality, we adopt an iterative process that minimizes the 3d separation within each of these components, finding N = 14 galaxies that are part of the first component (∆v = −340 +140 −50 km/s) and N = 25 galaxies that are part of the second component (∆v = +190 +120 −150 km/s).The median positions of these two groups are given by the black plus and cross in the left panel of Figure 3. Velocity offsets are calculated relative to the median spectroscopic redshift of the overdensity.

Stellar Population Modeling
Following the methodology outlined in Tacchella et al. (2022), we utilize the SED fitting code Prospector (v1.1.0;Johnson et al. 2021) to infer the stellar populations for the N = 81 objects that are part of our final  (Paxton et al. 2011(Paxton et al. , 2013(Paxton et al. , 2015(Paxton et al. , 2018)).We assume the MILES stellar spectral library (Falcón-Barroso et al. 2011;Vazdekis et al. 2015) and adopt a Chabrier (2003) initial mass function (IMF).Absorption by the intergalactic medium (IGM) is modeled after Madau (1995), where the overall scaling of the IGM attenuation curve is set to be a free parameter.Dust attenuation is modeled using a two-component dust attenuation model (Charlot & Fall 2000) with a flexible attenuation curve where the slope is tied to the strength of the ultraviolet (UV) bump (Kriek & Conroy 2013).Neb-ular emission (both from emission lines and continuum) is self-consistently modeled with the spectral synthesis code Cloudy (Byler et al. 2017).
To test the robustness of the inferred stellar populations and SFHs, we assume three different models for the SFH to see how our results depend on the assumed prior.Two of these models are non-parametric (one with the "continuity" prior, the other with the "bursty continuity" prior) while the third is parametric (with the shape of a delayed-tau function).For each of the non-parametric models, we assume that the SFH can be described by N SFR distinct time bins of constant starformation.The time bins are specified in units of lookback time and the number of distinct bins is fixed at N SFR = 6.The first two bins are fixed at 0−30 Myr and 30−100 Myr while the last bin is fixed between 0.15 t univ and t univ , where t univ is the age of the Universe at the galaxy's spectroscopic redshift, measured with respect to the formation redshift z form = 20.The rest of the bins are spaced equally in logarithmic time between 100 Myr and 0.85 t univ .Both the total stellar mass and the ratios of adjacent time bins are set to be free parameters.A summary of the parameters and priors associated with this Prospector model is presented in Table B1.We adopt the results from the non-parametric SFH with "continuity" prior as fiducial, since the Bayesian evidence does not strongly favor one model over another for the vast majority of the galaxies considered here.
Figure 4 shows the star-forming main sequence for the N = 81 objects at z = 5.2 − 5.5 that are part of the final spectroscopic sample identified in Section 2.3.The confirmed members of the field are given by the grey points while the confirmed members of the overdensity are given by the turquoise points.The reported stellar masses and SFRs are derived from the Prospector fits using the non-parametric SFH with "continuity" prior, where the SFRs are averaged over the last 100 Myr of lookback time.These stellar masses and SFRs are reported in Table C1, which gives a summary of the physical properties for our final spectroscopic sample.Based on these stellar masses and the mass-size relation reported in Shibuya et al. (2015), we find that the minimum separation between the sources in our sample is always larger than twice the effective radius of the sources in their sample.Therefore, each of the objects that is part of our final spectroscopic sample is likely an individual star-forming galaxy rather than individual starforming clumps within a much larger galaxy, although some are clearly merging systems in the final phase of coalescence.To be used as a point of comparison, the empirical star-forming main sequence at z = 5.4 derived by Popesso et al. (2023) is given by the solid black line.Additionally, the maximum allowed SFR assuming all of the stellar mass was formed in the last 100 Myr of lookback time is given by the solid magenta line.
We find that nearly all of the N = 81 objects at z = 5.2 − 5.5 that are part of our final spectroscopic sample agree with the empirically derived star-forming main sequence at z = 5.4 given by Popesso et al. (2023) within 1σ, despite our sample being biased as a result of our requirement to detect the Hα emission line at greater than 3σ.If we were to instead use Hα-based SFRs derived with the conversion factor from Kennicutt & Evans (2012), the objects within our final spectroscopic sample would be shifted upward by ∆ ∼ 0.5 dex.However, this is likely because the canonical hydrogen ionizing photon production efficiency (ξ ion ∝ L Hα /L UV ) used in Kennicutt & Evans ( 2012) is only ξ ion ∼ 10 25.1 erg −1 Hz, lower than measurements at z = 5 − 6 by ∆ξ ion ∼ 0.5 dex (e.g., Bouwens et al. 2016;Ning et al. 2023;Sun et al. 2023a).Furthermore, if we were to instead use UV-based SFRs derived with the conversion factor from Kennicutt & Evans (2012) alongside measurements of the rest-UV magnitudes derived as νL ν at λ rest = 1500 Å, it would not change the distribution of our final spectroscopic sample in Figure 4.This is ex-pected, since the Prospector-based SFRs largely rely on the available rest-UV photometry.
The most notable examples of galaxies disagreeing with the empirical relation include one well above the main sequence at high stellar mass (JADES−GS+53.13859−27.79025)and two somewhat below the main sequence at intermediate stellar mass (JADES−GS+53.06799−27.80816and JADES−GS+53.18328−27.77894)despite large uncertainties in the inferred SFRs.We note that the galaxy well above the main sequence at high stellar mass is an active galactic nucleus (AGN), originally identified as a broad line Hα emitter in Matthee et al. (2023) with a line width of 2200 ± 500 km/s.Combined with our relatively low SFR detection limit (∼ 2.1 M ⊙ /yr, see Section 2.3), the fact that our sample agrees with the empirical relation suggests that we are sampling the bulk of the star-forming population at these redshifts despite our selection criteria.However, we should mention that we are likely missing some amount of dusty star-forming galaxies (DSFGs), which historically have been used as tracers to identify protocluster candidates at these redshifts (for a review of environmental galaxy evolution, see Alberts & Noble 2022).
To test the impact of the assumed SFH, we compare the Prospector derived stellar masses and SFRs for the three different SFH models.We remind the reader that for Figure 4 and Table C1, the assumed SFH is non-parametric with the "continuity" prior.Compared to the non-parametric SFH with the "bursty continuity" prior, the derived stellar masses are a bit larger (∆ ≈ 0.1 dex) with large scatter (σ ≈ 0.7 dex) while the derived SFRs are broadly consistent (∆ ≈ 0.0 dex) with large scatter (σ ≈ 0.8 dex).Compared to the parametric SFH, both the derived stellar masses and SFRs are broadly consistent (∆ ≈ 0.0 dex) with large scatter (σ ≈ 0.8 dex and σ ≈ 1.0 dex, respectively).We find that the inferred Prospector parameters considered here do not depend strongly on the assumed SFH.
Figure 5 shows the comparison of SFRs in the two most recent time bins (corresponding to lookback times of 0−30 Myr and 30−100 Myr, respectively) for the N = 81 objects at z = 5.2 − 5.5.Once again, the confirmed members of the field are given by the grey points while the confirmed members of the overdensity are given by the turquoise points.A constant SFH is given by the solid black line, which assumes that the average SFRs in the two most recent time bins are equal.Values above this line represent a falling SFH over the last 100 Myr, while values below this line represent a rising SFH.The twin axes to the top and right are similar to those shown in Figure 4, where the median values for the confirmed  2023) is given by the solid black line.The maximum allowed SFR assuming all of the stellar mass was formed in the last 100 Myr of lookback time is given by the solid magenta line.
The median values of the inferred stellar masses and SFRs for the confirmed members of the field (overdensity) are given by the grey (turquoise) dashed lines in the twin axes to the top and right.The 1σ confidence interval for these parameters are given by the shaded regions in the same axes.The overdensity members appear to have larger inferred stellar masses and SFRs when compared to the field members.Additionally, nearly all of these objects agree with the empirically derived star-forming main sequence at z = 5.4 within 1σ.
members of the field (overdensity) are given by the grey (turquoise) dashed lines.By comparing members of the overdensity with members of the field in Figures 4 and 5, we can begin to explore the impact of environment on galaxy formation and evolution at z = 5.2 − 5.5.Table 1 gives a summary of percentiles for some of the Prospector inferred physical parameters for both the members of the field and members of the overdensity.In the star-forming main sequence and the summary of percentiles, overdensity members appear to have larger inferred stellar masses and SFRs when compared to the field members.We further compare the distributions of these inferred parameters by performing Kolmogorov-Smirnov (KS) and Anderson-Darling (AD) tests, which are two-sided tests for the null hypothesis that two independent samples are drawn from the same continuous distribution.The results of these tests indicate that the stellar masses and SFRs of the field and overdensity samples exhibit statistically significant differences at roughly the 2 − 3σ level (corresponding to p-values of 0.003 < p < 0.05).However, to fairly compare these two samples, the stellar mass distributions must be similar such that any observed differences can be attributed to external processes (e.g., environment) rather than internal processes (e.g., feedback).To create a mass-matched sample of field galaxies, we select the field galaxy that is closest in stellar mass to each overdensity galaxy, which produces a sample of galaxies with stellar masses that are consistent with being drawn from the same parent distribution as the overdensity members based on the KS-and AD- tests.The percentiles for the mass-matched members of the field also appear in Table 1.
When compared to the mass-matched field members, the overdensity members have similar inferred SFRs, consistent with being drawn from the same parent distribution based on the KS-and AD-tests.The ratio of SFRs in the two most recent time bins (corresponding to lookback times of 0 − 30 Myr and 30 − 100 Myr, respectively) also show an interesting trend, where overdensity members have SFHs consistent with constant while mass-matched field members have SFHs consistent with rising.Furthermore, KS-and AD-tests suggest that the SFHs of the mass-matched field and overdensity sample exhibit statistically significant differences at roughly the 2 − 3σ level (corresponding to p-values of 0.003 < p < 0.05).These results suggest that the physical processes associated with this extreme galaxy overdensity at z = 5.4 have induced earlier star formation and earlier stellar mass assembly relative to the field, although there are large uncertainties associated with all of these parameters, and the distributions of these parameters appear consistent within the associated uncertainties.

Dark Matter Halo Mass Estimates
To further understand the dynamical state of the overdensity identified in Section 3.1, the right panel of Figure 3 shows the phase-space diagram for the N = 39 confirmed members of the overdensity.Velocity offsets are calculated from the median redshift of the overdensity (see Figure 1) while spatial offsets are calculated from the median on-sky positions of the two groups that make up the bimodality (see the left panel of Figure 3).The median velocity offset of the overdensity is given by the grey dashed line while the median velocity offsets of the first and second groups are given by the green and blue dotted lines, respectively.For each of the groups, the 1σ confidence interval of the velocity offsets are given by the shaded regions (for group one, ∆v = −340 +140 −50 km/s; for group two, ∆v = +190 +120 −150 km/s).Following the methodology described in Long et al. (2020), we derive two different estimates of the dark matter halo mass for the overdensity at z = 5.4.For both of these methods, we use the stellar-to-halo abundance matching relation described in Behroozi et al. (2013) to convert stellar masses into halo masses.Uncertainties on our estimates are calculated by adopting the mean stellar-to-halo abundance matching relation from Behroozi et al. (2013) and propagating the uncertainties on the stellar masses for each individual galaxy.This is an underestimate of the true uncertainties since there is scatter in the stellar-to-halo abundance matching relation which we do not account for.
The first halo mass estimate is calculated by summing the halo masses for each individual galaxy that is part of the overdensity.This estimate assumes that (1) each galaxy is formed in a separate dark matter halo and (2) the individual galaxies are closer to virialization than the large-scale structure associated with the overdensity.We find that the first method gives a total halo mass of log 10 (M halo /M ⊙ ) = 12.8 ± 0.1.The second halo mass estimate is determined by summing the stellar masses for each individual galaxy that is part of the overdensity and converting to a halo mass.This estimate instead assumes that each galaxy is formed in the same dark matter halo.We find that the second method gives a total halo mass of log 10 (M halo /M ⊙ ) = 12.6 ± 0.3.
Weighing this kind of large-scale structure in the early Universe is challenging and requires a variety of assumptions.Typical methods for weighing galaxy clusters include (1) gravitational lensing, (2) the Sunyaev-Zel'dovich effect, and (3) using X-ray observations of the hot gas in the intracluster medium (ICM).However, current X-ray and sub-millimeter observations are not sensitive enough to properly weigh this extreme overdensity at z = 5.4.For this reason, we use the above dark matter halo mass estimates to infer a probable halo mass range of 12.6 ≲ log 10 (M halo /M ⊙ ) ≲ 12.8 for the overdensity identified in Section 3.1.
We note that none of these methods account for additional members of the overdensity that were not identified and included in the final spectroscopic sample, including objects that fall outside either the JADES or the FRESCO footprints.This is a non-negligible effect since the first component of the overdensity (see Figure 3) falls right at the edge of the JADES footprint (θ ≪ 1 ′ ).Additionally, since our final spectroscopic sample only includes galaxies with narrow photometric redshift constraints and secure Hα line detections, we are likely missing some additional subset of objects with relatively unconstrained photometric redshifts and/or low levels of star formation (e.g., DSFGs and/or obscured AGNs).There are zero known DSFGs at z = 5.4 in GOODS-S in the literature (e.g., Franco et al. 2018;Hatsukade et al. 2018;González-López et al. 2020;Gómez-Guijarro et al. 2022).However, wide and shallow surveys like GOODS-ALMA would miss DSFGs with infrared luminosities of L IR ≲ 3 × 10 11 L ⊙ at z = 5.4.Given that clusters induce earlier quenching for their constituent members, we cannot rule out the existence of a significant population of these kinds of objects.For these reasons, the halo mass range quoted above is likely an underestimate of the true halo mass associated with this extreme galaxy overdensity at z = 5.4.
Figure 6 shows the dark matter halo mass distribution for some representative groups and protoclusters at z < 8.The two halo mass estimates derived here in Section 3 are given by the turquoise triangles, where the triangles indicate that these are likely underestimates of the true halo mass.Groups and protoclusters at z < 6 from Li et al. (2022) are given in grayscale while protocluster candidates in the COSMOS field at z > 6 from Brinch et al. (2023) are given by the grey points for comparison, both selected based on photometric redshifts, with dark matter halo mass estimates that assume the same stellar-to-halo abundance matching relation used here.The magenta shaded region shows the expected halo mass evolution of a Coma-like cluster (Chiang et al. 2013) assuming a smooth evolution at z > 6.The black dashed line represents the typical threshold for shock stability assuming a spherical infall, below which the flows are predominantly cold and above which a shockheated ICM is expected (Dekel & Birnboim 2006).The black diagonal dashed line represents the typical threshold for penetrating cold gas flows.The overdensity identified in Section 3.1 is expected to eventually evolve Figure 6: The dark matter halo mass distribution for groups and protoclusters at z < 8.The two halo mass estimates derived in Section 3.3 are given by the turquoise triangles.For comparison, groups and protoclusters at z < 6 from Li et al. (2022) are given in grayscale while protoclusters at z > 6 from Brinch et al. (2023) are given by the grey points.The magenta shaded region shows the expected halo mass evolution of a Coma-like cluster (Chiang et al. 2013).The black horizontal dashed line represents the typical threshold for shock stability assuming a spherical infall, below which the flows are predominantly cold and above which a shock-heated ICM is expected (Dekel & Birnboim 2006).The black diagonal dashed line represents the typical threshold for penetrating cold gas flows.The overdensity at z = 5.4 is expected to evolve into a Coma-like cluster with log 10 (M halo /M ⊙ ) > 15 by z = 0.
into a Coma-like cluster with log 10 (M halo /M ⊙ ) > 15 by z = 0, rivaling some of the most massive galaxy clusters found in the local Universe (Ruel et al. 2014;Buddendiek et al. 2015).For a Coma-like cluster at z = 0, the effective radius at z = 5.4 would be R e ≈ 10 cMpc (Chiang et al. 2013).This value is much larger than the physical size of this extreme galaxy overdensity (see Figure 3).Thus, we conclude that the classification of this large-scale structure as a protocluster is justified.

Comparison with Previous Works
A number of previous works have found overdensities at similar redshifts to the one reported here.The earliest one identified at z > 5 is from Ouchi et al. (2005), who used narrow-band imaging with Subaru-Suprime Cam to identify galaxies with strong Lyα emission lines at z = 5.7 (∆z = 0.05).The distribution of sources was described as "clumpy" with one prominent overdensity significant at the 4.8σ level, and a second one at the 2.2σ level.Follow-up spectroscopy confirmed the narrow redshift range for these clumps (∆z = 0.03), each being about 2 ′ in diameter with a separation of about 9 ′ .The narrow redshift range combined with the small on-sky separations suggests that these two overdensities may be clumps within a large-scale structure analogous to the one reported here.
Similar structures have been reported more recently (e.g., Jiang et al. 2018;Chanchaiworawit et al. 2019;Harikane et al. 2019;Wang et al. 2021), with similar techniques in the first three cases, and based on overdensities in the sub-millimeter found with the South Pole Telescope in the fourth.Jiang et al. (2018) found two overdensities near that which was identified by Ouchi et al. (2005), one at z = 5.68 with a diameter of ∼ 15 ′ and a smaller one at z = 5.75.These both are presumably related to the overdensities reported by Ouchi et al. (2005), given their proximity in redshift and onsky.Chanchaiworawit et al. (2019) and Harikane et al. (2019) identified two similar structures at z = 6.5.Since the search method for these works starts with narrow-band imaging, the samples are from a narrow redshift range and the actual density of such overdensities on the sky must be significantly higher than implied by the existing detections.Thus, massive clusters of galaxies must be well on the way to formation by z = 6, which is something that we expect based on simulations (Chiang et al. 2013(Chiang et al. , 2017) ) and more recent observations (Laporte et al. 2022;Morishita et al. 2023).
One of the key breakthroughs illustrated here with respect to these previous works is the efficient spectroscopic confirmation of such a large sample of galaxies at z = 5.2 − 5.5, made possible by the powerful combination of deep imaging and WFSS provided by JWST/NIRCam.As also presented in Kashino et al. 2023, JWST/NIRCam WFSS enabled the discovery of three galaxy overdensities along the sight line of quasar J0100 + 2802 at z = 6.19, z quasar = 6.33, and z = 6.78 through the blind detections of [O III]emitting galaxies.It is also worth mentioning that a similar z = 5.2 galaxy overdensity was identified in the GOODS-N field with a probable halo mass range of 12.3 < log 10 (M halo /M ⊙ ) < 12.9 through Lyα and submillimeter spectroscopy (Walter et al. 2012;Calvi et al. 2021), which was also partially observed by FRESCO in February 2023 (see recent papers by Herard-Demanche et al. 2023;Sun et al. 2023b).Future observations with JWST will certainly find more structures similar to those highlighted here, allowing a more complete look at the progenitors of the most massive gravitationally bound structures in the local Universe: galaxy clusters.

SUMMARY & CONCLUSIONS
We have presented the discovery of an extreme galaxy overdensity at z = 5.4 in the GOODS-S field using data from the Near Infrared Camera (NIRCam) on JWST.These data consist of JWST/NIRCam imaging from the JWST Advanced Deep Extragalactic Survey (JADES) and JWST/NIRCam wide field slitless spectroscopy from the First Reionization Epoch Spectroscopic COmplete Survey (FRESCO).Our findings can be summarized as follows.
1. Galaxies were initially selected using HST+JWST photometry spanning λ = 0.4 − 5.0 µm.These data provide well-constrained photometric redshifts down to m ≈ 29 − 30 mag, particularly at z = 5.2 − 5.5, where Hα excess can be traced by comparing photometry in the F410M and F444W filters.Galaxies were subsequently selected using slitless spectroscopy over λ = 3.9 − 5.0 µm via a targeted emission line search for Hα around the best-fit photometric redshift.The final spectro-scopic sample of galaxies includes N = 81 objects at z = 5.2 − 5.5.
2. A Friends-of-Friends (FoF) algorithm was used to identify this extreme galaxy overdensity by iteratively looking for three-dimensional (3d) structural groupings within the final spectroscopic sample.One large-scale structure consisting of N = 39 galaxies was discovered, which is ∼ 10 times more dense in one-dimension (1d) and ∼ 12 times more dense in 3d than the N = 42 analogous field galaxies at z = 5.2 − 5.5.
3. The stellar populations for these N = 81 objects at z = 5.2 − 5.5 were inferred using the HST+JWST photometry spanning λ = 0.4 − 5.0 µm, the spectroscopic redshifts determined by the targeted line search, and the SED fitting code Prospector (Johnson et al. 2021).We constructed the starforming main sequence at z = 5.2 − 5.5 and found that nearly all the galaxies in our sample agree with the empirically derived star-forming main sequence at z = 5.4 derived by Popesso et al. (2023).Combined with our relatively low SFR detection limit, this suggests that we are sampling the bulk of the star-forming population at these redshifts despite our Hα selection criteria.By comparing members of the overdensity with a mass-matched sample of members of the field, we find evidence suggesting that environment has induced earlier star formation and earlier stellar mass assembly within the overdensity relative to the field.
4. Using two different methods, we estimated the total dark matter halo mass associated with this extreme galaxy overdensity at z = 5.4 to be within 12.6 ≲ log 10 (M halo /M ⊙ ) ≲ 12.8.As a result of our selection criteria, we are potentially missing objects that fall outside either the JADES or the FRESCO footprints, as well as some subset of objects with relatively unconstrained photometric redshifts and/or low levels of star formation.This means the total dark matter halo mass range quoted above is likely an underestimate of the true halo mass.This massive large-scale structure is expected to evolve into a Coma-like cluster with log 10 (M halo /M ⊙ ) > 15 by z = 0.
In this work, we have demonstrated the powerful combination of JWST/NIRCam imaging and slitless spectroscopy by efficiently confirming the redshifts for N = 81 galaxies at z = 5.2 − 5.5, inferring the physical properties of these galaxies, and assessing the largescale structure in which these galaxies reside.Follow-up spectroscopic observations using JWST and/or the Atacama Large Millimeter/sub-millimeter Array (ALMA) will (1) inform us about the chemical compositions of these (and similar) galaxies, (2) provide insight into the formation and evolution of extreme galaxy overdensities in the early Universe, and (3) constrain the total number of these kinds of large-scale structures immediately after the epoch of reionization (EoR; z > 6) when the Universe was less than a billion years old.

B. ASSUMED PHYSICAL MODEL FOR STELLAR POPULATION MODELING
Table B1 gives a summary of the parameter and priors assumed in the Prospector model that was used for the stellar population modeling described in Section 3.2.

C. PHYSICAL PROPERTIES OF FINAL SPECTROSCOPIC SAMPLE
Table C1 gives a summary of the physical properties for the 81 objects that are part of the final spectroscopic sample described in Section 2.3.

Figure 1 :
Figure1: The distribution of spectroscopic redshifts for the N = 81 objects at z = 5.2 − 5.5 that are part of the final spectroscopic sample identified in Section 2.3.As defined in Section 3.1, the grey histograms represent the N = 42 confirmed members of the field while the turquoise histograms represent the N = 39 confirmed members of the overdensity.The median redshift of the overdensity is given by the solid magenta line.Compared to the field members, the overdensity members appear much more clustered, representing a ∼ 10 times overdensity at z = 5.4.

Figure 2 :
Figure2: The on-sky distribution in angular units for the N = 81 objects at z = 5.2 − 5.5 that are part of the final spectroscopic sample.The grey points represent the N = 42 confirmed members of the field while the turquoise pluses represent the N = 39 confirmed members of the overdensity.The median position of the overdensity is given by the magenta point.The JADES deep (medium) footprint is illustrated by the dashed (solid) red line, the JEMS footprint by the solid orange line, and the FRESCO footprint by the solid yellow line.It is apparent that the overdensity falls near the edge of the JADES medium and FRESCO footprints, which means we cannot rule out the overdensity extending well beyond the region for which we currently have data.

Figure 3 :
Figure 3: Left panel: The on-sky distribution in physical units for the N = 39 confirmed members of the overdensity, color-coded by their spectroscopic redshift.The color-coded crosses (pluses) represent members of the first (second) component of the overdensity.The median position of the first (second) component is given by the black cross (plus).Right panel: The observed phase-space diagram for the N = 39 confirmed members of the overdensity.The median velocity offset of the first (second) component of the overdensity is given by the green (blue) dashed line.The 1σ confidence interval of the velocity offsets for each of the groups are given by the shaded regions.spectroscopicsample.Fits are performed on the rescaled Kron small photometry described in Section 2.1, while the redshift is fixed at the spectroscopic redshift determined in Section 2.2.Prospector uses a Bayesian inference framework and we choose to sample posterior distributions with the dynamic nested sampling code dynesty (v1.2.3; Speagle 2020).We use the Flexible Stellar Population Synthesis code (FSPS; Conroy et al. 2009; Conroy & Gunn 2010) via python-FSPS (Foreman-Mackey et al. 2014) with the Modules for Experiments in Stellar Astrophysics Isochrones and Stellar Tracks (MIST; Choi et al. 2016; Dotter 2016), which makes use of the Modules for Experiments in Stellar Astrophysics (MESA) stellar evolution package(Paxton et al. 2011(Paxton et al. , 2013(Paxton et al. , 2015(Paxton et al. , 2018)).We assume the MILES stellar spectral library(Falcón- Barroso et al. 2011;Vazdekis et al. 2015) and adopt aChabrier (2003) initial mass function (IMF).Absorption by the intergalactic medium (IGM) is modeled afterMadau (1995), where the overall scaling of the IGM attenuation curve is set to be a free parameter.Dust attenuation is modeled using a two-component dust attenuation model(Charlot & Fall 2000) with a flexible attenuation curve where the slope is tied to the strength of the ultraviolet (UV) bump(Kriek & Conroy 2013).Neb-

Figure 4 :
Figure4: The star-forming main sequence for the N = 81 objects at z = 5.2 − 5.5 that are part of the final spectroscopic sample.The stellar masses and SFRs reported here are derived from the Prospector fits.The SFRs are averaged over the last 100 Myr of lookback time.The grey points represent the N = 42 confirmed members of the field while the turquoise points represent the N = 39 confirmed members of the overdensity.The empirical star-forming main sequence at z = 5.4 derived byPopesso et al. (2023) is given by the solid black line.The maximum allowed SFR assuming all of the stellar mass was formed in the last 100 Myr of lookback time is given by the solid magenta line.The median values of the inferred stellar masses and SFRs for the confirmed members of the field (overdensity) are given by the grey (turquoise) dashed lines in the twin axes to the top and right.The 1σ confidence interval for these parameters are given by the shaded regions in the same axes.The overdensity members appear to have larger inferred stellar masses and SFRs when compared to the field members.Additionally, nearly all of these objects agree with the empirically derived star-forming main sequence at z = 5.4 within 1σ.

Figure 5 :
Figure 5: The comparison of SFRs in the two most recent time bins for the N = 81 objects at z = 5.2 − 5.5 that are part of the final spectroscopic sample.The stellar masses and SFRs reported here are derived from the Prospector fits.The grey points represent the N = 42 confirmed members of the field while the turquoise points represent the N = 39 confirmed members of the overdensity.A constant SFH is given by the solid black line.Values above (below) this line represent a falling (rising) SFH over the last 100 Myr of lookback time.The twin axes to the top and right are similar to those shown in Figure 4.The overdensity (field) members appear to have non-parametric SFHs consistent with a constant (rising) SFH with relatively small (large) scatter.

Figure A1 :
Figure A1: The NIRCam cutout images alongside the continuum subtracted 2d and 1d grism spectra of JADES-GS+53.08831-27.84042at z = 5.393.The upper-left panel shows the 1.2 ′′ × 1.2 ′′ F444W-F277W-F150W RGB thumbnail.The upper-right panel shows the extracted 2d spectrum around the Hα emission line detection, indicated by the solid red line.The lower-right panel instead shows the extracted 1d spectrum around the Hα emission line detection alongside the best-fit Gaussian profile given by the solid red line.The JADES ID and confirmed spectroscopic redshift are given in the lower-right panel for each galaxy.The complete figure set is available in the online journal.

Figure A2 :
Figure A2: The NIRCam cutout images alongside the continuum subtracted 2d and 1d grism spectra of JADES-GS+53.15686-27.86069at z = 5.327, just like in Figure A1.The complete figure set is available in the online journal.

Table 1 :
A summary of percentiles for some of the Prospector inferred physical parameters.

Table C1 :
A summary of the physical properties for the 81 objects in our final spectroscopic sample.