CEERS Key Paper. III. The Diversity of Galaxy Structure and Morphology at z = 3–9 with JWST

CEERS (S. L. Finkelstein et al. 2022, in preparation) is an Early Release Science program (proposal ID 1345) that will observe the Extended Groth Strip (EGS; Davis et al. 2007) extragalactic deep field (one of the five CANDELS fields; Grogin et al. 2011; Koekemoer et al. 2011) early in Cycle 1 with data made available to the public immediately. CEERS will obtain observations in several different modes with JWST, including a mosaic of 10 pointings with NIRCam, NIRSpec multi-object spectroscopic observations in parallel for six pointings, and six pointings with MIRI in parallel. The NIRCam imaging of CEERS will cover a total of 100 arcmin² with the F115W, F150W, F200W, F277W, F356W, F410M, and F444W filters down to a 5σ depth ranging from 28.8 to 29.7 (for a typical total exposure time of 2835 s filter^–1). The first set of CEERS observations were taken on 2022 June 21 in four pointings, hereafter referred to as CEERS1, CEERS2, CEERS3, and CEERS6.

We performed an initial reduction of the NIRCam images in all four pointings using version 1.7.2 of the JWST Calibration Pipeline ⁴⁸ with some custom modifications. We used the current set of NIRCam reference files, ⁴⁹ though we note that the majority were created preflight, including the flats. For details on the reduction steps, see Bagley et al. (2022). To summarize briefly, we applied detector-level corrections using stage 1 of the pipeline with default parameters and used custom scripts to remove 1/f noise, wisps, and snowballs from the count-rate maps. We then processed the cleaned count-rate maps with stage 2 of the pipeline and used a custom version of the TweakReg step to calibrate the astrometry. We then coadded the images using the drizzle algorithm with an inverse variance map weighting (Casertano et al. 2000; Fruchter & Hook 2002) in the Resample step in the pipeline. The final mosaics for each pointing in all filters have pixel scales of 003 pixel⁻¹. We then used a custom script to background subtract the images. These final background-subtracted images were used for the morphology measurements described in Section 3.

For the analysis in this paper, before updated catalogs incorporating JWST photometry are available, we use existing CANDELS v2 redshifts and stellar masses for the HST F160W-selected galaxies in the EGS field. Here we provide a brief summary of these measurements; for full details of the photometric redshift measurements and resulting catalogs, see Kodra et al. (2022).

The v2 photometric redshifts and stellar masses are based on the CANDELS EGS photometric catalog of Stefanon et al. (2017), which includes broadband TFIT (Laidler et al. 2007) photometry from the UV to IR imaging from both ground- and space-based telescopes. For this paper, we adopt the z_best column from the Kodra et al. (2022) catalog, which provides the overall best estimate of the redshift. This corresponds to the secure spectroscopic redshift if one is available or the mFDa4_z_weight photometric redshift otherwise. We use this photometric redshift value because it produces the most accurate confidence intervals (see Kodra et al. 2022 for further details).

We then use the Stefanon et al. (2017) photometry and the above estimated redshift to determine stellar masses for galaxies in the CANDELS EGS field. This was done using two different codes: Dense Basis ⁵⁰ and P12. Dense Basis (Iyer & Gawiser 2017; Iyer et al. 2019) is a python-based code that uses flexible stellar population synthesis to generate model spectra that correspond to a wide range of physically motivated nonparametric star formation histories, metallicities, and dust attenuation values. P12 (Pacifici et al. 2012) is a Fortran-based code that uses a Bayesian fitting algorithm and model spectral energy distributions (SEDs) generated using simple stellar population models and the reprocessing of these models using the photoionization code CLOUDY (Ferland et al. 1998). For more details on each of these two methods and detailed comparisons, see Pacifici et al. (2012), Iyer & Gawiser (2017), Iyer et al. (2019), and Pacifici et al. (2022).

We find that the stellar masses using these two methods agree well with one another, with the expected level of scatter (e.g., Pacifici et al. 2022). For this paper, we use the mean of the two measurements. We note here the caveats that these measurements are highly uncertain at the highest redshifts (z > 6) and the faintest magnitudes (F160W > 26) because they rely on the HST ACS and WFC3 photometry that does not trace the rest-frame optical light at these redshifts. Future work will improve on these measurements with the addition of JWST NIRCam fluxes.

In addition to the JWST CEERS images described above, we also use the existing CANDELS ACS and WFC3 images to compare the HST and JWST morphologies. Here we use the mosaics produced by the CEERS team ⁵¹ with updated astrometry tied to Gaia-EDR3 and a pixel scale of 003 pixel⁻¹.

For the sample analyzed in this paper, we select all galaxies with z_best > 3 within the CEERS1, CEERS2, CEERS3, and CEERS6 NIRCam short-wavelength (SW) or long-wavelength (LW) footprints. We cross-match these sources with those identified in NIRCam imaging using the Source Extractor segmentation map (see Section 3.1 for more details). We remove any spurious sources in the NIRCam imaging, such as those that result from the diffraction spikes of stars and those that are so close to the edge of the images that we cannot obtain reliable morphology measurements. This results in a total sample of 850 sources with a detection in any of the NIRCam filters. We note here that at the magnitude limit of the CANDELS WFC3 images (F160W < 27.6, 5σ), almost all of these galaxies have signal-to-noise ratios (S/Ns) in the NIRCam images high enough to enable morphological classifications, as discussed further in Section 3.

