COSMOS-Web: The Role of Galaxy Interactions and Disk Instabilities in Producing Starbursts at z < 4

For most of their lives, star-forming galaxies follow the main sequence (e.g., E. Daddi et al. 2007; D. Elbaz et al. 2007; K. G. Noeske et al. 2007), a tight relation between stellar mass and star formation rate (SFR) set by an equilibrium state between gas consumption, outflow, and inflow (R. Davé et al. 2012; S. J. Lilly et al. 2013; R. Feldmann 2015). The scatter of this relation (∼0.3 dex) is set by oscillations around that equilibrium, driven by constant adjustment of the galaxies' SFR to the available cold gas and replenishing gas inflows (e.g., N. Bouché et al. 2010; A. Dekel et al. 2013; S. Tacchella et al. 2016). Stronger oscillations are observed mostly in low-mass galaxies at low and local redshifts, while at higher redshifts, galaxies across all masses seem to show such a behavior (e.g., D. R. Weisz et al. 2012; N. Emami et al. 2019; A. L. Faisst et al. 2019). In addition to these oscillations, galaxies experience significant bursts of star formation, elevating their SFRs by factors of 10 and more above the main sequence (D. B. Sanders & I. F. Mirabel 1996). The fraction of starburst galaxies to the total population increases with redshift from ∼1% at z < 0.4 (N. Bergvall et al. 2016) to ∼2%–5% at z ∼ 0.5–1.0 (G. Rodighiero et al. 2011; L. Bisigello et al. 2018) and higher fractions at higher redshifts depending on stellar mass (M. T. Sargent et al. 2012; K. I. Caputi et al. 2017; L. Bisigello et al. 2018). Interactions between gas-rich galaxies lead to the inflow and compression of gas and are, therefore, usually the prerequisites for starbursts (e.g., D. B. Sanders & I. F. Mirabel 1996; T. J. Cox et al. 2008; R. Genzel et al. 2010; J. S. Kartaltepe et al. 2012; C.-L. Hung et al. 2013; P. F. Hopkins et al. 2018; A. Cibinel et al. 2019; J. Moreno et al. 2021; E. A. Shah et al. 2022). However, some simulations and observations suggest that starbursts and galaxy major mergers do not necessarily have to be causally connected (e.g., G. Rodighiero et al. 2011; S. Kaviraj et al. 2013; M. Sparre & V. Springel 2016; J. Fensch et al. 2017; S. Wilkinson et al. 2022; Y. A. Li et al. 2023; Z. Liu et al. 2024). In addition to starbursts, disk instabilities in gas-rich galaxies (A. Toomre 1964; F. Bournaud et al. 2007), which may be tightly related to increased gas inflows or minor mergers (N. Scoville et al. 2023), could increase the rate of star formation without external interactions, hence contributing to the starburst population (e.g., B. Wang & J. Silk 1994; A. Immeli et al. 2004; A. Dekel et al. 2009; K. Tadaki et al. 2018) and the buildup of bulge-dominated galaxies later on (e.g., R. Bouwens et al. 1999). The above shows that the origin of starbursts is, therefore, complex, and the definitive link between galaxy interactions, disk structure, increased star formation efficiency (SFE), and gas fraction out to high redshift has yet to be studied with robust measurements and comprehensive samples.

The works by D. Liu et al. (2019) and N. Scoville et al. (2023) are likely the most comprehensive submillimeter dust-based studies (in terms of sample size and redshift coverage) on the relation between gas and star formation properties of main-sequence and starburst galaxies to date. They combine robust SFRs from Herschel, gas fraction and SFE measurements from the Atacama Large Millimeter/submillimeter Array (ALMA), and other parameters (such as stellar masses) from comprehensive multiband photometry for >700 galaxies out to z = 6. Specifically, N. Scoville et al. (2023) found that starbursts at a given redshift and stellar mass are primarily triggered due to enhanced SFE rather than a high gas fraction. Similar results have been found in independent studies by J. D. Silverman et al. (2015) and L. J. Tacconi et al. (2018) using large-sample CO data (for a review, see L. J. Tacconi et al. 2020). On the other hand, an increase in the gas fraction is likely responsible for the maintained high specific SFRs at high redshift (e.g., L. J. Tacconi et al. 2010, 2018; R. Genzel et al. 2015; J. Freundlich et al. 2019; A. Gowardhan et al. 2019; D. Liu et al. 2019; M. Dessauges-Zavadsky et al. 2020).

In this work, we explore the population of starburst galaxies out to z = 4 in terms of their interactions, disk structure, SFE, and gas fractions. Specifically, we aim to study the role of galaxy–galaxy interactions and disk instabilities in producing starburst activity. To this end, we combine the N. Scoville et al. (2023) sample with high-resolution rest-frame optical imaging from the COSMOS-Web JWST/NIRCam imaging (C. M. Casey et al. 2023).

