Building the First Galaxies—Chapter 2. Starbursts Dominate the Star Formation Histories of 6 < z < 12 Galaxies

The HST and Beyond Committee (Dressler et al. 1996) chose a prime mission for the space telescope that would become the extraordinary JWST. The unanticipated reach of the "reborn" Hubble Space Telescope (HST) to galaxies with redshifts z = 2–3—only 2 billion yr after "the beginning"—held extraordinary promise for learning how the modern Universe actually began, when the first generations of stars collected into the first galaxies. Identifying this as the key moment in our origins—rather than, for example, the Big Bang—hinges on the idea that the miracle of this Universe and, more to the point, of life, is its incredible complexity. The chemical variety and abundance that made this possible could never have built to a critical level without multiple generations of heavy element–producing stars residing in these giant reservoirs we know as galaxies.

"The Hubble," as we now affectionately call it, took us to the trailhead of one of humanity's greatest journeys: traveling to this remote past to follow the story of our origins back to the present day.

Central to that mission is the simple notion that galaxies began in very different conditions—"cosmic environments"—at a time when dark matter was sculpting a terrain of mountains and valleys in a formerly smooth-as-glass Universe. There and then, "ordinary matter"—the baryons that we would be made of—were gathered and concentrated by gravity in a way that would fundamentally alter the very nature of the Universe. We had asked questions: how fast, how turbulent, how explosive, how dynamic were the forces that shaped these galaxies, so different from our times, when a slow and irreversible decline dooms this Universe to a grinding conclusion that we—fortunately—will never see.

The tremendous effort and challenge of building JWST is behind us. The unearthly, spectacular views we now have of the Universe—near and far—let us start that remarkable journey with a large, representative sample of faint, newborn galaxies, for a first look at how they survived their tumultuous births.

Where to start? A long-standing and popular approach to the study of galaxy formation and development has been to push to higher and higher redshift to compare the properties of ever-younger galaxies. JWST has already exceeded expectations with images and confirming spectroscopy of star-forming galaxies at z ∼ 13 (Curtis-Lake et al. 2023; Robertson et al. 2023), a scant 300 Myr after the Big Bang. Other studies identify candidates for even earlier galaxies; deeper searches and microlensing by rich clusters (Treu et al. 2022) promise to provide candidates as early as redshift z ∼ 18. However, by their very nature, such searches find extraordinary objects—the brightest sources at the earliest epochs. In order to assemble a representative sample of the first galaxies, a later starting point is required; early JWST results suggested z ∼ 12 as an epoch where dozens of galaxies could be identified in the JWST Advanced Deep Extragalactic Survey (JADES) GOODS-S field, with numbers increasing rapidly into the many hundreds for z ∼ 6. This motivated the choice of 6 < z < 12 for the present study of the star formation histories (SFHs). But for their number, these are the ancestors of the common galaxies we find today.

The paper is organized as follows. Section 2 briefly reviews the methodology—the program SEDz*—the details of which can also be found in Dressler et al. (2023, hereafter Paper I). Section 3 lays the foundation for this study, showing and explaining a new format for SEDz* results, with new data that highlight the four types of SFHs. Section 4 explains how a robust, representative sample of ∼900 galaxies was chosen from a catalog of GOODS-S sources with nine-band NIRCam fluxes (spectral energy distributions, SEDs) from JADES. Section 5 shows 72 examples of SFHs for four redshift ranges 6 < z < 12, as measured by SEDz*. Section 6 looks for changes in the mix of SFHs over time, and Section 7 presents a summary of the birth of stellar mass—when, where, and how much—from z ∼ 12 down to z ∼ 6. Section 8 investigates correlations of SFHs with environment, and Section 9 looks ahead to the implications of this study for future observational and theoretical programs. Appendix A shows histograms of SFH type with local density and nearest-neighbor distance, Appendix B shows examples of SFH templates used by SEDz*, and Appendix C discusses the impact of emission lines on the continuum-dominated SEDs of this study.

A standard cosmology with Ω_m = 0.3, Ω_Λ = 0.7, and H₀ = 70 km s⁻¹ Mpc⁻¹ is adopted throughout.

The buildup of stellar mass in a galaxy is expressed as its SFH, that is, the rate of star formation changing with time. Measuring the SFHs of galaxies has proven more than challenging, because while groups of young, massive stars are easy to recognize, populations of stars older than a billion years—though still young compared to the Universe—age in such a way that it is not possible to tell a population of old stars from another that is even older.

As part of a galaxy survey of 6 billion yr of cosmic history, Kelson et al. (2014) developed a tool for analyzing a galaxy's SED—basically, its rainbow of color. By measuring SEDs for galaxies observed ∼4 Gyr earlier in cosmic history, Kelson's analysis leveraged those observations to look back ∼2 Gyr further, to a time when the Universe was only around half its present age. ¹⁶ Because of that, the part of a galaxy's history that could be reliably measured was extended to when the Universe had half its present age, in this way revealing a substantial fraction of a galaxy's growth history. Applying this approach to a sample of galaxies observed at redshift z ∼ 0.4–0.8, Dressler et al. (2018) discovered that ∼20% of the galaxies were still growing rapidly, with increasing star formation rates (SFRs), in an era that was thought to be a time of decline for all galaxies as massive as the Milky Way.

The motivation for this study, then, emerged directly from that one. The ability to see only 1 billion yr (or at most 2) of SFH, a small part of a mature galaxy's lifetime, stretches to "longer than the age of the Universe" during the period when JWST's prime mission kicks in. Our job has been to gather sufficiently good SEDs, using JWST's near-infrared camera, NIRCam, to make accurate measurements of galaxy brightness at a series of colors—infrared light from 1–5 μm. Such data can fully constrain the ages of populations of stars that make up a galaxy, and lining up their stellar masses—epoch by epoch—constitutes an SFH of how the galaxy was built.

In this paper, we show the results of our second attempt (see Paper I) to transform the flux measurements of each galaxy—its SED—into SFHs, the buildup of stellar mass over the first billion years. Our subjects are a greatly expanded sample of 894 galaxies at redshifts 6 < z < 12 with nine-band near-infrared fluxes from JWST/NIRCam. We provide a short description of a program code—SEDz*—written expressly for the purpose of reconstructing the histories of galaxies from the rich information encoded in their SEDs. By choosing a sample of galaxies with ages less than a billion years, their fast-evolving stellar populations will be recorded for us to play back.

The development of SEDz* followed Kelson's maximum-likelihood Python code for analyzing galaxy SEDs from the Carnegie–Spitzer–IMACS Survey (Kelson et al. 2016) through the combination of stellar population templates. The program effectively isolated light from younger (≲1 Gyr) populations of main-sequence A stars and led to the discovery of "late bloomers"—the ≳20% of galaxies at z ∼ 0.5 that produced at least 50% of their stellar mass in the previous 1 Gyr, that is, rising SFRs instead of the falling SFRs that are conventionally described as predominant after z = 1. These late bloomers challenge theoretical models that tightly couple the growth of galaxies to that of their dark matter halos, because while it is possible to think of mechanisms that could delay star formation, it is not easy to imagine why some halos might form so much later than all the rest.

With this approach, the program SEDz* was developed to exclusively measure SFHs of the first galaxies. In Paper I, we described the challenges associated with measuring SFHs for stellar populations with ages of more than 1 Gyr, a deficiency turned on its head when the the population under study has an age of τ ≲ 1 Gyr. SEDz* takes advantage of this unique astrophysical opportunity that comes from the billion-year lifetimes of A stars (Dressler & Gunn 1983; Couch & Sharples 1987). Because they evolve rapidly over a Gyr, it is possible to derive SFHs from SEDs for A star–dominated populations, and vice versa. A stars are among the simplest main-sequence stars (Morgan & Keenan 1973), basically blackbody radiators with a relatively simple internal structure and an opacity produced by hydrogen absorption and thus free from the complications of chemical abundance and post-main-sequence evolution. These are critical and unique advantages for measuring SFHs during the first billion years.

The data input to SEDz* are SEDs—flux measurements in NIRCam's seven wide bands. ¹⁷ SEDz* uses a nonnegative least-squares (NNLS) engine (Lawson & Hanson 1995) and custom star formation templates (Robertson et al. 2010) ¹⁸ that are essentially a set of vectors that have a significant amount of "orthnomality," as can be seen in the plots in Appendix B, where a more complete discussion of templates can be found.

