The emergence of the Star Formation Main Sequence with redshift unfolded by JWST

P. Rinaldi Kapteyn Astronomical Institute, University of Groningen, P.O. Box 800, 9700AV Groningen, The Netherlands Steward Observatory, University of Arizona, 933 North Cherry Avenue, Tucson, AZ 85721, USA R. Navarro-Carrera Kapteyn Astronomical Institute, University of Groningen, P.O. Box 800, 9700AV Groningen, The Netherlands K. I. Caputi Kapteyn Astronomical Institute, University of Groningen, P.O. Box 800, 9700AV Groningen, The Netherlands Cosmic Dawn Center (DAWN), Copenhagen, Denmark E. Iani Kapteyn Astronomical Institute, University of Groningen, P.O. Box 800, 9700AV Groningen, The Netherlands G. Östlin Department of Astronomy, Stockholm University, Oscar Klein Centre, AlbaNova University Centre, 106 91 Stockholm, Sweden L. Colina Centro de Astrobiología (CAB), CSIC-INTA, Ctra. de Ajalvir km 4, Torrejón de Ardoz, E-28850, Madrid, Spain Cosmic Dawn Center (DAWN), Denmark S. Alberts Steward Observatory, University of Arizona, 933 North Cherry Avenue, Tucson, AZ 85721, USA J. Álvarez-Márquez Centro de Astrobiología (CAB), CSIC-INTA, Ctra. de Ajalvir km 4, Torrejón de Ardoz, E-28850, Madrid, Spain M. Annunziatella Centro de Astrobiología (CAB), CSIC-INTA, Ctra. de Ajalvir km 4, Torrejón de Ardoz, E-28850, Madrid, Spain INAF-Osservatorio Astronomico di Capodimonte, Via Moiariello 16, I-80131 Napoli, Italy L. Boogaard Max-Planck-Institut für Astronomie, Königstuhl 17, 69117 Heidelberg, Germany L. Costantin Centro de Astrobiología (CAB), CSIC-INTA, Ctra. de Ajalvir km 4, Torrejón de Ardoz, E-28850, Madrid, Spain J. Hjorth DARK, Niels Bohr Institute, University of Copenhagen, Jagtvej 128, 2200 Copenhagen, Denmark D. Langeroodi DARK, Niels Bohr Institute, University of Copenhagen, Jagtvej 128, 2200 Copenhagen, Denmark J. Melinder Department of Astronomy, Stockholm University, Oscar Klein Centre, AlbaNova University Centre, 106 91 Stockholm, Sweden T. Moutard Aix Marseille Univ, CNRS, CNES, LAM, Marseille, France European Space Agency (ESA), European Space Astronomy Centre (ESAC), Camino Bajo del Castillo s/n, 28692 Villanueva de la Canada, Madrid, Spain F. Walter Max-Planck-Institut für Astronomie, Königstuhl 17, 69117 Heidelberg, Germany
Abstract

We investigate the correlation between stellar mass (M) and star formation rate (SFR) across the stellar mass range log10(M/M)611subscript10subscriptMsubscriptMdirect-product611\log_{10}(\rm M_{\star}/M_{\odot})\approx 6-11roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≈ 6 - 11. We consider almost 50,000 star-forming galaxies at z37𝑧37z\approx 3-7italic_z ≈ 3 - 7, leveraging data from COSMOS/SMUVS, JADES/GOODS-SOUTH, and MIDIS/XDF. This is the first study spanning such a wide stellar mass range without relying on gravitational lensing effects. We locate our galaxies on the SFRMSFRsubscriptM\mathrm{SFR-M_{\star}}roman_SFR - roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT plane to assess how the location of galaxies in the star-formation main sequence (MS) and starburst (SB) region evolves with stellar mass and redshift. We find that the two star-forming modes tend to converge at log10(M/M)<7subscript10subscriptMsubscriptMdirect-product7\log_{10}(\rm M_{\star}/M_{\odot})<7roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) < 7, with all galaxies found in the SB mode. By dissecting our galaxy sample in stellar mass and redshift, we show that the emergence of the star-formation MS is stellar-mass dependent: while in galaxies with log10(M/M)>9subscript10subscriptMsubscriptMdirect-product9\log_{10}(\rm M_{\star}/M_{\odot})>9roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) > 9 the MS is already well in place at z=57𝑧57z=5-7italic_z = 5 - 7, for galaxies with log10(M/M)78subscript10subscriptMsubscriptMdirect-product78\log_{10}(\rm M_{\star}/M_{\odot})\approx 7-8roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≈ 7 - 8 it only becomes significant at z<4𝑧4z<4italic_z < 4. Overall, our results are in line with previous findings that the SB mode dominates amongst low stellar-mass galaxies. The earlier emergence of the MS for massive galaxies is consistent with galaxy downsizing.

Galaxies: formation, evolution, high-redshift, star formation, starburst, Epoch of Reionization
journal: ApJfacilities: HST, JWSTsoftware: Astropy (Astropy Collaboration et al., 2018), LePHARE (Arnouts & Ilbert, 2011), NumPy (Harris et al., 2020), pandas (pandas development team, 2020) Photutils (Bradley et al., 2021), SciPy (Virtanen et al., 2020) SExtractor (Bertin & Arnouts, 1996), TOPCAT (Taylor, 2005).

.

1 Introduction

In recent decades, galaxy surveys up to very high redshifts have significantly advanced our understanding of galaxy evolution, constraining galaxy physical properties such as stellar mass (MsubscriptM\mathrm{M_{\star}}roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT) and star formation rates (SFR) (Le Floc’h et al. 2005; Ellis et al. 2013; Oesch et al. 2014; Bouwens et al. 2015; Stefanon et al. 2019; Bowler et al. 2020; Bhatawdekar & Conselice 2021; Bouwens et al. 2021). These quantities have been extensively used to constrain the process of gas conversion into stars, i.e., the stellar mass assembly (Casey et al. 2012; L’Huillier et al. 2012; Bauer et al. 2013; Jackson et al. 2020). Statistical analysis of extensive galaxy samples has established a correlation between MsubscriptM\mathrm{M_{\star}}roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT and SFR in star-forming galaxies, revealing the “Main Sequence (MS) of star-forming galaxies”, and identified a passive cloud of galaxies with negligible star formation activity (Brinchmann et al. 2004; Noeske et al. 2007). These initial works triggered a vast amount of later papers studying galaxy evolution on the SFRMSFRsubscriptM\mathrm{SFR-M_{\star}}roman_SFR - roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT plane (e.g., Peng et al., 2010; Speagle et al., 2014; Salmon et al., 2015; Santini et al., 2017).

The existence of a star formation MS suggests that similar mechanisms could be responsible for growing low- and high-mass galaxies alike (Noeske et al., 2007). The MS galaxies grow continuously over a long time period from smooth gas accretion (e.g., Sánchez Almeida et al., 2014). The position of a galaxy on the SFRMSFRsubscriptM\mathrm{SFR-M_{\star}}roman_SFR - roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT plane has been proposed to be strictly correlated with its evolutionary stage (e.g., Tacchella et al., 2016), while the intrinsic scatter of the MS suggests some variety in the star formation histories (SFHs) for galaxies of a given stellar mass (e.g., Matthee & Schaye, 2019).

Recent works have shown that the normalization of the SFRMSFRsubscriptM\mathrm{SFR-M_{\star}}roman_SFR - roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT relation increases over cosmic time, particularly at z03𝑧03z\approx 0-3italic_z ≈ 0 - 3, reflecting higher gas accretion rate and, therefore, higher SFR in the past (e.g., Whitaker et al., 2012, 2014; Iyer et al., 2018; Popesso et al., 2023). This relation is typically modeled as a power-law, log10(SFR) = α𝛼\alphaitalic_αlog10(M) + β𝛽\betaitalic_β, with the slope (α𝛼\alphaitalic_α) ranging between 0.6 and 1.0 (see Speagle et al. 2014 for a more detailed study).

In addition, other works have analysed the presence of galaxies with significantly enhanced star formation activity at high redshifts, the so-called starbursts (SB), on the SFRMSFRsubscriptM\mathrm{SFR-M_{\star}}roman_SFR - roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT plane. Unfortunately, there is no absolute consensus in the literature of the starburst definition. Some studies simply considered that SB galaxies are those sources with a very high SFR (order of 10100Myr110100subscriptMdirect-productsuperscriptyr1\rm 10-100\;M_{\odot}\,yr^{-1}10 - 100 roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT roman_yr start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT; see Muxlow et al. 2006; Heckman 2006) and therefore located several σ𝜎\sigmaitalic_σ above the MS, where σ𝜎\sigmaitalic_σ denotes the observed scatter of the MS relation (e.g., Elbaz et al. 2011; Rodighiero et al. 2011; Schreiber et al. 2015; Lee et al. 2017; Orlitova 2020). One could adopt instead a criterion more in line with the textbook definition of starburst (Heckman, 2006), which refers to the galaxy birthrate parameter: a galaxy is in a SB episode when its ongoing SFR is much higher than its average past SFR. However, constraining the average past SFR of a galaxy is not trivial, as it assumes knowing its star formation history and age.

More recently, an alternative starburst definition for high-redshift galaxies has been recently proposed: SB are galaxies with high specific SFRs (sSFRs), with (empirically determined) log10(sSFR/yr1)>7.60subscriptlog10sSFRsuperscriptyr17.60\rm log_{10}(sSFR/yr^{-1})>-7.60roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_sSFR / roman_yr start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT ) > - 7.60 (Caputi et al., 2017, 2021). The inverse of the sSFR is the stellar mass doubling time, with the above limit corresponding to values <4×107yrabsent4superscript107yr<4\times 10^{7}\,\rm yr< 4 × 10 start_POSTSUPERSCRIPT 7 end_POSTSUPERSCRIPT roman_yr, which is consistent with the typical starburst timescales observed in the local Universe (e.g., Knapen & James, 2009).

In any case, the origin of the SB phenomenon is not completely understood. Current theories propose that starbursts may be driven by large-scale gravitational instabilities, which are influenced by stellar self-gravity (Inoue et al., 2016; Romeo & Fathi, 2016; Tadaki et al., 2018). Alternatively, starbursts might result from multiple star formation bursts, often triggered by galaxy mergers (Lamastra et al., 2013; Calabrò et al., 2019).