After summarizing the data (Section 2), we detail the morphological classification in Section 3. In Section 4, we present the results, and we conclude in Section 5. Throughout this work, we assume a ΛCDM cosmology with H₀ = 70 km s⁻¹ Mpc⁻¹, Ω_Λ = 0.7, and Ω_m = 0.3, and magnitudes are given in the AB system (J. B. Oke 1974). We use a G. Chabrier (2003) initial mass function for stellar masses and SFRs.

For this work, we use a sample of 387 galaxies out to z ∼ 4. This sample is a subsample of the 704 galaxies from N. Scoville et al. (2023; see more details in earlier works; N. Scoville et al. 2014, 2016, 2017) requiring existing JWST imaging. The selected galaxy sample is best for carrying out this study, as it includes

• 1.  
  robust measurements of total SFR from UV+far-IR Herschel and ALMA measurements,
• 2.  
  measured molecular gas masses from far-IR dust continuum observations with ALMA,
• 3.  
  measured stellar masses and accurate photometric redshifts from multiband photometry, and
• 4.  
  deep subkiloparsec-resolution JWST rest-frame optical and near-IR imaging data.

The galaxies reside in the COSMOS field (N. Scoville et al. 2007), which provides a wealth of ancillary data including X-ray and radio measurements, Hubble rest-frame UV and optical imaging (ACS/F814W; A. M. Koekemoer et al. 2007), and spectroscopic redshifts for more than half of the sample. In the following, we summarize the different data products, basic measurements, and sample properties.

2.1. Photometry and Imaging

2.2. Definitions and Sample Properties

We define the molecular gas fraction as f_gas = M_mol,gas/(M_* + M_mol,gas) and the SFE based on the total UV+IR SFR as SFE = M_mol,gas/SFR in units of Gyr⁻¹. The depletion time (t_depl) is the inverse of the SFE. To calculate Δ_MS (the logarithmic offset from the star-forming main sequence at a given stellar mass), we assume the main-sequence parameterization from N. Lee et al. (2015), which is based on Herschel measurements of individual galaxies and stacked samples on the COSMOS field. We emphasize that Herschel observations are crucial for deriving a robust Δ_MS of starburst galaxies, often dominated by dust-obscured star formation. The use of other parameterizations (e.g., J. S. Speagle et al. 2014; C. Schreiber et al. 2015) leads to similar Δ_MS and has a minimal impact on the final results of this work.

Figure 1 shows the dependence of SFE and f_gas on stellar mass and Δ_MS for our selected 387 galaxies in different redshift bins. The underlying hex bins show the full sample from N. Scoville et al. (2023). This shows that our selected subsample traces the full sample properties closely. In addition, note that Δ_MS correlates well with SFE and less with f_gas, as it has been found in various works including that of N. Scoville et al. (2014) and J. D. Silverman et al. (2015).

Zoom In Zoom Out Reset image size

There are mainly two selection biases to be noted. These are due to the requirement of Herschel and ALMA detections for the computation of accurate far-IR SFRs and gas masses. First, the high-redshift sample is biased toward high Δ_MS. Second, low stellar masses are biased toward high Δ_MS. While the first bias is not a large concern (our analysis will be done in separate redshift bins), the second bias may introduce unphysical relations with stellar mass. We therefore mitigate that bias by introducing mass-controlled samples for the following analysis by adopting two mass bins for each redshift range. We found that two stellar mass bins are sufficient to mitigate the bias, and more mass bins would reduce the sample and decrease the statistical robustness. However, due to the small sample at z > 3, we only adopt a single mass bin for the highest redshifts. We define the mass cuts in each of the redshift bins based on Figure 1 to remove any significant dependence between stellar mass and Δ_MS within the mass-binned subsamples. The adopted mass bins (in $\mathrm{log}(M/{M}_{\odot })$) are [10.1, 10.8] and [10.8, 11.7] for z ∼ 0.5, [10.4, 11.1] and [11.1, 11.8] for z ∼ 1.5, [10.7, 11.1] and [11.1, 11.7] for z ∼ 2.5, and [10.7, 11.2] for z ∼ 3.5 (all redshift bins with Δz = 1). Changing these mass bins in a reasonable range (±0.2 dex) does not affect the final conclusions of this work.

3.1. Visual Classification

We perform a visual classification of the 387 galaxies using all of the COSMOS-Web JWST NIRCam filters to minimize the effects of color-dependent morphology corrections. In the following, we define four morphological groups.