Figure 1 shows the redshift and mass distribution for this sample of objects. We note that 10 of these objects have existing spectroscopic redshifts from either the MOSDEF (Kriek et al. 2015) or DEEP2 (Newman et al. 2013) spectroscopic surveys; we use these spectroscopic redshifts for these 10 objects. Overall, this sample peaks at z ∼ 3 and has a long tail beyond z ∼ 6 out to z ∼ 9. There are a total of 40 sources in this sample with a photometric redshift of z > 6. We note here that the redshifts and stellar masses for these sources are uncertain, and some of these may turn out to be at lower redshift. The redshift and stellar mass estimates will be improved with the addition of JWST data to the SED modeling in the future.

Zoom In Zoom Out Reset image size

Galaxies were detected in the NIRCam images using Source Extractor ⁵² version 2.25.0 (Bertin & Arnouts 1996). The setup was optimized to detect the HST-selected z > 3 galaxies without over-deblending. We created empirical point-spread functions (PSFs) for each filter by stacking stars, and then the F115W, F150W, F200W, and F277W images were PSF-matched to the F356W image using the Python-based code PyPHER ⁵³ (Boucaud et al. 2016). In order to detect galaxies that may be very faint in some of the NIRCam filters, we used an inverse variance weighted combination of the PSF-matched F150W, F200W, F277W, and F356W images as the detection image. We first ran Source Extractor in "cold" mode to deblend nearby galaxies, then in "hot" mode to detect faint galaxies, following Stefanon et al. (2017). We then combined the "cold" and "hot" detections in order to keep all objects that were detected by at least one mode, using a factor of 2.5 to enlarge the cold isophotes. We visually inspected the segmentation map for the z > 3 sources and optimized the parameters to ensure that the sources were detected and adequately deblended from nearby neighbors without being shredded. We use the final segmentation map produced from this process for all of the measurements presented in this section.

Each of the 850 galaxies in our z > 3 sample were classified by three different people from among a total of 35 members of our team. We used the Zooniverse project builder ⁵⁴ to host images of each galaxy and designed a workflow of five tasks based on a modified version of the classification scheme of Kartaltepe et al. (2015). These five tasks ask classifiers to select options for the galaxy's main morphology class, interaction class, and various structural and quality flags, and finally to leave any specific comments about a particular object. In this paper, we focus on the first of these five tasks, the main morphology classification, which roughly corresponds to the typical Hubble type classification. The options for the main morphology class are (1) disk, (2) spheroid, (3) irregular/peculiar, (4) point source/unresolved, and (5) unclassifiable/junk. To reflect the overall complexity seen in high-redshift galaxies, these classes are not mutually exclusive, so a classifier can choose multiple options to best reflect the overall morphology of the galaxy. For example, a galaxy can have both a disk and a spheroid in the case where it is a disk galaxy with a bulge component. A galaxy can be both a disk and irregular if, for example, it is an asymmetric disk or a disk involved in an interaction. The exception is that if a galaxy is "unclassifiable" then it cannot also be one of the other classes. This level of complexity can make the interpretation of the various classes challenging, but the extra information provides an important level of nuance to the classifications.

Classifiers are presented with a collection of postage stamps for each galaxy being classified. These stamps include each of the JWST NIRCam filters (F115W, F150W, F200W, F277W, F356W, F410M, and F444W) at their native resolution, with an asinh scaling to bring out low surface brightness features, an RGB stamp made up of the filters that correspond to the rest-frame optical, a version of that stamp zoomed out by a factor of 2, the NIRCam Source Extractor segmentation map described above, three HST ACS/WFC3 filters (F814W, F125W, and F160W, also with an asinh scaling), and finally, an RGB stamp of these three HST filters and a similarly zoomed-out version. The stamps are scaled by the size of the galaxy as measured by Source Extractor, following Equations (2) and (3) of Häussler et al. (2007), with a minimum size of 100 × 100 pixels. An example set of stamps for one of the galaxies is shown in the Appendix.

The classifiers are asked to make a holistic judgment about the overall morphology of the galaxy, taking information across the full wavelength range into account. In a separate task, the classifiers can select flags to indicate that the morphology changes across the NIRCam filters or differs between JWST and HST images.

We perform parametric fits on the NIRCam images using both Galfit ⁵⁵ (Peng et al. 2002, 2010) and GalfitM ⁵⁶ (Häußler et al. 2013). Galfit is a least-squares fitting algorithm that finds the optimum Sérsic fit to a galaxy's light profile, and GalfitM is a modified version that uses images at different wavelengths to allow one to constrain the fit parameters that vary smoothly as a function of wavelength. The benefit of using GalfitM is that it fits all bands simultaneously and produces more physically consistent models. We performed fits using both codes to test for self-consistency, but since we use the rest-frame optical fit throughout this paper, we focus here on the GalfitM fits and describe the Galfit fits in the Appendix.