Earlier assessments that SB galaxies played a minor role in cosmic star-formation history were based on data from relatively massive galaxies with M1010M{}_{\star}\gtrsim 10^{10}\,\rm M_{\odot}start_FLOATSUBSCRIPT ⋆ end_FLOATSUBSCRIPT ≳ 10 start_POSTSUPERSCRIPT 10 end_POSTSUPERSCRIPT roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT up to z2𝑧2z\approx 2italic_z ≈ 2 (Rodighiero et al., 2011; Sargent et al., 2012; Lamastra et al., 2013). However, this view has now been reconsidered following the advent of deeper datasets that allow for the study of less massive galaxies.

In this regard, Caputi et al. (2017), by studying a sample of galaxies with M109Mgreater-than-or-equivalent-tosubscriptMsuperscript109subscriptMdirect-product\rm M_{\star}\gtrsim 10^{9}\,M_{\odot}roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT ≳ 10 start_POSTSUPERSCRIPT 9 end_POSTSUPERSCRIPT roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT at z45𝑧45z\approx 4-5italic_z ≈ 4 - 5, discovered a significant fraction of SB and found that star-forming galaxies displayed a bimodal distribution (SB/MS) on the SFR-MsubscriptM\rm M_{\star}roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT plane. Later on, Bisigello et al. (2018) analyzed a deep sample of star-forming galaxies at z03𝑧03z\approx 0-3italic_z ≈ 0 - 3 and determined that the fraction of SB galaxies increases both with redshift and towards lower stellar masses. Building on these findings, Rinaldi et al. (2022) further confirmed the existence of the SB/MS bimodality across a wider range of redshifts (z36.5𝑧36.5z\approx 3-6.5italic_z ≈ 3 - 6.5), showing that this bimodality is independent of the sSFR tracer used (e.g., UV or Hα𝛼\alphaitalic_α). These previous studies indicate that the SB population is much more significant than previously thought, underscoring the importance of studying low-mass galaxies to fully unveil the role of SBs in galaxy evolution.

In this study we make use of ultra-deep JWST data to extend the study of the star-forming galaxy distribution on the SFRMSFRsubscriptM\mathrm{SFR-M_{\star}}roman_SFR - roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT plane down to very low stellar masses. We complement our analysis with COSMOS/SMUVS data to trace the high-mass end of the SFRMsubscriptM-\rm M_{\star}- roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT plane, allowing us to cover in total five decades in stellar mass (log10(M/M)611subscript10subscriptMsubscriptMdirect-product611\rm\log_{10}(M_{\star}/M_{\odot})\approx 6-11roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≈ 6 - 11). Our main goal is to understand when the star-formation MS appeared in cosmic time for galaxies of different stellar masses.

The structure of the paper is as follows. Section §2 provides an overview of our sample, comprising over 50,000 star-forming galaxies, and Section §3 details our methodology. In Section §4, we present our findings. Specifically, in Section §4.2, we show the SFRMSFRsubscriptM\rm SFR-M_{\star}roman_SFR - roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT plane across four redshift bins: z2.83.2𝑧2.83.2z\approx 2.8-3.2italic_z ≈ 2.8 - 3.2, z3.23.9𝑧3.23.9z\approx 3.2-3.9italic_z ≈ 3.2 - 3.9, z3.95𝑧3.95z\approx 3.9-5italic_z ≈ 3.9 - 5, and z57𝑧57z\approx 5-7italic_z ≈ 5 - 7 and validate the SB/MS bimodality at these cosmic epochs. Subsequently, we analyze the sSFR distribution. In Section §4.3, we study the emergence of the MS as a function of stellar mass. Finally, in Section §4.4, we analyze the emergence of the MS as a function of both cosmic time and stellar mass. Section §5 summarizes our key findings.

Throughout this paper, we consider a cosmology with H0=70kms1Mpc1subscript𝐻070kmsuperscripts1superscriptMpc1H_{0}=70\;\rm km\;s^{-1}\;Mpc^{-1}italic_H start_POSTSUBSCRIPT 0 end_POSTSUBSCRIPT = 70 roman_km roman_s start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT roman_Mpc start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT, ΩM=0.3subscriptΩ𝑀0.3\Omega_{M}=0.3roman_Ω start_POSTSUBSCRIPT italic_M end_POSTSUBSCRIPT = 0.3, and ΩΛ=0.7subscriptΩΛ0.7\Omega_{\Lambda}=0.7roman_Ω start_POSTSUBSCRIPT roman_Λ end_POSTSUBSCRIPT = 0.7. All magnitudes are total and refer to the AB system (Oke & Gunn, 1983). A Chabrier (2003) initial mass function (IMF) is assumed (0.1–100 M). We define SB as galaxies with sSFR>107.60yr1sSFRsuperscript107.60superscriptyr1\rm sSFR>10^{-7.60}\,yr^{-1}roman_sSFR > 10 start_POSTSUPERSCRIPT - 7.60 end_POSTSUPERSCRIPT roman_yr start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT (Caputi et al., 2017, 2021).

2 Dataset

In this work, we made use of data from two fields: the JWST Advanced Deep Extragalactic Survey (JADES) in the GOODS-SOUTH region (67.7 arcmin2) and the COSMOS/SMUVS survey (0.66 deg2).

2.1 JADES/GOODS-SOUTH

For this field, we considered the complete dataset from JADES/GOODS-SOUTH data release 2 (DR2), covering a total of 67.7 arcmin2. This dataset encompasses the Hubble eXtreme Deep Field (XDF), which includes ultra-deep MIRI data at 5.6 μ𝜇\muitalic_μm from MIRI Deep Imaging Survey (MIDIS; Östlin et al. 2024, in prep.).

2.1.1 MIRI data

In this work, we made use of the MIDIS/F560W ultra-deep observations carried out in December 2022 over Hubble XDF. These data will be described in Östlin et al. (2024, in prep.) and consist of 40absent40\approx 40≈ 40 hours on source taken in the HUDF, which allowed us to reach 28.6 mag (5σ𝜎\sigmaitalic_σ) for point-like sources measured in an r=0.arcsecond\farcsstart_ID start_POSTFIX SUPERSCRIPTOP . ′ ′ end_POSTFIX end_ID23 circular aperture (the radius being chosen to ensure a 70absent70\approx 70≈ 70% encircled energy). For a more comprehensive view on the MIRI data reduction, refer to Section 2.1.2 in Rinaldi et al. (2023a).

2.1.2 NIRCam data

We also considered the NIRCam imaging taken by the JWST Advanced Deep Extragalactic Survey, JADES (Eisenstein et al., 2023a), Data Release 2 (DR2; Eisenstein et al. 2023b), which includes observations from the JWST Extragalactic Medium-band Survey (JEMS; Williams et al. 2023) and the First Reionization Epoch Spectroscopically Complete Observations (FRESCO, Oesch et al. 2023). This dataset provides a total of 14 bands from 0.9 to 4.8 μ𝜇\muitalic_μm (6 at short-wavelength, SW, and 8 at long-wavelength, LW), with 5σ𝜎\sigmaitalic_σ depths ranging from 30.5 to 30.9 mag (measured in a 0.2″  radius circular aperture). We remark that JADES is one of the deepest NIRCam surveys on the sky, only matched in depth (in some bands) by the MIDIS/NIRCam-parallel project (Pérez-González et al., 2023) and The Next Generation Deep Extragalactic Exploratory Public Near-Infrared Slitless Survey, NGDEEP (Bagley et al., 2024), which means that we can have access to very low-mass galaxies that, prior JWST’s launch, were accessible only thanks to lensed fields (e.g., Santini et al. 2017; Rinaldi et al. 2022).

2.1.3 HST data

We obtained all the HST images over GOODS-SOUTH from the Hubble Legacy Field (HLF/GOODS-SOUTH). The HLF/GOODS-SOUTH provides 13 HST bands covering a wide range of wavelengths (0.2--1.6μ𝜇\muitalic_μm), from the UV (WFC3/ UVIS F225W, F275W, and F336W filters), optical (ACS/ WFC F435W, F606W, F775W, F814W, and F850LP filters), to near-infrared (WFC3/IR F098M, F105W, F125W, F140W and F160W filters). In this work, we only made use of the deepest ones (i.e., F435W, F606W, F775W, F814W, F850LP, F105W, F125W, F140W, F160W). We refer the reader to Whitaker et al. (2019) for a more detailed description of these observations111The HLF/GOODS-SOUTH) imaging is available at https://archive.stsci.edu/prepds/hlf/;.

2.2 COSMOS/SMUVS

As a complement, we also considered deep imaging data from the COSMOS field (Scoville et al., 2007). By design, these data allow us to explore a very different region in parameter space, as they cover a much wider area than the ones probed by JADES/GOODS-SOUTH (about 35×\approx 35\times≈ 35 × larger), but are about three magnitudes shallower. Therefore, these blank fields are useful to probe the high-mass end, which is not properly probed by smaller fields.

We leveraged the Spitzer Matching survey of the Ultra-VISTA ultra-deep stripes (SMUVS; Ashby et al. 2018), which utilized Spitzer’s IRAC at 3.6 and 4.5 μ𝜇\muitalic_μm across 0.66deg20.66superscriptdeg20.66\;\rm deg^{2}0.66 roman_deg start_POSTSUPERSCRIPT 2 end_POSTSUPERSCRIPT of the COSMOS field, integrating \approx 25 h/pointing and achieving 80% completeness at 25.5absent25.5\approx 25.5≈ 25.5 mag (Deshmukh et al. 2018). This survey overlaps with the deepest near-IR and optical data from UltraVISTA (McCracken et al. 2012) and Subaru (Taniguchi et al. 2007), respectively.

Our analysis employs the SMUVS galaxy catalogue, initially compiled by Deshmukh et al. (2018) and later updated by van Mierlo et al. (2022), encompassing 300,000absent300000\approx 300,000≈ 300 , 000 Spitzer sources with 28-band photometry from U𝑈Uitalic_U band to 4.5 μ𝜇\muitalic_μm. SED fitting was performed using LePHARE, utilizing synthetic templates from Bruzual & Charlot (2003) and a reddening law by Calzetti et al. (2000), as detailed in Deshmukh et al. (2018).

3 Photometry and SED fitting