For this study, SEDz* divides redshifts 6 < z < 12 into seven integer steps, 12, 11...6, lasting 42, 50, 61, 79, 101, 134, and 186 Myr, respectively. SEDz* operates with two sets of SED templates, one with 10 Myr bursts (unresolved after subsequent star formation) and another characterized by constant star formation (CSF). The program builds up an SED by adding stellar population templates (starting at z = 12 and working down) as needed to improve the fit and evolving them forward—adding up subsequent populations to improve the fit, as measured by χ². The epoch of observation (OE) is chosen as the lowest χ², and the SFH is read off as the stellar mass contributed by each scaled template that, in combination, make the best fit. This can include the addition of a CSF template that signals CSF at OE. SEDz* can combine a CSF template with a final burst template to expand its fitting possibilities, as shown in Appendix B.

SEDz* requires no "priors" in the important sense that every trial solution, as lower-redshift templates are added one by one, is fully independent of previous ones. This allows for modeling free of any particular functional form (e.g.,nonparametric) and allowing for a range of galaxy histories. However, that range is bounded; only 33 templates are used (28 bursts, z = 12 to z = 6, evolving, and 5 CSF, z = 10, 9, 8, 7, 6—at the OE). Since flux ratios are fixed for each template, choosing a template constrains all seven bands. There is no fitting of individual bands, so an acceptable solution is not guaranteed; finding one also indicates that the SED can be well described by the templates. ¹⁹

Paper I describes the working of the code in considerable detail, so we do not repeat it here. Also discussed in Paper I are tests of SEDz*, including its ability to reproduce the SEDs of synthesized galaxies in a simulated deep field, from the NIRCam "Data Challenge" (Williams+2018). The Appendix of Paper I shows how test SFHs generated by combining the stellar population templates were recovered by SEDz*, and how the distinction between bursts and extended SFHs is robust.

As in Paper I, we neglect the potential impact of dust. We note, as before, that galaxies in our sample seem to be well described by SEDs with little or no dust, consistent with the results of several papers showing that these initial JWST-selected samples are uniformly blue (Nanayakkara et al. 2023) and fairly dust-free. Our library of stellar population templates includes sets with extinctions of A_v = 0.4 and 1.0. We show in Appendix Figure 13 that extinction at the level of A_V = 0.4 is ruled out; in fact, even A_V = 0.1 would produce a reddened SED not found in the 894-galaxy sample. While we did find 27 other examples of steeply rising SEDs in the current work, their identification as galaxies at high redshift is not possible with SEDz* because of its reliance on the A-star model. ²⁰ We conclude from this that our sample is by construction a "no-extinction" sample, and that it does describe more than 95% of galaxies found by deep NIRCam imaging in this field.

However, we do speculate here that the absence of rising SFHs in this and our previous study could indicate that this phase of galaxy building is largely hidden by dust, and that the SFHs we find began with 10–20 Myr of intense, dust-enshrouded star formation, then cleared by its explosive feedback.

In Paper I, we introduced and explained the SEDz* code and applied it to data from the GLASS JWST/ERO program (Treu et al. 2022) from parallel imaging with NIRCam of NIRISS grism spectroscopy of the cluster Abell 2744 (Merlin et al. 2022). Some challenges in processing and calibrating one of the first deep-imaging programs—a crucial rationale for the ERO program—limited the targets of Paper I to only 24 galaxies that were judged suitable for a first measurement of SFHs.

To that point, SEDz* had only been tested on simulated data of the NIRCam deep-imaging program (the Data Challenge; (Williams et al. 2018) and on simulated SFHs produced using SEDz* itself, an admittedly easier test to pass. Valid questions had been raised about how different the SEDs of the earliest galaxies might be from their descendants, particularly because the nature of stellar populations at these early times was largely unknown. However, the application of SEDz* to the first such data produced encouraging, and surprising, results, in the sense that the program was able to reproduce 24 complex nine-band SEDs with the code's limited library of stellar population templates. ²¹ In other words, SFHs with recognizable characteristics, bursts—single and multiple, with contiguous "triples," and continuous "long SFHs," with reasonable masses in the range of 10⁸–10⁹ M_⊙—fit all 24 SEDs within the errors of the photometry. Considering the limited number of stellar population templates available to SEDz* and their unique shapes, obtaining good fits from the "get-go" was surprising and remarkable.

With the comparatively exquisite photometry from NIRCam imaging of GOODS-S for JADES in late 2022, deeper imaging data have led to an ∼900-galaxy sample at redshifts 6 < z < 12, and our analysis of these new data confirms the basic conclusions of Paper I: good reproduction of observed SEDs and four SFH "types" and, in particular, well establishing the unexpected prominence of starburst SFHs over longer, steadier runs of star formation. Moreover, we now have sufficient data to begin to examine the dependence of these SFHs with redshift, mass, environment, and large-scale structure.

We begin by revisiting the variation in SFH type. Figure 1 introduces a new data format for SEDz* output and gives examples of the four types of SFHs found in Paper I and now in this paper, demonstrating both how the code derives SFHs and that it does so with considerable fidelity. The left-hand box displays the observed SED—fluxes in each of the nine bands with 1σ error bars. The NNLS solution SEDz* found by combining stellar population templates is the magenta band—the quartile range of 21 iterations, each a random perturbation of the SED by its errors. The χ² of the fit is inset in the upper left; a prominent minimum in χ² defines the observed redshift, or OE. The solution and operation of SEDz* are recorded in the right box, which shows the stellar mass added at each of seven epochs (integer redshifts z = 12, 11...6). This is how the SFH is calculated, by scaling and combining stellar population templates to make the best NNLS fit to the observed SED. The four plots of Figure 1 show this procedure through animation; by "clicking" on a plot, the steps SEDz* takes to derive the best-fit SFH are displayed in sequence.

In Figure 1, this "best fit" for the upper left SED is a starburst, a single epoch of star formation observed at z = 7—in this case, a combination of a burst and CSF. As shown in the top-left panel of Figure 13 (Appendix B), the slope of CSF alone is more level and a burst alone much steeper than what is observed. SEDz* has determined that a ratio of 2:1 (more CSF than burst) is an excellent fit to the observed SED. That is, this two-component stellar population model from present-epoch stars in our Galaxy "works." While it might seem remarkable that present-epoch stellar populations provide an almost perfect fit, in one sense this is completely unremarkable: the A stars that dominate the light of this starburst at the OE are the least-complicated stars along the main sequence, a fully convective core and a fully radiative envelope with opacity from hydrogen ions—no metals required.

The stellar mass calculated for star formation at z = 7 is ∼3 × 10⁸ M_⊙. The green arrow signifies the CSF contribution, with ongoing star formation at OE. ²² The CSF template accounts for emission at OE, detected as the modest elevation of F356W and F444W from the stellar continuum level (F410M) due to [O iii] emission. The CSF template only serves at OE, because stars formed in a prior epoch are, by definition, old. In fact, although a burst or CSF is the way SEDz* models the growth of stellar mass through the epochs, these are indistinguishable by the next epoch, becoming an "old"' population, regardless of how the mass was spread over the interval. It is also important to remember that, despite representing continuing CSF at OE, most of the light from the CSF template is from young A stars; the O and B stars, and any nonstellar continuum and emission are gone, except for the last 10–20 Myr.

This starburst example begs the question, could another SFH—for example, a combination of bursts—also reproduce the upper left SED? The upper right panel shows that the answer is "no" by means of an SED that is perfectly fit with a combination of two bursts: one of ∼5 × 10⁸ M_⊙ at z ≈ 12 and an ∼2 × 10⁸ M_⊙ CSF+burst ∼500 Myr later (OE z ≈ 6). The integrated mass of this SFH is the red-encircled black dot at ∼7 × 10⁸ M_⊙. Epochs without star formation are marked as lower limits of 10⁶ M_⊙. In this study, at the depth of the JADES data for the GOODS-S field, detections as low as 10⁷ M_⊙ have been made but are incomplete below ∼5 × 10⁷ M_⊙(limited by the "noise" from other epochs of star formation), and severely so below 2 × 10⁷ M_⊙. Nevertheless, it is clear that multiple epochs of star formation cannot combine to produce the SED on the left, nor can a single epoch of star formation imitate a double.

Perhaps the most revealing feature of this double-burst example is how clearly it shows SEDz* combining two widely separated epochs of star formation to match the observed SED. The two colored curves below the SED are A-star templates that correspond to the age of the two bursts. The z = 12 burst in orange is that of an "old" stellar population (∼500 Myr before OE), while the purple curve has the strong ultraviolet flux of a young population (∼100 Myr), demonstrating the ability of SEDz* to resolve the 6 < z < 12 epoch into old and young populations. This time it is the burst at z = 6 that dominates CSF, by 2:1, which accounts for the steeper slope of the purple template, as shown in the right plot in Figure 13 (Appendix B). It is the striking differences between the two templates that together reproduce the observed SED that shows how accurately SEDz* can decompose an SED into separate stellar populations, and how well an SFH can be matched with a comparatively small number of A-star templates. ²³

Before moving on to examples of longer SFHs, we want to clarify why these two examples are called "starbursts," given that each represents a few hundred million M_⊙ in approximately 100 Myr. An averaged SFR of only a few M_⊙ yr⁻¹ is not usually considered a "burst." However, we have no information on the history of star formation (SFR) over these individual epochs; we are not measuring SFRs but rather the total stellar mass associated with that epoch. What makes the upper left example a "burst" is that the z = 7 epoch contains only young A stars; that is the template that fits the SED. No older A stars—stars born in previous epochs—are detected. As we discuss in Section 9, we think it more likely that the stellar mass we find within a single epoch was born in ∼10–20 Myr early within the epoch, with SFRs in the "tens" M_⊙ yr⁻¹), and a falling SFR (or even a shutdown) in the 50–100 Myr after. Such "more-likely" SFHs would account for the A star–dominated SED and a corresponding absence of older A stars, τ ≳ 150 Myr. And that is a "burst," by any definition.