GalfitM fits were performed using the IDL program Galapagos-2 from the MegaMorph Project ⁵⁷ (Bamford et al. 2011; Häußler et al. 2013, 2022; Vika et al. 2013). Galapagos-2 is a wrapper that enables GalfitM to be run over larger survey images. We used the Source Extractor setup described above with Galapagos-2. As input, we used all seven NIRCam filters (F115W, F150W, F200W, F277W, F356W, F410M, and F444W) and used the NIRCam WHT images produced by the JWST pipeline to create rms images to be used as the input sigma image. We used the F200W Source Extractor catalog for initial guesses and used the NIRCam empirical PSFs. In addition to the final output catalog with the Sérsic fit parameters, Galapagos-2 also outputs the original stamp, the GalfitM model, and the residual image for each galaxy in each filter. Out of the 850 z > 3 galaxies in our sample, 37 (4%) were flagged because GalfitM reached one of the constraint limits in one of the filters.

We measure nonparametric morphologies using the Python package Statmorph ⁵⁸ (Rodriguez-Gomez et al. 2019). For each NIRCam filter, we create 100 × 100 pixel cutouts of the 850 galaxies in our z > 3 sample to use as input to Statmorph, along with a cutout of the segmentation map, the empirical PSF, and the gain.

Statmorph measures a wide range of morphology statistics commonly used in astrophysics. The ones that we use for the analysis in this paper are concentration (C), asymmetry (A), and clumpiness/smoothness (S; Bershady et al.2000; Conselice et al. 2000; Conselice 2003); the Gini coefficient (G) and the second moment of the region of the galaxy containing 20% of the total flux (M₂₀; Abraham et al. 2003; Lotz et al. 2004); the Gini–M₂₀ bulge and merger statistics (Rodriguez-Gomez et al. 2019); the S/N per pixel; and quality flags.

As for the parametric fits, we use the fits for the filter corresponding to the rest-frame optical emission at the redshift of the galaxy: F277W for galaxies at z = 3.0–4.0, F356W for galaxies at z = 4.0–4.5, and F444W for galaxies at z > 4.5. Of the 850 galaxies fit, 81% have a S/N per pixel of >2.5 in the corresponding rest-frame optical filter; below this value, the fit results may not be reliable (Lotz et al. 2006). We compare these commonly used measures of galaxy morphology to our visual classifications in Section 4.3.

For each of the 850 galaxies in our z > 3 sample, we assign a galaxy a given visual classification if two out of three people select a given option as the main morphological class. There is only one object in our sample for which all three classifiers disagree, meaning one selected only "disk," one selected only "spheroid," and one selected only "irregular." This object is therefore not included in any of the figures presented here. As noted above, since the main morphological classes are not mutually exclusive, various combinations are possible. Throughout this paper, we break things down into the following nonexclusive morphological groups (highlighted in Figure 2).

• 1.  
  Galaxies with disks. The disk category contains galaxies classified as disk-only (without a spheroid or irregular component), disk+spheroid (a galaxy with both a disk and spheroid component; a separate structural flag indicates whether the disk or the bulge is dominant), disk+irregular (a disk galaxy with irregularities such as asymmetries, a warp, or disturbance by a nearby companion), and disk+spheroid+irregular (a disk galaxy with a spheroid component that also has some irregularities; note that these are fairly rare). When we refer to "all disks," we are referring to the sum of the galaxies in all of these categories.
• 2.  
  Galaxies with spheroids. The spheroid category contains galaxies classified as spheroid-only (without a disk or irregular component), spheroid+disk (same as disk+spheroid above), spheroid+irregular (a spheroid galaxy with irregularities such as asymmetries or surrounding low surface brightness features), and spheroid+disk+irregular (same as disk+spheroid+irregular above). When we refer to "all spheroids," we are referring to the sum of the galaxies in all of these categories.
• 3.  
  Galaxies with irregular features. The irregular category contains galaxies classified as irregular-only (no discernible disk or spheroid component), irregular+disk (same as disk+irregular above), irregular+spheroid (same as spheroid+irregular above), and irregular+disk+spheroid (same as disk+spheroid+irregular above). When we refer to "all irregulars," we are referring to the sum of the galaxies in all of these categories. Note that the irregular category may include merging or interacting systems, as well as galaxies that are irregular for other reasons, such as clumpy star formation. Mergers and interactions themselves will be discussed in a future paper (C. Rose, et al. 2022, in preparation).

Zoom In Zoom Out Reset image size

Our sample contains the full range of morphological types across all redshift and stellar masses. Over the entire redshift range, only 16 and 18 galaxies are classified as point source/unresolved or unclassifiable, respectively. Figure 3 shows the fraction of the total number of galaxies that each morphological class makes up as a function of redshift. For a fair comparison across redshifts, we limit this to the 666 galaxies with stellar masses greater than 10⁹ M_⊙, since the galaxies with lower stellar masses are only present at the low-redshift end of our sample (see Figure 1). We emphasize that this represents a mass-selected sample but not a mass-complete sample, as there are likely to be many galaxies identified by JWST in this mass range that are undetected by CANDELS HST imaging.