• Noninteracting. Isolated galaxies of various kinds, such as disk or compact galaxies, that are not in an interacting state.
• Disk galaxies. These are extended galaxies (to be distinguished from compact galaxies), preferentially with a semiminor-to-semimajor axis ratio of less than 0.8 and a disk structure. The disk can either be face-on or edge-on.
• Clumpy disks. This is a subset of the disk galaxies category, showing more than one star-forming clump. The clumps are pronounced in the bluer bands but visible throughout redder wavelengths. Due to limitations in resolution and confusion with interacting systems, we have not classified clumpy disks at redshifts z > 3 in this work.
• Interacting. Interacting or irregular galaxy systems are identified by multiple nuclei, irregular structures, or tidal features. This category should include mainly galaxy systems in premerger and postmerger phases.

It is important to note that we must make use of all available NIRCam filters for a reliable morphological classification. One example to showcase the importance of multiple photometric bands for a morphological classification is discussed in Z. Liu et al. (2024). As shown in their Figure 1, this starburst galaxy could be classified as a clear merger in JWST filters F115W and F150W. However, at longer rest-frame wavelengths (specifically F227W), this source shows a clear disk structure—a finding that is supported by kinematic measurement of the CO gas.

For consistency across redshifts, we therefore visually classified interacting galaxies primarily at rest-frame optical wavelengths (0.8–1 μm) and redder bands where available. This choice is motivated by the fact that those wavelengths trace the stellar mass of galaxies and hence are more ideal to discriminate between interacting systems and galaxies with a significant disk disk component as shown in the case by Z. Liu et al. (2024). We use F150W for z < 1.0, F277W for 1 < z < 3, and F444W for z > 3 but also consider redder bands where available. In addition, red, green, and blue (RGB) color images (in addition to redshifts where available) are used to identify companions (mergers) and reject background/foreground galaxies. We note that the PSF varies significantly between those bands (from FWHM 0$\mathop{.}\limits^{\unicode{x02033}}$06 to 0$\mathop{.}\limits^{\unicode{x02033}}$16 for F115W to F444W; M.-Y. Zhuang & Y. Shen 2024), which may be a concern for high redshifts. We tested our classification by carefully comparing PSF-matched images to rule out significant biases in our visual classification due to such PSF differences. Although we made use of all available bands in the morphological classification, we additionally tested the dependence of our morphological classification on single bands. We find, as expected, that using filters at bluer wavelengths breaks up the galaxies in clumps (see also below) and biases the morphological classification against disk galaxies. On the other hand, the use of filters redder than F227W for z > 1 galaxies does not change the classification. In addition, ∼5% of galaxies are covered by the PRIMER-COSMOS observations, which include four more filters in addition to COSMOS-Web. The morphological classification was also tested against these filters (specifically F356W, providing an intermediate band between F277W and F444W at a better PSF resolution than F444W), and we do not find biases using redder bands for the classification of interacting systems.

On the other hand, we use bluer bands (F115W and F150W) for the identification of clumpy disk galaxies. This is motivated by the fact that such clumps are generally star-forming and rest-UV-bright (and at the same time, the PSF of bluer bands is smaller). However, there may be a bias in clump identification in high-redshift galaxies caused by their smaller physical sizes and the large PSF at the same rest-frame filters. Furthermore, a differentiation between clumpy disks and interacting galaxies without kinematic information from spectroscopy is nearly impossible in these cases (G. C. Jones et al. 2021). We therefore do not attempt to identify clumpy disk galaxies at z > 3. In addition to the above method, we use residual images for identifying clumps. These are generated in two ways: (i) by subtracting a smooth profile fit (Sérsic) performed using statmorph (V. Rodriguez-Gomez et al. 2019) from the galaxy image and (ii) by directly subtracting a long-wavelength from a short-wavelength image (e.g., F150W–F444W), which highlights blue rest-UV structure and removes the smooth ISM. We classify a galaxy as a clumpy disk if there are more than three bright (S/N > 3) clumps detected. We use both methods to support the visual classification of clumpy disks from the individual or RGB images.

Finally, we compare the above visual classification with more statistical classifications such as the G–M₂₀ (Gini versus M₂₀) or C–A (concentration versus asymmetry) diagnostics. We find that these methods are most efficient in selecting interacting galaxies; however, they have a significantly lower purity and completeness for selecting clumpy disks (see the Appendix).

In the following, we denote the fraction of interacting galaxies in our sample with f_int. Note that, as mentioned above, without spectroscopic confirmation of the detailed kinematics, we are not able to identify close-pair mergers. The interacting group therefore mostly focuses on the pre- and postmerger phases. The fraction of clumpy disk galaxies compared to the total disk galaxy population in a given selection bin is denoted as f_clumpy. We derive the uncertainties on these fractions using the formalism of (Bayesian) binomial confidence intervals suggested by E. Cameron (2011) and shown in code form in their Appendix A.