The triple-epoch SFH shown in the lower left panel is common in the sample, particularly at z ∼ 6–8. This example is particularly revealing because of the flat (purple) SFH at z = 6 from a CSF template that also produces [O iii] and Hα emission, a burst that contributes the peak flux at λ ∼ 1 μm, and an ∼200 Myr old aging burst from the z ∼ 8 (green) burst template (but indistinguishable as a burst after z = 8) that produces the rising flux beyond ∼3 μm. In addition to "contiguous," the three epochs of star formation are close in mass at 2 × 10⁸ M_⊙, a near-constant level of SFH that might continue to z = 5 and beyond. We note that a small but significant fraction of these "triples" could be fit with only two contiguous epochs of star formation, but none of them can be fit by a single epoch of SFH. And, at the redshift most are seen, two or three epochs amount to 300–400 Myr in duration, long enough to host multiple generations of star formation.

This leads naturally into the example in the lower right, where four epochs of the possible six have substantial star formation that is both continuous and, like the triple SFH, a coherent decline in added stellar mass. The most important thing to recognize about these SFHs is that, again, they cannot be reproduced by the single or a pair of bursts of star formation like in the upper left and right or by the shorter contiguous episodes that are shown in the lower left panel (see also Appendices B and C of Paper I). A long continuous SFH is defined here as four epochs or more of star formation (1) with less than a factor of 2 uncertainty in stellar mass and (2) free of two-epoch gaps. The example shown here, and the more that follow in Figures 2–5 , are representative of 8% of our 6 < z < 10 sample. It is important to reiterate that what matters for these longer, continuous SFHs is that they exist, rather than a dissection into gaps, spikes, dips, and wiggles. The data for this study are insufficient for SEDz* to deliver such detail, and indeed, this might not even be possible with higher-resolution photometric filters or spectra due to limitations in NNLS fitting. For example, a single-epoch gap in one of the continuous histories is likely a consequence of nonnegative solutions only: at any given redshift, adding star formation at that epoch may not improve the SED fit and might even degrade it. For this reason, the longer SFHs we find here are indicative of the general, not the detailed, behavior of star formation in the galaxy.

The SEDs and solutions shown Figure 1 are exemplary; they display the very good fits of a majority of SEDs to single-epoch templates (bursts) and the ability of combination templates to fit a range of SED shapes. That the fraction of cases that fail to find a good solution is small is, we think, strong evidence that the principle behind SEDz* is valid: A star–dominated stellar populations can be used to recover accurate SFHs for galaxies in the τ < 1 Gyr period of cosmic history.

In the examples to follow, we will see that SFHs of significant duration are all found to be declining or nearly flat. What is almost absent from the SFHs shown here is a stretch of rising star formation. As we explained in Paper I, it is possible to attribute this to the difficulty of detecting a declining but "current" older population against a young one that was rising in flux but is now fading—a simple selection effect. However, a more intriguing possibility is that the flat/declining phases of nearly all first-galaxy SFHs—from bursts to long SFHs—were preceded by a strong, rapid burst of dusty star formation lasting ≲50 Myr. Among the galaxies that were excluded from our sample we find 27 with steep, red SEDs, but SEDz* cannot be used to determine their redshifts, so their connection to the unobscured galaxies in our sample remains unknown.

NIRCam's exceptional performance and diversity of modes (Rieke et al. 2023) have supported Guaranteed Time Observations that cover a wide variety of science programs, including exoplanets, star formation, and our Galaxy with its neighbors. However, the largest component is devoted to the study of the early Universe (Eisenstein et al. 2023)—the era of galaxy birth. The program described in this paper comes from deep-field imaging that is beginning to answer long-standing questions about how the first generation of stars—collected into galaxies that served as reservoirs for the buildup of heavy chemical elements—fundamentally changed the evolution of our Universe.

This study uses nine-band NIRCam imaging of an ∼25 arcmin² field of the JADES GOODS-S survey. The images have been intensively processed for calibration and to remove instrument signatures, combine dithered exposures, and characterize noise. These data were used to produce catalogs of objects (here v0.7.2) that identify stars and galaxies, deblend overlapping images, and generate a number of different radial extractions of photometric measurements in nine filters, NIRCam bands F090W, F115W, F150W, F200W, F277W, F335M, F356W, F410M, and F444W. The catalog used here reached a depth of ∼30.0 AB mag and contained 24,350 extended sources. Our study selected galaxies with a signal-to-noise ratio (S/N) > 5σ (F200W and F277W flux) within a 4 pixel diameter circular aperture, appropriate for the small size of z > 5 galaxies. Details of data quality, data reduction, and the production of photometric catalogs can be found in Tacchella et al. (2023) and Rieke et al. (2023). To establish our sample, SEDz* was run on the complete GOODS-S v0.7.2 catalog, with results sorted into four different redshift ranges: 6 < z < 7, 7 < z < 8, 8 < z < 10, and 10 < z < 13. To be precise, the actual ranges were shifted down by 0.25 in z, for example, 5.75 < z < 6.75, in order to center on the epoch redshift (to match the templates), in this case, z = 6. From the comparison of SEDz* redshifts with the "known redshifts" of the synthetic Data Challenge, SEDz* was found to have an accuracy of Δσ_z ≈ 0.15 for z > 6 galaxies. This better-than-expected performance meant that interpolation between the templates used by SEDz* was justified, and in fact, not interpolating could add to systematic errors. For this reason, boundaries were set as follows: z1: 5.75 < z < 6.75, z2: 6.75 < z < 7.75, z3: 7.75 <  z < 9.75, and z4: 9.75 < z < 12.75. Thus, when a galaxy's redshift was found to be within the range 5.75 < z < 6.25, the z = 6 template was used, but when 6.25 < z < 6.75, a 50/50 interpolation between the z = 6 and 7 templates was made (and so on for the other samples). In effect, this is equivalent to saying the χ² for z = 6 and 7 are the same (within errors), so the program averages the SED solutions (resulting SFHs) for those two redshifts. ²⁴ Eventually, data for the separate redshift ranges would be combined into a catalog spanning the full 5.75 < z < 12.75 range, but the previous step allowed an investigation of performance over the redshift range that was helpful.