Zoom In Zoom Out Reset image size

Overall, 56% of the galaxies above this mass cut at z > 3 have a visually identifiable disk component, dropping from ∼60% at z = 3–4, to ∼45% at z ∼ 5, to ∼30% at z > 6. Again, we note that the photometric redshifts at z > 6 are uncertain, so some fraction of these sources may actually be at lower redshift; we caution the reader when interpreting the results at z > 6. The disk-only and disk+irregular groups each make up ∼20% (and slightly less at z > 6), while the disk+spheroid group makes up ∼10%, and the disk+spheroid+irregular group makes up <5%. Of the galaxies at z > 3, 38% have a visually identifiable spheroidal component, decreasing from 42% to 26% between z = 3 and 4.5, then varying between ∼30% and 40% beyond z = 4.5. This is largely driven by the similar decrease then increase in the spheroid-only group. Part of this apparent trend at higher redshifts may be due to small number statistics, and part may be due to a number of selection effects. For example, there is a possibility that we miss fainter extended features in some of these systems at high redshift. It is also possible that a larger fraction of galaxies at higher redshift are small enough to be at the resolution limit of NIRCam, given the expected size evolution of galaxies, and are therefore more round and compact in appearance.

Of the galaxies at z > 3, 43% have irregular features and this fraction remains roughly constant across the full redshift range due in part to the fraction of disk+irregular galaxies being roughly constant at 20% and then decreasing, while the fraction of irregular-only galaxies is roughly 10%–15% and then increases to 20% by z = 4.5. Note that the total fractions of objects that are all disks, all spheroids, or all irregulars do not add up to one due to the overlapping objects in each of these classes.

Finally, we note that the fraction of point sources and unclassifiable objects remains at 0%–2% across most of the redshift range. At z > 6, 13% of galaxies are unresolved and 8% are unclassifiable, corresponding to five and three individual galaxies, respectively, in this redshift bin. We remind the reader that the above percentages correspond to galaxies that were bright enough to be detectable with HST CANDELS imaging and may not be representative of the overall galaxy sample detectable by JWST at these redshifts.

One of the major advances that JWST NIRCam imaging brings to morphological analyses of galaxies is that the broad wavelength coverage enables us to probe the rest-frame optical morphologies of galaxies across a wide redshift range. As described in Section 3.3, we used GalfitM to perform multiwavelength parametric fits across all of the NIRCam filters. For a fair comparison of these parameters at different wavelengths, we select the NIRCam filter closest to the rest-frame optical at different redshifts. Throughout this section, we use the F277W filter for galaxies at 3.0 < z < 4.0, F356W for galaxies at 4.0 < z < 4.5, and F444W for galaxies at z > 4.5. We compare these measurements to the visual classifications in Figures 4 and 5.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Overall, the distribution of Sérsic indices (where n = 0.5 corresponds to a Gaussian profile, n = 1 to an exponential profile, and n = 4 to a de Vaucouleurs profile) tracks with the expectations from the visual classifications. Disk galaxies with no apparent spheroid or irregular features (disk-only) peak at low Sérsic indices with a long tail out to higher values, as expected ($\langle n\rangle ={1.16}_{-0.40}^{+0.88}$, where the error bar denotes the 16th–84th percentile range of the distribution). Galaxies that are pure spheroids (spheroid-only) have a much broader distribution and peak at higher n ($\langle n\rangle ={2.46}_{-1.41}^{+2.32}$), as has been noted at lower redshift based on HST imaging (e.g., Kartaltepe et al. 2015; Vika et al. 2015). Galaxies with both a disk and spheroidal component (disk+spheroid) peak at intermediate values ($\langle n\rangle ={2.39}_{-0.98}^{+1.06}$). This illustrates that a cut at a fixed Sérsic index would not cleanly select disk- or spheroidal-dominated galaxies. For example, a dividing line of n = 2 would identify 71% of the visually identified disks and only 45% of the visually identified spheroids. However, it is worth noting that a fraction of the objects visually identified as spheroids with low n might have extended low surface brightness disks that are difficult to pick out by eye.

Irregular galaxies with no disk or spheroid component peak at very low n, with a substantial fraction at n < 1 and a long tail out to higher values ($\langle n\rangle ={1.19}_{-0.82}^{+1.51}$). Disk galaxies with irregular features (disk+irregular) peak closer to n = 1 with a narrower distribution that more closely resembles that of the disk-only galaxies. Likewise, the distribution of Sérsic indices for the spheroidal galaxies with irregular features (spheroid+irregular) closely resembles that of the spheroid-only group. A visual inspection of the models and residuals for the irregular galaxy population reveals that, unsurprisingly, irregular features are not well fit by a Sérsic profile. For disks and spheroidal galaxies with irregular features, the model fits the disk/spheroidal component well and leaves behind features in the residuals, while the irregular-only population is not well fit at all. We caution the reader against overinterpreting Sérsic indices for irregular galaxies and using Sérsic indices to select disk galaxies without first checking the images (and residuals) for irregular features.