Figure 2 shows examples of our classification scheme in four redshift bins. The redshifts and galaxy types are indicated on each cutout, as well as f_gas and SFE. Note the clear distinction of disk galaxies with and without clumps.

Zoom In Zoom Out Reset image size

4.1. The Impact of Galaxy–Galaxy Interactions on Starburst Activity

We first compare the fraction of interacting systems to the offset from the main sequence, the gas fraction, and the SFE.

The top panel in Figure 3 shows the fraction of interacting galaxies as a function of Δ_MS in two different mass bins per redshift. The fraction clearly increases off the main sequence (indicated by the gray area) toward higher Δ_MS and reaches 20%–40% at z < 3 and 80% at z = 3.5. However, we note that, due to the lack of near-IR observations at these redshifts, the merger classification may not be as robust as at lower redshifts (for example, due to dust obscuration for starburst galaxies). This can lead to a classification bias (e.g., due to clumpy disks; see discussion in Section 3.1) that may lead to an overprediction of the fraction of merging galaxies at z = 3.5.

Zoom In Zoom Out Reset image size

Note that here we focus on the change in fractions across different properties and do not attempt to compare their absolute values to the literature. This is because relative changes are more reliable than the absolute values of f_int, which depend on sample selection and subjective visual selection thresholds and therefore are difficult to compare. The two mass bins show very similar behaviors. Keeping this in mind, an increase of the fraction of interacting systems toward higher redshifts is found, which is in agreement with a global increase in the merger fraction studied in other works (e.g., L. A. M. Tasca et al. 2014; M. Romano et al. 2021). Overall, this is consistent with the theory of galaxy–galaxy interactions playing some role in inducing starburst activity. Here we show that this may be the case out to z ∼ 4. We also note that this is in line with and directly related to the increased fraction of interacting galaxies found at higher infrared luminosities (see the study at z < 1.5 by C.-L. Hung et al. 2013).

It is suggested that the gas fraction is weakly correlated with starburst activity (e.g., N. Scoville et al. 2023), and we therefore would not expect a significant correlation between f_int and f_gas. The middle panels of Figure 3 show that this is indeed the case by displaying f_int as a function of the gas fraction of galaxies relative to the population mean (i.e., for a given stellar mass and redshift). No significant correlation between those quantities is seen in either stellar mass bin.

Finally, the bottom panels in Figure 3 show f_int as a function of SFE (relative to the mean of the population at a given redshift and stellar mass). Up to z = 2, we find a clear trend of galaxies with a higher-than-average SFE residing more frequently in interacting systems, specifically in the lower mass bin. At higher redshifts, this trend is less significant. As discussed in Section 5, theoretical works indicate that the SFE increase may be less correlated with interactions prior to cosmic noon. Furthermore, as argued later, disk fragmentation in gas-rich environments at high redshifts could play a more significant role in increasing the galaxies' SFE.

Overall, we found that the fraction of interacting galaxies is increased in the starburst regime out to z ∼ 4 and correlates with SFE but not gas fraction (at least out to z = 2). Galaxy–galaxy interactions therefore represent a veritable way to push galaxies into the starburst region at redshifts beyond cosmic noon. However, our analysis also shows that interacting galaxies only make up at most 40% of the z < 3 starburst galaxies in our sample. This implies that noninteracting systems contribute significantly to the starburst population. The increase of star formation (efficiency) through disk instabilities could provide another avenue for galaxies to reach the starburst regime (see Section 1). In this case, we would expect a higher fraction of disk galaxies with pronounced star-forming clumps in the starburst regime. This is studied in the next section.

4.2. The Impact of Disk Instabilities on Starburst Activity

A significant fraction (50%–80%) of galaxies in the starburst regime in our sample are noninteracting (i.e., isolated) disk galaxies. These may have been interacting in the past; however, taking the current morphological evidence at face value shows that they are currently not in a merging state. A possible way to reach the starburst regime is through an increased SFE due to the instability in gas-rich disks. The idea (see, for example, A. Dekel et al. 2009; A. B. Romeo & K. Fathi 2016) here is that gas-rich streams increase the gas density of the disk, which then becomes unstable and starts to fragment into clumps (Toomre instability; A. Toomre 1964). These clumps can contribute to several percent of the total disk mass, and star formation is maintained in dense subclumps over timescales of several hundred Myr. Steady inbound gas streams on the disk maintain the instability of the disk and replenish gas over several Gyr. Eventually, the clumps might migrate toward the center due to dynamical friction and may form spheroid-dominated galaxies later on. Recent measurements with JWST using multiple broad and medium bands suggest that the majority of clumps are less than ∼200 Myr old, suggesting a similarly long survival time (e.g., A. Claeyssens et al. 2023).