By fitting combinations of stellar population templates, SEDz* found maximum-likelihood fits to SEDs that yielded redshifts in these four redshift intervals. A well-recognized problem in deriving z > 4 redshifts from SEDs is a degeneracy with z ∼ 2 galaxies, where NIRCam's range of 1–5 μm translates into a rest-frame coverage of ∼0.3–1.7 μm, typically covering a Balmer break over an otherwise flat SED. When very faint galaxies are the targets, these are easily mistaken for z ∼ 6–8 galaxies with a Lyman break. We used three different methods to mitigate the problem. First, we wrote code in SEDz* that compared the shape of the χ² curve when two minima were found, one for z ∼ 2 and another for the higher redshift. During the Data Challenge tests, we found that using the overall slope of the run of χ² with redshift, the depth and width of the minima, and the color of the SEDz* itself, removed roughly 50% of cases. We also used archived HST data to find at least two 2.5σ detections in three visible bands of the WFC3 imaging, F606W, F775W, and F814W. Flux below the Lyman break is the best rejection method, but we have found that, for galaxies this faint in the near-infrared, only about one-third of low-redshift galaxies are detected in these bands, even with the deepest HST imaging available. The third check was to use EAZY "photo-z" redshifts for the cataloged objects, part of the JADES team data processing for internal and eventual community use (Hainline et al. 2020).

The procedure was to run SEDz* for each redshift interval using only (1) the rejection of low-z objects by detection of visible flux and (2) the internal SEDz* z < 4 rejection. This produced the four subsets z1, z2, z3, and z4 with 759, 374, 277, and 61 galaxies, respecitvely. Taking the SEDz* redshift (hereafter z_SED) as the adopted redshift was required, because that is the value for which the SFH is derived, but comparing z_SEDwith the EAZY-derived redshift, z_a, and finding them consistent was taken as the next level of "qualification." This was termed the "gold sample" and amounted to 446, 183, 115, and 15 "confirmed" objects. This left 313, 191, 162, and 46 objects with z_SEDunconfirmed by z_a. Their SEDs were inspected, one by one, to decide if the Lyman-break range was sufficiently well defined to suggest a low redshift for the galaxy. If so, the seven wideband images (readily accessed through a "FitsMap" viewer from the JADES image-processing team), in particular, the F090W, F115W, and F150W images, were inspected in order to evaluate visual evidence for a Lyman break. From these ∼712 inspections, 225, 123, 103, and 38 galaxies were rejected. The remaining 88, 68, 59, and 8 galaxies were added to the "gold sample" based on judgment that at least half of these were at the higher redshift found by SEDz*. The count for the four redshift slices was then z1—534, z2—251, z3—174, and z4—23, a total of 982 galaxies.

Subsequent to the first release of this paper, we "pruned" the data sample further. An additional 72 galaxies were deleted from z1, 18 were deleted from z2 and 2 were moved from z2 to z3. About half had been removed after inspection for proximity to diffraction spikes, and the remaining half was split between marginal "S/N ≈ 5" cases and poor overall fits of the model to the data. The final numbers in this study, z1—462, z2—233, z3—176, and z4—23 (unchanged)–a total of 894 galaxies.

To expand on the introduction of the four SFH types we showed in Figure 1, we show 18 examples for each of four redshift ranges in Figures 2–5. In each, there is a row-by-row progression, starting from the top, from a single-burst population, through stochastic (multiple bursts), to contiguous epochs of star formation, and finally longer, continuous over four or more epochs.

SEDz* plots like the ones that follow, and their associated data, will be available for the full 894-galaxy sample at https://obs.carnegiescience.edu/SEDz-star/SFHs.

Beginning with Figure 2, we see three examples of a z = 6 starburst, each with a mass of ∼2 × 10⁸ M_⊙. As in the discussion for the burst SFH in Figure 1, a single epoch of star formation—both occurring and observed at z ∼ 6—is a very good fit to each observed SED. ²⁵ For all three, the solution is a combination of a burst and continuing CSF with less than a 50% contribution to the flux (denoted by the small green arrow). (A confirmation of the redshift found by SEDz* is the [O iii] emission in the top left case, detected in the medium-band filter F335M.) Again, the excellent fit found for such cases validates that the stellar population templates are correct for the task. Single bursts are the most common SFHs in our study—473 cases, 53% of the full sample.

Zoom In Zoom Out Reset image size

We also identify the next two examples (second row, left and center, hereafter 2l and 2c) as single bursts. Their two adjacent epochs of star formation may be unresolved, that is, one event that is best fit by using consecutive templates. In both 2l and 2c, the burst at z = 7 would have been very blue, but one epoch later—when these galaxies are observed—their contribution rises to the red, causing the elevation in each spectra. It seems reasonable that the star formation episode of τ ∼ 100 Myr could produce this SEDz* result in the relatively long "integer-redshift" epochs of z = 6 and 7. Under this interpretation of a single event, we include later examples where the two masses differ by as much as an order of magnitude, implying a sharply rising or falling episode of star formation.

The next examples, which exhibit well-separated bursts, are additional demonstrations of the ability of SEDz* to reveal dominant early star formation in galaxies observed ∼500 Myr later, giving confidence in its identification of SFHs that extend over the 6 < z < 12 range of our study. Comparing the stellar population templates shown in Appendix B, it is clear that these three examples show histories that cannot be reproduced by any single stellar population. Figure 2 shows three examples (2r, 3l, and 3c) of strong, well-separated bursts of star formation. We call these "stochastic" SFHs because it seems unlikely that these separate bursts can be "tied together," not even if there was weaker, undetected star formation in between (which would be at least an order of magnitude smaller in mass). The following SFH (3r) is interpreted as either three or four epochs of star formation starting at z = 10. All four are graphic illustrations of very young, blue stellar populations whose red flux has been "boosted" by older populations at z = 10–12. For all four of these stochastic SFHs, redshifts are again confirmed by [O iii] emission in the F335M band. Stochastic star formation, by our definition, is almost as common as the single bursts (in total epochs of star formation), adding 165 cases (19%) of two or more bursts. Together, "bursts" and "stochastic" make up 72% of the SFHs found here; clearly, this is a strong and, we believe, unexpected dominant mode of galaxy building in its beginnings.

However, the bottom three rows of Figure 2 remind us that a significant fraction of early galaxies are undergoing more orderly, gradual growth. In examples 4l, 4c, 4r, and 5l, three stellar populations are required to fit the observed SEDs. SFH3 star formation is usually "contiguous" but sometimes with a one-epoch gap, likely due to "noise" or the asymmetry in NNLS solutions—no negative star formation. (See Figure 3, the gap in 4c, but note the z = 9 star formation 4l that shows the opposite—a detection that is dubious.) SFH3 describes 21% of the full sample (25% of the z = 6 sample). All three have final star formation at z = 6 with added CSF components, hence the very flat purple CSF contributions (see Appendix B). Again, it is also clear that no amount of z = 6 or 7 star formation can produce the "rising to the red" in all but 4r, with its uncertain star formation at z = 8 and noisy SED. Although they look short, these last three epochs add up to half a billion years, or half of the time since the Big Bang, so these are a substantial departure from what seems to be the dominant SFH mode—≲100 Myr starbursts.

The final five examples, 5c, 5r, and row 6, are longer, more continuous histories (SFH4) that stretch over the full redshift range, with four to six epochs of star formation. Although two have a two-epoch gap that could have landed them in the "stochastic" category, for all, there is an orderly history of CSF or declining star formation. These are pronounced examples of young stellar populations with a history of very early star formation that accounts for their strongly "rising-to-the-red" SEDs. All five show z = 11 and 12 star formation; however, it is likely that this period of ∼100 Myr is not "resolved," so the SED could be reproduced by only one, with the combined mass. Although they represent only 3% of the z = 6 population and only 69 galaxies in the 894-galaxy sample (8%), the persistence of SFH4 over the full time range explored here suggests a different environment, one where galaxies can evolve more slowly and relatively undisturbed, for example, experiencing only minor mergers.

From this first set, we see the prevalent signature of early (z ≳ 10) star formation: an SED rising to the red, often accompanied by substantial star formation at OE. This is best seen by the different levels for the medium bands F335M and F410M, though sometimes only from the broad bands F356W and F444W. This raises the question of whether strong emission lines could be boosting the far red fluxes and mimicking early star formation. We show in Appendix C that this is not the case, the main point being that the CSF templates include active star formation, so that the level of the broad bands cannot be "raised" by emission—it is already included. ²⁶ The medium bands F335M and F410M (open circles; not used in the NNLS fit) provide good reference points for the presence of emission and the level in the continuum, because emission-line strengths were not calculated for the narrower bands. Appendix C provides more examples and discussion.