The middle panels of Figure 4 show the distribution of sizes (effective radii, R_e ) measured by GalfitM for each of the morphological types. Galaxies with disks and irregular features generally have larger sizes than those with spheroids. For example, galaxies with disks only have a median effective radius of ${1.19}_{-0.48}^{+0.61}$ kpc, galaxies with irregulars only have $\langle {R}_{e}\rangle ={1.32}_{-0.57}^{+0.79}$ kpc, while spheroid-only galaxies have $\langle {R}_{e}\rangle ={0.54}_{-0.32}^{+0.42}$ kpc. These trends are seen more clearly in Figure 5. The disk+irregular and disk+spheroid+irregular groups have size distributions that more closely match the distribution for the irregular-only group, while the spheroid+irregular group has a smaller median size ($\langle {R}_{e}\rangle ={0.74}_{-0.34}^{+0.81}$ kpc).

The distribution of the axis ratios is shown in the right panel of Figure 4 and in Figure 5 and offers another way to compare our visual morphologies to a quantitative measurement. A population of disks with exponential profiles and random orientations is expected to have a relatively flat distribution of axis ratios that falls off at low values, while triaxial ellipsoids are expected to have a distribution that is peaked at higher values, b/a ∼ 0.6 (e.g., Elmegreen et al. 2005; Ravindranath et al. 2006; Padilla & Strauss 2008; Law et al. 2012; Robertson et al. 2022). The mean values we see for the different morphological groups follow this trend. The spheroid-only group has the largest median axis ratio (${0.64}_{-0.18}^{+0.17}$), while the disk+irregular group has the smallest (0.36${}_{-0.11}^{+0.21}$).

Figure 6 shows the axis ratio as a function of effective radius split into several redshift bins. In each redshift bin from z = 3 to 6, the spheroid-only galaxies have the smallest median effective radius and largest median axis ratios, suggestive of a population of true triaxial ellipsoids. Overall, we see the general trend of small galaxies being rounder, for all morphology types, as seen by Padilla & Strauss (2008) and Zhang et al. (2019), in each of these redshift bins. We cannot draw any conclusions at z > 6 due to the small sample size and the previously mentioned uncertainties.

Zoom In Zoom Out Reset image size

As described in Section 3.4, we used statmorph to measure nonparametric image statistics for all of the HST-selected z > 3 galaxies in our sample across all NIRCam filters. We use the same filters corresponding to the rest-frame optical emission for each galaxy as we did for the above parametric comparison. In total, 81% of the galaxies in our sample have a reliable fit from statmorph; Figures 7 and 8 highlight two of the commonly used methods to separate galaxies into the standard Hubble types and identify mergers based on these image statistics (e.g., Bershady et al. 2000; Conselice 2003; Lotz et al. 2004, 2008a).

Zoom In Zoom Out Reset image size

The top panels of Figure 7 show the location of each galaxy on the asymmetry–concentration plane, with the classic lines used to mark the boundaries between disk galaxies, elliptical galaxies, intermediate galaxies, and mergers for moderate-redshift HST images (Bershady et al. 2000). While these boundaries do not cleanly separate z > 3 galaxies into different types relative to their visual classifications, a few trends can be seen. On average, galaxies with a spheroid have a higher concentration (C) than those with disks and irregulars. Similarly, irregular galaxies have a higher asymmetry value (A), on average. Very few galaxies lie above the classic demarcation for mergers (Conselice 2003), and those that do span the full range of visual morphologies, albeit with a higher fraction of irregulars. Figure 8 shows the distribution of the same galaxies on this plane but color coded by the median Sérsic index for each bin. This distribution highlights the correlation of the concentration value with the Sérsic index for the sample.

Zoom In Zoom Out Reset image size

The bottom panels of Figure 7 show the location of each galaxy on the Gini–M₂₀ plane. The lines mark the boundaries between disk galaxies, ellipticals/S0s, and mergers based on nearby galaxies and adjusted for galaxies at z ∼ 1 (Lotz et al. 2004, 2008a). While there is no discernible difference between z > 3 disks and spheroids with this diagnostic (similar to what has been seen at lower redshift and in simulations; e.g., Lotz et al. 2008b; Kartaltepe et al. 2010; Pearson et al. 2019), irregular galaxies have higher Gini and M₂₀ values, on average. Of galaxies with irregular features, 28% lie above the merger line, whereas only 17% of disks and 15% of spheroids do. The right panel of Figure 8 shows the distribution of the same galaxies on this plane but color coded by the median Sérsic index for each bin. Galaxies that occupy the ellipticals/S0s portion of this plane have a higher Sérsic index, on average. A significant fraction of the galaxies in the merger region of the plane also have a higher average Sérsic index.