The occurrence of UV-bright star-forming clumps in galaxies is ubiquitous at high redshifts (about 60% of z = 1–3 galaxies contain bright UV clumps; L. L. Cowie et al. 1995; C. J. Conselice et al. 2004; B. G. Elmegreen et al. 2013; Y. Guo et al. 2015; E. Soto et al. 2017; A. Zanella et al. 2019). Several studies support the formation of these clumps through in situ physical processes based on differences in the evolution of the clump fraction and minor/major mergers; the disk nature of their host galaxies and, similarly, the kinematic properties of ordered disk rotation with high velocity dispersion; the distribution of clumps with scale height for edge-on galaxies; and the stellar mass function being similar to local star clusters and H ii regions (see T. Shibuya et al. 2016; B. G. Elmegreen et al. 2017; M. Dessauges-Zavadsky & A. Adamo 2018; M. Girard et al. 2020). A recent study of lensed clumpy galaxies at z = 1 (M. Dessauges-Zavadsky et al. 2019, 2023) corroborates the picture in which clumps are formed in situ through disk instabilities.

Along these lines, B. Wang & J. Silk (1994) present a simple recipe to link star formation to disk instability. In this simple theoretical model, the SFE (defined here as SFE/M_gas) is proportional to

$$\begin{eqnarray}&&{\rm{SFE}}\equiv \frac{{\rm{SFR}}}{{M}_{{\rm{gas}}}}\propto \frac{{(1-{Q}^{2})}^{1/2}}{Q},\end{eqnarray} \tag{ 1 }$$

where $Q=\frac{\kappa {v}_{{\rm{gas}}}}{\pi G{{\rm{\Sigma }}}_{{\rm{gas}}}}$ is the Toomre disk instability parameter, with κ the epicyclic frequency, v_g the radial cloud velocity dispersion, ${{\rm{\Sigma }}}_{{\rm{gas}}}=\frac{{M}_{{\rm{gas}}}}{\pi {r}_{{\rm{disk}}}^{2}}$ the gas surface density, and G the gravitational constant. A similar expression can be derived from the analytical model presented in A. Dekel et al. (2009; see their Equation (47)). Equation (1) is valid for an unstable disk, i.e., Q < 1. According to this model, the more unstable the disk, the more stars per surface are formed. Consequently, SFE rises as 1/Q as Q → 0 (see Equation (1)). We emphasize the simplicity of this model (see, for example, the discussion in A. B. Romeo & J. Wiegert 2011; A. B. Romeo & N. Falstad 2013; S. E. Meidt et al. 2023) as stars can also be formed if Q > 1; however, Q ≪ 1 may be a condition of exceptional star formation. Note that by simply applying empirical correlations, we get

$$\begin{eqnarray}&&Q\propto {{\rm{\Sigma }}}_{{\rm{gas}}}^{-1}\propto \frac{{r}_{{\rm{disk}}}^{2}}{{M}_{{\rm{gas}}}}\propto \frac{{M}_{* }^{0.4}}{\mu {M}_{* }}\propto {M}_{* }^{-0.24},\end{eqnarray} \tag{ 2 }$$

where we have used $\mu =\frac{{f}_{{\rm{gas}}}}{1-{f}_{{\rm{gas}}}}\propto {M}_{* }^{-0.36}$ (L. J. Tacconi et al. 2018) and ${r}_{{\rm{disk}}}\propto {M}_{* }^{\alpha }$ with α ∼ 0.2 (e.g., L. Yang et al. 2021). We therefore would expect only a weak stellar mass dependence. In addition, we note that the dispersion Q ∝ v_g may be weakly positively correlated with stellar mass (H. Übler et al. 2019), making the dependence of Q on stellar mass even weaker. If disk instabilities are at work to push disk galaxies into the starburst regime, we would expect an increase of galaxies with a number of dense star-forming clumps embedded in their disks for higher Δ_MS. That fraction would not significantly depend on stellar mass.

The top panel of Figure 4 shows the fraction of clumpy disk galaxies (with respect to the isolated disk galaxy population) as a function of Δ_MS for the same stellar mass bins as in Figure 3 and three redshift bins. We indeed observe an increase in the clumpy disk fraction toward the starburst regime. In the middle and bottom panels of Figure 4, we provide a comparison of f_clumpy with gas fraction and SFE normalized to the population mean per bin. We find that the fraction of clumpy disks increases with both f_gas and SFE; however, the correlation with the latter is more pronounced. As expected, we do not find significant differences in that behavior for the different stellar mass bins (see Equation (2)), except for the comparison between f_clumpy and f_gas, for which we do not find a significant correlation for the lower mass bin. We find some weak trend of f_clumpy with redshifts (as shown by the different panels), which is in agreement with studies of larger galaxy samples (Y. Guo et al. 2015).

Zoom In Zoom Out Reset image size