For further examples of SFHs, Figure 3 shows for z ∼ 7 a virtual replay of the z ∼ 6 population: five individual bursts (four single, one "twin"); four well-separated "stochastic" events; three more triples; three noisy SFH3/SFH4 histories; three long SFHs, none with a gap of more than one epoch; and two cases with star formation in all six epochs. (Some have large error bars that remind us that these are "representative" SFHs, not to be considered faithful "epoch by epoch.)" As in Figure 2 but not explicitly called out are very high S/N examples; in particular, examples 1l, 2l, 2r, 3r, and 4c are well-defined SED solutions that match the fluxes of the data with remarkable fidelity. Considering the relatively crude redshift resolution of the stellar population templates, the agreement of model and data is convincing evidence that the modeling works and the templates are fully "descriptive"; three or four templates are usually all that is needed, and that is good, because that is all that is available.

Zoom In Zoom Out Reset image size

Extending the redshift coverage up to z = 8–10 in Figure 4 shows that the classification into the four types of SFHs continues to earlier times. All are represented in similar frequency to those at lower observed redshift. The Lyman break is, of course, pronounced, which helps with constraining the SEDz* fiting. Most of the the bursts (single and double) in the top three rows rise sharply to the (rest-frame) ultraviolet, signaling very young populations, but conversely, the bottom two rows show the steepest rise to the red in the sample, in this case signaling strongly declining rates of star formation since z = 12.

Zoom In Zoom Out Reset image size

Finally, in Figure 5, we see 18 of the 23 galaxy sample of the highest redshifts. The Lyman breaks are clifflike, and the SEDs are all blue within this sample that covers only 150 Myr, only ∼200 Myr since z = 20 when the first "modern" galaxies were likely born. Rapid evolution of the Universe at this time, especially the strong growth of dark matter halos and rapidly decreasing density of large-scale structure, suggests star formation that might be itself changing in character or composition. Yet, remarkably, the stellar populations we observe are all matched by the stars that we have all around us today, some 13 billion yr later.

Zoom In Zoom Out Reset image size

Does the population of SFHs itself "evolve"? Figure 6 plots the SFH type for 894 different galaxies as measured at their OE, color-coded to show burst SFHs as red and orange (single and multiple) and contiguous "triples" and continuous "long SFHs," as green and blue. To first order, we see two things: (1) the "final" masses of all 894 galaxies in our SFH sample are largely confined to 10⁸–10¹⁰ M_⊙, and (2) all four types are represented over the full redshift range covered in this study. That there is some "accumulated" stellar mass below 10⁸ M_⊙ but none below 5 × 10⁷ M_⊙ is a simple detection limit: the fluxes would generally fall below the S/N > 5 limit we have chosen. On the other hand, the steep falloff in numbers above M > 3 × 10⁹ M_⊙ is probably a reflection of astrophysics, for example, stellar feedback suppression in the environment of rapid star formation in the compressed volumes of these subkiloparsec-sized sources. Powerful feedback is probably expected in the case of the bursts, but it is less obvious why the longer histories, with lifetimes of hundreds of Myr, would be subject to the same limitations, and yet the distribution here suggests, interestingly, that they are. This, and the clear result that most of the stellar mass in this formative time is made in relatively large starbursts, should place strong constraints on numerical simulations of galaxy growth.

Zoom In Zoom Out Reset image size

The "evolution" of the proportions of SFHs over time is harder to parameterize, but it is clear from the fact that all types show up over the diagram that there is no strong evolution of the population of SFHs. Figure 7 shows a crude graph of the frequency of the four SFH types for the sample in Figure 6. With the exception of the 70% fraction of bursts at the highest redshifts, ²⁷ the variation in the rest of the plot is factors of 2 and not much in the way of trends.

Zoom In Zoom Out Reset image size

We should add that there is some nonstatistical "bunchiness" in Figure 6, such as the "swells" of SFH3 "contiguous" at z ∼ 6.5 and SFH4 "long SFHs" at z ∼ 7.5 and the "valley" of single bursts at z ≈ 6.5, but we believe that these are more likely SEDz* "preferences" associated with interpolating between the epochs, rather than a real effect. A denser sampling of epochs, with twice the time resolution for the templates, would likely resolve this, but for now we focus on broad trends.

The lack of strong change over time (for this early period) provides little in the way of clues about the nature or "causes" of the different histories. This suggests that dependencies based on space, rather than time, might provide more insight. We investigate this in Section 8.

Returning to our theme, we want to explore how the results obtained in this study can help in understanding how galaxies began to grow and build up the essential elements for planets and life. Specifically, how has the progress of star formation proceeded in the first billion years of cosmic history, birthing and growing new galaxies in a environment more gas-rich and more turbulent than we easily imagine? At the same time, these young galaxies are likely fed continuously by a smooth inflow of gas, incorporating less massive galaxies and adding their moderately metal-enriched stellar populations. All this continues as young galaxies are vulnerable to violent major mergers and the huge energy released in massive starbursts and possible large-scale changes through black hole formation and growth. It is not likely that such questions will be answered solely through observations, but rather—as always—theory will be required to explore the physics of each of these elements. Numerical simulations should benefit greatly from these kinds of data, replacing previous and various speculations, where the manifestations of star formation that set the course are reproduced and, we hope, understood.

What we have found in this study already confirms the dynamism of the epoch where galaxies achieved masses of 10⁸–10¹⁰ M_⊙. We find large contributions to the growth of stellar populations by bursts unlike any we see today, strong enough to produce a stellar mass of 10⁸–10⁹ M_⊙ in an episode lasting only ≲100 Myr—little more than a dynamical time—and, apparently, strong enough to ward off further star formation for more than ∼500 Myr, and maybe even a Gyr. And yet, we also see common cases of multiple bursts over which the total stellar mass can reach well over 10¹⁰ M_⊙. What is the difference, then, between these and systems making similar amounts of stellar mass but over the same long period of time where bursting galaxies seem to go dormant?

A good way to appreciate the power of these data can be seen in Figures 8 and 9, where we plot the mass buildup over 6 < z < 12 from these different modes of star formation. These plots use the SEDz* SFHs for our 894-galaxy sample to graph the onset of star formation and, in the cases of longer SFHs, its subsequent addition of stellar mass. Here we have "connected the dots"—from when star formation began to the last epoch where it is detected—to show stellar mass buildup. (Color coding is the same as for previous figures.)

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Figure 8 shows this for (left) single bursts and bursty "stochastic" histories, SFH1 and SFH2, and (right) longer "contiguous" and "continuous" histories, SFH3 and SFH4. In an effort to provide guidance to numerical simulation modelers trying to answer such questions, we now express the data we have described here in terms of the growth of galaxies with such different histories in mass-buildup diagrams. The plot is for mass versus time, as in Figure 6, but now with tracks that connect the first and last epochs of recorded star formation within the 6 < z < 12 era. In Figure 8, SFH1 (a single burst) shows up as a dot, with a handle marking its epoch, while SFH2 ("stochastic") in orange, appears as a shallow rise, since the two bursts are usually well separated. It is also easy to see that, while the single bursts are the most common SFH type, multiple bursts are both larger in mass to begin with (we find mostly declining SFHs) and growing with subsequent bursts, such that they add the most mass over most of this era. A critical point, though, is that the mass from these burst-dominant galaxies is growing through the appearance of new bursts. The most frequent number of bursts for SFH2 galaxies is just two—a single added burst rather than several. Stellar mass is primarily growing by adding bursting objects, single and multiple, not by many smaller bursts in each galaxy. ²⁸ The same diagram for the three epochs "contiguous" (SFH3) and four epochs or more "continuous" (SFH4) histories shows more tilt in its tracks; this is most obvious in the three-epoch tracks that become more and more dominant from redshift z = 9 down.

Putting it all together in Figure 9, we combine the burst and longer histories, summing and integrating to learn how the stellar mass of this collection of 894 galaxies grew from z = 12 to z = 6. We see the contribution of single bursts grow strongly, while the average mass for SFH1—the lowest of the four types—does not increase over the full redshift range. By z ∼ 6, the single-burst, stochastic, and continuous histories have each added ∼10¹¹ M_⊙ to this volume of space. The SFH4 galaxies have contributed a factor of 3 less mass, but most of this is the result of their lower frequency of 8%, so this means that all four types contribute a similar amount of stellar mass from 6 < z < 12. Finally, we plot the integrated stellar mass for this special epoch, which reaches 4 × 10¹¹ M_⊙ in this volume at z = 6 and is growing at a rate of ∼500 M_⊙ Myr⁻¹.

We hope that both the rates and manner of star formation in these youngest of galaxies will provide the first meaningful constraints for numerical modeling studies of the evolution of the Universe at the end of the first billion years of cosmic history.

SFHs are known to have a strong spatial variance, in the sense that different kinds of galaxies (for example, ellipticals as opposed to spirals) are found to dominate in different environments (Dressler 1980), and the different SFHs of these have clearly differed greatly. Recognizing this, we looked at the spatial distributions to see if galaxy "environment" could be connected to the four types of SFHs. Figure 10 shows the distribution on the sky of the four SFHs (again with the same color coding). There is obvious large-scale structure in the z1, z2, andz3 maps, especially in the lowest redshift z ∼ 6 map, where the contrast between large voids and substantial clustering is strong. This is probably both a result of the growing clustering with epoch and because our much larger sample makes any contrast more discernible. Still, there is a strong impression of substantial large-scale structure over the period 6 < z < 10—covering most of the epoch of reionization. A map of only the 23 galaxies in the z4 sample is not useful, of course, but we have made another view of the z = 11–12 Universe that we discuss below.

Zoom In Zoom Out Reset image size

It is not surprising that there are no visual spatial distribution differences in the four SFHs, but it is reasonable to expect that correlations of SFH types with local density or nearest-neighbor distance might provide some insight into whether the environments of these galaxies influence their SFHs. In Figure 11 (Appendix A), we show histograms of local density and nearest-neighbor distance for the z ∼ 6, z ∼ 7, z ∼ 8, 9 redshift ranges (a.k.a., z1, z2, z3). The left three panels show density ranges that vary from a few to tens of galaxies arcmin⁻². As with Figure 6, the most notable feature of these diagrams are ups and downs likely associated with density fluctuations, with no clean separation by type. Perhaps there is a slight preference for the longer SFHs to be in denser regions or closer to their neighbors, but nothing clear enough to be helpful.

However, in the sky maps of Figure 10, the fourth panel may hold an important clue about the role of environment. Here we have done something new, based on the ability of SEDz* to identify prior star formation from galaxies observed at later epochs—galaxy SFHs. We made this map by choosing galaxies for which substantial star formation has been detected at z = 11 or 12, ²⁹ and in doing so, we selected 201 galaxies of our sample of 894 that had their first epoch of star formation at that time; we plot in the fourth panel what a sky with galaxies forming stars at z = 11–12 might have looked like.

Remarkably, we see a relatively tight galaxy group of 30 members in an area of less than ∼2 arcmin⁻². Eleven of those, 37%, are SFH4—long, continuous SFHs, only 8% of the full sample in this study. We claim that this detection of many long SFHs "spatially," while not statistically conclusive, is a unique data point in the search for an environmental dependence of SFH type. Furthermore, this distribution further confirms that SEDz* actually works; there is no other explanation for how these two quantities—position on the sky and SFH—could be well correlated. If the z = 11–12 star formation SEDz* detected were bogus, these galaxies should spread randomly across the field. ³⁰ The clustering, as well as what appears to be large-scale structure for the full sample of galaxies, confirms these "earliest" SFHs. In addition, above and below the circled area, the density remains high and appears to contain a higher fraction of SFH4 galaxies, compared to the much larger sample of the left side of the figure.

Subsequently, we noticed that this this higher galaxy density region also appears in the upper right panel, 6.75 < z < 7.75. In fact, this galaxy group had already been independently discovered by members of the JADES team (Endsley et al. 2023) at an observed redshift of z ∼ 7.5, further confirming their long SFHs. Even more members are observed at the later epoch.

This fortunate "feature" of the GOODS-S field called to mind that, for the original morphology–density relationship, Dressler (1980) and Postman & Geller (1984) found no gradient in morphological type beyond the rich-cluster boundaries, that is, no slowly changing populations beyond the effective radius of the cluster. This and other evidence convinced Dressler (see Section IV-C) that the origin of the morphology–density relation did not lie in the late-time evolution of the cluster or its galaxies. Rather, it pointed to an environmental dependence arising at the epoch of galaxy birth: the different morphological types were foretold by their birth environment, a kind of early nurture that is in fact nature.

Soon there will be more deep-field imaging like GOODS-S, and more groups will be found at—incredibly—z ∼ 12. This result, if confirmed by many other cases, might suggest a rather obvious conclusion about the different SFHs we have found here: burst histories are most common in the equivalent of the lower-density "field" of the modern Universe, probably the result of stochastic merger events that reflect the sparser environment. In this picture, long SHFs are destined for the richer, denser environments of the future. They will be galaxies that, unlike the "bursts," were built up in a more orderly series of accretion events and minor mergers.

Observations of the early Universe with the incomparable JWST have already had a profound effect on our ideas about how the first galaxies were born and how they grew. In scarcely more than a year, hundreds of studies have feasted on galaxy samples that are large, deep, various, and multiplexed. They have been probed with diverse modes of cameras and spectrographs that are unprecedented in a space telescope, offering factors of hundreds to thousands greater sensitivity in the near-to-mid-IR: the early Universe is viewed anew.

Our study is but one of multitudes focused on the rise of stars and galaxies, the elements that redefined the Universe at the very epoch we now explore, for the first time. So far, this wealth of JWST data tells a new story: that the birth of galaxies appears to be a volatile affair. Specifically, a central theme has been that "starbursts" are prevalent in this first billion years, or at the very least, that they play a critical role in understanding what see. Before comparing the present work to a small sample of these papers, it is important to reiterate why—for a majority of galaxies in our GOODS-S field—we identify their SFHs as "bursts." Our study has particular significance because it was the first to trace the buildup of the stellar mass in the first galaxies (Dressler et al. 2023), and because it is the only one to rely solely on the light of main-sequence A stars. The mass-to-light ratios of these stars are thus accurately known and independent of chemical composition. Furthermore, a stellar population with an age of between 100 and 1000 Myr is not a dusty one; even while O and B stars are still forming in H ii regions, A stars have already migrated into near-dustless environments.

Our study of 894 newborn galaxies at 6 < z < 12 has shown that their star formation was predominately in bursts; we have called these histories "starbursts" (single) and "stochastic" (multiple). Because of the methodology of this study and the limitations of SEDs with seven broadband fluxes, we actually know little about the nature of these bursts that we find dominating the SFHs of the first galaxies. For these, we observe primarily, almost exclusively, the light of A stars whose ages range between ∼50 and 150 Myr. It is likely that the mass of stars born during the full epoch was much less than that of those born in the first ≲50 Myr of explosive star formation that preceded it. In other words, most of the observed A stars came from a shorter period that started early in the epoch. By this reasoning, we suggest that SFRs of 10–30 M_⊙ yr⁻¹ characterized the dust-enshrouded phase we have imagined, and it is during this phase that the ∼50–120 Myr old A stars that define our SEDs were born. We can estimate the SFRs of the CSF at OE and imagine them to be typical of the SFR over the whole epoch; that might only be 1–10 M_⊙ yr⁻¹ (see, for example, Emami et al. 2019; Faisst et al. 2019;  Rezaee et al. 2023), but we will be hard-pressed to make an SED measurement, or any other we can think of, that can recover that history of the SFR with a resolution of tens of megayears for any given object. ³¹

It is more accurate, then, to say that our study has concluded that starbursts dominate not because we see signatures of 10⁸–10⁹ M_⊙ of stars forming over an epoch but because we see a population of A stars in that epoch that are limited to ages of less than ≲150 Myr. Most likely, they were born in a much shorter period early in that epoch. What defines our result that "starbursts dominate" is that we see no A stars older than, say, 150 Myr, that is, from a prior epoch. That is why we call these SFHs bursts, not because we are seeing the burst (although we may be, in some cases), but because the only stars we do see are young.

We think that the principal benefit of our observations of the prominent role of starbursts will be to help inform theory and numerical simulation models about the growth of the baryonic component of the Universe that winds up as galaxies. The prevalence of bursts of star formation in the first billion years should be influential in guiding theoretical work to understand galaxy growth in a dynamic environment. It also seems that there are environmental factors that might send galaxies down one path rather than another that could help constrain numerical simulations.

Early results from JWST have focused on a large population of bright galaxies, particularly at z > 10, e.g., Donnan et al. (2023), Harikane et al. (2024), and Finkelstein et al. (2023), suggesting that the luminosity function of the first galaxies evolved less slowly with redshift. ³² Observations of galaxies that are brighter than expected, more numerous than expected, or evolving less rapidly than expected have promoted "burstiness" as a way to reconcile these findings (e.g., Looser et al. 2023) with previous models (e.g., Wilkins et al. 2022). However, those galaxies are higher redshift and/or rarer (and more massive?) than those studied here, so our sample adds little to that discussion.

On the other hand, Sun et al. (2023a, 2023b) have in particular suggested that lower-mass galaxies, with an abundance and mass predicted with pre-JWST models, have been elevated in luminosity by bursts, specifically that their light-to-mass ratios have risen through substantial bursts of star formation. While again, those samples are not comparable our own, we can by analogy question this explanation. Our study finds a high level of "burstiness" through measurements of stellar mass—the mass of A stars from the SEDs. Likewise, the z > 10 "bright" examples may be cases of galaxies forming more stellar mass in bursts.

In a study that also seeks to measure a change in SFHs over 6 < z < 12, Ciesla et al. (2023) use a very different method than our own to extract such information from SEDs. Their analysis correlates SED shape with position on the star-forming main sequence to characterize whether a galaxy's SFH is "stochastic," like our definition—bursty over multiple epochs—or "secular"—longer and steadier, like our SFH3 and SFH4 types. They posit a smooth transition from mostly stochastic to mostly secular at around z ∼ 9. Although our sample is relatively small beyond z > 9, Figure 6 shows no clear transition at z > 9, with a similar fraction of SFH1+SFH2 and SFH3+SFH4 down to z ∼ 7. There seems to be a strong upturn in both SFH1 and SFH3 below z ∼ 7 (see Figure 7 with the ratio between "stochastic," SFH1, and "secular," SFH3, basically unchanged). Again, our methodology is more straightforward—SEDz* actually derives SFHs from ages of A-star populations, but this important issue requires further study.

About this important matter, what seems clear is that most of the youngest galaxies do not grow steadily in a calm, peaceful environment. Rather, their journey to what we today regard as "galaxy-sized" may be more chaotic, or even explosive. Perhaps these are gas-rich mergers dominating in lower-density regions, manifesting in "one- or two-event" growth spurts in the first billion years, reaching ∼10⁹ M_⊙(what we now call galaxy-sized) when this first phase completes. Importantly, however, by z ∼ 2—"cosmic noon"—these galaxies should grow by an order of magnitude in mass to reach a halfway point for L*. It will be important, and challenging, to relate that later, more easily observed Universe to what we have witnessed for the earliest galaxies. These starburst galaxies at 6 < z < 12 seem neither poised for another burst nor prepared to settle into steadier SFHs. The simple fact that most of our 894 galaxies do not show star formation over many epochs precludes the notion that these bursty objects are "picking up" again before z ∼ 6. Perhaps the "contiguous" SFH3 galaxies, common in our z = 6 sample, are destined to become those L* galaxies. We note that almost all show the signatures of continuing star formation at z = 6, so perhaps they will grow into the most common galaxies of today, just beginning in earnest their journey to maturity. It is undeniable that a transition from "stochastic" to "secular" (to use Ciesla et al.'s term) did occur, so a focused study on 4 < z < 6 galaxies, by whatever means their SFHs can be characterized, seems a priority, and imperative for Chapter 3.

Finally, pursuing our suggestion that the starbursts we have observed were the consequence of intense, extreme-density, dust-enshrouded bursts of star formation, we hope an effort can be made to match up this population of growing galaxies to a coeval population of heavily dust-obscured protogalaxies that will come from JWST/MIRI and perhaps ALMA. That our sample seems effectively dust-free strongly suggests that there is a "just-before" phase of tens of millions of years in which explosive star formation ignited in gas-rich protogalaxies, leading directly to the objects we have been studying. Forming stars rapidly in these subkiloparsec volumes is a formidable challenge to our considerable knowledge of how stars form. In particular, the feedback from the supernovae of such starbursts should blow things apart.

We have taken note, with great interest, of work suggesting that such feedback can be suppressed if, at that time, supernovae of mass M ≳ 10 M_⊙ collapsed directly into black holes. Renzini et al. (2022) have suggested that this extraordinary explanation is nearly "mandatory" to explain the multiple stellar populations observed for many globular clusters, traditionally thought to be limited to a single generation by supernova-driven feedback. It appears that the same challenge applies toward our understanding of the earliest galaxies. We find first "bursts" of star formation with stellar masses ≳10⁸ M_⊙—a scaled-up version of the largest globular clusters—similarly in need of a physics miracle (see Dekel et al. 2023) to avoid destruction before fulfillment.

Travelers on JWST, the ultimate spacetime ship, we seek answers to questions first asked a century ago. Wonderfully, we are fortunate to be grasping this "once-in-a-species" opportunity to know our origins.

A.D. gratefully acknowledges the support of the NIRCam team for the opportunity to contribute to the NIRCam science program on the earliest galaxies and for prior funding from JWST/NIRCam contract NAS5-02015. The DC2 "deep-field" simulation was a Herculean task that provided the essential guidance for the program described here and for the NIRCam GTO program in general.

This work is based on observations made with the NASA/ESA/CSA James Webb Space Telescope, using imaging data processed by the NIRCam project team in collaboration with JADES, a extragalactic project team consisting of members of the NIRCam and NIRSpec teams. D.S., K.M., M.R., and A.D. have been supported by JWST/NIRCam contract to the University of Arizona NAS5-02015.

N.B. acknowledges support from the Cosmic Dawn Center (DAWN), funded by the Danish National Research Foundation under grant No. 140. K.B. acknowledges partial support by the Australian Research Council Centre of Excellence for All Sky Astrophysics in 3 Dimensions (ASTRO 3D) through project No. CE170100013. A.B. acknowledges funding from the "FirstGalaxies" Advanced Grant from the European Research Council (ERC) under the European Union's Horizon 2020 research and innovation program (grant agreement No. 789056). S.C. acknowledges support by European Union HE ERC starting grant No. 101040227—WINGS. D.E. is supported as a Simons Investigator and by (JWST/NIRCam) contract to the University of Arizona NAS5-02015. R.H. acknowledges funding for this research from the Johns Hopkins University, Institute for Data Intensive Engineering and Science (IDIES).

We use data from the JADES High Level Science Product at MAST (Rieke et al. 2023).

The histograms of Figure 11 show distributions of the four SFH types with local density calculated using the area of the 10 nearest-neighbor galaxies and distance to the single nearest-neighbor galaxy. The absence of any obvious trends between SFH type and environment seems at odds with results for environmental dependencies later in cosmic history, for example, the obviously difference in the SFHs of rapidly growing elliptical galaxies compared to the drawn-out SFHs of spirals like the Milky Way. However, from another perspective, the lack of environmental dependence is reminiscent of studies of the environment around rich clusters in the early 1980s (Dressler 1980; Postman & Geller 1984) that found that a steep dependence of galaxy morphology with local density within the effective radius of a cluster did not continue to the lower-density "field" beyond. In this study, Figure 10 shows a hint of the same behavior; in the region of higher density, a z ∼ 7 rich group, the long-SFH types, are strongly represented, a correlation that is not expressed in lower-density surroundings.

Zoom In Zoom Out Reset image size

In this section, we show samples of the stellar population templates used by SEDz* to characterize the SFHs of 5 < z < 12 galaxies. Figure 12 plots the fluxes of 10 stellar population templates with a 10 Myr burst of star formation (at 1 M_⊙ yr⁻¹ = 10⁷ M_⊙) at the start of epoch z = 12, "observed" to evolve at epochs z = 11...3—later epochs without further star formation. For this study, SEDz* uses seven templates, for bursts starting at z = 11, 10...6. The principal feature of this plot is that the templates are largely nonconformal, that is, not a conformal set of curves scaled by some parameter or set of parameters. These "vectors" describing stellar populations are different enough—sufficiently orthonormal—that a least-squares combination of them is substantially "resolved" from any other combination. This property allows SEDz* to "calculate" the history of a stellar population, essentially, by vector algebra: finding the coefficients of the vector sum that best represent the observed SED. ³³ What makes this particular application of the method potent is the nonconformal character of SEDs for stellar populations of ages τ < 1 Gyr—the templates covering the early Universe for z = 12 to z = 5 whose light is dominated by main-sequence A stars. The figure shows why, as has been known for half a century, only finding the ages of stellar populations with stars older than 2 Gyr is, in practice, impossible; note how the templates z = 5, 4, and 3 are becoming a simple scaling of a single shape, as the Universe reaches an age of 2 Gyr at z ∼ 3.

Zoom In Zoom Out Reset image size

The signature of a burst of star formation is a very blue SED at that epoch, but for the subsequent epochs, the history of star formation within that epoch is unresolved. Therefore, SEDz* accumulates the sum of bursts as the stellar masses of each; this is indistinguishable from continuous star formation over the prior epochs.

A final point of note is that SEDz* works because it is strongly constrained by the shape of each template, which means that variation of the coefficients in the maximum-likelihood solution cannot either make or break the fit. If these templates, made from present-epoch stars in our Galaxy, were not representative of stellar populations at z > 5, this attempt to reproduce observed SEDs would fail badly. Quite the opposite is true.

The top-right panel of Figure 13 shows templates for six templates, from z = 10 to z = 5, corresponding to the flux resulting from 1.0 M_⊙ yr⁻¹ of CSF over that particular epoch. Unlike bursts, there is no evolution of stellar population over subsequent epochs, because ongoing star formation can only be recorded in the SED from the OE. Prior epochs of CSF are indistinguishable from bursts of the same mass. With the time resolution offered by broadband SEDs, no additional information is available; there is no signature to distinguish CSF from a more complex behavior over the ∼100 Myr duration of each epoch.

Zoom In Zoom Out Reset image size

Clearly, the distinguishing feature of the six CSF templates is that they are all flat—as conformal as it gets—in comparison to the burst templates. The modulation that is apparent comes from the Balmer break—moving from ∼2 μm at z ∼ 6 to ∼3.5 μm at z ∼ 10, and from the jaggedness of the SED from 3 to 5 μm—due to [O iii] and Hα emission lines. This flatness, when combined with bursts of previous epochs, is responsible for much of the "character" of the long SFHs. The right plot of Figure 13 shows that combining a burst early in an epoch with CSF at that same epoch produces signatures that are found in many z > 5 SEDs. The ratios of 3:1 to 1:2 for the burst/CSF flux, shown in this case for z ∼ 7, are apparent in hundreds of the SEDs in our sample. It is worth remembering that this combination is just equivalent to ongoing star formation that is declining, rather than constant, over the epoch.

The bottom plot shows the attenuation in the SED expected for different values of A_V , as described in Section 2.

In general, the influence of emission lines in photometric studies using broad bands, like this one, are only significant if the line-to-continuum ratio is very large. For example, a moderate-resolution spectrograph (with resolutions of hundreds) is effective at detecting even weak lines only when the continuum level is low, and that is true only for young stellar populations, τ ≤ 20 Myr. In contrast, broad photometric bands like the ones used in this study select against such populations, because the continuum flux is weak (the objects are faint).

Because our detections of faint galaxies rely on the sensitivity to small fluxes, the equivalent width of an emission line must be enormous, that is, a low continuum flux. Since the "integer epochs" of this study cover ∼100 Myr of cosmic history, any object selected through broadband photometry will by necessity require a large burst of star formation compared to the stellar mass generated over the epoch. Thus, the contribution to the flux of the youngest populations, through emission lines, should be modest.

Figure 14 verifies that by selecting cases of relatively strong star emission lines (the right-side examples) and comparing to SEDs on the left with little evidence of star formation. Our SEDs are made up of seven broad bands (filled circles) and the two medium bands, F335M and F410M (open circles; not used in the SEDz* fitting). In particular, these and the F356W and F444W broad bands are sampling the [O iii] and Hα lines over the redshift range 6 < z < 9.

Zoom In Zoom Out Reset image size

The medium bands often provide good evidence of emission because their width is ∼40% of the broad bands, so they can also—in cases that exhibit emission—establish the continuum level at that color, sometimes considerably below that of the broad bands. ³⁴ The close agreement in flux of both broad and narrow bands in the top and bottom left examples shows that there is little or no detected emission in these cases, which is common in our 894-galaxy sample. Note how the levels of F356W and F444W are close. The middle left SED shows moderate emission in F444W and detection of the continuum in F335M, and here F444W appears to be elevated above F356W by [O iii] emission. The case in the top right is a stronger example: strong detected emission in F356W and F335W, also from [O iii], at this lower redshift.

When compared to the flat blue (CSF) templates on the left and purple on the right (excluding middle right), these SEDs all show a substantial rise in flux over the F356–F444 region. What is the contribution of emission lines to this rise?

For most of our sample, the answer is little to none. This is because the CSF templates include the emission-line fluxes for the appropriate stellar population for both the broad and narrow bands. As explained in Stark et al. (2013), the SEDs of these templates include hydrogen lines (Robertson et al. 2010; based on Osterbrock tabulations and case B recombination), metallic lines (Z = 0.2 Z_⊙, likely appropriate for our high-z galaxies), and continuum radiation, calculated and described by Stark et al. (2013). Of course, the strength of these lines should vary from object to object, and perhaps systematically from z ∼ 0 to z > 5, but examination of dozens of our sample suggests that such variations, though present, are smaller than the effect of "line versus no-line." However, the existence of variation, as well as not knowing the redshift well enough to place the emission lines accurately, prevents us from making detailed arguments about whether these low-z line ratios are in fact a good match for very young galaxies.

Knowing that the emission lines are appropriately included in the templates we use, is it easier to understand the examples we show here. Since emission is included, the F356W and F444W fluxes of the top and bottom left SEDs show no emission; indeed, if their "elevation" from a flat SED were due to emission, the narrow bands would show higher, not at the same level. The middle left example does show elevation of F444W compared to its more typical closeness to the level of F356W, but the SEDz* model passes through both points because emission is included (this also indicates that the present-epoch templates are appropriate). The same goes for the top right example with strong emission in F335M. In this case, and most others, the continuum level of F410M, compared to the elevated value of F356W, comes in at the proper ratio of 2.5:1.

In the middle and bottom right examples, we again see the SEDz* model passing through F356W and F444W, and the emission—here [O iii] and Hα —is included. The continuum below is sampled in each by both medium bands. Here we also see that the continuum level is provided by a single, older, z = 11 population, which is boosting the level of F356W and F444 to give the SED its distinctive shape, one that is inconsistent with any single population. For the bottom right example, the continuum is matched by z = 11–12 and z = 8 flux to complete the fit. The only way emission could contribute significantly would be if [O iii] and Hα emission were substantially stronger in these regions of star formation in z > 5 galaxies than in the galaxies of today.

To look for galaxies with stronger-than-present-day emission, we examined SEDs to find cases where the points indicating emission had fluxes systematically higher than the values of the five bands F090W, F115W, F150W, F200W, and F277W. There were few among the 690 cases of the z1 and z2 samples, z < 7.75, but 50 of the 128 galaxies in the redshift range 7.75 < z <9.75 showed conspicuous [O iii] emission in the F444 band. Among these, 28 had SEDs that rose steadily through F356W, which SEDz* attributed to earlier star formation, z ≥ 10. For the remaining 22 the lack of such a clear trend suggested that these were, in fact, cases of [O iii] emission much stronger than represented by the templates. For such cases the "earlier" stellar mass found by SEDz* was incorrect. By examining the z ≥ 10 contribution, we determined that for such galaxies total stellar mass had been overestimated by factors of 1.5 to 3.5. This is important for future SEDz studies. At the same time, it does not alter the conclusions of this study, in fact, it moves some SFH2 types (an early and late burst) back to the SFH1 category, reaffirming our principal result of the prominence of starbursts in the early universe.

Since these strong emission cases are a small part of the study, this result, added to the cases with little or no emission, confirms that "red rise" we find for ∼100 galaxies is flux from older stellar populations, in other words, SEDz* star formation histories confirmed.