We find that galaxies detected by both HST and JWST in the z > 3 universe have a wide diversity of morphologies. Overall, ∼60% of these galaxies have disks (including those that also have spheroids and/or irregular features) at z = 3–4, and this fraction has an apparent downward trend with increasing redshift. Other early JWST studies have identified candidate disk galaxies at these redshifts (e.g., Ferreira et al. 2022a, 2022b; Robertson et al. 2022) and find similar fractions. Galaxies with spheroids make up ∼40% over this redshift range, with some variations in the higher-redshift bins that are likely related to the small numbers in these bins overall or the difficulty in identifying low surface brightness features at these redshifts. Galaxies with a pure spheroid, i.e., without discernible disks or irregular features, make up ∼20% across the full redshift range, roughly consistent with the findings of Ferreira et al. (2022a). While the fraction of all galaxies with irregular features is roughly constant at all redshifts (∼40%–50%), the fraction of galaxies that are purely irregular (i.e., those that have no discernible disk or spheroidal features) increases from ∼12% at z = 3–3.5 to ∼20% at z > 4.5. This fraction is lower than that reported by Ferreira et al. (2022a, 2022b). Slight differences among these early studies likely arise due to the different classification schemes and mass ranges used.

The distribution of axis ratios and sizes presented in Section 4.2 and Figures 5 and 6 suggests that our z > 3 sample indeed contains a mix of true disks and spheroids. The distribution of axis ratios for the galaxies classified as spheroid-only is consistent with that expected from a population of triaxial ellipsoids, while the relatively flat distribution and lower median for the disk galaxies is expected for a population of disks with exponential profiles and random orientations (e.g., Padilla & Strauss 2008). The axis ratio distribution for the spheroid-only group peaks at b/a > 0.6 across the entire redshift range of our sample, while the median for the disk-only group remains at b/a ∼ 0.4. The disk+spheroid and spheroid+irregular populations have axis ratio and size distributions that are intermediate. Previous theoretical work has found that galaxy shapes have evolved over time from prolate to oblate as they transition from having dark matter–dominated interiors to baryonic matter interiors following a compaction event (e.g., Ceverino et al. 2015; Tomassetti et al. 2016) and that this transition happens earlier for more massive galaxies (e.g., Zhang et al. 2019).

It is worth noting that selection effects may be partially responsible for the axis ratio distributions observed for these objects. For example, it has been shown that the distribution of axis ratios has a strong dependence on the mass and luminosity of the galaxy population (e.g., Padilla & Strauss 2008; Zhang et al. 2019). Galaxy orientation also plays an important role in this distribution, as face-on disks may be more difficult to detect than edge-on disks at the magnitude limit (e.g., Elmegreen et al. 2005), and the presence of dust can impact the measured axis ratios (Padilla & Strauss 2008). The size of the current sample does not allow binning by mass, luminosity, or finer morphology groupings; however, future work with larger sample sizes will allow greater exploration of this parameter space.

To summarize, we see evidence for galaxies with established disks and spheroidal morphologies across the full redshift range of our sample. We emphasize that the fractions quoted here are apparent fractions only and that several observational effects likely play a role in these measurements (as discussed in Section 5.2). Further work is needed to quantify our ability to pick out disk features and resolve spheroidal galaxies in JWST images at varying image depths in order to quantify the true fraction of galaxies with disks and spheroids at these redshifts. Likewise, larger samples, particularly at z > 6, will be needed to truly establish when the first disks began to form, when disks grew their bulges, and when spheroids emerged.

Figure 9 shows the morphological fractions as a function of redshift based on CANDELS HST imaging and using the visual classifications of Kartaltepe et al. (2015) for all z > 3 galaxies in all five CANDELS fields (1375 galaxies in total). The HST classifications were limited to galaxies with F150W < 24.5 because fainter galaxies could not be reliably classified, so only 59 galaxies out of the 850 in our sample are bright enough to make that cut. A comparison of the JWST morphological fractions with those 59 specific galaxies is shown in Figure 10. Based on the HST imaging alone, a smaller fraction of galaxies at z = 3.0–4.5 have disks (∼40%), and a larger fraction are pure spheroids (∼20% at z = 3.0–5.0). The fraction of galaxies that are only irregular is small and drops with redshift, from ∼5% at z > 6. The fraction of galaxies that are unclassifiable rises sharply, from ∼5% at z = 3.5, to ∼35% at z = 5.5, to ∼80% at z > 7. Likewise, ∼30% are unresolved at z = 5–7. Among the 59 galaxies with both HST and JWST classifications, a higher fraction are classified as spheroids and a lower fraction are classified as disks or irregulars with HST than with JWST. At z > 6, all of the objects were unclassifiable or unresolved with HST.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

The large difference seen between the HST and JWST morphologies at these redshifts is expected and due to the difference in depth and wavelength coverage. A total of 488 galaxies were flagged by at least one classifier as having a different morphology in the JWST images compared to the HST images (159 galaxies were flagged by two out of three classifiers; see Figure 11 for examples). A significant number of galaxies with disks were previously identified as spheroids because of their compact central morphologies, with low surface brightness disk features that only became visible with deeper imaging (see, for example, Conselice et al. 2011; Mortlock et al. 2013; Kartaltepe et al. 2015). This suggests that some fraction of the spheroidal galaxies observed with JWST, particularly those that are faint and/or at higher redshifts, possibly have unobserved disks as well. It is not likely that these disks would previously have been identified as irregular, except for some at the low-redshift end, as these irregular features are also too faint to be easily identified at these redshifts with HST. At the low-redshift end (z = 3–3.5), some disks may have been classified as irregular if the HST data only picked up the brighter star-forming clumps rather than the underlying disk structure.