Figure 5 shows the above results more quantitatively by indicating how much f_gas and SFE differ (expressed in σ values) between disk galaxies on the main sequence and in the starburst regime (defined as >10× above the main sequence). For this test, we mass-match the main-sequence and starburst disk sample to ${\rm{log}}(M/{M}_{\odot })=10.8\pm 0.3$ to remove potential dependence on stellar mass. Averaged over redshifts up to z = 2, the difference between f_gas of main-sequence and starburst galaxies is <3σ, while they differ by >4σ in SFE. This suggests that disk galaxies appear in the starburst population mainly due to an increased SFE (achieved through disk instabilities) and not an increase in the gas fraction. An increase in f_gas may contribute at higher redshifts as the difference in SFE and f_gas becomes comparable between two types of galaxy populations.

Zoom In Zoom Out Reset image size

In this work, we have collected a sample of galaxies to investigate the properties of starbursts out to z ∼ 4 in terms of their morphology, gas fraction, and SFEs. Our sample is well fit for this investigation, as all galaxies have robust total (UV+far-IR) SFRs, gas mass measurements, and high-resolution optical-to-near-IR imaging from JWST. In summary, these are our two main findings.

• 1.  
  We find a significant enhancement of the fraction of interacting galaxies in the starburst regime out to z ∼ 4. This fraction correlates with SFE (but less with f_gas) out to z ∼ 2, but no significant correlation is found in the range 2 < z < 4.
• 2.  
  A large fraction of starburst galaxies in our sample are noninteracting disks (50%–80% depending on redshift). We find that the fraction of clumpy disk galaxies with respect to the disk population is enhanced in the starburst regime. This increase correlates more significantly with SFE than f_gas.

We note that with our visual classification, we are most sensitive to major to equal-mass mergers. This is suggested by the flux ratios of the merging components in the red filters, which are tracing their stellar mass. As mentioned earlier, kinematic spectroscopic measurements and deeper imaging data would likely be necessary to identify minor mergers robustly. The interacting systems identified here therefore do not contribute to the more continuous accretion of mass and gas as minor mergers would do (e.g., N. Scoville et al. 2023).

The former of these findings suggests that galaxy–galaxy interactions are a way to push galaxies into the starburst regime by increasing the efficiency of star formation in the galaxies. This has been observed at lower redshifts and is suggested by simulations (D. B. Sanders & I. F. Mirabel 1996; R. Genzel et al. 2010; J. S. Kartaltepe et al. 2012; C.-L. Hung et al. 2013; P. F. Hopkins et al. 2018; A. Cibinel et al. 2019; F. Renaud et al. 2019; J. Moreno et al. 2021; E. A. Shah et al. 2022; S. Wilkinson et al. 2022; H. He et al. 2023). Specifically, we showed that the fraction of interacting galaxies increases with increasing Δ_MS out to z ∼ 4. However, we also find that the interacting fraction does correlate with SFE only out to z = 2, while the correlation becomes weaker in the earlier Universe and lower stellar masses. This is also indicated by some simulation results (e.g., J. Fensch et al. 2017; D. R. Patton et al. 2020) and could suggest a weaker correlation between interaction and SFE increase, or that such a correlation occurs on a shorter timescale, before the cosmic noon. It also suggests that galaxy–galaxy interactions may not be the prime avenue to cause starbursts at z > 2—instead (as shown in Section 4.2), disk instabilities could be a more common path as f_clumpy clearly correlates with Δ_MS out to at least z = 4.

We note that even in the case of a strong one-to-one correlation between f_int and SFE, this signal could be diluted, leading to no significant correlation between these two quantities. The main reason for this is the different timescales probed here. Specifically, we would expect an increase in star formation (and hence SFE, assuming that little changes in the gas mass) triggered through the compression of gas upon close interaction within a few hundred Myr (J. M. Lotz et al. 2008; R. Teyssier et al. 2010; E. Cenci et al. 2024). Star formation traced by optical emission line indicators (e.g., Hα) may provide the closest link between SFE increase and close interactions. In fact, simulations show that infrared SFR tracers can lag behind as far-IR light is emitted mostly after coalescence when sufficient amounts of metals and dust have been produced (e.g., J. M. Lotz et al. 2008). Therefore, while a strong increase in SFE can be observed in close mergers, the late- and early-stage mergers (averaging over many hundreds of Myr) selected in our classification may not directly trace an SFE increase (see also R. Teyssier et al. 2010; E. A. Shah et al. 2022). Thus, averaged over the whole sample of interacting galaxies, the correlation between f_int and SFE can be significantly diluted due to the mismatch in timescales and the short time span of increased SFE (or likewise SFR) during interactions.