In this section, we provide an overview of the photometry and spectral energy distribution (SED) fitting performed for our sources in JADES/GOODS-SOUTH. A comprehensive explanation of the photometry and SED fitting for sources in JADES/GOODS-SOUTH will be presented in Navarro-Carrera et al. (in prep.). For the SMUVS sources, detailed methodologies are referenced in Deshmukh et al. (2018) and van Mierlo et al. (2022).

3.1 Photometry

We adopted SExtractor (Bertin & Arnouts, 1996) to perform the source detection and extract the photometry across all HST and JWST images (with a pixel scale of 30 mas for both HST and NIRCam images and 60 mas for MIRI/F560W). We adopted a detection strategy similar to the one adopted in Rinaldi et al. (2023a, b) and Navarro-Carrera et al. (2024). The detection has been performed by using a super stack image combining F277W, F335M, F356W, F410M, and F444W.

We used SExtractor by considering a hot-mode configuration, following Galametz et al. (2013) which has been proved being well suited for detecting very faint sources. Particularly, we adopted a detection threshold of 2σ2𝜎2\sigma2 italic_σ and a minimum number of contiguous pixels of 9.

We recovered more than 85%percent8585\%85 % of the JADES DR2 official catalogue. Interestingly, when anti-crossmatching with the official JADES DR2 catalogue (Eisenstein et al. 2023b), we found a substantial number of real sources (visually inspected) which have not been reported in the official JADES DR2 catalogue (due to the different detection and deblending strategies, especially close to very bright and extended sources).

We built up our photometric catalogue following the same approach as the one adopted in Rinaldi et al. (2023a, b) and Navarro-Carrera et al. (2024). Briefly, we combined aperture-corrected photometry, adopting circular apertures (i.e., MAG_APER) of 0.5″  diameter, and Kron apertures (i.e., MAG_AUTOKron 1980). We chose a circular-aperture flux over a Kron flux when the sources were fainter than a given magnitude. In this case, as we were dealing with very deep images, we adopted mag=lim27\rm{}_{lim}=27start_FLOATSUBSCRIPT roman_lim end_FLOATSUBSCRIPT = 27 as our faint limit for the Kron aperture. The aperture correction for JWST images has been obtained by using the software WebbPSF (Perrin et al. 2014). When SExtractor, for a given source, was not able to recover any flux, we estimated an upper limit (3σ3𝜎3\sigma3 italic_σ): we placed random apertures (0.250.250.250.25″ṙadius) in a square box (20202020×\times× 20202020″) around each source and measured the local background, assigning -1 as an error. For those sources, instead, with no photometric coverage, we simply put -99. All fluxes have been corrected for Galactic extinction by using dustmaps (Green 2018). Finally, for each source, we imposed a minimum error of 0.05 mag to account for photometric calibration uncertainties and the errors being underestimated by SExtractor (Sonnett et al., 2013).

We double-checked our photometry with the one from the official JADES DR2 catalogue and found good agreement. See Navarro-Carrera et al. (in prep.) for more details.

3.2 SED fitting

Once we built up our photometric catalogue, we performed SED fitting for all our sources. We adopted a two-step process. We first derived photometric redshifts by using Eazy (Brammer et al., 2008) and, then, derived the stellar properties with LePHARE (Arnouts & Ilbert, 2011).

The Eazy SED fitting has been performed by adopting a linear combination of different templates. Particularly, we adopted the same templates as the ones used in Eisenstein et al. (2023b); Hainline et al. (2023). This approach allowed us to reach a superlative outlier fraction (8%absentpercent8\approx 8\%≈ 8 %), considering a large amount of spectroscopic redshifts (>1000absent1000>1000> 1000) in GOODS-SOUTH from different programs (e.g., JADES/NIRSpec MSA, 3D-HST, MUSE, and CANDELS/GOODS-SOUTH; Brammer et al. 2012; Guo et al. 2013; Bacon et al. 2023; Bunker et al. 2023).

Once we derived the photometric redshifts from Eazy and established the goodness of our results, we inferred the stellar properties of our sources with LePHARE by fixing the redshifts to the photometric estimate by Eazy and adopting the same approach as the one used in Rinaldi et al. (2023a, b). Briefly, we considered the stellar population synthesis (SPS) models proposed by Bruzual & Charlot (2003), based on the Chabrier IMF (Chabrier, 2003). We made use of two different star formation histories (SFHs): a standard exponentially declining SFH (known as “τ𝜏\tauitalic_τ-model”) and an instantaneous burst adopting a simple stellar population (SSP) model. We opted for two distinct metallicity values, a solar metallicity (Z = 0.02) and a fifth of solar metallicity (Z = 0.2Z = 0.004). We considered the Calzetti et al. (2000) reddening law in combination with Leitherer et al. (2002). In particular, we adopted the following color excess values: 0E(BV)1.50𝐸𝐵𝑉1.50\leq E(B-V)\leq 1.50 ≤ italic_E ( italic_B - italic_V ) ≤ 1.5, with a step of 0.1 mag.

A more detailed description of the adopted methodology for the SED fitting will be presented in Navarro-Carrera et al. (in prep.).

Refer to caption
Figure 1: The SFRMSFRsubscriptM\mathrm{SFR}-\mathrm{M_{\star}}roman_SFR - roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT plane, showcasing all sources (JADES/GOODS-SOUTH + COSMOS/SMUVS) analyzed in this study, divided in redshift bins as indicated. The pale blue region marks the lower envelope for SB galaxies, based on the criteria from Caputi et al. (2017, 2021). Fits for the MS and SB are derived from Rinaldi et al. (2022). The gray shaded area represents the SFR threshold derived from the 2σ𝜎\sigmaitalic_σ detection of the JADES images used in this work. The vertical dashed line in each panel refers to the stellar mass completeness (75%) of JADES sample at each redshift. The error bar showed in gray (upper right panel) indicate the median uncertainties on M and SFR. White contours are also presented to show the bimodality between MS and SB.

4 Results

We considered our independent determinations of MsubscriptM\rm M_{\star}roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT (from SED fitting) and SFRs (from rest-frame UV dust-corrected luminosities) to locate our sample of galaxies (JADES/GOODS-SOUTH) on the SFRMSFRsubscriptM\rm SFR-M_{\star}roman_SFR - roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT plane at z37𝑧37z\approx 3-7italic_z ≈ 3 - 7 complemented with the one from COSMOS/SMUVS, as shown in Figure 1. Additionally, we investigated the sSFR distribution of these galaxies, illustrated in Figure 2 and Figure 3. Finally, we inspected the evolution of the MS and SB percentages as a function of redshift at different stellar mass regimes (Figure 4).

4.1 Star Formation Rate estimates from the UV

We derived the SFRs for our sample independently of their SED fitting by considering their rest-frame UV luminosities (Lνsubscript𝐿𝜈L_{\nu}italic_L start_POSTSUBSCRIPT italic_ν end_POSTSUBSCRIPT). To do so, we estimated Lνsubscript𝐿𝜈L_{\nu}italic_L start_POSTSUBSCRIPT italic_ν end_POSTSUBSCRIPT at a reference wavelength λrestsubscript𝜆𝑟𝑒𝑠𝑡\lambda_{rest}italic_λ start_POSTSUBSCRIPT italic_r italic_e italic_s italic_t end_POSTSUBSCRIPT = 1500 Å from the photometry of every galaxy at the filter with the closest effective wavelength to λobs=λrest×(1+z)subscript𝜆𝑜𝑏𝑠subscript𝜆𝑟𝑒𝑠𝑡1𝑧\lambda_{obs}=\lambda_{rest}\times(1+z)italic_λ start_POSTSUBSCRIPT italic_o italic_b italic_s end_POSTSUBSCRIPT = italic_λ start_POSTSUBSCRIPT italic_r italic_e italic_s italic_t end_POSTSUBSCRIPT × ( 1 + italic_z ), where z𝑧zitalic_z is the photometric redshift of that galaxy. We corrected the UV fluxes for dust extinction following the Calzetti et al. (2000) reddening law in order to recover the intrinsic UV fluxes. In order to do so, we adopted E(BV)𝐸𝐵𝑉E(B-V)italic_E ( italic_B - italic_V ) values from the SED fitting analysis. Then we converted them into a monochromatic luminosity (Lνsubscript𝐿𝜈L_{\nu}italic_L start_POSTSUBSCRIPT italic_ν end_POSTSUBSCRIPT). Finally we obtained SFRUV using the prescription given by Kennicutt (1998):

SFR(Myr1)=1.4×1028Lν(ergs1Hz1).SFRsubscriptMdirect-productsuperscriptyr11.4superscript1028subscriptL𝜈ergsuperscripts1superscriptHz1\mathrm{SFR(M_{\odot}\;yr^{-1})}=1.4\times 10^{-28}\;\mathrm{L_{\nu}}\mathrm{(% erg}\;\mathrm{s^{-1}Hz^{-1}).}roman_SFR ( roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT roman_yr start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT ) = 1.4 × 10 start_POSTSUPERSCRIPT - 28 end_POSTSUPERSCRIPT roman_L start_POSTSUBSCRIPT italic_ν end_POSTSUBSCRIPT ( roman_erg roman_s start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT roman_Hz start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT ) . (1)

The conversion formula (Eq. 1) has an intrinsic scatter of 0.3 dex. Therefore, we propagated that error into our uncertainty on the SFR. Furthermore, Eq. 1 is based on a Salpeter IMF (Salpeter, 1955), while, in this work, we used a Chabrier one. Therefore, we scaled our SFRs from a Salpeter IMF to a Chabrier by considering a conversion factor (i.e., 1.55 as reported in Madau & Dickinson 2014).

4.2 The bimodality between MS and SB galaxies

Over recent decades, numerous studies have identified that star-forming galaxies predominantly lie along the so-called “Main-Sequence of star-forming galaxies” (Brinchmann et al., 2004; Elbaz et al., 2007; Speagle et al., 2014; Whitaker et al., 2014; Salmon et al., 2015). Only a subset of these sources was classified as starburst galaxies, ranging from about 2% (Rodighiero et al., 2011) to about 30% (Caputi et al., 2017), depending on the stellar mass cut.