Zoom In Zoom Out Reset image size

To explore the impact of observed wavelength on the classifications, classifiers were also asked to flag objects for which the morphology changes (i.e., they would have selected different main morphology classes) between the NIRCam SW filters and the LW filters. At least one classifier chose this flag for 190 galaxies, and two out of three chose it for 37 (see Figure 11 for examples). Note that some of the differences seen across the different filters may be due in part to the increased resolution of the SW bands. This flag was rarely chosen, suggesting that the depth of the JWST images is the primary driver in the morphological differences observed between JWST and HST.

We compare the results of our surface brightness profile fits and nonparametric fits with the results from mock images and catalogs based on several different simulations in Figure 12 in three different redshift bins.

Zoom In Zoom Out Reset image size

First, we use a mock galaxy catalog based on the Santa Cruz semianalytic model (SAM) and publicly available as part of the CEERS simulated data release SDR3. ⁵⁹ The CEERS mock galaxy catalog is an augmented version of the EGS light cone presented by Yung et al. (2022), which spans 782 arcmin² and contains galaxies with −16 ≳ M_UV ≳ −22 in the range 0 < z ≲ 10. The physical properties of the galaxies are modeled with the physics-based Santa Cruz SAM (Somerville et al. 2015, 2021; Yung et al. 2019). The sizes of the disk components of the galaxies are computed based on the ansatz that the specific angular momentum of the halo gas is equal to that of the dark matter halo and that it is conserved during disk formation (Mo et al. 1998; Somerville et al. 2008).

We also use the publicly available ⁶⁰ mock images and derived morphological catalogs of Costantin et al. (2022b) and Rose et al. (2022), which use the IllustrisTNG cosmological simulation. ⁶¹ The IllustrisTNG project (Marinacci et al. 2018; Naiman et al. 2018; Nelson et al. 2018; Pillepich et al. 2018; Springel et al. 2018) is a series of large cosmological magnetohydrodynamical simulations of galaxy formation and is an update to the original Illustris-1 simulation (Vogelsberger et al. 2014). It consists of three different runs that span a range of cosmological volumes and resolutions: TNG50, TNG100, and TNG300.

Costantin et al. (2022b) produced synthetic images of TNG50 galaxies and used the Multi Instrument Ramp Generator (MIRaGe ⁶² ) to simulate raw NIRCam images for the CEERS depth and filter combination. These images were then reduced using the JWST pipeline. Morphological measurements were made using statmorph. Rose et al. (2022) produced noiseless synthetic images with TNG100 galaxies using the public visualization API (Nelson et al. 2019) in each of the CEERS filters. These images were then convolved with the model PSF for each filter using WebbPSF (Oschmann et al. 2014). Poisson and background noise estimated from the JWST exposure time calculator (Pontoppidan et al. 2016) were added to create mock images at the CEERS depth. Parametric models were fit using Galapagos-2 and GalfitM, while nonparametric fits were performed using statmorph.

Figure 12 compares the distribution of the Sérsic indices and sizes of galaxies from the SAM, TNG50, and TNG100 to the distribution measured from CEERS galaxies (Section 4.2 and Figure 4). The overall distributions from the SAM have very similar peaks with a narrower distribution, which holds for all three redshift bins. The Sérsic index for both TNG50 and TNG100 peak at lower values than the CEERS galaxies and have narrower distributions at all redshifts. At z = 3–4, TNG50 galaxies have larger sizes than TNG100 galaxies and even larger than both the SAM galaxies and the observed CEERS galaxies. At z > 4, the distributions match more closely. At all redshifts, the simulations do not contain the smaller (lower R_e ), more compact (larger n) galaxies that we observe with JWST CEERS imaging.

Figure 12 also compares the measured axis ratio and asymmetry value for the TNG50 and TNG100 galaxies to the distribution from CEERS. In all three redshift bins, the axis ratios of the TNG50 and TNG100 galaxies match each other well but peak at higher b/a (∼0.6) and fall off more sharply at lower values than the observed CEERS galaxies. At z = 3–4, the asymmetry distributions for TNG50, TNG100, and CEERS are well matched, but the TNG50 and TNG100 distributions shift toward lower (more negative) values at higher redshift. Negative asymmetry values are unphysical and typically result from low-S/N sources, where the source is very close to the background level that is being subtracted when making the asymmetry measurement.

Overall, the agreement between our measurements for the z > 3 JWST CEERS galaxies and the various simulations is encouraging. The differences seen (for example, the difference in axis ratio and the lack of small compact galaxies in the simulations) are worthy of a more in-depth look in order to determine if there are selection effects impacting the results or there is an actual physical difference between the galaxies in these simulations and those in the real observed universe.