The significant fraction of noninteracting disk galaxies in our sample defined as starbursts indicates another avenue other than galaxy–galaxy interactions at work to offset galaxies from the main sequence (see starburst rotators found at z = 2; e.g., N. M. Förster Schreiber et al. 2011). Other studies (e.g., S. Wilkinson et al. 2022) have found similar fractions, and in this work, we have studied the properties of these galaxies. As shown by Equation (1) and also by recent simulation work (E. Cenci et al. 2024), disk instabilities can cause an increase in SFE, which can push galaxies above the main sequence. We expect that this effect may be largely independent of stellar mass Equation (2). The enhanced fraction of clumpy disk galaxies in the starburst population of isolated disk galaxies in our sample and the coupling between f_clumpy and SFE are direct observational indications of disk instabilities in action.

It is interesting that while the fraction of clumpy disks shows a significant correlation with SFE, the fraction of interacting galaxies does not. This could suggest that the SFE increase due to disk instabilities is maintained over longer timescales (compared to merger-driven enhancement as discussed above); thus, it is statistically easier to observe a correlation. The factor of ∼3 higher fraction of galaxies with clumps at high Δ_MS compared to interacting galaxies is an additional indication of a longer sustained increase of SFE through disk instabilities. As mentioned above, recent JWST observations suggest that the clump survival time is likely less than 200 Myr. If the SFE increase is indeed maintained over longer timescales, this would mean that the gas in the disk needs to be replenished on similar timescales as the clump lifetime in order to maintain disk fragmentation.

More detailed analyses, including kinematic measurements with JWST or ALMA, are necessary to study the effect of interactions on starburst galaxies in more detail. Specifically, Hα observations would provide an important constraint on the more instantaneous SFE. At higher redshifts, a statistical study of star-forming clumps in starburst disk galaxies would reveal their properties, specifically their ages and star formation densities, similarly to what has been achieved at z = 2 by multiple integral field unit programs (e.g., N. M. Förster Schreiber & S. Wuyts 2020).

Further characterization of star-forming clumps in larger galaxy samples will additionally become important in understanding the role of disk instabilities in growing the most massive galaxies at the highest redshifts. A possible model to explain the high SFEs needed to produce the most massive galaxies during the Epoch of Reionization (e.g., Y. Harikane et al. 2023; C. M. Casey et al. 2024; M. Xiao et al. 2024) are localized regions of essentially feedback-free star formation (e.g., A. Dekel et al. 2023; Z. Li et al. 2024). This model specifically works in low-metallicity regions that are self-shielded from significant external stellar winds. In such a regime, stars are efficiently formed before the feedback introduced by supernovae starts (for more details, see A. Dekel et al. 2023). A way to form these dense pockets of star formation could be through violent disk instabilities in gas-rich galaxies fed by pristine low-metallicity gas. Interestingly, such a scenario may have been recently observed in a lensed z = 6 galaxy hosting 15 star-forming clumps (S. Fujimoto et al. 2024). The high gas surface density of the clumps (Σ_gas ∼ 10³⁻⁵ M_⊙ pc⁻²) suggests an overall SFE (assuming star formation takes place in the clumps) that is 6–9× ​​​​​​higher (H. Fukushima & H. Yajima 2021) than for typical main-sequence galaxies with densities Σ_gas < 10³ M_⊙ pc⁻², assuming star formation takes place in their extended disks (e.g., L. J. Tacconi et al. 2018; M. Dessauges-Zavadsky et al. 2020). The contribution of the clumps to the integrated star formation is expected to be significant (>50%) in this case, in contrast to other studies finding contributions in the low tens of percent in galaxies after the cosmic noon (J. A. Hodge et al. 2016; A. Cibinel et al. 2017; D. Elbaz et al. 2018). The integrated f_gas of that galaxy is ∼0.75, thus not significantly higher than main-sequence galaxies at a similar redshift (0.4–0.8; M. Dessauges-Zavadsky et al. 2020). If disk instabilities play an important role in increasing the SFE (as shown by our work and recently by S. Fujimoto et al. 2024), the most massive galaxies in the early Universe (although they look compact) may show similar structures on resolved scales of ∼100 pc (unless other formation channels are present, such as monolithic collapse; e.g., O. J. Eggen et al. 1962). Such scales can currently only be probed by gravitational lensing; however, later, this may be possible with diffraction-limited observations using the next generation of large telescopes such as the Thirty Meter Telescope or the European Extremely Large Telescope.

The JWST data presented in this article were obtained from the Mikulski Archive for Space Telescopes (MAST) at the Space Telescope Science Institute. The specific observations analyzed can be accessed via doi:10.17909/ph8h-qf05.