However, these prior studies were limited by their observational depth. In this work, we leveraged deep observations from HST and JWST in GOODS-SOUTH, with the latter reaching 5σ5𝜎5\sigma5 italic_σ depth of about 30.5 mag (Eisenstein et al. 2023b). These ultra-deep observations enabled us, for the first time, to explore a region of the parameter space that, previously, was only accessible through the study of lensed fields (i.e., by exploiting the gravitational lensing effect). The latter, prior JWST’s launch, was the only way to probe very low-mass galaxies (e.g., Santini et al., 2017; Atek et al., 2018; Bhatawdekar & Conselice, 2021; Rinaldi et al., 2022). Notably, in this study, we expanded our dataset with data from SMUVS (Deshmukh et al., 2018), which allowed us to extend our study toward the high-mass end of the SFRMSFRsubscriptM\rm SFR-M_{\star}roman_SFR - roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT plane. Overall, our statistical sample consists of more than 50,000absent50000\approx 50,000≈ 50 , 000 galaxies, giving us the opportunity, for the first time, to study the SFRMSFRsubscriptM\rm SFR-M_{\star}roman_SFR - roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT plane over five decades in stellar mass across cosmic time (z37𝑧37z\approx 3-7italic_z ≈ 3 - 7) in blank fields.

Following the criteria established in Caputi et al. (2017, 2021), our analysis reveals that approximately 41%percent4141\%41 % of the galaxies in our sample is located along the MS (i.e., they have log10(sSFR/yr1)<8.05subscriptlog10sSFRsuperscriptyr18.05\rm log_{10}(sSFR/yr^{-1})<-8.05roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_sSFR / roman_yr start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT ) < - 8.05), while about 48%percent4848\%48 % of the sources is placed within the SB cloud (i.e., they have log10(sSFR/yr1)>7.60subscriptlog10sSFRsuperscriptyr17.60\rm log_{10}(sSFR/yr^{-1})>-7.60roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_sSFR / roman_yr start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT ) > - 7.60). The residual 11%percent1111\%11 % of our sample is located in the so-called “Star Formation Valley” (SFV). The SFV region in the SFRMSFRsubscriptM\rm SFR-M_{\star}roman_SFR - roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT plane is indicative of galaxies either transitioning from the MS to SB, experiencing a rejuvenation effect (Rosani et al., 2020; Zhang et al., 2023), or moving from SB back to MS, thereby returning to a steady state after a burst of intense star formation, probably triggered by disk instabilities or (minor and/or major) mergers.

Interestingly, when we apply the same mass cut to galaxies as used by Rodighiero et al. (2011), we find that SB galaxies constitute only 2%percent22\%2 % at z2.8greater-than-or-equivalent-to𝑧2.8z\gtrsim 2.8italic_z ≳ 2.8. This is in agreement with the findings reported by Rodighiero et al. (2011) at z2𝑧2z\approx 2italic_z ≈ 2, regardless of which sSFR cut is adopted.

All in all, this comparison suggests that focusing solely on massive galaxies (i.e., log10(M/M)10greater-than-or-equivalent-tosubscriptlog10subscriptMsubscriptMdirect-product10\rm log_{10}(M_{\star}/M_{\odot})\gtrsim 10roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≳ 10) might lead to a significant underestimation of the role of SB galaxies in the Cosmic Star Formation History, since the starburst phase in high mass galaxies is typically triggered by significant events such as major mergers, as highlighted by studies like Pearson et al. (2019) and Renaud et al. (2022).

Notably, until the present, the SB cloud has been defined only down to log10(M/M)9subscriptlog10subscriptMsubscriptMdirect-product9\rm log_{10}(M_{\star}/M_{\odot})\approx 9roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≈ 9 (Caputi et al., 2017, 2021). In Figure 1, we report an extrapolation toward lower stellar masses (log10(M/M)=6subscriptlog10subscriptMsubscriptMdirect-product6\rm log_{10}(M_{\star}/M_{\odot})=6roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) = 6) of the SB lower envelope (pale blue shade). Interestingly, from Figure 1, we can observe that this extrapolation nicely follows the separation in two clouds between MS and SB galaxies. Indeed, one should note that this SB lower envelope corresponds to stellar-mass doubling times of approximately 4×107yr4superscript107yr\rm 4\times 10^{7}\;yr4 × 10 start_POSTSUPERSCRIPT 7 end_POSTSUPERSCRIPT roman_yr, aligning with the typical timescales of local SB episodes as proposed in Knapen & James (2009). Therefore, classifying all galaxies with log10(sSFR/yr1)>7.60subscriptlog10sSFRsuperscriptyr17.60\rm log_{10}(sSFR/yr^{-1})>-7.60roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_sSFR / roman_yr start_POSTSUPERSCRIPT - 1 end_POSTSUPERSCRIPT ) > - 7.60 as SB galaxies is consistent across all stellar mass ranges, as demonstrated in Figure 1.

We then divided our galaxy sample into four distinct redshift bins, following the same approach as the one used in Rinaldi et al. (2022). Specifically, 27%percent2727\%27 % of our sample is located within the redshift range z2.83.2𝑧2.83.2z\approx 2.8-3.2italic_z ≈ 2.8 - 3.2, 32%percent3232\%32 % falls in the z3.23.9𝑧3.23.9z\approx 3.2-3.9italic_z ≈ 3.2 - 3.9 range, 24%percent2424\%24 % are in the z3.95𝑧3.95z\approx 3.9-5italic_z ≈ 3.9 - 5 interval, and the remaining 17%percent1717\%17 % is located within z57𝑧57z\approx 5-7italic_z ≈ 5 - 7 range. The above division in redshift bins has been adopted in order to probe the same amount of time (400absent400\approx 400≈ 400 Myr) in each redshift bin. Notably, we retrieve the bimodality between MS and SB in each redshift bin. This finding further supports the existence of this bimodality, as previously suggested in Caputi et al. (2017, 2021); Rinaldi et al. (2022).

As a sanity check, we also inspected all sources in our sample that show a photometric excess ascribable to the Hα𝛼\alphaitalic_α emission line by using the same approach as the one adopted in Rinaldi et al. (2023a); Caputi et al. (2023). We first estimated their SFRs(Hα𝛼\alphaitalic_α) and, then, their corresponding sSFR(Hα𝛼\alphaitalic_α). Finally, we looked at the sSFR distribution, confirming that the SB/MS bimodality is evident and does not depend on which SFR tracer it has been adopted, as it has been already shown in Rinaldi et al. (2022). A complete analysis of the Hα𝛼\alphaitalic_α emission of those sources will be presented in Navarro-Carrera et al. (in prep.).

Finally, by looking at Figure 1, we observe that, with some scatter, our sample consistently remains above a certain SFR threshold. The gray shaded area in each panel, indeed, represents the lowest SFRUV that we can probe at each redshift bin, using the JADES/GOODS-SOUTH data’s 2σ2𝜎2{\sigma}2 italic_σ detection limit. The vertical dashed line indicates where we achieve completeness in terms of stellar mass (75%percent7575\%75 % completeness).

Interestingly, below log10(M/M)<7subscriptlog10subscriptMsubscriptMdirect-product7\rm log_{10}(M_{\star}/M_{\odot})<7roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) < 7, there is no clear distinction between MS and SB galaxies; instead, galaxies appear to mostly populate the SB cloud. This trend could be explained either by an increased burstiness of star formation that becomes important as we progressively go to lower stellar masses (Atek et al. 2022, Navarro-Carrera et al., in prep), or by the current limitations in observational depth.

Refer to caption
Figure 2: The sSFR distribution of the entire sample (JADES/GOODS-SOUTH + COSMOS/SMUVS) divided in four distinct stellar mass bins. The entire plane is colour coded following the regions derived by Caputi et al. (2017): the star-formation MS for sSFR >>> 10-8.05 yr-1, the Starburst cloud for sSFR >>> 10-7.60 yr-1, and the Star Formation Valley for 10-8.05 yr-1 \leq sSFR \leq 10-7.60 yr-1. The sSFR distribution are color coded by age, as derived by LePHARE. To consider the different areas covered by JADES/GOODS-SOUTH (67.7 arcmin2) and COSMOS/SMUVS (0.66 deg2), we normalized the JADES/GOODS-SOUTH counts to match the COSMOS/SMUVS survey area, which is approximately 35 times larger than that of JADES/GOODS-SOUTH.

4.3 The emergence of the MS across the stellar mass range

Since our sample spans five decades in MsubscriptM\rm M_{\star}roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT within the SFRMSFRsubscriptM\rm SFR-M_{\star}roman_SFR - roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT plane, we investigated the evolution of the SB/MS bimodality as a function of MsubscriptM\rm M_{\star}roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT to identify whether there is a specific stellar mass regime where the “Main Sequence of star-forming galaxies” takes place.

Our approach involved dividing our sample into the following stellar mass bins: log10(M/M)7subscriptlog10subscriptMsubscriptMdirect-product7\rm log_{10}(M_{\star}/M_{\odot})\leq 7roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≤ 7, log10(M/M)78subscriptlog10subscriptMsubscriptMdirect-product78\rm log_{10}(M_{\star}/M_{\odot})\approx 7-8roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≈ 7 - 8, log10(M/M)89subscriptlog10subscriptMsubscriptMdirect-product89\rm log_{10}(M_{\star}/M_{\odot})\approx 8-9roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≈ 8 - 9, and log10(M/M)>9subscriptlog10subscriptMsubscriptMdirect-product9\rm log_{10}(M_{\star}/M_{\odot})>9roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) > 9.

Specifically, we analyzed the sSFR distribution within each stellar mass bin. To address the significant differences in survey area between our datasets, we normalized the JADES/GOODS-SOUTH (67.7 arcmin2) counts to match the COSMOS/SMUVS survey area (0.66 deg2), which is approximately 35 times larger than that covered by JADES/GOODS-SOUTH.

From Figure 2, we note that in the lowest mass regime (log10(M/M)7subscriptlog10subscriptMsubscriptMdirect-product7\rm log_{10}(M_{\star}/M_{\odot})\leq 7roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≤ 7), MS galaxies are almost absent, constituting less than 3%percent33\%3 % of the total population at those stellar masses. Interestingly, we find the same result when looking at the sample of Hα𝛼\alphaitalic_α-emitters at z2.8greater-than-or-equivalent-to𝑧2.8z\gtrsim 2.8italic_z ≳ 2.8 (Navarro-Carrera et al., in prep). By looking at their stellar ages222The ages for our galaxies directly come from LePHARE and they are purely based on the formation time as given by Bruzual & Charlot (2003) models (i.e., the models that we used to perform the SED fitting); (as indicated by the color bar on the right in Figure 2), we note that MS galaxies, albeit rare, preferentially show an older population compared to the starbursts in the same stellar mass bin. This suggests that MS galaxies may have already reached a steady state (e.g., Wang et al. 2019), with no significant feedback mechanisms triggering new star formation episodes (e.g., Renaud et al. 2019) and, therefore, young and massive stars are not entirely dominating their light.