We used Galfit to compute parametric fits on the F150W and F200W images. As initial guesses, we use the source location, magnitude, size, position angle, and axis ratios from the Source Extractor catalogs and segmentation maps. We use cutouts of the error array (ERR extension) produced by the JWST pipeline as the input sigma images for each source, which includes Poisson noise from the sources themselves, as well as the usual instrument noise. As input PSFs, we create empirical PSFs for each filter from stacked stars from all four CEERS pointings. For each galaxy, the Kron radius measured by Source Extractor was used to scale the size of the cutout used as input to Galfit. All galaxies in the cutout within 3 mag of the primary source were simultaneously fit, down to a magnitude limit of 27, with all other sources masked. Based on our testing, we find that sources fainter than this limit were not reliably fit. We assign each fit a quality flag. A flag of 0 indicates a good fit; a flag of 1 indicates that the fit is suspect, meaning the resulting Galfit magnitude differs substantially from the input Source Extractor magnitudes; a flag of 2 indicates a poor fit, where one or more parameters reached a constraint limit; a flag of 3 indicates that the fit failed to find a solution; and a flag of 4 indicates that the source was not fit because it was either an artifact or located too close to the edge.

Of the 850 galaxies in the z > 3 sample, 74% have a Galfit flag of 0, 13% have a flag of 1, 8% have a flag of 2, 2% have a flag of 3, and none have a flag of 4 (these would already have been removed by our initial sample selection). For the comparisons discussed below, we use all of the galaxies with a flag of 0 or 1, representing 87% of the total sample.

As a check on the quality of our fits, and for consistency, we compared the Galfit and GalfitM fits for the F200W filter and visually inspected the model and residuals. Overall, we see a high level of agreement with no significant offsets for the Sérsic index, size, axis ratio, and magnitudes between the two measurements for the sources that do not reach a constraint limit. Throughout the paper, we use the GalfitM measurements for the filter closest to the rest-frame optical at the redshift of the galaxy. The median and 16th and 84th percentile range for the Sérsic index (n), effective radius (R_e ), and axis ratio (b/a) for each morphological type are given in Table 1 and used throughout the text.

Table 1. Median Properties for Each Morphology Class

Morphology Group	No.	Redshift	log(Mass (M⊙))	na	Re (kpc) a	b/aa
All disks	467	3.84 ± 0.87	9.31 ± 0.63	${1.41}_{-0.65}^{+1.38}$	${1.19}_{-0.51}^{+0.75}$	${0.40}_{-0.14}^{+0.22}$
Disk-only	192	3.86 ± 0.90	9.28 ± 0.47	${1.16}_{-0.40}^{+0.88}$	${1.19}_{-0.48}^{+0.61}$	${0.40}_{-0.13}^{+0.21}$
Disk+spheroid	88	3.86 ± 0.98	9.29 ± 1.12	${2.38}_{-0.98}^{+1.06}$	${0.83}_{-0.33}^{+0.73}$	${0.51}_{-0.16}^{+0.17}$
Disk+irregular	155	3.82 ± 0.82	9.32 ± 0.41	${1.12}_{-0.51}^{+1.10}$	${1.33}_{-0.42}^{+0.78}$	${0.36}_{-0.11}^{+0.21}$
Disk+spheroid+irregular	32	3.81 ± 0.60	9.49 ± 0.48	${2.63}_{-1.14}^{+1.63}$	${1.33}_{-0.45}^{+1.37}$	${0.39}_{-0.10}^{+0.18}$
All spheroids	323	3.94 ± 1.04	9.30 ± 0.71	${2.48}_{-1.34}^{+2.12}$	${0.72}_{-0.36}^{+0.66}$	${0.56}_{-0.20}^{+0.21}$
Spheroid-only	156	4.03 ± 1.11	9.29 ± 0.47	${2.46}_{-1.41}^{+2.32}$	${0.54}_{-0.32}^{+0.42}$	${0.64}_{-0.18}^{+0.17}$
Spheroid+irregular	47	3.90 ± 1.12	9.21 ± 0.39	${2.92}_{-1.96}^{+2.30}$	${0.72}_{-0.34}^{+0.81}$	${0.55}_{-0.20}^{+0.18}$
All irregulars	376	3.90 ± 0.93	9.29 ± 0.43	${1.37}_{-0.80}^{+1.73}$	${1.28}_{-0.57}^{+0.84}$	${0.40}_{-0.13}^{+0.21}$
Irregular-only	142	4.02 ± 1.02	9.24 ± 0.44	${1.12}_{-0.82}^{+1.51}$	${1.32}_{-0.57}^{+0.79}$	${0.42}_{-0.13}^{+0.17}$
Point source/unresolved	16	5.27 ± 1.82	9.66 ± 0.65	⋯	⋯	⋯
Unclassifiable	18	4.57 ± 1.54	9.47 ± 0.84	⋯	⋯	⋯

Note.

^a Sérsic index (n), effective radius (R_e in kpc), and axis ratio (b/a) as measured in the NIRCam filter that most closely represents the rest-frame optical for the redshift of the galaxy: F277W for 3.0 < z < 4.0, F356W for 4.0 < z < 4.5, and F444W for z > 4.5.

Download table as:  ASCIITypeset image