Abstract
During the epoch of reionization, the first galaxies were enshrouded in pristine neutral gas, with one of the brightest emission lines in star-forming galaxies, Lyman α (Lyα), expected to remain undetected until the Universe became ionized. Providing an explanation for the surprising detection of Lyα in these early galaxies is a major challenge for extragalactic studies. Recent James Webb Space Telescope observations have reignited the debate about whether residence in an overdensity of galaxies is a sufficient and necessary condition for Lyα to escape. Here, we take unique advantage of both high-resolution and high-sensitivity images from the James Webb Space Telescope Near Infrared Camera to show that all galaxies in a sample of Lyα emitters with redshift >7 have close companions. We exploit on-the-fly radiative-transfer magnetohydrodynamical simulations with cosmic ray feedback to show that galaxies with frequent mergers have very bursty star formation histories that drives episodes of high intrinsic Lyα emission and facilitates the escape of Lyα photons along channels cleared of neutral gas. We conclude that the rapid buildup of stellar mass through mergers presents a compelling solution to the long-standing puzzle of the detection of Lyα emission deep in the epoch of reionization.
Similar content being viewed by others
Main
Young, vigorously star-forming galaxies have been identified in the very early Universe1,2,3. These galaxies should be excellent sources of Lyman-α emission (Lyα; wavelength (λ) = 1,215.67 Å)—the intrinsically brightest emission line4—which stems from the recombination of hydrogen that has been ionized by their young stellar populations. However, deep in the epoch of reionization, galaxies are expected to be exceptionally gas-rich, such that their stellar nurseries are enshrouded in copious amounts of neutral hydrogen, which leads to extreme damped absorption of Lyα5. Furthermore, the intergalactic medium (IGM) is increasingly neutral as we probe to higher redshift6,7, and this neutral gas is expected to resonantly scatter Lyα emission. Hence, as a result of ‘local’ attenuation by a gas-rich interstellar medium (ISM) and scattering by a neutral IGM, Lyα emission should be detectable only towards the end of the reionization era, about one billion years after the Big Bang4,8,9,10.
Although the decreasing observability of Lyα emission with increasing redshift has been repeatedly observed11,12,13, this picture has been challenged by the occasional, surprising detection of Lyα emission in several galaxies deep in the reionization era14,15. It has been suggested that Lyα can escape through the neutral IGM if the galaxies reside in sufficiently large ionized bubbles embedded in the neutral IGM, driven either by active galactic nuclei (AGN)16,17,18 or by an enhanced radiation field produced by an overdensity of associated objects19,20,21,22,23,24,25. However, a solution to the escape of Lyα emission through the ISM and circumgalactic medium of what are expected to be very gas-rich galaxies remains elusive. Before the advent of the James Webb Space Telescope (JWST), the sensitivity and resolution of imaging instruments meant that studies of Lyα emitters (LAEs) at high redshift were spatially unresolved. Hence, it was not possible to probe the physical processes that could explain the escape of Lyα emission from the ISM.
To address this crucial issue, we study a sample of nine galaxies that have been spectroscopically confirmed with the detection of Lyα emission at redshift (z) > 7 with high-resolution spectrographs and that have been observed with publicly available JWST Near Infrared Camera (NIRCam)26 imaging (the properties of which are reported in Table 1). These galaxies fall in the GOODS-North, GOODS-South, EGS and COSMOS fields and have been observed as part of five different programmes: Public Release IMaging for Extragalactic Research (PRIMER, primary investigator (PI) Dunlop), the First Reionization Epoch Spectroscopic Complete Survey (FRESCO, PI Oesch)27, the Cosmic Evolution Early Release Science survey (CEERS, PI Finkelstein)28, the JWST Advanced Deep Extragalactic Survey (JADES, PI Eisenstein)29 and Director’s Discretionary Time (DDT) programme 4426 (providing NIRSpec integral field unit (IFU) observations of GN-z11, PI Maiolino)30. Six of these nine Lyα-emitting galaxies are known to lie within overdensities21,23,24,25,31, three of which also probably host accreting black holes16,17,18. Moreover, it has been shown that the remaining galaxies are incapable of alone blowing a large enough ionized bubble to facilitate the escape of Lyα through the neutral IGM, indicating that an underlying population of faint galaxies must surround these LAEs25.
However, the presence of spectroscopically confirmed galaxies within known ionized bubbles that do not show Lyα emission13,23,25 indicates that there must be further local processes at play driving Lyα emission deep into the epoch of reionization. With this in mind, we use the unparalleled sensitivity and resolution of JWST/NIRCam to accurately determine the properties of these nine LAEs with equivalent width greater than 10 Å, at z > 7.
Remarkably, NIRCam images of these galaxies show the presence of several components of all these LAEs, as shown in Fig. 1. Given the tight redshift constraints provided by strong [O iii]4959,5007 emission falling into medium-band filters, we estimate the probability of these being unrelated, high-redshift objects and find this is always less than 14% for each of our systems and is often just 1% (see Methods for further details). We estimate the redshift of each candidate companion by carefully extracting the spectral energy distribution (SED) of both the main component and its companions in Extended Data Figs. 1 and 2. We then fit each SED by assuming several parametric star formation histories (SFHs: burst, constant, delayed, exponential), allowing a large range for all parameters32. In all cases, we deduce that the photometric redshifts of the companions are consistent with that of their spectroscopically confirmed main component.
To spectroscopically confirm their identity as companions, we exploit existing X-shooter/Very Large Telescope33, MOSFIRE/Keck34,35, NIRSpec on-board JWST36 and NIRCam wide-field slitless spectroscopy (WFSS) observations37 of these galaxies. For three systems in our sample, we have no spectroscopic information on their companions given the limitations of seeing at ground-based telescopes and the small aperture size of JWST/NIRSpec observations. For the remaining six systems in our sample, we are able to spectroscopically confirm the companion galaxy as it falls within the field of view (FoV) of the observations of the main LAE. For two companions in our sample, we see tentative (∼3σ) evidence of Lyα emission in Extended Data Fig. 3 (one of which has recently been confirmed by NIRSpec IFU observations38). For a further two companions, we observe the [O iii]5007 emission line at >8σ (as seen in Extended Data Fig. 4, with fluxes reported in Table 2), and NIRSpec IFU observations spectroscopically confirm the companions of the final two systems (refs. 30,39 and Carniani et al., manuscript in preparation). Therefore, for all systems for which we have resolved spectroscopic observations of the companion (∼67%), we spectroscopically confirm its nature as a physical companion (see Methods for further details of the spectroscopic analysis).
To confirm that the presence of a close companion is the primary factor governing the visibility of Lyα, we determine the fraction of companions seen in a mass-matched sample of z > 7 galaxies with high-resolution spectroscopic data for which Lyα is not detected13. We find that 43% of these galaxies have photometric-candidate companions within 5″ of the central galaxy, which is consistent with the companion fraction determined by a more comprehensive study (Puskás et al., manuscript in preparation). The lower fraction of close companions among samples that are not selected for Lyα emission is evidence that the 100% rate of companions for our sample of LAEs is atypical of z > 7 galaxies.
To reinforce the observational evidence supporting the idea that continuing interactions drive an increase in the detectability of Lyα emission during the epoch of reionization, we explore comparable galaxy mergers using the Azahar simulation suite (Martin-Alvarez et al., manuscript in preparation). This simulation suite was performed before we obtained the NIRCam data, but remarkably, we find many simulated galaxies that match the photometry and spatial geometry of our objects. Azahar is a cosmological, high-resolution, zoom-in simulation that uses a magnetohydrodynamical solver together with state-of-the-art cosmic ray feedback (see the pathfinder Pandora project40 and Methods). Most importantly for this work, Azahar also features on-the-fly radiative transfer41,42,43 capable of self-consistently reproducing reionization while fully modelling the ISM of galaxies (with a maximum spatial resolution of 10 pc) and, crucially, resolving the propagation and escape of ionizing radiation on the ISM scales44,45. We post-process our simulations with the publicly available RAdiative SCattering in Astrophysical Simulations (RASCAS) code46 to account for both the production and resonant scattering of Lyα photons as well as scattering or absorption of Lyα photons by dust47.
As a result of a substantial overdensity of galaxies in our simulated volume, the main progenitor undergoes repeated mergers with other galaxies brought in by the cosmic filaments, often involving several companions or mergers in rapid succession as has been observed in ref. 48. To highlight one such merger occurring at z ≈ 7.3, Fig. 2b–e show different projections of Azahar. This particular interaction features three merging galaxies that will constitute the main progenitor of the spiral galaxy formed by z ≈ 1. These are very good analogues to the EGSY8p68 observations shown on the left (modelling of the observations with GALFIT in Extended Data Fig. 5 confirms the presence of three components), with very similar individual galaxy sizes, mutual distances and merger configuration. Specifically, the three simulated galaxies have stellar masses of M*,1 = 2 × 108 M⊙, M*,2 = 3 × 108 M⊙ and M*,3 = 8 × 107 M⊙ at this redshift, which are comparable to the stellar masses of the EGSY8p68 galaxies (see Table 1). Furthermore, the simulated star formation rates (SFRs) vary between 1 and 10 M⊙ yr−1 (shown in Fig. 3) and are consistent with the SFRs of the EGSY8p68 galaxies (7.62 M⊙ yr−1, 2.22 M⊙ yr−1, 4.26 M⊙ yr−1). We also consider the total Lyα luminosity of the simulated system and find that although it varies, its average value is ∼1043 erg s−1 (shown in Fig. 3), in agreement with the value for EGSY8p68 reported in Table 1. Finally, we note that simulated volume hosts a substantial overdensity of galaxies, which is consistent with the local environments of the observed LAEs in our sample (as discussed earlier), including EGSY8p6821. We thus conclude that the simulated system matches all key observed photometric and spectroscopic properties of EGSY8p68 (see also Methods).
The observations shown in Fig. 2a clearly indicate the superior resolution and sensitivity of JWST (top) compared to the Hubble Space Telescope (HST, bottom). The tricomponent nature of EGSY8p68 is entirely unresolved in existing HST observations but clearly identifiable in the F150W NIRCam imaging. Fig. 2b shows the large-scale view (150 kpc across) of the three filaments encompassing the merging system, where large amounts of H i gas (blue) along the filaments are feeding the star formation in these systems, resulting in substantial NIRCam F150W emission (yellow-white). This active star formation is driving ionized hydrogen bubbles (grey) away from the filaments and into the low-density regions. Fig. 2c shows the simulated JWST-like 150W observation (with dust extinction modelled as an absorption screen along the line of sight (LOS)), where we select the LOS to illustrate the resemblance with EGSY8p68. Note that the simulated merging system is displayed at a slightly earlier stage of the merger with respect to EGSY8p68, with our simulated galaxies approaching the probable physical separation of the observed galaxies at z ≈ 7. Interestingly, the mock F150W emission shows compact galactic cores analogous to these JWST observations, with diffuse emission (blending with the background) emerging from the extended, low-density stellar disks. Fig. 2d shows a synthetic colour-composite observation using the JWST filters at the full resolution of our simulation. Although the two companion galaxies are more diffuse, with stars actively forming in their disks, the main galaxy has a much more compact core due to the preceding rapid growth through repeated mergers. Fig. 2e is a colour composite of various simulated properties relevant for Lyα detectability. The emission from LyC ionizing photons (purple) and He ii ionizing photons (orange) is distributed on a scale of a few (physical) kiloparsecs. Importantly, there is a clear separation of the stellar emission and the H i gas during the merger event, with the ionizing radiation escaping from star-forming regions.
To provide a quantitative confirmation and in-depth physical understanding of this effect, we explore in Fig. 3 the time evolution of the three merging galaxies in detail. In Fig. 3a, the colours encode the same quantities as in Fig. 2e, and Fig. 3b shows the intrinsic Lyα emission (blue) and resonantly scattered Lyα emission (orange) produced by post-processing our simulations with RASCAS. The first column shows a larger-scale view of the main progenitor (galaxy 2) and its approaching companions at z ≈ 8.1. The second column shows a close-up view of the main progenitor at z ≈ 8.1 undergoing a minor merger, with a complex topology of neutral gas and channels through which Lyα photons are escaping. The third column shows the continuing merger of our galaxies 1, 2 and 3 at z ≈ 7.3 (with a different orientation with respect to Fig. 2), where bursty star formation drives ionized channels through the ISM, with a substantial amount of gas also tidally stripped from the galaxies. The fourth column shows a post-merger view of the newly formed disk-dominated galaxy at z ≈ 6.1. Here, the H i is confined to the merger remnant galaxy, and the radiation can only escape through the channels opened by collimated, cosmic-ray-driven biconical outflows. Focusing now on the Lyα radiation, all three galaxies have substantial intrinsic Lyα emission due to their continuing star formation. However, at this particular viewing angle, because of absorption by dust, Lyα scattered emission is only present around merging galaxy 2 at z ≈ 8.1 and around galaxies 1 and 3 at z ≈ 7.3. At z ≈ 7.3, the remaining emission, not absorbed by dust, is very diffuse on large scales and hence not observable.
To understand what drives the observable Lyα emission, Fig. 3c shows the SFR evolution for our three merging galaxies, together with the SFR history of the merger remnant at z ≈ 6.1. SFRs of all three galaxies undergo repeated cycles of bursts, driven by numerous mergers (‘major’ mergers (filled dots) and ‘minor’ mergers (empty dots)) and a high SFR ‘plateau’ during the final merger episode from z ≈ 7 to 6. This is mirrored in the evolution of the intrinsic Lyα emission of each galaxy as well as the total emission integrated across our three simulated galaxies, as shown in Fig. 3d, where there is a clear correspondence between the SFR peaks and peaks of intrinsic Lyα. Our simulations hence demonstrate that high SFR triggered by gas-rich mergers results in brighter intrinsic Lyα emission, which enhances the probability of these systems being detected.
To constrain this more robustly, we select an observable Lyα flux threshold of 10−18erg s−1 cm−2 consistent with our observations and compute the ratio of scattered to intrinsic Lyα emission over the spatial extent of our mock observations of the galaxies. This is shown in Fig. 3e, where the dark green curve denotes the median of the ratio of scattered to intrinsic Lyα flux and the shaded bands (quartiles and minimum-maximum values) show its distribution computed over 12 different LOSs. As the fraction of escaping Lyα emission has several high peaks, we conclude that the Lyα emission should be detectable for several of the redshifts studied here. Interestingly, during close galactic encounters (as shown by the distance trajectories of three main galaxies shown in Fig. 3e by coloured lines) and mergers, substantial gas tidal stripping and bursty star formation feedback lead to the opening of more low-column density sightlines, hence facilitating the escape of Lyα radiation. This effect can be clearly seen when comparing the most important merger instances to the peaks (high-tails) in the distribution of the Lyα escape fraction. Our findings about higher intrinsic Lyα emission and enhanced Lyα escape fraction in mergers are robust to the assumed dust model (Methods and Extended Data Fig. 6). It is worth noting that inclusion of AGN feedback in our simulations would only reinforce our findings, as expected black hole fuelling during mergers would power AGN49.
The rightmost panels of Fig. 3c–e show various physical relations for the observations (black points) and simulations (pink ‘violin’ symbols, denoting the distribution of Lyα escape fraction over 12 LOSs) using the same flux filter as above, demonstrating a very good correspondence between our simulations and our observed LAEs. Interestingly, for both observations and simulations we do not find clear correlations between the Lyα escape fraction and any other galaxy properties (distance, SFR, stellar mass), further corroborating our finding that ‘favourable’ LOSs with evacuated neutral hydrogen channels lead to directional Lyα photon leakage.
In this Article, we introduce a new interpretation to explain the unexpected detection of Lyα in z ≥ 7 galaxies in an epoch where the IGM is mostly neutral and gas-rich early galaxies are heavily ‘locally’ Lyα damped. Combining the high-resolution and high-sensitivity of JWST data with state-of-the-art radiative-transfer on-the-fly magnetohydrodynamics simulations, we demonstrate that three ingredients are key to making Lyα emission detectable from our sample of galaxies deep in the epoch of reionization: galactic mergers driving high intrinsic Lyα emission in the host galaxy, a ‘favourable’ LOS cleared of local neutral hydrogen in the host galaxy by tidal interactions with companions and by star formation feedback and a sufficiently large ionized bubble facilitating the escape of Lyα emission through the IGM.
Methods
Data
The emergence of JWST with its substantial technical advancements over both ground- and space-based telescopes has the potential to provide new hints as to the driving forces behind the observed Lyα emission discussed in this Article. We therefore decided to analyse all spectroscopically confirmed LAEs at z > 7 with JWST/NIRCam imaging. This search originally showed 11 such LAEs; however, further constraints on the SED of these galaxies from NIRCam imaging shows that the Lyα emission line detected in two of these LAEs is in fact not Lyα (discussed in more detail below). Therefore, the sample becomes nine LAEs (seen in Table 1): eight have companions that are identifiable in publicly available NIRCam imaging, and GN-z11 has spectroscopically confirmed companions that do not show an ultraviolet (UV) continuum30,39. These galaxies are found in the GOODS-North, GOODS-South, EGS and COSMOS fields, which have been observed with NIRCam by FRESCO (PI Oesch, ID 1895), JADES (PI Eisenstein, ID 1180), CEERS (PI Finkelstein, ID 1345) and PRIMER (PI Dunlop, ID 1837), respectively.
The PRIMER survey covers ∼150 arcmin2 of the COSMOS field with the F090W, F115W, F150W, F200W, F277W, F356W, F410M and F444W filters all at a 5σ limiting depth of AB magnitude (mAB) > 28. The CEERS survey took imaging of 100 arcmin2 of the EGS field in the F115W, F150W, F200W, F277W, F356W, F410M and F444W filters, reaching a 5σ limiting depth of mAB > 29 in the short-wavelength filters and mAB > 28.4 in the long-wavelength (LW) filters50. The June 2023 JADES data release surveys a ∼36 arcmin2 area of the GOODS-S field using the F090W, F115W, F150W, F200W, F277W, F335M, F356W, F410M and F444W filters, reaching 5σ limiting depths of mAB ≈ 30 (ref. 51).
The FRESCO survey’s primary aim is to obtain NIRCam LW Grism spectra over the F444W filter for galaxies across GOODS-S and GOODS-N. However, the simultaneous imaging in the short-wavelength channel provides imaging in both the F182M and F210M filters, and direct imaging obtained in the F444W filter also provides us with a total of three filters for our galaxies in GOODS-S/GOODS-N. FRESCO’s LW WFSS brings substantial advantages as the F444W filter covers the wavelength range of 3.9–5 μm, therefore observing both [O iii]4959,5007 and Hβ at z ≈ 7–9. These emission lines ordinarily contaminate imaging in the F444W filter, so constraints on the fluxes of these emission lines allow us to reduce the uncertainty in galaxy properties such as stellar mass, SFR and age. We use these constraints on the emission lines in F444W to create a mock F430M filter. Furthermore, when we fit the SED (discussed further below), we check the fluxes of the [O iii]4959,5007 and Hβ emission lines in the best-fit SED model and confirm that these are consistent with the measured fluxes from the WFSS spectrum.
Data reduction
Throughout our analysis, we use the final reduced data products created by the PRIMER team for the COSMOS field, the publicly available reduced data products from the CEERS team for the EGS field, the publicly available reduced data products from the JADES team for the deep HST region of GOODS-S and our own reduction of the GOODS-S FRESCO imaging for regions outside of the publicly available JADES imaging. To produce the GOODS-S images, we take all available NIRCam imaging data of GOODS-S available from the FRESCO survey and reduce them using v.1.8.5 of the JWST pipeline. We first download the uncalibrated files from the MAST archive and follow the steps of the JWST pipeline: (Stage 1) detector-level corrections and ramp fitting, (Stage 2) instrument-level and observing-mode corrections producing fully calibrated exposures and (Stage 3) producing the final mosaic. We take more care after Stage 1 to ensure the removal of horizontal and vertical striping that can be present in the rate files.
We extract the photometry of our objects using SEXTRACTOR52 on point-spread function (PSF)-matched images, extracting all objects with more than 4 px (DETECT_MINAREA 4) detected at more than 1.5σ (DETECT_THRESH and ANALYSIS_THRESH 1.5). We used large deblending parameters to allow efficient separation of each component (DEBLEND_NTHRESH 64 and DEBLEND_MINCOUNT 0.000001); however, the small separation between the components of the EGSY8p68 and JADES-GS-z7-LA systems makes the extraction of each component’s photometry impossible with SEXTRACTOR. For JADES-GS-z7-LA, we use the photometry of the LAE and its companion reported in ref. 53, and for EGSY8p68, we model the shape of each component using GALFIT54 assuming a Sersic profile and NIRCam PSFs, simulated using WEBBPSF. We first run GALFIT on one of the NIRCam images in which all the components are clearly resolved (F115W), and we then fix the resultant position, Sersic index, axis ratio and angle for each component when extracting the photometry for the remaining images. Extended Data Fig. 5 shows the modelling of the components for this system.
As discussed earlier, it is desirable to spectroscopically confirm as many companion galaxies as possible. The use of the Multi-Shutter Array (MSA) with NIRSpec to observe three of the systems (EGSY8p68, CEERS-1027 and JADES-GS+53.15682-27.76716) means that given the small shutter size, the companion is not coincident with the shutter position. Instead, we use all existing ground-based observations of our sample of LAEs, four of which have been observed by MOSFIRE (z7-13433, ID N199, PI Jung; EGSY8p68, ID C228M, PI Zitrin; GSDY, ID U069, PI Treu; z7-GSD-3811, ID C182M, PI Scoville) and one with X-shooter (COSY, ID 097.A-0043, PI Ellis).
We reduce the spectroscopic data using the standard MOSFIRE data-reduction pipeline and EsoReflex. These pipelines perform wavelength calibration, sky subtraction and flat fielding as well as further standard data-reduction steps to produce two-dimensional (2D) spectra. The resultant 2D spectra span the length of the original slit used for the observations and have spectral and spatial resolutions of R ≈ 3,380, 0.1798″per px and R ≈ 8,900, 0.158″per px for MOSFIRE and X-shooter, respectively, but are limited by the seeing of each observation, which can be as poor as ∼0.9″.
Two of our sample are located in the GOODS-S field (z7-GSD-3811 and GSDY). This field was the target of WFSS observations by the NIRCam LW Grism as part of the FRESCO observing programme. FRESCO’s LW WFSS covers the wavelength range of 3.9–5 μm using the F444W filter, therefore observing both the [O iii]4959,5007 and Hβ emission lines for z ≈ 7–9 galaxies, which are expected to be bright in high-redshift galaxies. These lines have previously been observed in z > 7 galaxies using the publicly available FRESCO data27,55. Moreover, the substantially superior spatial resolution of NIRCam LW WFSS (60 Mas allows us to resolve emission lines where MOSFIRE and X-shooter would fail to do so.
We reduce the Grism data following the steps described in ref. 56, including flat fielding, background subtraction, 1/f noise subtraction and World Coordinate System assignment. We also perform a careful astrometric analysis to avoid any offsets between the LW direct imaging in F444W and the LW channel Grism spectra in F444W. This allows us to extract the 2D Grism spectra associated with the astrometric position of any source; these spectra are then stacked, and the one-dimensional spectrum at the position of the source with an aperture of 3 px is extracted. We also perform a continuum subtraction from the one-dimensional spectrum to remove any contaminant continuum, and we focus on measuring [O iii]4959,5007 and Hβ emission line fluxes.
Redshift determination
We initially take care to ensure that for galaxies that have only Lyα emission detected in their spectra, the emission line detected is indeed Lyα. There are initially five LAEs at z > 7 in the literature, imaged by NIRCam, for which Lyα is the only observed emission line (z7-GSD-3811, z8-19326, z7-13433, z7-20237 and GSDY). To confirm that this emission line is Lyα, we use the SED-fitting code BAGPIPES32 on the available NIRCam data to produce an improved photometric redshift of each LAE. One of the key diagnostics in this fitting becomes the bright [O iii]4959,5007 emission lines expected from high-redshift objects, which fall into either the medium-band filter F410M (for z ≈ 6.8–7.6) or F430M (for z ≈ 7.3–7.8) or a lower flux consistent with continuum emission in this filter but a boost in the F444W filter (for z ≈ 6.8–9). This photometric-redshift fitting results in a redshift that is consistent with the redshift of the emission lines observed to be Lyα in all but two of these galaxies (z8-19326 and z7-20237). The photometric redshift, constrained by a boost in F410M (z ≈ 7.1) for z8-19326, becomes inconsistent with the emission line detected being Lyα, and therefore this galaxy is not included in our sample. Boosts in both F356W and F444W strongly indicate that z7-20237 is in fact at z ≈ 6. The original HST photometry reported in ref. 31 also slightly favoured a z ≈ 6 solution, making the observed emission line probably C iv. The removal of these two galaxies as contaminants underlines the challenges in determining galaxy redshifts with just one emission line and wide-band HST filters pre-JWST. The JWST offers an important advancement in this area, given NIRSpec’s sensitivity, allowing the detection of several lines, but also given the frequency of programmes that exploit the NIRCam medium-band filters around ∼4 μm. All the galaxies in our final sample show a clear [O iii]4959,5007 excess in NIRCam imaging. We can therefore be confident that the emission line detected is indeed Lyα.
Following the previous step, to confirm that all components of each system are at a similar redshift, we again use BAGPIPES to estimate their photometric redshifts. We then compare these with the spectroscopic redshift of the main component. This fitting confirms that all the companion galaxies preferentially have high photometric redshifts, all of which are consistent with the spectroscopic redshift of the massive galaxy. In combination with the prior of the companion being closely located to a confirmed high-redshift galaxy, these become strong-photometric candidate companion galaxies. We also use all available ground- and space-based spectroscopic observations that cover the companions to spectroscopically confirm the companion galaxy where possible. Below we discuss the available redshift constraints on our sample.
EGSY8p68, CEERS-1027 and JADES-GS+53.15682-27.76716
These galaxies constitute the subsample of LAEs that have photometric candidate companions with no spectroscopic confirmation. This lack of confirmation is due to a combination of factors. The observations of CEERS-1027 and JADES-GS+53.15682-27.76716 were made with the NIRSpec MSA, and the positioning of the shutter misses the companion; as such, spectroscopic constraints cannot be placed on these companions. EGSY8p68 has been observed by both the NIRSpec MSA and the ground-based telescope MOSFIRE/Keck. The MOSFIRE observations of the system cover the companion galaxies, but as a result of seeing the observations (∼0.7″) and the small separation between the components of the system (∼0.1″), resolving the components was not possible. Follow-up observations with the NIRSpec MSA show EGSY8p68 is probably an AGN, given the presence of high-ionization emission lines and a broadening of the Hβ emission line17. Although one of the companion galaxies lies in the shutter of this observation, the available spatial resolution and strength of the emission lines from EGSY8p68 makes it challenging to disentangle the emission from each component. However, the reported NIRSpec MSA flux is seven times fainter than that observed in the MOSFIRE spectrum57. Given that the MSA shutter size is notably smaller than the MOSFIRE slit, one solution to this discrepancy in flux is that some of the Lyα emission is extended outside of the MSA slit. We note that the component EGSY8p68-C sits outside of the MSA slit position; as such, this component may also host Lyα emission at the same redshift as the main component. Moreover, ref. 17 claims that the system is consistent with being involved in a major merger, and recent NIRSpec IFU observations of EGSY8p68 spectroscopically confirm that they are physical companions (Carniani et al., manuscript in preparation). However, because these results are not yet published and the data are not yet public, we continue our analysis of the redshift constraints assuming that we do not have this spectroscopic confirmation.
Despite the lack of spectroscopic constraints on this subsample, the similarities in the SEDs of the main and companion galaxies are distinct. The galaxies in each system show similar colours, Lyman breaks that are consistent and an emission-line boost, driven by [O iii]4959,5007 and Hβ emission, in the same filters, as can be seen in Extended Data Figs. 1 and 2. This results in a photometric redshift for each of the companions that clearly favours a z > 7 solution, and given the emission line feature in the F444W filter, the photometric redshift is further constrained to z = 7.6–9. Although the redshifts of these components are very confidently constrained to within this Δz = 1.4, further constraining their redshift is challenging without spectroscopy. However, we use the UV luminosity function from ref. 58 to estimate the expected number of galaxies within the volume to which the [O iii]4959,5007 emission redshift constraint constrains the companion galaxies. We note that the photometric redshift from BAGPIPES (reported in Table 1) is much more tightly constrained than that provided by the [O iii]4959,5007 boost alone, but we use this wider redshift constraint (Δz = 1.4) to provide an upper limit on the expected number of galaxies. Given that we have searched for companions in a 3″ × 3″ region, we use this as the area for this estimate and the [O iii]4959,5007 redshift constraints to produce a volume. We then integrate the z = 7 or z = 8 (dependent on the spectroscopic redshift of the main component) UV luminosity function down to the 5σ magnitude limit of the NIRCam images. We estimate that the expected number of galaxies in this FoV that are consistent with Hβ and [O iii]4959,5007 producing the excess F444W flux and with observed UV magnitude is \(0.0{9}_{-0.04}^{+0.16}\) for JADES-GS+53.15682-27.76716, \(0.00{7}_{-0.002}^{+0.009}\) for CEERS-1027 and \(0.00{3}_{-0.001}^{+0.006}\) for EGSY8p68. It is therefore very unlikely that another galaxy that is not associated with the spectroscopically confirmed LAE would be present in this volume.
JADES-GS-z7-LA
JADES-GS-z7-LA is an extremely high equivalent width LAE (400 Å) discovered in JADES NIRSpec observations of galaxies in the GOODS-S field53. The presence of a photometric candidate companion galaxy to this LAE has already been reported53, but the companion lies outside of the NIRSpec MSA shutter and hence is not spectroscopically confirmed. We perform the same analysis as above, noting that the [O iii]4959,5007 emission is in the F410M filter, limiting the redshift to between z = 6.8 and z = 7.6, and find the expected number of galaxies to be \(0.1{3}_{-0.05}^{+0.14}\). Therefore, once again, it is unlikely that the observed companion is at a redshift that makes it unrelated to the spectroscopically confirmed LAE. Moreover, we note the observed excess in the F115W filter that may be driven by Lyα emission in both the main and companion galaxies53, supporting the conclusion that these two galaxies are at the same redshift.
COSY and z7-13433
COSY, z7-13433 and their companions show excess fluxes in the F410M filter driven by [O iii]4959,5007 emission. SED-fitting of their companions strongly prefers a z > 7 solution, and the presence of [O iii]4959,5007 emission confines the redshift to z = 6.8–7.6. We again perform the same analysis as above, estimating the expected number of galaxies in the FoV. The expected number of galaxies is \(0.0{3}_{-0.01}^{+0.03}\) for COSY and \(0.0{4}_{-0.01}^{+0.03}\) for z7-13433, and thus we conclude that it is very likely that the central and companion galaxies are associated with each other.
Moreover, the companions of both COSY and z7-13433 are within the slits used for the observations with X-shooter and MOSFIRE, respectively. The two components of each system are shown in Extended Data Fig. 3, where we can identify both the positive (white) and negative (black) traces of the objects caused by the ABBA nodding pattern of the observations. We find that for z7-13433, two emission lines are offset by ∼3 px, which is consistent with the ∼0.5″ offset between the two components in the slit. The Lyα emission coincident with the position of the companion is detected at 2.8σ and corresponds to z = 7.479, and the Lyα flux associated with the main galaxy is at 4σ corresponding to z = 7.482. Therefore, these two components are offset by ∼125 km s−1. In the case of COSY, we expect a similar spatial offset between the two components, but the ∼0.7″ seen during the observations makes distinguishing between the two components very challenging. However, there is a second, previously unidentified, component to the Lyα emission of COSY at ∼0.9898 μm. We attribute this 2.7σ emission line to COSY-B and note a velocity offset of ∼280 km s−1 between COSY-A and COSY-B. Although these Lyα detections are clearly tentative, the fact that they are coincident with the expected position in the 2D spectra of strong-photometric candidate galaxies, together with their proximity to confirmed LAEs, boosts their likelihood of being a true detection.
Finally, we note that recent NIRSpec IFU observations of COSY spectroscopically confirm its companion galaxy38. However, as with EGSY8p68, because the data is not yet public, we continue our analysis without this spectroscopic information.
z7-GSD-3811 and GSDY
Both galaxies in this sample are in the FoV of the FRESCO NIRCam WFSS survey. Their final 2D WFSS spectra (Extended Data Fig. 4) show [O iii]4959,5007 emission lines for z7-GSD-3811, GSDY and their companions. For z7-GSD-3811, both components of the [O iii]4959,5007 emission are present, and there is an ∼3σ faint detection of Hβ. Moreover, [O iii]4959,5007 emission is also detected in the companion galaxy. This emission is offset from the central galaxy by 170 ± 70 km s−1. The 2D spectrum of GSDY shows a clear detection of [O iii]5007 in both the central and companion galaxies, and there is a very faint indication of [O iii]4959 and Hβ in the companion galaxy, but at a less than 3σ confidence for both. The emissions detected in the central and companion galaxies are offset by just ∼26 km s−1, but given the low spectral resolution of the NIRCam WFSS, this is smaller than the ∼70 km s−1 uncertainty on this measurement.
We discuss in the main text the use of Hβ to estimate the intrinsic Lyα flux and hence the escape fraction of Lyα photons; however, it can also be used as a probe of the instantaneous SFR. We assume a Hα to Hβ ratio of 2.85 and Case B recombination at Te = 10,000 K and Ne = 104 cm−3, and following ref. 59, we estimate the instantaneous SFR reported in Table 2. Given the lack of constraints on the reddening of these galaxies, we do not dust-correct the Hβ flux, and as such, these measured instantaneous SFRs are lower bounds.
The observation of several spatially resolved components of each emission line is taken as spectroscopic confirmation of the companion galaxies. Moreover, the velocity offset observed in all of these galaxies is indicative of either a LOS separation or, more likely, a difference in the LOS velocities of the two objects as they interact with each other. We therefore consider this further evidence that the systems are indeed interacting, as opposed to being two components of the same stable system.
GN-z11
Recent observations of GN-z11 with the NIRSpec IFU have shown the presence of three regions emitting either He ii (GN-z11-C)30 or Lyα emissions (GN-z11-B and a more tentative candidate)39. Although none of these regions show evidence of a UV continuum in deep NIRCam imaging, a tentative indication of Civλλ1548,1551 emission in GN-z11-B and analysis of the emission in GN-z11-C30 indicate that they are indeed being driven by a UV-faint stellar population in these regions. As such, we conclude that GN-z11 does have spectroscopically confirmed companion galaxies. Moreover, the presence of several LAEs in the local vicinity of GN-z11 indicates that peculiar velocity is not playing a vital role in the escape of Lyα through the neutral IGM (see ref. 39 for a detailed analysis of the Lyα emission present in this system). The density of objects in this field also makes GN-z11 probably a protocluster core39.
SED fitting
We then use BAGPIPES with the redshift fixed to that of the spectroscopic redshift of the system, using different SFHs (constant, delayed, burst and burst+constant) to fit the SED of all central and companion galaxies, as seen in Extended Data Figs. 1 and 2. The best-fit model is then obtained by taking the SFH that leads to the smallest Bayesian information criterion, as described in ref. 60. We report the galaxy properties returned from this method in Table 1.
Close-companion fraction of non-LAEs
To ascertain the companion fraction of a mass-matched sample of z > 7 galaxies, we take 30 galaxies with stellar masses from \(7.4 < \log [{{{{{M}}}}}_{* }({{{{{M}}}}}_{\odot })] < 9.3\) from refs. 23,61,62,63. We find that 47% of these galaxies have photometric-candidate companions within 5″ of the central galaxy. However, many of these galaxies have not been spectroscopically followed up and as such, we do not have constraints on their Lyα emission.
Unfortunately, examples of spectroscopically confirmed z > 7 non-Lyα-emitting galaxies that have been observed by high-resolution spectrographs are rare in the literature. As such, it is challenging to probe the companion fraction for known non-LAE galaxies. Moreover, the challenges of slit/aperture spectroscopy mean it is incredibly rare to have spectroscopic information on the local environment of such galaxies. Thus, although a galaxy may be confirmed as non-LAE, there is no information about the Lyα emission of any companion that may or may not be present.
The recent JADES data release of NIRSpec and NIRCam observations in GOODS-S includes seven z > 7 galaxies, observed using the NIRSpec high-resolution G140M grating, that show no Lyα emission13. This does not provide any Lyα information on any would-be companion galaxies; however, the companion fraction for this sample, ∼43%, is again consistent with both our larger sample and that of Puskás et al. (manuscript in preparation). We therefore conclude that the 100% companion fraction of our LAE sample is clearly atypical compared to the companion fractions of samples not selected for Lyα emission.
Close companions are necessary but not sufficient for Lyα emission
We note the recently reported z = 9.3 merging system64 with an absence of detected Lyα emission. The reported specific SFR of this galaxy (approximately−8) is consistent with the high specific SFR seen in our sample of merging LAEs, and therefore, one may naively posit that this is a contradiction to our results. We argue that non-detection of Lyα emission from a merging system is entirely reasonable and in fact our simulation modelling predicts this. Although we see boosts in the SFR and hence intrinsic Lyα emission from merging systems, the escape of Lyα emission is only possible when the escape fraction of Lyα emission out of both the host galaxy and the surrounding IGM is non-zero. We explore the effect of viewing angle on the observed Lyα emission escaping our simulated merging systems and find that although most viewing angles facilitate the escape of Lyα, several viewing angles exist for which the only Lyα emission escaping the galaxy is diffuse and therefore unobservable. Moreover, without the presence of a large ionized bubble, the neutral IGM will not facilitate the escape and observation of Lyα emission. Ultimately, we conclude that close companions are necessary but not sufficient for the observation of Lyα from galaxies deep in the epoch of reionization. A combination of a large-scale ionized bubble and a preferable LOS, stripped of neutral hydrogen by the interaction with the companion, are also required.
Azahar simulations
The new Azahar simulations are a suite of several models spanning various combinations of canonical hydrodynamics, magnetic fields, radiative transfer and cosmic ray physics with the aim of understanding their effects on galaxy formation as well as their complex interplay. The simulations are generated using the magnetohydrodynamical code Ramses65, which uses a constrained transport method for divergenceless evolution of the magnetic field66. The main target of study for Azahar is a massive spiral galaxy with Mhalo (z = 1) ≈ 2 × 1012 M⊙ in a relatively large zoom-in region, about 8 cMpc across along its largest axis. Azahar has a maximum spatial resolution in cells with a full cell width of Δx ≈ 20 pc (or, equivalently, an approximate cell half-size or radius of 10 pc), and refinement is triggered whenever its size is larger than a quarter of its Jeans length or it contains a total mass larger than 8 mDM, where mDM ≈ 4.5 × 105 M⊙ is the mass of dark-matter particles in the zoom region. The stellar mass particle is m* ≈ 4 × 104 M⊙. Azahar follows the set of physics presented by the pathfinder Pandora40. We provide here a brief summary, and we refer the reader to that reference and Martin-Alvarez et al. (manuscript in preparation) for further details.
In this work, we investigate one of the most complete simulations in the Azahar suite, the RTnsCRiMHD model. In this Article, we refer to the RTnsCRiMHD simulation simply as Azahar. Radiative cooling is modelled both above and below 104 K (refs. 67,68). The model adopts a magneto-thermo-turbulent prescription for star formation69,70, mechanical supernova (SN) feedback44 and astrophysically seeded magnetic fields through magnetized SN feedback71. This feedback injects 1% of the SN energy as magnetic energy and intends to reproduce the approximate magnetization of SN remnants. This particular simulation of Azahar also includes cosmic rays modelled as an energy density allowed to anisotropically diffuse, evolved with an implicit solver72,73. In Azahar, cosmic rays are exclusively sourced by SN feedback, with each event injecting 10% of its total energy as cosmic rays. We assume a constant diffusion coefficient for cosmic rays of κCR = 3 × 1028 cm2 s−1 and do not account for cosmic ray streaming in this model. This model also includes on-the-fly radiative transfer41, with a configuration similar to that of the SPHINX simulations42, featuring three energy bins for radiation spanning the ionization energy intervals 13.6–24.59 eV (H i ionization), 24.59–54.42 eV (HeI ionization) and above 54.42–eV (He ii ionization). Finally, recent work has shown that when expanding simulation models to account for cosmic ray feedback physics at comparable resolutions, the escape fractions of LyC photons from galaxies are lower than in their absence74.
Therefore, Azahar features a realistic (yet computationally taxing) ISM model that considerably approaches the resolution required to converge in the propagation and escape of ionizing photons from galaxies44. We note that although close, our resolution still falls short of that regime and that even more sophisticated ISM models75 may be important for a complete understanding of photon propagation in the ISM. Opportunely, a considerable part of our results relies on gas being ejected from galaxies during merger events, which implies that our estimate of LAE detectability during these events will be more resilient to these caveats.
Simulated ionized bubble evolution in the neutral IGM
In the high-resolution patch of the Universe probed by our simulation (≲8 Mpc on a side), centred on the main progenitor of our z ≈ 1 galaxy, the first ionized bubbles emerge at z ≈ 20. By z ≈ 15, several ionized bubbles exist that are a few (physical) kpc in radius and rapidly begin merging as redshift approaches z ≈ 10. This drastically increases the mean free path of ionizing photons, with the coalesced bubble reaching scales of 100 kpc. The resulting main bubble continues to grow, reaching a radius ∼0.5 Mpc at z ≈ 7.3.
For intrinsically bright LAEs, as studied here, a ≲0.5 Mpc local ionized bubble may be sufficient for these sources to be detectable76. We further note that because of continuing mergers, our simulated galaxies have considerable relative peculiar velocities of up to 230 km s−1. Sufficiently large velocities relative to the intervening neutral IGM gas may reduce damping wing suppression for particular LOSs, hence requiring smaller local ionized bubbles to facilitate the detection of Lyα. Our simulation also does not consider star formation outside of the high-resolution zoom-in region. Hence, our simulated ionized bubbles could be even larger: for example, as a result of the presence of a neighbouring large-scale overdensity or large-scale clustering77 on spatial scales beyond our zoom-in region.
Nevertheless, knowing from our observations that the LAEs sit in sufficiently large ionized bubbles, here we primarily focus on the probability of Lyα photons escaping the host galaxy and its intermediate surroundings. We find that a large amount of gas remains neutral within the local ionized bubble, actively feeding the galaxies along the cosmic web.
Galaxy selection and measurements in the simulations
To select galaxies in our simulation, we use the Halomaker software78 to detect and characterize dark-matter halos. We identify the three progenitor galaxies of interest (as well as the galaxy merging with the main progenitor at z ≈ 8.1) and follow them through time by tracking their innermost stellar particles. Their centres are determined using a shrinking-spheres algorithm applied to their stellar component79, and they are assigned their corresponding halo as obtained by Halomaker. To select the system for study, we broadly reviewed the evolution of all galaxies with stellar masses M* > 107 M⊙ in the redshift interval z ∈ (9.0, 6.0). The most promising candidate, the main progenitor of the main Azahar galaxy, was selected for further investigation.
To assign measurements to individual galaxies, we measure values within their galactic region, defined by the radius rgal < 0.2 rrvir. During the merging stages, we use r = min(rgal, 0.45Dij), where Dij is the distance between the i and j progenitors, down to a minimum distance of 1.5 kpc.
We have briefly reviewed the observability of lower-mass pairs in the simulated domain. Although these present behaviour during mergers similar to that of our main studied system, their observability and that of other galaxies with similar masses remains low. This is because their low stellar masses generate low luminosities. Their low stellar masses also make them more vulnerable to disruption during powerful star formation bursts, which increases their escape fractions43. Despite this, most isolated systems present quiet star formation histories and, in the absence of mergers, have low observabilities throughout their evolution.
RASCAS post-processing
We post-process the simulation with the publicly available, massively parallel code RASCAS46 for modelling the Lyα emission and resonant scattering in our simulations. RASCAS accounts for two sources for the Lyα emission: recombination and collisional excitation. For recombination, RASCAS adopts the case B recombination coefficient from ref. 80 and the fraction of recombinations producing Lyα photons from ref. 81. For collisional excitation, RASCAS uses the fitting function for the collisional excitation rate from level 1s to 2p from ref. 82. We cast NMC = 106 Monte Carlo photon packets to sample the real Lyα photon distributions from recombination and collisional excitation, respectively. We sample NMC = 107 photon packets for the images along the LOS in Extended Data Fig. 3.
After each scattering event, the Lyα photon changes its frequency and direction according to a phase function as implemented in RASCAS. RASCAS adopts the phase function in refs. 83,84 for the scattering of Lyα photons around the line centre and Rayleigh scattering for Lyα photons in the line wing. At high H i column density, Lyα photons will scatter many times locally until they shift in frequency enough that they can have a long mean free path. To reduce the associated computation cost, RASCAS adopts a core-skipping mechanism85 to transit the photons to line wings while avoiding local scattering in space. RASCAS also implements the recoil effect and the transition due to deuterium, with an abundance of D/H = 3 × 10−5. The dust distribution is modelled according to the dust model described in ref. 47. Dust can either scatter or absorb Lyα photons. The probability of scattering is given by the dust albedo adust = 0.32, following ref. 86. The dust scattering with the Henyey–Greenstein phase function87 and the asymmetry parameter is set to g = 0.73, following ref. 86. We generate synthetic images and spectra with a peeling algorithm described in refs. 49,88,89 along 12 LOSs uniformly sampled using the healpix algorithm90.
Finally, we note that the dust distribution has a large effect on the Lyα radiative transfer as well as the [O iii]4959,5007 and Hβ emission used to infer the intrinsic Lyα emission from observations. In Extended Data Fig. 6, we show how reducing the amount of dust changes the observed Lyα images. We see that the dust is able to completely mask the Lyα emission of some individual galaxies as resonant scattering of H i increases the probability of Lyα photons encountering dust grains. The dust can suppress the Hα emission by up to 20% around the centre of the three merging galaxies, if we assume intrinsic Lyα can be converted directly to intrinsic Hα emission. We therefore conclude that inclusion of dust modelling, as done here, is fundamental for understanding observed Lyα emission. However, because the channels are cleared of local neutral hydrogen, our results about the boosted intrinsic Lyα emission in star formation bursts and enhanced Lyα escape fraction in mergers are robust to the assumed dust model.
We further consider the Lyα profiles escaping along the LOS shown in the upper panel of Extended Data Fig. 3 at z = 7.3 after applying a reasonable IGM damping wing from ref. 91. We note that the shape of the Lyα profile is very complex and sensitive to both the LOS92 and IGM attenuation93. As such, further analysis of this is needed to constrain the expected Lyα profiles from simulations requiring a larger simulated galaxy sample with careful IGM attenuation modelling as well as a larger observational sample to compare to, which is beyond the scope of this work. However, as an immediate consistency test, we take the Lyα profile from our simulated system at z = 7.3 after IGM attenuation91 and compare the velocity offset and full-width at half-maximum (FWHM) to that reported in our observational comparison, EGSY8p68. We find that the observed Lyα profile17 is well matched with the simulated profile’s velocity offset of ∼200 km s−1 and FWHM of ∼200 km s−1 and shows a clear asymmetric line shape as observed.
Data availability
The JWST data used in this analysis is publicly available from the STScI MAST Archive.
Code availability
The WEBBPSF tool can be found at https://www.stsci.edu/jwst/science-planning/proposal-planning-toolbox/psf-simulation-tool. The standard MOSFIRE data-reduction pipeline can be found at https://keck-datareductionpipelines.github.io/MosfireDRP/and EsoReflex can be found at https://www.eso.org/sci/software/esoreflex/.
Change history
01 February 2024
A Correction to this paper has been published: https://doi.org/10.1038/s41550-024-02210-1
References
Curtis-Lake, E. et al. Spectroscopic confirmation of four metal-poor galaxies at z = 10.3–13.2. Nat. Astron. 7, 622–632 (2023).
Arrabal Haro, P. et al. Confirmation and refutation of very luminous galaxies in the early Universe. Nature 622, 707–711 (2023).
Williams, H. et al. A magnified compact galaxy at redshift 9.51 with strong nebular emission lines. Science 380, 416–420 (2023).
Dijkstra, M. Lyα Emitting galaxies as a probe of reionisation. Publ. Astron. Soc. Aust. 31, e040 (2014).
Heintz, K. E. et al. Extreme damped Lyman-α absorption in young star-forming galaxies at z = 9–11. Preprint at arxiv.org/abs/2306.00647 (2023).
Robertson, B. E., Ellis, R. S., Furlanetto, S. R. & Dunlop, J. S. Cosmic reionization and early star-forming galaxies: a joint analysis of new constraints from Planck and the Hubble Space Telescope. Astrophys. J. 802, L19 (2015).
De Barros, S. et al. VLT/FORS2 view at z 6: Lyman-α emitter fraction and galaxy physical properties at the edge of the epoch of cosmic reionization. Astron. Astrophys. 608, A123 (2017).
Endsley, R. et al. The REBELS ALMA Survey: efficient Ly-α transmission of UV-bright z = 7 galaxies from large velocity offsets and broad line widths. Mon. Not. R. Astron. Soc. 517, 5642–5659 (2022).
Bolan, P. et al. Inferring the intergalactic medium neutral fraction at z 6-8 with low-luminosity Lyman break galaxies. Mon. Not. R. Astron. Soc. 517, 3263–3274 (2022).
Hsiao, T. Y.-Y. et al. JWST NIRSpec spectroscopy of the triply-lensed z = 10.17 galaxy MACS0647—JD. Preprint at arxiv.org/abs/2305.03042 (2023).
Stark, D. P. et al. Lyα and C III] emission in z = 7-9 galaxies: accelerated reionization around luminous star-forming systems? Mon. Not. R. Astron. Soc. 464, 469–479 (2017).
De Barros, S. et al. The GREATS H β + [O III] luminosity function and galaxy properties at z ∼ 8: walking the way of JWST. Mon. Not. R. Astron. Soc. 489, 2355–2366 (2019).
Jones, G. C. et al. JADES: The emergence and evolution of Ly-alpha emission & constraints on the IGM neutral fraction. Preprint at arxiv.org/abs/2306.02471 (2023).
Bunker, A. J. et al. JADES NIRSpec Spectroscopy of GN-z11: Lyman-α emission and possible enhanced nitrogen abundance in a z = 10.60 luminous galaxy. Astron. Astrophys. 677, A88 (2023).
Jung, I. et al. CEERS: Diversity of Lyman-alpha emitters during the epoch of reionization. Preprint at arxiv.org/abs/2304.05385 (2023).
Laporte, N. et al. A spectroscopic search for AGNactivity in the reionization era. Astrophys. J. 851, 40 (2017).
Larson, R. L. et al. A CEERS discovery of an accreting supermassive black hole 570 Myr after the Big Bang: identifying a progenitor of massive z > 6 quasars. Astrophys. J. 953, L29 (2023).
Maiolino, R. et al. A small and vigorous black hole in the early Universe. Preprint at arxiv.org/abs/2305.12492 (2023).
Castellano, M. et al. First observational support for overlapping reionized bubbles generated by a galaxy overdensity. Astrophys. J. 818, L3 (2016).
Tilvi, V. et al. Onset of cosmic reionization: evidence of an ionized bubble merely 680 Myr after the Big Bang. Astrophys. J. 891, L10 (2020).
Leonova, E. et al. The prevalence of galaxy overdensities around UV-luminous Lyman-α emitters in the epoch of reionization. Mon. Not. R. Astron. Soc. 515, 5790–5801 (2022).
Morishita, T. et al. Early Results from GLASS-JWST. XIV. A spectroscopically confirmed protocluster 650 million years after the Big Bang. Astrophys. J. 947, L24 (2023).
Tang, M. et al. JWST/NIRSpec spectroscopy of z = 7-9 star-forming galaxies with CEERS: new insight into bright Lyα emitters in ionized bubbles. Mon. Not. R. Astron. Soc. 526, 1657 1686 (2023).
Tacchella, S. et al. JADES Imaging of GN-z11: revealing the morphology and environment of a luminous galaxy 430 Myr after the Big Bang. Astrophys. J. 952, 74 (2023).
Witstok, J. et al. Inside the bubble: exploring the environments of reionisation-era Lyman-α emitting galaxies with JADES and FRESCO. Preprint at arxiv.org/abs/2306.04627 (2023).
Rieke, M. J. et al. Performance of NIRCam on JWST in flight. Publ. Astron. Soc. Pac. 135, 028001 (2023).
Oesch, P. A. et al. The JWST FRESCO survey: legacy NIRCam/grism spectroscopy and imaging in the two GOODS fields. Mon. Not. R. Astron. Soc. 525, 2864–2874 (2023).
Bagley, M. B. et al. CEERS epoch 1 NIRCam imaging: reduction methods and simulations enabling early JWST science results. Astrophys. J. 946, L12 (2023).
Rieke, M. J. et al. JADES initial data release for the Hubble ultra deep field: revealing the faint infrared sky with deep JWST NIRCam imaging. Astrophys. J. Suppl. Ser. 269, 16 (2023).
Maiolino, R. et al. JWST-JADES. Possible population III signatures at z = 10.6 in the halo of GN-z11. Preprint at arxiv.org/abs/2306.00953 (2023).
Jung, I. et al. New z > 7 Lyman-alpha emitters in EGS: evidence of an extended ionized structure at z ∼ 7.7. Preprint at arxiv.org/abs/2212.09850 (2022).
Carnall, A. C., McLure, R. J., Dunlop, J. S. & Davé, R. Inferring the star formation histories of massive quiescent galaxies with BAGPIPES: evidence for multiple quenching mechanisms. Mon. Not. R. Astron. Soc. 480, 4379–4401 (2018).
Vernet, J. et al. X-shooter, the new wide band intermediate resolution spectrograph at the ESO very large telescope. Astron. Astrophys. 536, A105 (2011).
McLean, I. S. et al. Design and development of MOSFIRE: the multi-object spectrometer for infrared exploration at the Keck Observatory. In Proc.SPIE, Ground-based and Airborne Instrumentation for Astronomy III, 77351E (SPIE, 2010).
McLean, I. S. et al. MOSFIRE, the multi-object spectrometer for infra-red exploration at the Keck Observatory. In Proc. SPIE, Ground-based and Airborne Instrumentation for Astronomy IV, 84460J (SPIE, 2012).
Jakobsen, P. et al. The Near-Infrared Spectrograph (NIRSpec) on the James Webb Space Telescope. I. Overview of the instrument and its capabilities. Astron. Astrophys. 661, A80 (2022).
Rieke, M. J., Kelly, D. & Horner, S. Overview of James Webb Space Telescope and NIRCam’s role. In Proc. SPIE, Cryogenic Optical Systems and Instruments XI, 590401 (SPIE, 2005).
Übler, H. et al. GA-NIFS: JWST discovers an offset AGN 740 million years after the Big Bang. Preprint at arxiv.org/abs/2312.03589 (2023).
Scholtz, J. et al. GN-z11: The environment of an AGN at z = 10.603. Preprint at arxiv.org/abs/2306.09142 (2023).
Martin-Alvarez, S. et al. The Pandora project – I. The impact of radiation, magnetic fields, and cosmic rays on the baryonic and dark matter properties of dwarf galaxies. Mon. Not. R. Astron. Soc. 525, 3806–3830 (2023).
Rosdahl, J. & Teyssier, R. A scheme for radiation pressure and photon diffusion with the M1 closure in RAMSES-RT. Mon. Not. R. Astron. Soc. 449, 4380–4403 (2015).
Rosdahl, J. et al. The SPHINX cosmological simulations of the first billion years: the impact of binary stars on reionization. Mon. Not. R. Astron. Soc. 479, 994–1016 (2018).
Rosdahl, J. et al. LyC escape from SPHINX galaxies in the epoch of reionization. Mon. Not. R. Astron. Soc. 515, 2386–2414 (2022).
Kimm, T. & Cen, R. Escape fraction of ionizing photons during reionization: effects due to supernova feedback and runaway OB stars. Astrophys. J. 788, 121 (2014).
Garel, T. et al. Ly α as a tracer of cosmic reionization in the SPHINX radiation-hydrodynamics cosmological simulation. Mon. Not. R. Astron. Soc. 504, 1902–1926 (2021).
Michel-Dansac, L. et al. RASCAS: radiation scattering in astrophysical simulations. Astron. Astrophys. 635, A154 (2020).
Laursen, P., Sommer-Larsen, J. & Andersen, A. C. Lyα radiative transfer with dust: escape fractions from simulated high-redshift galaxies. Astrophys. J. 704, 1640–1656 (2009).
Asada, Y. et al. JWST catches the assembly of a z ∼ 5 ultra-low-mass galaxy. Mon. Not. R. Astron. Soc. 523, L40–L45 (2023).
Costa, T. et al. AGN-driven outflows and the formation of Lyα nebulae around high-z quasars. Mon. Not. R. Astron. Soc. 517, 1767–1790 (2022).
Finkelstein, S. L. et al. CEERS key paper. I. An early look into the first 500 Myr of galaxy formation with JWST. Astrophys. J. 946, L13 (2023).
Eisenstein, D. J. et al. Overview of the JWST Advanced Deep Extragalactic Survey (JADES). Preprint at arxiv.org/abs/2306.02465 (2023).
Bertin, E. & Arnouts, S. SExtractor: software for source extraction. Astron. Astrophys. Suppl. Ser. 117, 393–404 (1996).
Saxena, A. et al. JADES: Discovery of extremely high equivalent width Lyman-α emission from a faint galaxy within an ionized bubble at z = 7.3. Astron. Astrophys. 678, A68 (2023).
Peng, C. Y., Ho, L. C., Impey, C. D. & Rix, H.-W. Detailed structural decomposition of galaxy images. Astron. J. 124, 266–293 (2002).
Laporte, N., Ellis, R. S., Witten, C. E. C. & Roberts-Borsani, G. Resolving ambiguities in the inferred star formation histories of intense [O III] emitters in the reionization era. Mon. Not. R. Astron. Soc. 523, 3018–3024 (2023).
Sun, F. et al. First sample of Hα+[O III]λ5007 line emitters at z > 6 through JWST/NIRCam slitless spectroscopy: physical properties and line-luminosity functions. Astrophys. J. 953, 53 (2023).
Zitrin, A. et al. Lymanα emission from a luminous z = 8.68 galaxy: implications for galaxies as tracers of cosmic reionization. Astrophys. J. 810, L12 (2015).
Bouwens, R. J. et al. UV luminosity functions at redshifts z ∼ 4 to z ∼ 10: 10,000 galaxies from HST legacy fields. Astrophys. J. 803, 34 (2015).
Kennicutt, J. & Robert, C. Star formation in galaxies along the Hubble Sequence. Annu. Rev. Astron. Astrophys. 36, 189–232 (1998).
Laporte, N. et al. Probing cosmic dawn: ages and star formation histories of candidate z ≥ 9 galaxies. Mon. Not. R. Astron. Soc. 505, 3336–3346 (2021).
Laporte, N. et al. A lensed protocluster candidate at z = 7.66 identified in JWST observations of the galaxy cluster SMACS0723 − 7327. Astron. Astrophys. 667, L3 (2022).
Harikane, Y. et al. A comprehensive study of galaxies at z 9-16 found in the early JWST data: ultraviolet luminosity functions and cosmic star formation history at the pre-reionization epoch. Astrophys. J. Suppl. Ser. 265, 5 (2023).
Bouwens, R. et al. UV luminosity density results at z > 8 from the first JWST/NIRCam fields: limitations of early data sets and the need for spectroscopy. Mon. Not. R. Astron. Soc. 523, 1009–1035 (2023).
Boyett, K. et al. A massive interacting galaxy 525 million years after the Big Bang. Preprint at arxiv.org/abs/2303.00306 (2023).
Teyssier, R. Cosmological hydrodynamics with adaptive mesh refinement. A new high resolution code called RAMSES. Astron. Astrophys. 385, 337–364 (2002).
Teyssier, R., Fromang, S. & Dormy, E. Kinematic dynamos using constrained transport with high order Godunov schemes and adaptive mesh refinement. J. Comput. Phys. 218, 44–67 (2006).
Ferland, G. J. et al. CLOUDY 90: numerical simulation of plasmas and their spectra. Publ. Astron. Soc. Pac. 110, 761–778 (1998).
Rosen, A. & Bregman, J. N. Global models of the interstellar medium in disk galaxies. Astrophys. J. 440, 634 (1995).
Kimm, T. et al. Feedback-regulated star formation and escape of LyC photons from mini-haloes during reionization. Mon. Not. R. Astron. Soc. 466, 4826–4846 (2017).
Martin-Alvarez, S., Slyz, A., Devriendt, J. & Gómez-Guijarro, C. How primordial magnetic fields shrink galaxies. Mon. Not. R. Astron. Soc. 495, 4475–4495 (2020).
Martin-Alvarez, S., Katz, H., Sijacki, D., Devriendt, J. & Slyz, A. Unravelling the origin of magnetic fields in galaxies. Mon. Not. R. Astron. Soc. 504, 2517–2534 (2021).
Dubois, Y. & Commerçon, B. An implicit scheme for solving the anisotropic diffusion of heat and cosmic rays in the RAMSES code. Astron. Astrophys. 585, A138 (2016).
Dubois, Y., Commerçon, B., Marcowith, A. & Brahimi, L. Shock-accelerated cosmic rays and streaming instability in the adaptive mesh refinement code Ramses. Astron. Astrophys. 631, A121 (2019).
Farcy, M., Rosdahl, J., Dubois, Y., Blaizot, J. & Martin-Alvarez, S. Radiation-magnetohydrodynamics simulations of cosmic ray feedback in disc galaxies. Mon. Not. R. Astron. Soc. 513, 5000–5019 (2022).
Katz, H. et al. PRISM: A Non-equilibrium, multiphase interstellar medium model for radiation hydrodynamics simulations of galaxies. Preprint at arxiv.org/abs/2211.04626 (2022).
Hayes, M. J. & Scarlata, C. On the sizes of ionized bubbles around galaxies during the reionization epoch. The spectral shapes of the Lyα emission from galaxies. Astrophys. J. 954, L14 (2023).
Weinberger, L. H., Haehnelt, M. G. & Kulkarni, G. Modelling the observed luminosity function and clustering evolution of Ly-α emitters: growing evidence for late reionization. Mon. Not. R. Astron. Soc. 485, 1350–1366 (2019).
Tweed, D., Devriendt, J., Blaizot, J., Colombi, S. & Slyz, A. Building merger trees from cosmological N-body simulations. Towards improving galaxy formation models using subhaloes. Astron. Astrophys. 506, 647–660 (2009).
Power, C. et al. The inner structure of ΛCDM haloes - I. A numerical convergence study. Mon. Not. R. Astron. Soc. 338, 14–34 (2003).
Hui, L. & Gnedin, N. Y. Equation of state of the photoionized intergalactic medium. Mon. Not. R. Astron. Soc. 292, 27–42 (1997).
Cantalupo, S., Porciani, C. & Lilly, S. J. Mapping neutral hydrogen during reionization with the Lyα emission from quasar ionization fronts. Astrophys. J. 672, 48–58 (2008).
Goerdt, T. et al. Gravity-driven Lyα blobs from cold streams into galaxies. Mon. Not. R. Astron. Soc. 407, 613–631 (2010).
Hamilton, D. R. On directional correlation of successive quanta. Phys. Rev. 58, 122–131 (1940).
Dijkstra, M. & Loeb, A. Lyα-driven outflows around star-forming galaxies. Mon. Not. R. Astron. Soc. 391, 457–466 (2008).
Smith, A., Safranek-Shrader, C., Bromm, V. & Milosavljević, M. The Lyman α signature of the first galaxies. Mon. Not. R. Astron. Soc. 449, 4336–4362 (2015).
Li, A. & Draine, B. T. Infrared emission from interstellar dust. II. The diffuse interstellar medium. Astrophys. J. 554, 778–802 (2001).
Henyey, L. G. & Greenstein, J. L. Diffuse radiation in the Galaxy. Astrophys. J. 93, 70–83 (1941).
Yusef-Zadeh, F., Morris, M. & White, R. L. Bipolar reflection nebulae : Monte Carlo simulations. Astrophys. J. 278, 186–194 (1984).
Wood, K. & Reynolds, R. J. A model for the scattered light contribution and polarization of the diffuse Hα galactic background. Astrophys. J. 525, 799–807 (1999).
Górski, K. M. et al. HEALPix: A framework for high-resolution discretization and fast analysis of data distributed on the sphere. Astrophys. J. 622, 759–771 (2005).
Keating, L. C. et al. JWST observations of galaxy damping wings during reionization interpreted with cosmological simulations. Preprint at arxiv.org/abs/2308.05800 (2023).
Blaizot, J. et al. Simulating the diversity of shapes of the Lyman-α line. Mon. Not. R. Astron. Soc. 523, 3749–3772 (2023).
Witten, C. E. C., Laporte, N. & Katz, H. Evidence for a low Lyman continuum escape fraction in three massive, ultraviolet-bright galaxies at z > 7. Astrophys. J. 944, 61 (2023).
Roberts-Borsani, G. W. et al. z ≳ 7 Galaxies with red Spitzer/IRAC [3.6]-[4.5] colors in the full CANDELS data set: the brightest-known galaxies at z ∼ 7-9 and a probable spectroscopic confirmation at z = 7.48. Astrophys. J. 823, 143 (2016).
Song, M. et al. Keck/MOSFIRE spectroscopy of z = 7-8 galaxies: Lyα emission from a galaxy at z = 7.66. Astrophys. J. 826, 113 (2016).
Roberts-Borsani, G. et al. Nature and nurture? Comparing Lyα detections in UV-bright and fainter [O III]+Hβ emitters at z ∼ 8 with Keck/MOSFIRE. Astrophys. J. 948, 54 (2023).
Acknowledgements
We acknowledge the work of the PRIMER core team in obtaining and reducing the NIRCam data for the COSMOS field used in this work (PI Dunlop, ID 1837). We also acknowledge the work of the CEERS team in obtaining and reducing the NIRCam data of the EGS field used in this work (PI Finkelstein, ID 1345). We acknowledge the FRESCO team’s work in obtaining observations of the GOODS-S and GOODS-N fields (PI Oesch, ID 1895). We also thank J. Rosdahl for support with the generation of radiative transfer RAMSES simulations and T. Garel for support with the use of RASCAS. We thank M. Tang, D. Stark and H. Übler for providing further information on the properties of some of our samples. Support for this work was provided by NASA through grant no. JWST-GO-01837 awarded by the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS 5-26555. This research has used the Keck Observatory Archive, which is operated by the W. M. Keck Observatory and the NASA Exoplanet Science Institute under contract with the National Aeronautics and Space Administration and observations collected at the European Southern Observatory. This work used the DiRAC@Durham facility managed by the Institute for Computational Cosmology on behalf of the Science and Technology Facilities Council (STFC) DiRAC HPC Facility (www.dirac.ac.uk). The equipment was funded by BEIS capital funding by means of STFC capital grant nos. ST/P002293/1, ST/R002371/1 and ST/S002502/1, Durham University, and STFC operations grant no. ST/R000832/1. DiRAC is part of the National e-Infrastructure. C.W. thanks the Science and Technology Facilities Council for a PhD studentship funded by UKRI grant no. 2602262. N.L. acknowledges support from the Kavli foundation. D.S., M.G.H and J.S.D. acknowledge STFC support. P.S. acknowledges INAF mini grant no. 2022, ‘The evolution of passive galaxies through cosmic time’. R.S.E. acknowledges financial support from European Research Council advanced grant no. FP7/669253. P.G.P.-G. acknowledges support from Spanish Ministerio de Ciencia e Innovación MCIN/AEI/10.13039/501100011033 through grant no. PGC2018-093499-B-I00. D.P. acknowledges support by the Huo Family Foundation through a P.C. Ho PhD studentship. R.M., W.B. and W.M. acknowledge support by the STFC and European Research Council advanced grant no. 695671, ‘QUENCH’. R.M. also acknowledges funding from a research professorship from the Royal Society. This research was supported in part by the National Science Foundation under grant no. NSF PHY-1748958.
Author information
Authors and Affiliations
Contributions
C.W. reduced and analysed NIRCam imaging and WFSS, MOSFIRE and X-shooter data. N.L. extracted the photometry of each candidate and performed the SED-fitting analysis. S.M.-A. developed the Azahar simulations, performed their basic analysis and extracted the Lyα luminosity for all the merging systems. Y.Y. ran the RASCAS simulations and showed the influence of LOSs. D.S. and M.G.H. coordinated the simulations part of this paper. C.W., N.L., S.M.-A., D.S., Y.Y. and M.G.H. wrote the paper and developed the main interpretation of the results. R.M., G.R.-B., R.S.E., W.M., W.M.B., D.P., C.S. and H.K. were key in interpreting the observational results. J.S.D., R.S.E., N.A.G., G.I., A.M.K., D.M., P.G.P.-G. and P.S. were involved in obtaining and reducing the data for the PRIMER programme used in this study. All authors discussed the results and commented on the manuscript.
Corresponding author
Ethics declarations
Competing interests
The authors declare no competing interests.
Peer review
Peer review information
Nature Astronomy thanks the anonymous reviewers for their contribution to the peer review of this work.
Additional information
Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Extended data
Extended Data Fig. 1 SEDs of the sample.
Best SED-fits for all systems with strong Lyα. The yellow dots are the observed photometry, the grey line shows the best fit. The error bars indicate the uncertainty on the NIRCam photometry. The physical parameters of the best fit are also indicated (see text for details).
Extended Data Fig. 2 SEDs of the sample.
Same as Extended Figure 1.
Extended Data Fig. 3 Ground-based spectroscopy of.
The X-shooter (above) and MOSFIRE (below) Gaussian-smoothed 2D spectra of COSY and z7-13433, respectively. These cover the wavelength range of the previously spectroscopically confirmed Lyα emission of the system. White circular apertures indicate the proposed two components of the Lyα emission that are coincident with the expected positions of the main and companion galaxies in the slit. The positive trace of the emission is indicated by white pixels, the negative trace by black pixels, this effect of positive and negative traces is caused due to the ABBA nodding procedure employed by both instruments.
Extended Data Fig. 4 NIRCam Grism spectroscopy.
NIRCam WFSS Grism spectra of z7-GSD-3811 (above) and GSDY (below). Within each panel, the continuum subtracted 2D spectrum (above) and 1D spectrum (below) of the central galaxy are shown. The upper left sub-panel shows the F444W image of the system where two components can be identified. The horizontal dashed lines indicate the position of the central galaxy in both the image and 2D spectrum. The three dashed lines in the 1D spectrum denote the expected positions of the [OIII] doublet (central and right dashed lines) and Hβ (left dashed line) emission lines given the object’s spectroscopic redshift, while the grey shaded region denotes the one-sigma noise level.
Extended Data Fig. 5 Modelling of morphology.
Example of GALFIT modelling for the EGSY8p68 system for which the separation between each component is too small for the photometry to be extracted by SEXTRACTOR.
Extended Data Fig. 6 Effects of dust opacity on Lyα emission.
Dust opacities superimposed on synthetic Lyα images at z = 7.3. The top row shows the intrinsic Lyα emission, the middle row shows the scattered Lyα emission generated with a default dust model, and the bottom row shows the scattered Lyα emission with practically no dust (dust opacity reduced by a factor of 106). The three different columns show the galaxy from three different directions. Red and orange iso-contours correspond to a dust optical depth of 0.2 at the Lyα wavelength and Hα wavelength, correspondingly. Note that the dust does not act as a screen for the Lyα emission because of the resonant nature of Lyα scattering.
Rights and permissions
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article’s Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article’s Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.
About this article
Cite this article
Witten, C., Laporte, N., Martin-Alvarez, S. et al. Deciphering Lyman-α emission deep into the epoch of reionization. Nat Astron 8, 384–396 (2024). https://doi.org/10.1038/s41550-023-02179-3
Received:
Accepted:
Published:
Issue Date:
DOI: https://doi.org/10.1038/s41550-023-02179-3