Lastly, by looking at the HAEs that fall in this stellar mass regime, we find that galaxies show a very high Hα𝛼\alphaitalic_α equivalent width (EW) and are preferentially located in the starburst cloud, which confirms that these sources are going trough a violent episode of star formation. On the contrary, the low-mass HAEs that show a low EW(Hα𝛼\alphaitalic_α) preferentially lie along the MS and represent a very minor fraction of the entire population of HAEs (1%less-than-or-similar-toabsentpercent1\lesssim 1\%≲ 1 %) (Navarro-Carrera et al. in prep).

However, we caution the reader that the fact that galaxies preferentially appear to be in a SB phase at lower stellar masses (log10(M/M)7subscriptlog10subscriptMsubscriptMdirect-product7\rm log_{10}(M_{\star}/M_{\odot})\leq 7roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≤ 7) could be explained either by an increase in the burstiness of star formation, which becomes more significant as stellar mass decreases (Atek et al. 2022, Navarro-Carrera et al. in prep), or by the limitations of current observational depth.

Moving to higher stellar mass bins, MS galaxies (i.e., those sources with log10(sSFR/yr)<8.05subscriptlog10(sSFR/yr)8.05\text{log}_{10}\text{(sSFR/yr)}<-8.05log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT (sSFR/yr) < - 8.05) begin to emerge marking the onset of the “Main Sequence of star-forming galaxies”, as we observe for galaxies with stellar mass log10(M/M)78subscriptlog10subscriptMsubscriptMdirect-product78\rm log_{10}(M_{\star}/M_{\odot})\approx 7-8roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≈ 7 - 8, and the SB/MS bimodality becomes more prominent, especially at log10(M/M)89subscriptlog10subscriptMsubscriptMdirect-product89\rm log_{10}(M_{\star}/M_{\odot})\approx 8-9roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≈ 8 - 9.

As we move to the highest stellar mass bin (M109Mgreater-than-or-equivalent-tosubscriptMsuperscript109subscriptMdirect-product\rm M_{\star}\gtrsim 10^{9}\,M_{\odot}roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT ≳ 10 start_POSTSUPERSCRIPT 9 end_POSTSUPERSCRIPT roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT), MS galaxies dominate the entire sample, now accounting for approximately 64%percent6464\%64 % at those stellar masses, which favours the picture pointed out in Rodighiero et al. (2011) where SB galaxies becomes increasingly rare as we move toward higher mass (especially at log10(M/M)10greater-than-or-equivalent-tosubscriptlog10subscriptMsubscriptMdirect-product10\rm log_{10}(M_{\star}/M_{\odot})\gtrsim 10roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≳ 10). On the other hand, the HAEs that fall in this stellar mass regime are characterized by having a very low EW(Hα𝛼\alphaitalic_α), due to the strong anti-correlation between MsubscriptM\rm M_{\star}roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT and EW(Hα𝛼\alphaitalic_α) (already pointed out in Atek et al. (2022)), indicating a lack of intense star formation activity (Navarro-Carrera et al. in prep). The latter may explain why SB galaxies are less common at log10(M/M)>9subscriptlog10subscriptMsubscriptMdirect-product9\rm log_{10}(M_{\star}/M_{\odot})>9roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) > 9, now comprising only 26%percent2626\%26 % of the population, with the remaining 10%absentpercent10\approx 10\%≈ 10 % in the SFV region. Surprisingly, the fraction of SB galaxies primarily accumulates around log10(M/M)9subscriptlog10subscriptMsubscriptMdirect-product9\rm log_{10}(M_{\star}/M_{\odot})\approx 9roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≈ 9 and strongly diminishes at higher stellar masses.

This result, indeed, indicates that the starburst phenomenon, at very high stellar masses (e.g., M1010Mgreater-than-or-equivalent-tosubscriptMsuperscript1010subscriptMdirect-product\rm M_{\star}\gtrsim 10^{10}\,M_{\odot}roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT ≳ 10 start_POSTSUPERSCRIPT 10 end_POSTSUPERSCRIPT roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT), likely requires galaxies to undergo either a major merger event (although the merger fraction is still not well constrained at high redshift), which could trigger a violent episode of star formation activity (e.g., Cibinel et al. 2019; Renaud et al. 2022), or violent disk instabilities, gas accretion, and/or interactions (flybys and minor mergers) (e.g., Bournaud et al. 2012; Dannerbauer et al. 2017; Ho et al. 2019; Pan et al. 2019).

Interestingly, the results are independent of the assumed SFH. In addition to the initial assumptions discussed in Section §3.2, we explored also a delayed SFH with LePHARE. Despite a different SFH, the trends illustrated in Figure 2 remain overall unchanged. Finally, by following Iani et al. (2024) (see their Appendix C), we investigated whether our results are influenced by the introduction of a correction factor (k) to the UV-derived SFR adopted in this study. Specifically, we analyzed the theoretical evolution of the SFR(Hα𝛼\alphaitalic_α)/SFR(UV) ratio as a function of time on a log-log scale, considering the prescriptions of Kennicutt (1998). This was done for the two metallicities used in our work, namely solar and sub-solar. We employed the same SFH models adopted for the SED fitting of our sources, including single burst and τ𝜏\tauitalic_τ-models with τ𝜏\tauitalic_τ values of 0.001, 0.01, 0.03, 1, 2, 3, 5, 8, 10, and 15 Gyr, along with a constant SFH model. These models incorporate the assumptions from the BPASS models (Eldridge & Stanway, 2022), tailored for a Chabrier IMF, a cutoff mass of 100 MsubscriptMdirect-product\rm M_{\odot}roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT, and excluding binary stars. By correcting the SFRs for galaxies younger than 15 Myr—where the UV estimates of SFR are known to be significantly underestimated (e.g., Kennicutt 1998; Calzetti 2013)—we observed that the overall sSFR distribution remains unaffected.

4.4 The emergence of the MS over cosmic time

Refer to caption
Figure 3: The sSFR distribution of the entire sample (JADES/GOODS-SOUTH + COSMOS/SMUVS) divided, this time, in four distinct stellar mass bins and four redshift bins. Each column refers to a specific redshift bin, while each row refers to a specific stellar mass bin. All 16 panels are color coded following Caputi et al. (2017): the star-formation MS for sSFR >>> 10-8.05 yr-1, the Starburst cloud for sSFR >>> 10-7.60 yr-1, and the Star Formation Valley for 10-8.05 yr-1 \leq sSFR \leq 10-7.60 yr-1. Also in this case, as we did in Figure 2, to consider the different areas covered by JADES/GOODS-SOUTH (67.7 arcmin2) and COSMOS/SMUVS (0.66 deg2), we normalized the JADES/GOODS-SOUTH counts to match the COSMOS/SMUVS survey area, which is approximately 35 times larger than that of JADES/GOODS-SOUTH.

We then inspected the emergence of the MS galaxies as a function of cosmic time. For this purpose, as we already assumed in Figure 1, we divided our sample into four different redshift bins, each spanning approximately 400Myr400Myr\rm 400\;Myr400 roman_Myr: z2.83.2𝑧2.83.2z\approx 2.8-3.2italic_z ≈ 2.8 - 3.2, z3.23.9𝑧3.23.9z\approx 3.2-3.9italic_z ≈ 3.2 - 3.9, z3.95𝑧3.95z\approx 3.9-5italic_z ≈ 3.9 - 5, and z57𝑧57z\approx 5-7italic_z ≈ 5 - 7.

In Figure 3, we illustrate the sSFR distribution as a function of redshift, segmented into four stellar mass bins: we show 16 panels in total, with each row corresponding to a specific stellar mass bin and each column to a specific redshift bin.

As evident from Figure 3, there is a clear evolution of the SB/MS bimodality as a function of cosmic time. Going from the highest redshift bin (z57𝑧57z\approx 5-7italic_z ≈ 5 - 7) to the lowest one (z2.83.2𝑧2.83.2z\approx 2.8-3.2italic_z ≈ 2.8 - 3.2), the SB/MS bimodality evolves remarkably, showing that the onset of the MS over cosmic time clearly depends on the stellar mass regime considered.

Refer to caption
Figure 4: The evolution of the MS, SB, and SFV percentages with cosmic time in four stellar mass bins. Each panel has color bands corresponding to the four redshift bins analyzed in this study Also in this case, to consider the different areas covered by JADES/GOODS-SOUTH (67.7 arcmin2) and COSMOS/SMUVS (0.66 deg2), we normalized the JADES/GOODS-SOUTH counts to match the COSMOS/SMUVS survey area, which is approximately 35 times larger than that of JADES/GOODS-SOUTH.

Interestingly, as we move toward lower redshifts, SB galaxies become increasingly less common. Concurrently, we observe the onset of the “Main Sequence of star-forming galaxies”: for galaxies M108Mgreater-than-or-equivalent-tosubscriptMsuperscript108subscriptMdirect-product\rm M_{\star}\gtrsim 10^{8}\,M_{\odot}roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT ≳ 10 start_POSTSUPERSCRIPT 8 end_POSTSUPERSCRIPT roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT, the MS is already in place at z57𝑧57z\approx 5-7italic_z ≈ 5 - 7 (more prominently for M109Mgreater-than-or-equivalent-tosubscriptMsuperscript109subscriptMdirect-product\rm M_{\star}\gtrsim 10^{9}\,M_{\odot}roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT ≳ 10 start_POSTSUPERSCRIPT 9 end_POSTSUPERSCRIPT roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT). Notably, MS galaxies are nearly absent in the lowest stellar mass regime, a trend consistent across all redshift bins.

This evolution over cosmic time suggests that at lower redshifts, the SB phenomenon becomes less common as galaxies have already assembled most of their M. Significant events, such as major mergers, are thus required to trigger intense star formation episodes in massive galaxies. However, such events become increasingly rare as we move toward lower redshifts (z3less-than-or-similar-to𝑧3z\lesssim 3italic_z ≲ 3; e.g., Ventou et al. 2019).