We thank the anonymous referee for the suggestions, which substantially improved this manuscript. Support for this work was provided by NASA grants JWST-GO-01727 and HST-AR15802 awarded by the Space Telescope Science Institute, operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS 5-26555. C.M.C. thanks the National Science Foundation for support through grants AST-1814034 and AST-2009577 as well as the University of Texas at Austin College of Natural Sciences for support. C.M.C. also acknowledges support from the Research Corporation for Science Advancement from a 2019 Cottrell Scholar Award sponsored by IF/THEN, an initiative of Lyda Hill Philanthropies. The French part of the COSMOS team is partly supported by the Centre National d'Etudes Spatiales (CNES). O.I. acknowledges the funding of the French Agence Nationale de la Recherche for the project iMAGE (grant ANR-22-CE31-0007). This work was made possible by utilizing the CANDIDE cluster at the Institut d'Astrophysique de Paris. The cluster was funded through grants from the PNCG, CNES, DIM-ACAV, the Euclid Consortium, and the Danish National Research Foundation Cosmic Dawn Center (DNRF140). It is maintained by Stephane Rouberol. S.G. acknowledges financial support from the Villum Young Investigator grants 37440 and 13160 and the Cosmic Dawn Center (DAWN), funded by the Danish National Research Foundation (DNRF) under grant DNRF140.

Facilities: ALMA - Atacama Large Millimeter Array, Herschel - European Space Agency's Herschel space observatory, JWST - James Webb Space Telescope.

Software: astropy (Astropy Collaboration et al. 2013, 2018); PhotUtils (L. Bradley et al. 2024); SExtractor (E. Bertin & S. Arnouts 1996); statmorph (V. Rodriguez-Gomez et al. 2019); WebbPSF (M. D. Perrin et al. 2014).

For this work, we performed a visual classification of our galaxies in mainly disks, clumpy disks, and interactive systems. Here, we compare this visual classification with a statistical morphological classification. Specifically, the Gini (G) versus M₂₀ (e.g., J. M. Lotz et al. 2004) or concentration (C) versus asymmetry (A) (e.g., M. A. Bershady et al. 2000) diagnostics are some of the most-used classifiers of interacting and noninteracting galaxy systems.

We use statmorph (V. Rodriguez-Gomez et al. 2019) to measure these nonparametric morphological quantities. We first create a mask for each source to remove contaminants. A mask for each filter was generated using the detect_sources function provided by the PhotUtils (L. Bradley et al. 2024) Python package and requiring a minimum of 20 connected pixels and a threshold of >2.5σ above background. The masks are then carefully combined to a "master mask" to optimally enclose the flux of the source and remove contaminants. We then run statmorph on all galaxies for all their available bands using the masks as created above.

The left panel of Figure A1 shows the resulting G–M₂₀ diagram for the galaxies at rest-frame ∼1 μm. A commonly used differentiation between merging and nonmerging galaxies (e.g., J. M. Lotz et al. 2008) is indicated as a dashed line. It can be seen that 60%–70% of visually classified interacting systems are recovered. This number can be improved slightly by optimizing the statistical selection criteria. However, while the completeness is high, the purity is not. Specifically, as shown in the figure, the clumpy disk galaxies occupy a similar space as the interacting systems. On the other hand, some interacting systems would not be classified as mergers (e.g., panel (B1) in Figure A1). We notice that these are usually interacting galaxies identified visually by their tidal tail structures. These may be mergers in the final approach. We also investigated the C–A classification diagram (M. A. Bershady et al. 2000). However, we found that the separation between interacting and noninteracting is less clear; therefore, we do not go into further details here.

Zoom In Zoom Out Reset image size

Panels (A1)–(A3) of Figure A1 show an example of a clumpy disk galaxy. It would not be identified as such using the G–M₂₀ statistical classification as it is too similar to an isolated disk galaxy. Panel (A1) shows the F150W image. Panel (A2) shows the same image but with a subtracted Sérsic model fit by statmorph, revealing the disk structure and clumps in the residual. Panel (A3) shows the F444W image subtracted from the PSF-matched F150W image. The PSF matching was achieved with the PhotUtils function create_matching_kernel() using the WebbPSF (M. D. Perrin et al. 2014) theoretical PSFs of the two bands. The star-forming clumps are clearly identifiable in the Sérsic-subtracted residual image and the color image (as blue clumps). This is an example of the technique used for visual classification and identification of clumpy disk galaxies. A more in-depth study of these clumps and their identification will be provided in a follow-up paper (B. Kalita et al. 2025, in preparation).

All in all, this analysis shows that our visual classification of interacting systems is robust. Furthermore, we believe that a statistical morphological identification would not be as reliable as our visual classification by members of our team, given the difficulties in selecting clumpy disk galaxies. This is especially true for such a relatively small sample of galaxies. A more sophisticated treatment of this selection using machine learning algorithms is therefore also beyond the scope of this work.