This result is in line with what one should expect by looking at the Cosmic Star Formation Rate Density over comsic time (e.g., Lilly et al., 1996; Madau et al., 1996; Madau & Dickinson, 2014), where the SFR density (ρSFRsubscript𝜌SFR\rm\rho_{SFR}italic_ρ start_POSTSUBSCRIPT roman_SFR end_POSTSUBSCRIPT) starts declining at z2less-than-or-similar-to𝑧2z\lesssim 2italic_z ≲ 2. In particular, this result is also in agreement with what has been presented in Bisigello et al. (2018), where they clearly observe that at z3less-than-or-similar-to𝑧3z\lesssim 3italic_z ≲ 3 (see their Figure 8) the dominant mode of forming stars is the one that regulates the galaxies along the Main-Sequence, while SB galaxies are almost absent. Finally, in Figure 3, we also observe that within a fixed redshift bin, the SB/MS bimodality evolves as a function of stellar mass , as already shown in Figure 2.

The emergence of the Main Sequence of star-forming galaxies can be further observed in Figure 4, where we show the percentage evolution of MS, SB, and SFV galaxies across redshift and stellar mass bins. To estimate uncertainties in the percentages of MS, SB, and SFV galaxies, we employed Monte Carlo (MC) simulations following the approach used in Rinaldi et al. (2022). We generated 1000 mock catalogs from our initial dataset (JADES/GOODS-SOUTH + COSMOS/SMUVS), perturbing M and SFR within their error bars for each simulation run. Subsequently, these mock galaxies were divided according to the same four stellar mass and redshift bins analyzed in Figure 4, adopting the criteria for MS, SB, and SFV classifications. Then, for each bin, we built up a distribution of the percentages for MS, SB, and SFV galaxies. We then determined the 1σ𝜎\sigmaitalic_σ uncertainty as the half-distance between the 16th and 84th percentiles of each distribution. The errors were found to be up to a maximum of 5-6%.

Consistent with our analysis in Figures 2 and Figure 3, MS galaxies are almost absent at log10(M/M)7subscript10subscriptMsubscriptMdirect-product7\rm\log_{10}(M_{\star}/M_{\odot})\leq 7roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≤ 7. This may be due to either to the fact that very low-mass galaxies generally build up their MsubscriptM\rm M_{\star}roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT through more violent episodes of star formation (e.g., Asada et al. 2024), or due to observational constraints that prevent the detection of such low-mass MS galaxies with correspondingly low SFR (below the threshold imposed by the current depth of our observations).

As we consider higher stellar masses (log10(M/M)>7subscript10subscriptMsubscriptMdirect-product7\rm\log_{10}(M_{\star}/M_{\odot})>7roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) > 7), the Main Sequence starts taking place, becoming the dominant mode of forming stars at log10(M/M)>9subscript10subscriptMsubscriptMdirect-product9\rm\log_{10}(M_{\star}/M_{\odot})>9roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) > 9 at z37𝑧37z\approx 3-7italic_z ≈ 3 - 7. Figure 4 highlights again that the evolution of a galaxy along the MS is strictly correlated with MsubscriptM\rm M_{\star}roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT. The onset of the MS becomes evident well before the Cosmic Noon (z2𝑧2z\approx 2italic_z ≈ 2), particularly for galaxies with M>109MsubscriptMsuperscript109subscriptMdirect-product\rm M_{\star}>10^{9}\,M_{\odot}roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT > 10 start_POSTSUPERSCRIPT 9 end_POSTSUPERSCRIPT roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT, for which the MS is already in place by z57𝑧57z\approx 5-7italic_z ≈ 5 - 7. This result suggests that galaxies can start evolving through secular processes well before the peak of Cosmic Star Formation, as suggested both in the past (e.g., Smit et al., 2016) and more recently (Langeroodi & Hjorth, 2024). The transition from bursty to secular star formation is especially evident in higher-mass galaxies, which generally achieve a steady state of star formation at earlier epochs. This emphasizes that as galaxies build up their MsubscriptM\rm M_{\star}roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT, they increasingly favor secular evolution, because significant enhancements in their star formation activity typically require events like (minor and/or major) mergers or disk instabilities (e.g., Bournaud et al. 2012; Dannerbauer et al. 2017; Cibinel et al. 2019; Ho et al. 2019; Rodríguez Montero et al. 2019; Pan et al. 2019; Renaud et al. 2022). This implies that more massive galaxies are less likely to enter a SB phase without such significant events.

5 Summary and Conclusions

In this paper, we investigated the emergence of the Main Sequence of star-forming galaxies with stellar mass and redshift. For this purpose, we analyzed over 50,000 galaxies at z37𝑧37z\approx 3-7italic_z ≈ 3 - 7 on the SFRMSFRsubscriptM\rm SFR-M_{\star}roman_SFR - roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT plane, spanning a wide stellar mass range log10(M/M)611subscriptlog10MsubscriptMdirect-product611\rm log_{10}(M/M_{\odot})\approx 6-11roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≈ 6 - 11. This study has been made possible thanks to the joint analysis of ultra-deep JWST data in GOODS-SOUTH (FRESCO, JADES, JEMS, and MIDIS), as well as shallower COSMOS/SMUVS data over an area 35×\approx 35\times≈ 35 × larger, which allowed us to populate the high-mass end of the SFRMSFRsubscriptM\rm SFR-M_{\star}roman_SFR - roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT plane. This represents the first study that allows to reach very low stellar masses (log10(M/M)<7subscript10subscriptMsubscriptMdirect-product7\rm\log_{10}(M_{\star}/M_{\odot})<7roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) < 7) in blank fields without relying on gravitational lensing effects. We divided our sample into four redshift bins: z2.83.2𝑧2.83.2z\approx 2.8-3.2italic_z ≈ 2.8 - 3.2, z3.23.9𝑧3.23.9z\approx 3.2-3.9italic_z ≈ 3.2 - 3.9, z3.95𝑧3.95z\approx 3.9-5italic_z ≈ 3.9 - 5, and z57𝑧57z\approx 5-7italic_z ≈ 5 - 7. This division allowed us to probe similar amounts of time in each redshift bin (400absent400\approx 400≈ 400 Myr).

Our key findings are summarized as follows:

  • In agreement with previous results (Caputi et al. 2017; Rinaldi et al. 2022), we find a bimodality on the SFRMSFRsubscriptM\rm SFR-M_{\star}roman_SFR - roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT plane, such that the vast majority of star-forming galaxies lie either on the MS or SB cloud. This pattern is clearly observed for all galaxies with log10(M/M)>7subscript10subscriptMsubscriptMdirect-product7\mathrm{\log_{10}(M_{\star}/M_{\odot})>7}roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) > 7 across all redshifts (Figure 1).

  • Instead, at stellar masses log10(M/M)<7subscript10subscriptMsubscriptMdirect-product7\rm\log_{10}(M_{\star}/M_{\odot})<7roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) < 7, the two star-formation modes appear to converge in the SFRMSFRsubscriptM\rm SFR-M_{\star}roman_SFR - roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT plane (across all redshift bins), such that all low stellar-mass galaxies lie in the SB zone (Figure 2). Although our sample progressively loses completeness at such low stellar masses, we note that, at z2.93.2𝑧2.93.2z\approx 2.9-3.2italic_z ≈ 2.9 - 3.2, it is still 50% complete down to log10(M/M)6.3subscript10subscriptMsubscriptMdirect-product6.3\rm\log_{10}(M_{\star}/M_{\odot})\approx 6.3roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≈ 6.3, and 75% complete down to log10(M/M)6.65subscript10subscriptMsubscriptMdirect-product6.65\rm\log_{10}(M_{\star}/M_{\odot})\approx 6.65roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≈ 6.65. Therefore, this is not merely a selection effect. Actually, this is expected on physical grounds. If the star formation is driven by the stochastic collapse of giant molecular clouds, a low stellar-mass galaxy will by nature be bursty. For example, a single minor gas accretion event could lead to the formation of 107Msuperscript107subscriptMdirect-product10^{7}\,\rm M_{\odot}10 start_POSTSUPERSCRIPT 7 end_POSTSUPERSCRIPT roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT in 10 Myr, doubling the stellar mass of a low-mass galaxy. Instead, for a high stellar-mass galaxy, a coordinated global SF event will be required required to produce a SB (e.g., Gerola et al., 1980; Östlin et al., 2001; Atek et al., 2022).

  • As the stellar mass increases, the MS of star-forming galaxies becomes more prominent and dominant, especially for log10(M/M)>9subscriptlog10subscriptMsubscriptMdirect-product9\rm log_{10}(M_{\star}/M_{\odot})>9roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) > 9. Notably, at higher stellar masses (log10(M/M)8greater-than-or-equivalent-tosubscriptlog10subscriptMsubscriptMdirect-product8\rm log_{10}(M_{\star}/M_{\odot})\gtrsim 8roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≳ 8), the MS is already established by z4greater-than-or-equivalent-to𝑧4z\gtrsim 4italic_z ≳ 4. More importantly, for galaxies with log10(M/M)9greater-than-or-equivalent-tosubscriptlog10subscriptMsubscriptMdirect-product9\rm log_{10}(M_{\star}/M_{\odot})\gtrsim 9roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≳ 9, the MS is already in place at z57𝑧57z\approx 5-7italic_z ≈ 5 - 7, aligning with the concept of galaxy downsizing (Sparre et al. 2015; Franco et al. 2020). In contrast, galaxies with lower stellar masses (log10(M/M)8less-than-or-similar-tosubscriptlog10subscriptMsubscriptMdirect-product8\rm log_{10}(M_{\star}/M_{\odot})\lesssim 8roman_log start_POSTSUBSCRIPT 10 end_POSTSUBSCRIPT ( roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT / roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT ) ≲ 8) at these redshifts (z47𝑧47z\approx 4-7italic_z ≈ 4 - 7) are not yet on the MS and remain in the starburst phase (Figure 3). Thus, the emergence of the MS for star-forming galaxies at different redshifts heavily depends on MsubscriptM\rm M_{\star}roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT (Figure 4), suggesting that galaxies with higher stellar mass can achieve a steady state well before the Cosmic Noon (z2𝑧2z\approx 2italic_z ≈ 2; e.g., Smit et al. 2016; Langeroodi & Hjorth 2024).

In conclusion, current ultra-deep JWST observations indicate that very low-mass galaxies (107Mless-than-or-similar-toabsentsuperscript107subscriptMdirect-product\rm\lesssim 10^{7}\;M_{\odot}≲ 10 start_POSTSUPERSCRIPT 7 end_POSTSUPERSCRIPT roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT) predominantly experience bursty star formation, with only a few rare cases that can be classified as MS, typically showing signs of an evolved population. More importantly, our findings confirm that the emergence of the “Main Sequence of star-forming galaxies” in the SFRMSFRsubscriptM\rm SFR-M_{\star}roman_SFR - roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT plane is heavily dependent on MsubscriptM\rm M_{\star}roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT. The emergence of the “Main Sequence of star-forming galaxies” varies across cosmic time, with more massive galaxies (M109Mgreater-than-or-equivalent-tosubscriptMsuperscript109subscriptMdirect-product\rm M_{\star}\gtrsim 10^{9}\,M_{\odot}roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT ≳ 10 start_POSTSUPERSCRIPT 9 end_POSTSUPERSCRIPT roman_M start_POSTSUBSCRIPT ⊙ end_POSTSUBSCRIPT) establishing a more regulated, secular evolution already at z57𝑧57z\approx 5-7italic_z ≈ 5 - 7. This highlights the critical role of MsubscriptM\rm M_{\star}roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT in determining galaxy evolution.

Deeper and wider JWST observations will be instrumental in further constraining the emergence of the MS galaxies within the SFRMSFRsubscriptM\rm SFR-M_{\star}roman_SFR - roman_M start_POSTSUBSCRIPT ⋆ end_POSTSUBSCRIPT plane at different cosmic times.

The authors thank Nadav Peleg Brochstein for the useful discussion and tests on the effects of different SFHs on the stellar masses. This work is based on observations made with the NASA/ESA/CSA James Webb Space Telescope. The data were obtained from the Mikulski Archive for Space Telescopes at the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS 5-03127 for JWST. These observations are associated with programs GO #1963, GO #1895 and GTO #1283. The authors acknowledge the team led by coPIs C. Williams, M. Maseda and S. Tacchella, and PI P. Oesch, for developing their respective observing programs with a zero-exclusive-access period. Also based on observations made with the NASA/ESA Hubble Space Telescope obtained from the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS 5–26555. The specific observations analyzed can be accessed via:. (catalog DOI: 10.17909/gdyc-7g80, 10.17909/fsc4-dt61, 10.17909/fsc4-dt61, 10.17909/T91019, 10.17909/1rq3-8048, 10.17909/z2gw-mk31) K.I.C. and R.N.C. acknowledge funding from the Dutch Research Council (NWO) through the award of the Vici Grant VI.C.212.036. K.I.C. and E.I. acknowledge funding from the Netherlands Research School for Astronomy (NOVA). S.A. acknowledges support from the JWST Mid-Infrared Instrument (MIRI) Science Team Lead, grant 80NSSC18K0555, from NASA Goddard Space Flight Center to the University of Arizona. This work was supported by research grants (VIL16599, VIL54489) from VILLUM FONDEN. L.C. acknowledges support by grants PIB2021-127718NB-100 and PID2022-139567NB-I00 from the Spanish Ministry of Science and Innovation/State Agency of Research MCIN/AEI/10.13039/501100011033 and by “ERDF A way of making Europe”. J.A.-M., and L.C., acknowledge support by grant PIB2021-127718NB-100 from the Spanish Ministry of Science and Innovation/State Agency of Research MCIN/AEI/10.13039/501100011033 and by “ERDF A way of making Europe”. L.C. thanks the support from the Cosmic Dawn Center received during visits to DAWN as international associate .

References

  • Arnouts & Ilbert (2011) Arnouts, S., & Ilbert, O. 2011, LePHARE: Photometric Analysis for Redshift Estimate, , , ascl:1108.009
  • Asada et al. (2024) Asada, Y., Sawicki, M., Abraham, R., et al. 2024, MNRAS, 527, 11372
  • Ashby et al. (2018) Ashby, M. L. N., Caputi, K. I., Cowley, W., et al. 2018, ApJS, 237, 39
  • Astropy Collaboration et al. (2018) Astropy Collaboration, Price-Whelan, A. M., Sipőcz, B. M., et al. 2018, AJ, 156, 123
  • Atek et al. (2022) Atek, H., Furtak, L. J., Oesch, P., et al. 2022, MNRAS, 511, 4464
  • Atek et al. (2018) Atek, H., Richard, J., Kneib, J.-P., & Schaerer, D. 2018, MNRAS, 479, 5184
  • Bacon et al. (2023) Bacon, R., Brinchmann, J., Conseil, S., et al. 2023, A&A, 670, A4
  • Bagley et al. (2024) Bagley, M. B., Pirzkal, N., Finkelstein, S. L., et al. 2024, ApJ, 965, L6
  • Bauer et al. (2013) Bauer, A. E., Hopkins, A. M., Gunawardhana, M., et al. 2013, MNRAS, 434, 209
  • Bertin & Arnouts (1996) Bertin, E., & Arnouts, S. 1996, A&AS, 117, 393
  • Bhatawdekar & Conselice (2021) Bhatawdekar, R., & Conselice, C. J. 2021, ApJ, 909, 144
  • Bisigello et al. (2018) Bisigello, L., Caputi, K. I., Grogin, N., & Koekemoer, A. 2018, A&A, 609, A82
  • Bournaud et al. (2012) Bournaud, F., Juneau, S., Le Floc’h, E., et al. 2012, ApJ, 757, 81
  • Bouwens et al. (2015) Bouwens, R. J., Illingworth, G. D., Oesch, P. A., et al. 2015, ApJ, 803, 34
  • Bouwens et al. (2021) Bouwens, R. J., Oesch, P. A., Stefanon, M., et al. 2021, AJ, 162, 47
  • Bowler et al. (2020) Bowler, R. A. A., Jarvis, M. J., Dunlop, J. S., et al. 2020, MNRAS, 493, 2059
  • Bradley et al. (2021) Bradley, L., Sipőcz, B., Robitaille, T., et al. 2021, astropy/photutils: 1.0.2, v1.0.2, Zenodo, doi:10.5281/zenodo.4453725
  • Brammer et al. (2008) Brammer, G. B., van Dokkum, P. G., & Coppi, P. 2008, ApJ, 686, 1503
  • Brammer et al. (2012) Brammer, G. B., van Dokkum, P. G., Franx, M., et al. 2012, ApJS, 200, 13
  • Brinchmann et al. (2004) Brinchmann, J., Charlot, S., White, S. D. M., et al. 2004, MNRAS, 351, 1151
  • Bruzual & Charlot (2003) Bruzual, G., & Charlot, S. 2003, MNRAS, 344, 1000
  • Bunker et al. (2023) Bunker, A. J., Saxena, A., Cameron, A. J., et al. 2023, A&A, 677, A88
  • Calabrò et al. (2019) Calabrò, A., Daddi, E., Fensch, J., et al. 2019, A&A, 632, A98
  • Calzetti (2013) Calzetti, D. 2013, in Secular Evolution of Galaxies, ed. J. Falcón-Barroso & J. H. Knapen, 419
  • Calzetti et al. (2000) Calzetti, D., Armus, L., Bohlin, R. C., et al. 2000, ApJ, 533, 682
  • Caputi et al. (2017) Caputi, K. I., Deshmukh, S., Ashby, M. L. N., et al. 2017, ApJ, 849, 45
  • Caputi et al. (2021) Caputi, K. I., Caminha, G. B., Fujimoto, S., et al. 2021, ApJ, 908, 146
  • Caputi et al. (2023) Caputi, K. I., Rinaldi, P., Iani, E., et al. 2023, arXiv e-prints, arXiv:2311.12691
  • Casey et al. (2012) Casey, C. M., Berta, S., Béthermin, M., et al. 2012, ApJ, 761, 140
  • Chabrier (2003) Chabrier, G. 2003, PASP, 115, 763
  • Cibinel et al. (2019) Cibinel, A., Daddi, E., Sargent, M. T., et al. 2019, MNRAS, 485, 5631
  • Dannerbauer et al. (2017) Dannerbauer, H., Lehnert, M. D., Emonts, B., et al. 2017, A&A, 608, A48
  • Deshmukh et al. (2018) Deshmukh, S., Caputi, K. I., Ashby, M. L. N., et al. 2018, ApJ, 864, 166
  • Eisenstein et al. (2023a) Eisenstein, D. J., Willott, C., Alberts, S., et al. 2023a, arXiv e-prints, arXiv:2306.02465
  • Eisenstein et al. (2023b) Eisenstein, D. J., Johnson, B. D., Robertson, B., et al. 2023b, arXiv e-prints, arXiv:2310.12340
  • Elbaz et al. (2007) Elbaz, D., Daddi, E., Le Borgne, D., et al. 2007, A&A, 468, 33
  • Elbaz et al. (2011) Elbaz, D., Dickinson, M., Hwang, H. S., et al. 2011, A&A, 533, A119
  • Eldridge & Stanway (2022) Eldridge, J. J., & Stanway, E. R. 2022, ARA&A, 60, 455
  • Ellis et al. (2013) Ellis, R. S., McLure, R. J., Dunlop, J. S., et al. 2013, ApJ, 763, L7
  • Franco et al. (2020) Franco, M., Elbaz, D., Zhou, L., et al. 2020, A&A, 643, A30
  • Galametz et al. (2013) Galametz, A., Grazian, A., Fontana, A., et al. 2013, ApJS, 206, 10
  • Gerola et al. (1980) Gerola, H., Seiden, P. E., & Schulman, L. S. 1980, ApJ, 242, 517
  • Green (2018) Green, G. 2018, The Journal of Open Source Software, 3, 695
  • Guo et al. (2013) Guo, Y., Ferguson, H. C., Giavalisco, M., et al. 2013, ApJS, 207, 24
  • Hainline et al. (2023) Hainline, K., Robertson, B., Tacchella, S., et al. 2023, in American Astronomical Society Meeting Abstracts, Vol. 55, American Astronomical Society Meeting Abstracts, 212.02
  • Harris et al. (2020) Harris, C. R., Millman, K. J., van der Walt, S. J., et al. 2020, Nature, 585, 357. https://doi.org/10.1038/s41586-020-2649-2
  • Heckman (2006) Heckman, T. 2006, in Encyclopedia of Astronomy and Astrophysics (Murdin)
  • Ho et al. (2019) Ho, S. H., Martin, C. L., & Turner, M. L. 2019, ApJ, 875, 54
  • Iani et al. (2024) Iani, E., Caputi, K. I., Rinaldi, P., et al. 2024, ApJ, 963, 97
  • Inoue et al. (2016) Inoue, S., Dekel, A., Mandelker, N., et al. 2016, MNRAS, 456, 2052
  • Iyer et al. (2018) Iyer, K., Gawiser, E., Davé, R., et al. 2018, ApJ, 866, 120
  • Jackson et al. (2020) Jackson, T. M., Pasquali, A., Pacifici, C., et al. 2020, MNRAS, 497, 4262
  • Kennicutt (1998) Kennicutt, Robert C., J. 1998, ARA&A, 36, 189
  • Knapen & James (2009) Knapen, J. H., & James, P. A. 2009, ApJ, 698, 1437
  • Kron (1980) Kron, R. G. 1980, ApJS, 43, 305
  • Lamastra et al. (2013) Lamastra, A., Menci, N., Fiore, F., & Santini, P. 2013, A&A, 552, A44
  • Langeroodi & Hjorth (2024) Langeroodi, D., & Hjorth, J. 2024, arXiv e-prints, arXiv:2404.13045
  • Le Floc’h et al. (2005) Le Floc’h, E., Papovich, C., Dole, H., et al. 2005, ApJ, 632, 169
  • Lee et al. (2017) Lee, N., Sheth, K., Scott, K. S., et al. 2017, MNRAS, 471, 2124
  • Leitherer et al. (2002) Leitherer, C., Li, I. H., Calzetti, D., & Heckman, T. M. 2002, ApJS, 140, 303
  • L’Huillier et al. (2012) L’Huillier, B., Combes, F., & Semelin, B. 2012, A&A, 544, A68
  • Lilly et al. (1996) Lilly, S. J., Le Fevre, O., Hammer, F., & Crampton, D. 1996, ApJ, 460, L1
  • Madau & Dickinson (2014) Madau, P., & Dickinson, M. 2014, ARA&A, 52, 415
  • Madau et al. (1996) Madau, P., Ferguson, H. C., Dickinson, M. E., et al. 1996, MNRAS, 283, 1388
  • Matthee & Schaye (2019) Matthee, J., & Schaye, J. 2019, MNRAS, 484, 915
  • McCracken et al. (2012) McCracken, H. J., Milvang-Jensen, B., Dunlop, J., et al. 2012, A&A, 544, A156
  • Muxlow et al. (2006) Muxlow, T., Beswick, R. J., Richards, A. M. S., & Thrall, H. J. 2006, in Proceedings of the 8th European VLBI Network Symposium, ed. W. Baan, R. Bachiller, R. Booth, P. Charlot, P. Diamond, M. Garrett, X. Hong, J. Jonas, A. Kus, F. Mantovani, A. Marecki, H. Olofsson, W. Schlueter, M. Tornikoski, N. Wang, & A. Zensus, 31
  • Navarro-Carrera et al. (2024) Navarro-Carrera, R., Rinaldi, P., Caputi, K. I., et al. 2024, ApJ, 961, 207
  • Noeske et al. (2007) Noeske, K. G., Weiner, B. J., Faber, S. M., et al. 2007, ApJ, 660, L43
  • Oesch et al. (2014) Oesch, P. A., Bouwens, R. J., Illingworth, G. D., et al. 2014, ApJ, 786, 108
  • Oesch et al. (2023) Oesch, P. A., Brammer, G., Naidu, R. P., et al. 2023, MNRAS, 525, 2864
  • Oke & Gunn (1983) Oke, J. B., & Gunn, J. E. 1983, ApJ, 266, 713
  • Orlitova (2020) Orlitova, I. 2020, arXiv e-prints, arXiv:2012.12378
  • Östlin et al. (2001) Östlin, G., Amram, P., Bergvall, N., et al. 2001, A&A, 374, 800
  • Pan et al. (2019) Pan, H.-A., Lin, L., Hsieh, B.-C., et al. 2019, ApJ, 881, 119
  • pandas development team (2020) pandas development team, T. 2020, pandas-dev/pandas: Pandas, vlatest, Zenodo, doi:10.5281/zenodo.3509134. https://doi.org/10.5281/zenodo.3509134
  • Pearson et al. (2019) Pearson, W. J., Wang, L., Alpaslan, M., et al. 2019, A&A, 631, A51
  • Peng et al. (2010) Peng, Y.-j., Lilly, S. J., Kovač, K., et al. 2010, ApJ, 721, 193
  • Pérez-González et al. (2023) Pérez-González, P. G., Costantin, L., Langeroodi, D., et al. 2023, ApJ, 951, L1
  • Perrin et al. (2014) Perrin, M. D., Sivaramakrishnan, A., Lajoie, C.-P., et al. 2014, in Society of Photo-Optical Instrumentation Engineers (SPIE) Conference Series, Vol. 9143, Space Telescopes and Instrumentation 2014: Optical, Infrared, and Millimeter Wave, ed. J. Oschmann, Jacobus M., M. Clampin, G. G. Fazio, & H. A. MacEwen, 91433X
  • Popesso et al. (2023) Popesso, P., Concas, A., Cresci, G., et al. 2023, MNRAS, 519, 1526
  • Renaud et al. (2019) Renaud, F., Bournaud, F., Agertz, O., et al. 2019, A&A, 625, A65
  • Renaud et al. (2022) Renaud, F., Segovia Otero, Á., & Agertz, O. 2022, MNRAS, 516, 4922
  • Rinaldi et al. (2022) Rinaldi, P., Caputi, K. I., van Mierlo, S. E., et al. 2022, ApJ, 930, 128
  • Rinaldi et al. (2023a) Rinaldi, P., Caputi, K. I., Costantin, L., et al. 2023a, ApJ, 952, 143
  • Rinaldi et al. (2023b) Rinaldi, P., Caputi, K. I., Iani, E., et al. 2023b, arXiv e-prints, arXiv:2309.15671
  • Rodighiero et al. (2011) Rodighiero, G., Daddi, E., Baronchelli, I., et al. 2011, ApJ, 739, L40
  • Rodríguez Montero et al. (2019) Rodríguez Montero, F., Davé, R., Wild, V., Anglés-Alcázar, D., & Narayanan, D. 2019, MNRAS, 490, 2139
  • Romeo & Fathi (2016) Romeo, A. B., & Fathi, K. 2016, MNRAS, 460, 2360
  • Rosani et al. (2020) Rosani, G., Caminha, G. B., Caputi, K. I., & Deshmukh, S. 2020, A&A, 633, A159
  • Salmon et al. (2015) Salmon, B., Papovich, C., Finkelstein, S. L., et al. 2015, ApJ, 799, 183
  • Salpeter (1955) Salpeter, E. E. 1955, ApJ, 121, 161
  • Sánchez Almeida et al. (2014) Sánchez Almeida, J., Elmegreen, B. G., Muñoz-Tuñón, C., & Elmegreen, D. M. 2014, A&A Rev., 22, 71
  • Santini et al. (2017) Santini, P., Fontana, A., Castellano, M., et al. 2017, ApJ, 847, 76
  • Sargent et al. (2012) Sargent, M. T., Béthermin, M., Daddi, E., & Elbaz, D. 2012, ApJ, 747, L31
  • Schreiber et al. (2015) Schreiber, C., Pannella, M., Elbaz, D., et al. 2015, A&A, 575, A74
  • Scoville et al. (2007) Scoville, N., Aussel, H., Brusa, M., et al. 2007, ApJS, 172, 1
  • Smit et al. (2016) Smit, R., Bouwens, R. J., Labbé, I., et al. 2016, ApJ, 833, 254
  • Sonnett et al. (2013) Sonnett, S., Meech, K., Jedicke, R., et al. 2013, PASP, 125, 456
  • Sparre et al. (2015) Sparre, M., Hayward, C. C., Springel, V., et al. 2015, MNRAS, 447, 3548
  • Speagle et al. (2014) Speagle, J. S., Steinhardt, C. L., Capak, P. L., & Silverman, J. D. 2014, ApJS, 214, 15
  • Stefanon et al. (2019) Stefanon, M., Labbé, I., Bouwens, R. J., et al. 2019, ApJ, 883, 99
  • Tacchella et al. (2016) Tacchella, S., Dekel, A., Carollo, C. M., et al. 2016, MNRAS, 457, 2790
  • Tadaki et al. (2018) Tadaki, K., Iono, D., Yun, M. S., et al. 2018, Nature, 560, 613
  • Taniguchi et al. (2007) Taniguchi, Y., Scoville, N., Murayama, T., et al. 2007, ApJS, 172, 9
  • Taylor (2005) Taylor, M. B. 2005, in Astronomical Society of the Pacific Conference Series, Vol. 347, Astronomical Data Analysis Software and Systems XIV, ed. P. Shopbell, M. Britton, & R. Ebert, 29
  • van Mierlo et al. (2022) van Mierlo, S. E., Caputi, K. I., Ashby, M., et al. 2022, A&A, 666, A200
  • Ventou et al. (2019) Ventou, E., Contini, T., Bouché, N., et al. 2019, A&A, 631, A87
  • Virtanen et al. (2020) Virtanen, P., Gommers, R., Oliphant, T. E., et al. 2020, Nature Methods, 17, 261
  • Wang et al. (2019) Wang, E., Lilly, S. J., Pezzulli, G., & Matthee, J. 2019, ApJ, 877, 132
  • Whitaker et al. (2012) Whitaker, K. E., van Dokkum, P. G., Brammer, G., & Franx, M. 2012, ApJ, 754, L29
  • Whitaker et al. (2014) Whitaker, K. E., Franx, M., Leja, J., et al. 2014, ApJ, 795, 104
  • Whitaker et al. (2019) Whitaker, K. E., Ashas, M., Illingworth, G., et al. 2019, ApJS, 244, 16
  • Williams et al. (2023) Williams, C. C., Tacchella, S., Maseda, M. V., et al. 2023, ApJS, 268, 64
  • Zhang et al. (2023) Zhang, J., Li, Y., Leja, J., et al. 2023, ApJ, 952, 6