Cosmic Evolution of Gas and Star Formation^*

Galaxy evolution is driven by three processes: the conversion of interstellar gas into stars, the accretion of intergalactic gas to replenish the galactic interstellar gas reservoirs, and the interactions and merging of galaxies. The latter redistributes the accreted gas, promotes starburst (SB) activity, fuels active galactic nuclei (AGNs), and can transform the morphology from disks (rotation dominated) to ellipsoidal systems. The interstellar gas plays a defining role in these processes. The gas, being dissipative in its dynamics, will concentrate in the nucleus where it may fuel a central SB or AGN.

The gas supply and its evolution at high redshifts is only loosely constrained (see the reviews by Solomon & Vanden Bout 2005; Carilli & Walter 2013; Tacconi et al. 2020). At present, only ∼200 galaxies at z >1 have been observed in the CO lines, and most are not in the CO (1−0) line, which is a well-calibrated mass tracer of cold molecular gas. To properly understand the early evolution requires significantly larger samples of galaxies. In particular, one must probe the dependence of star formation rates (SFRs) and star formation efficiency (SFE) on the multiple independent variables:

• 1.  
  the variation with redshift or cosmic time;
• 2.  
  the dependence on galaxy stellar mass (M_*); and
• 3.  
  the differences between the galaxies with "normal" star formation (SF) activity on the main sequence (MS) and SB galaxies.

At each redshift, the majority of star-forming galaxies populate a relatively narrow locus in SFR versus M_*, but this so-called MS of star-forming galaxies migrates dramatically to higher SFRs at higher redshift (z). At z ≃ 2, the SB galaxies, with SFRs elevated to 2–100 times higher levels than the MS galaxies, constitute only 5% of the total population of star-forming galaxies (Rodighiero et al. 2011). However, their contribution to the total SF is 8%–14% (Sargent et al. 2012) at z = 2 and likely larger at higher z. Constraining these evolutionary dependences requires both large samples and consistently accurate estimates of the redshifts, gas contents, stellar masses, and SFRs, including both the SF sampled in the optical/UV and the dust-embedded SF probed by far-infrared observations.

The galaxies on the MS at high z would also be classified as SBs if they were at low redshift since they have 3–20 times shorter gas depletion times (M_gas/SFR) than typical spiral galaxies at low z. One might ask whether their shorter gas depletion times are due to nonlinear processes within more gas-rich galaxies (e.g., due to cloud–cloud collisions) or external mechanisms such as galactic merging.

In the work presented here, we develop and analyze a sample of 708 galaxies in the COSMOS survey field (Scoville et al. 2007) within all the Atacama Large Millimeter/submillimeter Array (ALMA) archive pointings (as of 2021 June) in Bands 6 and 7 (1.3 mm and 850 μm). The ALMA data are used to estimate the gas masses of each galaxy from the observed Rayleigh–Jeans (RJ) dust continuum fluxes (Scoville et al. 2016)—a technique calibrated here (see Figure 18). This technique is likely more reliable than the excited-state CO emission lines that are commonly used.

The major questions we address are as follows:

• 1.  
  How do the gas contents depend on the stellar mass of the galaxies?
• 2.  
  How do these gas contents evolve with cosmic time, down to the present, where they constitute typically less than 5%–10% of the galactic mass?
• 3.  
  In the SB phase, is the prodigious SF activity driven by increased gas supply or increased efficiency for converting the existing gas into stars?
• 4.  
  How does the efficiency of star formation vary with cosmic epoch, stellar mass of the galaxy, and whether the galaxy is on the MS or undergoing an SB?
• 5.  
  How rapidly is the gas being depleted? The depletion timescale is characterized by the ratio M_gas/SFR, but to date this has not been measured in broad samples of galaxies owing to the paucity of reliable interstellar gas measurements at high redshift.
• 6.  
  At present there are virtually no observational constraints on the accretion rates needed to maintain the SF at high redshifts, only theoretical predictions. Here we provide estimates on the accretion rates needed to maintain the observed gas contents.

1.1. Modifications from Our Analysis in Earlier Work

The galaxy sample used here has ALMA continuum observations in Band 6 (240 GHz) and/or Band 7 (345 GHz); they are all within the COSMOS survey field and thus have uniform quality and deep ancillary data (Weaver et al. 2022). The ALMA pointings are noncontiguous, but their fields of view (FOVs), totaling 102.9 arcmin², include 708 galaxies with measured far-infrared fluxes from the Spitzer and Herschel Observatories. This sample also has calibrations that are uniform across the full sample without the need for zero-point corrections. Their redshift and stellar mass distributions are shown in Figure 1. All of the Herschel sources within the ALMA pointings are detected by ALMA. There are, of course, sampling biases in the different areas of the parameter space, but the multivariable fitting described below should continuously join those areas.

Zoom In Zoom Out Reset image size

The most uncertain parameter analyzed here is the SFR of the galaxies. The photometric redshift fitting yields estimates for the SFRs derived from fitting the observed optical/UV spectra with SB models plus an overlying extinction. In the case of the gas-rich galaxies in our sample, such estimates typically involve large corrections, by factors of 5–100, for the dust extinction correction. We will refer to this as SFR_opt/UV. An alternative is to remove the extinction correction from SFR_opt/UV, multiplying by the factor 100^−A_nuv/5, to yield an observed SFR_obs.opt/UV. In either case, one adds this to the SFR_IR estimated from the luminosity absorbed by dust which is reradiated in the far-infrared. We refer to these alternatives as SFR(ext cor opt/uv+IR) and SFR (obs. opt/uv+IR).

In Figure 2 these alternatives for the total SFR (observed + obscured) are plotted as a function of the extinction-corrected SFR from the photometric redshift fitting. The left panel clearly shows lower dispersion in the vertical direction; on the other hand, it is implicitly double counting the obscured SFR. In the following work, we have consistently used the left panel SFR (SFR(ext cor opt/uv+IR)) for the SFRs. However, in Table 1, we provide the fitting results for both formulations of the SFRs.

The dusty SF activity is, in virtually all cases, dominant (5–10 times) over the unobscured SF probed in the optical/UV. This is clearly demonstrated in Figure 2. The use of SFRs based on optical/UV data alone, with extinction corrections derived from the optical/UV data, can result in an order-of-magnitude underestimation of the total SFR, even for galaxies close to the MS. ²⁴

Zoom In Zoom Out Reset image size

2.1. Distributions of M_gas and SFE

The ALMA fluxes were translated into gas masses using Equation (A3). These gas masses are shown in Figure 3 for the entire sample of 708 detected sources. As clearly seen in the redshift color-coding of the points, there must be strong dependence of the gas masses on z. Yet this cannot be the entire explanation of the scatter, since the high-redshift galaxies with a large range of gas masses exhibit similarly high SFRs, implying that the rate of SF per unit gas mass (the SFE) is also varying (see Figure 4, right panel).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Figures 4–6 show the trends in the gas contents and SFE with cosmic time, stellar mass, and whether the galaxies are on the MS or have elevated SFRs and are in the SB region. However, these plots do not adequately take account of the joint, simultaneous dependences on the independent parameters (z, M_*, and MS vs. SB galaxies). (The latter is parameterized by SB below.) We fit for all of these dependencies simultaneously in Section 2.2 below.

Zoom In Zoom Out Reset image size

2.2. Simultaneous Fitting for the Joint Dependences of M_gas, SFR, and SFE

Simultaneous fits were done using two techniques: Markov Chain Monte Carlo (MCMC, MLINMIX_ERR; Kelly 2007), and Levenberg–Marquardt (LM) least-squares fitting (lm_fit.pro in IDL). For each fitting, the terms for which power-law coefficients were obtained included the following:

• 1.  
  SFR_MS(z), representing the temporal evolution of the MS as a function of z;
• 2.  
  SFR_MS(M_*), encapsulating the mass-dependent shape of the MS (shown in Figure 7 by the z = 0 MS curve);
• 3.  
  SB to quantify dependences for SB galaxies above the MS—specifically, SB is equal to the ratio of a galaxy's total SFR (=SFR_opt/UV + SFR_IR) to that of a galaxy on the MS at the same redshift and stellar mass; and
• 4.  
  M_*, to capture any simple dependence on the stellar mass of the galaxies.

Two of these terms represent variations relative to the MS: (1) with respect to cosmic age (SFR_MS(z)), and (2) with respect to SFR as a function of stellar mass along the MS (SFR_MS(M_*)). As noted above, SB measures the total SFR relative to the MS. The use of separate functions for the evolution in time and mass relative to the MS and elevation above the MS allows one to probe the evolutionary dependence and the difference between the SB population and the MS galaxies. The steepness of the power-law coefficients for each term enables judgment of the relative importance to changes in the gas contents and the SFE.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Figure 7 shows the MS loci in the M_* and SFR plane for redshifts z = 0–5. The shape of the MS as a function of mass is taken from the z = 1.2 MS of Lee et al. (2015), and the evolution as a function of cosmic age (i.e., z) is from Speagle et al. (2014) fit # 49, normalized to a fiducial mass of log M_* = 10.5 M_⊙. The same mass dependence of the MS is used for all z.

Both the MCMC and LM solutions gave similar coefficients within the uncertainties (2%–10% of the values of the coefficients). The primary advantage of the MCMC fitting is that it fully probes the parameter space, while the LM fitting takes account of uncertainties in the dependent variables and is much faster. The fitting results are shown in Figure 8 and listed in Table 1; the coefficients listed in Table 1 are from the LM fit. Additional relations for the gas depletion time and gas mass fraction are easily derived from those equations, so they are not listed here. The last relation in the table for the net gas accretion rate is derived in Section 3.

Table 1.  M_gas, SFR, and SFE Fitting

	Equation No.
Using SFR = SFR(ext cor opt/uv+IR)
${M}_{\mathrm{gas}}=2.68\times {10}^{9}\,{M}_{\odot }\times {\mathrm{SFR}}_{\mathrm{MS}}{\left(z\right)}^{0.70}\times {\mathrm{SFR}}_{\mathrm{MS}}{\left({M}_{* }\right)}^{-0.56}\times {\left({M}_{* \ 10}\right)}^{0.45}\times {\left(\mathrm{SB}\right)}^{0.29}$	1
$\mathrm{SFE}=0.25\times {\mathrm{SFR}}_{\mathrm{MS}}{\left(z\right)}^{0.33}\times {\mathrm{SFR}}_{\mathrm{MS}}{({M}_{* })}^{1.81}\times {\left({M}_{* \ 10}\right)}^{-0.60}\times {\left(\mathrm{SB}\right)}^{0.71}$	2
(M⊙ yr−1 per 109 M⊙ of gas)
$\mathrm{SFR}=0.66\,{M}_{\odot }\,{\mathrm{yr}}^{-1}\times {\mathrm{SFR}}_{\mathrm{MS}}{(z)}^{1.03}\times {\mathrm{SFR}}_{\mathrm{MS}}{({M}_{* })}^{1.25}\times {\left({M}_{* \ 10}\right)}^{-0.15}\times {\left(\mathrm{SB}\right)}^{1.0}$	3
$\mathrm{ACC}=0.22\,{M}_{\odot }\,{\mathrm{yr}}^{-1}\times {\mathrm{SFR}}_{\mathrm{MS}}{(z)}^{1.70}\times \left[{{M}_{* \ 10}}^{0.62}-0.31{\left({M}_{* \ 10}\right)}^{1.13}\right]$	4
Using SFR = SFR(obs. opt/uv+IR)
${M}_{{\rm{g}}{\rm{a}}{\rm{s}}}=4.34\times {10}^{9}\,{M}_{\odot }\times {{\rm{S}}{\rm{F}}{\rm{R}}}_{{\rm{M}}{\rm{S}}}{\left(z\right)}^{0.59}\times {{\rm{S}}{\rm{F}}{\rm{R}}}_{{\rm{M}}{\rm{S}}}{\left({M}_{\ast }\right)}^{-0.54}\times {\left({M}_{\ast 10}\right)}^{0.42}\times {\left({\rm{S}}{\rm{B}}\right)}^{0.53}$	1-alt
${\rm{S}}{\rm{F}}{\rm{E}}=0.20\times {\rm{S}}{\rm{F}}{{\rm{R}}}_{{\rm{M}}{\rm{S}}}{(z)}^{0.52}\times {\rm{S}}{\rm{F}}{{\rm{R}}}_{{\rm{M}}{\rm{S}}}{({M}_{\ast })}^{1.83}\times {({M}_{\ast 10})}^{-0.56}\times {({\rm{S}}{\rm{B}})}^{0.49}({M}_{\odot }\,{\rm{y}}{{\rm{r}}}^{-1}\,\text{per\unicode{x000A0}}{10}^{9}\,{M}_{\odot }\,\text{of gas})$	2-alt

Note. M_{* 10} is M_stellar in units of 10¹⁰ M_⊙. The fitting was done using both MCMC and LM techniques as detailed in the text, yielding similar exponents for the fit terms (within ∼5%). Equations (1)–(A2) are from the LM fitting. Equation (A2) is obtained from the product of the first two equations. The bottom two equations were derived using the alternative formulation for the SFRs. SFR​MS​(z) is the MS evolution as a function of z from Speagle et al. (2014); SFR​MS​(M_*) is the dependence of the MS as a function of M_* from Lee et al. (2015); and SB is the ratio of the galaxy's SFR to that if it were on the MS.

Download table as:  ASCIITypeset image

Equations (1) and (2) in Table 1 quantify major results regarding the gas contents and SFE in high-redshift galaxies relative to those in low-redshift galaxies and SBs versus MS galaxies:

• 1.  
  The gas masses clearly increase at higher z, varying as SFR_MS(z)^0.70 (i.e., 70% of overall temporal dependence of the MS SFR). The remaining 30% of the increase in MS SFR is due to increased SFE at higher redshifts (based on the exponents of the SFR_MS(z) terms in Equations (1) and (2)).
• 2.  
  The reverse is true for the SB galaxies above the MS: the dominant driver for the activity is an increased SFE—the gas contents increase only as the 0.29 power of SB, while the SFE increases as the 0.71 power of SB. Thus, the SB activity is largely driven by increased efficiency in forming stars per unit mass of gas and only ∼30% due to increased gas masses.
• 3.  
  From the terms involving M_*, it is also apparent that the higher stellar mass galaxies have higher gas contents, but not in proportion to M_* (i.e., the gas mass fraction decreases in the highest-mass galaxies).

The first conclusion clearly implies that the SFE must increase at high redshift (as discussed below). The second conclusion indicates that the galaxies above the MS have higher gas contents, but not in proportion to their elevated SFRs. This confirms that these galaxies are undergoing bursts of activity, rather than long-term elevated SFRs. The galaxies with higher stellar mass likely use up their fuel at earlier epochs and have lower specific accretion rates (see Section 3) than the low-mass galaxies. This is a new aspect of "downsizing" in the cosmic evolution of galaxies. Perhaps many of the high-mass galaxies undergo a last fatal SB that rapidly exhausts their gas supplies.

2.3. Increased SFR versus Higher SFE

In fitting for the SFR dependencies, we intend to clearly distinguish the obvious intuition that when there is more gas there will be both more SF and a higher efficiency for converting the gas to stars. Thus, in solving for the SFE we impose a fixed, linear dependence of the SFR on M_gas. We are then effectively fitting for the SFEs (SFR/M_gas) for star formation per unit gas mass as a function of z, SB, and M_*. This isolates the SFE variation with redshift, SB, and M_* from the variation in M_gas with the same three parameters.

The lack of strong dependence of the SFE on galaxy mass (${M}_{* }^{0.6}$ in Equation (2)) is reasonable. If star-forming gas at high redshift is in self-gravitating GMCs as at low z, the very local gravity (environment) near the GMCs may influence the internal SF, but the GMC will not "care" that it is in the distant potential well of a more or less massive galaxy.

Figure 9 shows the relative evolutionary dependencies of the SFRs, the gas masses, and the SFE per unit mass of gas normalized to unity at z = 0. The fundamental conclusion is that the elevated rates of SF activity both at high redshift and above the MS are due to both increased gas contents and increased efficiencies for converting the gas to stars.

Zoom In Zoom Out Reset image size

Using the MS continuity principle (Appendix B) and the gas contents obtained from Equation (1) for the MS (i.e., SB ≡ SFR/SFR_MS = 1), we can now derive the net accretion rates of MS galaxies required to maintain the MS evolutionary tracks (see Figure 10). This accretion may be purely gaseous or via minor mergers. This analysis is based on the simple logic that we now have estimates of the gas contents that the MS galaxies should have at each redshift as they evolve to lower redshift. And if through their SF activity the galaxies use up their gas at too great a rate and arrive at a later MS curve with a lower M_gas than Equation (1), then the difference must, on average, be made up by accretion. Along each evolutionary track (dashed lines in Figure 10), the rate balance must be given by

$$\begin{eqnarray}&&\begin{array}{l}\displaystyle \frac{{{dM}}_{\mathrm{gas}}}{{dt}}={\dot{M}}_{\mathrm{acc}}-(1-{f}_{\mathrm{mass}\,\mathrm{return}})\times \mathrm{SFR},\end{array}\end{eqnarray} \tag{ 1 }$$

assuming that major merging events are rare. Parameter f_{mass return} is the fraction of stellar mass returned to the gas through stellar mass loss, taken to be 0.3 here (Leitner & Kravtsov 2011). Since these paths are following the galaxies in a Lagrangian fashion, the time derivatives of a mass component M must be taken along the evolutionary track and

$$\begin{eqnarray}\begin{array}{rcl}\displaystyle \frac{{dM}}{{dt}} & = & \displaystyle \frac{{dM}}{{dz}}\displaystyle \frac{{dz}}{{dt}}+\displaystyle \frac{{dM}}{{{dM}}_{* }}\displaystyle \frac{{{dM}}_{* }}{{dt}}\\ & = & \displaystyle \frac{{dM}}{{dz}}\displaystyle \frac{{dz}}{{dt}}+\mathrm{SFR}\displaystyle \frac{{dM}}{{{dM}}_{* }}.\end{array}\end{eqnarray} \tag{ 2 }$$

Figure 11 shows the required gas accretion rates (contours) from Equation (4) in Table 1. The background colors indicate the SFRs (left panel) and the gas contents (right panel). The required rates are ∼100 M_⊙ yr⁻¹ at z ∼2. The combination of two separate stellar mass terms in Equation (4) is required to match the curvatures shown in Figure 11; they do not have an obvious physical justification. The first term dominates at low mass and the second at higher masses.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Two important points to emphasize are as follows: (1) these accretion rates should be viewed as net rates (that is, the accretion from the halos and satellite galaxies minus any outflow rate from SF or AGN feedback), and (2) these rates refer only to the MS galaxies where evolutionary continuity is a valid assumption.

The derived accretion rates are required in order to maintain the SF in the early universe galaxies. Even though the existing gas contents are enormous compared to present-day galaxies, the observed SFRs will deplete this gas within ∼5 ×10⁸ yr; this is short compared to MS evolutionary timescales. The large accretion rates are comparable to the SFRs. The higher SFEs deduced for the high-z galaxies and for galaxies above the MS may be dynamically driven by the infalling gas and galactic merging. These processes will shock-compress the galaxy disk gas, since the induced velocities are likely larger than the internal supersonic turbulence within the clouds.

It is worth noting that although one might think that the accretion rates could have been readily obtained simply from the evolution of the MS SFRs, this is not the case. One needs the mapping of M_gas and its change with time in order to estimate the first term on the right-hand side of Equation (1). Figure 11 shows the relative evolution of each of the major rate functions over cosmic time and stellar mass. Comparing the proximity of the curves as a function of redshift, one can see modest differential change in the accretion rates and SFRs as a function of redshift and stellar mass.

One might ask, what are reasonable accretion rates to adopt for the SB galaxies above the MS, since they are not necessarily obeying the continuity assumption? Here there appear to be two possibilities: either the accretion rate is similar to the MS galaxy at the same stellar mass, or if the elevation above the MS was a consequence of galactic merging, one might assume a rate equal to twice that of an individual galaxy with half the stellar mass. The latter assumption would imply $\sqrt{2}$ higher accretion rate, thus being consistent with the higher gas masses of the galaxies above the MS. In this case, the higher SFRs will be maintained longer than the simple depletion time it takes to reduce the preexisting gas mass back down to that of an MS galaxy. The same $\sqrt{2}$ factor of increase in the gas mass will arise from the merging of the preexisting gas masses of two galaxies of half the observed mass. This follows from the dependence of M_gas on stellar mass, varying only as ${M}_{* }^{0.30}$, rather than linearly. Thus, the notion that the SB galaxies are the result of galaxies merging is favored.

Using the mass functions of star-forming and passive galaxies (Ilbert et al. 2013), we estimate the total cosmic mass density of gas as a function of redshift using Equation (1). (This is the equivalent of the Lilly–Madau plot for the SFR density as a function of redshift.) We do this for the redshift range z = 0–4 and M_* = 10¹⁰–10¹² M_⊙, a modest extrapolation of the ranges covered in the data presented here. Figure 12 shows the derived cosmic mass densities of stars (star-forming and passive galaxies) and gas as a function of redshift. We applied Equation (1) only to the star-forming galaxies and did not include any contribution from the passive galaxy population; to include the SB population, we multiplied the gas mass of the normal star-forming population by a factor of 1.1. If the galaxy distribution is integrated down to stellar mass equal to 10⁹ M_⊙, the stellar mass (and presumably the gas masses) is increased by 10%–20% (Ilbert et al. 2013).

Zoom In Zoom Out Reset image size

The evolution of the gas mass density shown in Figure 12 is similar in magnitude to the theoretical predictions based on semianalytic models by Obreschkow et al. (2009), Lagos et al. (2011), and Sargent et al. (2013) (see Figure 12 in Carilli & Walter 2013). However, all of their estimations exhibit a more constant density at z > 1. The empirically based, prescriptive predictions of Popping et al. (2015) exhibit closer agreement with the evolution found by us; they predict a peak in the gas at z ≃ 1.8 and a falloff at higher and lower redshift. (All of those previous estimations have much larger uncertainties.)

The gas mass fractions computed for galaxies with M_* = 10¹⁰–10¹² M_⊙ are shown in Figure 13. The gas is dominant over the stellar mass down to z ≃ 1.5. At z = 3–4 the gas mass fractions get up to ∼80% when averaged over the galaxy population. Thus, gas contents that peak at z ≃ 2 are likely responsible for the peak in SF and AGNs at that epoch, since the latter are dependent on the gas for fueling. At z = 4 down to 2, the buildup in the gas density is almost identical to that of the cosmic SFR density. (Note that the gas density point at z ∼0.3 is uncertain since it relies on extrapolation of Equation (1) to low M_* and low z, where there exist relatively few galaxies in our sample.)

Zoom In Zoom Out Reset image size

All of the foregoing has focused on the global properties of galaxies at high redshift. The internal distributions of gas, star formation, and stars are critical to developing a physical understanding of the nature of the gas clouds and the star formation within them. To inspire this discussion, a visible-wavelength, multiband image of M74 is provided in Figure 14. It shows the dramatic organization of the star cluster formation along the spiral arms. The thinness of such galactic disks is due to the dissipative nature of the star-formation GMCs, which damp the gas motions perpendicular to the disk. In high-z galaxies, which are accreting gas at high rates and have large energy and momentum input from both the accretion and the active star formation, the disk is also likely to be more irregular in structure (e.g., Forster Schreiber et al. 2011).

Zoom In Zoom Out Reset image size

In nearby galaxies, essentially all of the SF takes place within the GMCs, not in the H i. In the Milky Way disk the GMCs have extraordinary properties—some of which are yet to be understood. The half-mass point for the GMC population is at 2 × 10⁵ M_⊙, diameter 40 pc, and mean density n_H2 = 300 cm⁻³. Although the gas thermal temperature is ∼10–20 K (sound speed ∼0.1–0.2 km s⁻¹), the molecular emission-line widths are typically 5 km s⁻¹, indicating Mach ≳ 10 supersonic turbulence. The magnitude of the turbulence is marginally contained by the self-gravity of the molecular gas given the sizes, masses, and velocity dispersions. Hence, they are self-gravitating and do not fly apart within the cloud-crossing timescale of a few megayears (see Scoville & Sanders 1987 for a summary of GMC properties). The source of the turbulent pressure support is not well constrained. Although we routinely visualize the clouds as spherical for ease of analysis, that is certainly not the case—they are often filamentary, indicating that they are only marginally gravitationally bound. Both the filamentary structure and the supersonic turbulent support may reflect the importance of internal magnetic fields.

In the central 300 pc of the Galaxy, there are a number of molecular clouds with more extreme properties—10⁶ M_⊙ and velocity dispersions 10–20 km s⁻¹. We have found in the high-z galaxies total molecular gas masses 10–100 times that of the Galaxy. We expect that the gas clouds there will be much more massive than those in the local galaxies—this intuition guided our modeling of the IR radiative transfer in clouds of 10⁶–10⁸ M_⊙ (see Appendix A).

For how long do the GMCs and their H₂ survive? One might think that, given the marginal self-gravity of the GMCs and the energy releases from SFR, the GMCs last no more than a few cloud-crossing times (i.e., 10⁶−10⁷ yr). However, there is a very simple argument suggesting that they last more than 10⁸ yr. In the inner disks of local star-forming galaxies, it is always the case that the molecular gas dominates the atomic gas in mass. Within each annulus around the center of the galaxy there must be conservation of the mass flux from H₂ to H i+H ii and vice versa (since the SF consumption is relatively small for each orbit):

$$\begin{eqnarray}&&\begin{array}{l}\langle {\tau }_{{\rm{H}}2}\rangle =\langle {\tau }_{{\rm{H}}\,{\rm\small{I}}}\rangle \times \displaystyle \frac{{M}_{{\rm{H}}2}}{{M}_{{\rm{H}}\,{\rm\small{I}}}},\end{array}\end{eqnarray} \tag{ 3 }$$

where 〈τ〉 is the mass-weighted average timescale of each phase and H ii is neglected since it is generally a minor mass component. Since the typical timescale for the H i is at least ∼10⁸ yr to get to the next spiral arm and the H₂/H i mass ratio is typically greater than 10:1 in the interior regions of the galactic disks, the lifetime of H₂ should be 10⁸–10⁹ yr. An alternative way to phrase this point is that one cannot have the dominant mass component confined to a narrow range of azimuth in the spiral arms, since the azimuthal velocities are different by less than a factor 2 between arm and interarm regions. The lifetime estimate above is of course the mean molecular lifetime—it is not the lifetime of the cloud structures, which may well break up into smaller molecular clouds upon leaving the spiral arms. (This issue is discussed in more quantitative detail in Koda et al. 2016.)

Why is SF concentrated in the spiral arms? If the H₂ clouds have a long life span and persist into the interarm regions, why is the visible SF so beautifully concentrated in the spiral arms of M74 and other local galaxies (see Figure 14)? Despite the fact that the GMCs are self-gravitating, they are not individually on the verge of collapse to form star clusters owing to the supersonic turbulent pressure support, corresponding to a few kilometers per second. However, within the spiral arms, the galactic orbits converge and the cloud–cloud velocity dispersions are increased; cloud–cloud collisions can then compress the internal motions. The molecular gas is very dissipative of the shock energy, leading to collapse of large masses, and in some cases precipitating massive cluster formation. In high-redshift gas-rich galaxies there may not be such well-organized spiral structure, but there will be a much higher rate of cloud collisions, simply due to the larger number density of cloud structures and the large dispersive motions. A similar scenario probably occurs in the nuclear regions of the local ultraluminous IR galaxies, which have much higher SFEs (the role of GMC collisions in triggering the formation of star clusters is discussed in Scoville & Hersh 1979; Scoville & Sanders 1987; Tan 2000; Fukui et al. 2013, 2015; Wu et al.2017).

Why is the overall cycling time so long for galaxies like the Milky Way? Taking the molecular gas content of the Milky Way and dividing by the overall SFR of a few M_⊙ yr^–1 yields a very long cycling time of ∼2 Gyr. Why is this mean SFE so low? Here one should recall the internal structure of the GMCs. As mentioned above, the gas is in filamentary structures, and when two such clouds collide, the filaments in each are unlikely to be aligned—thus, only a small fraction of the gas content is likely to be compressed into dissipative collapse. A good example of this is provided by the two Orion GMCs associated with Orion A and B (see Lombardi et al. 2014). These elongated clouds have their massive star clusters (M42 and NGC 2024) in the nearest regions of their respective clouds, where they might have collided ∼10 Myr in the past. Very possibly, the low Galactic efficiency for forming stars is due to the internal filamentary structure of the molecular clouds. The filamentary structure is also likely to reduce the effectiveness of energetic SF activity in disrupting the clouds—thus extending their lifetimes.

In summary of the discussion here, we have pointed out that in galaxies where the star-forming molecular gas is abundant, the molecular gas and perhaps the clouds are likely to be long-lived, and a significant fraction of the SF is likely to be associated with and triggered by compressive collisions of the clouds. Both of these phenomena are likely even more true in the early universe, where the gas densities are much higher and the cloud motions are more disordered.

The variations of gas masses, accretion, and their relation to star formation have been explored with the most extensive sample yet of high-redshift galaxies. The deduced estimates are "consistent" with most existing studies using the CO lines (Daddi et al. 2010; Genzel et al. 2010, 2015; Tacconi et al. 2010, 2013; Ivison et al. 2011; Riechers et al. 2011; Magdis et al. 2012; Carilli & Walter 2013; Saintonge et al. 2013; Bolatto et al. 2015). However, the sample of galaxies used here is vastly more extensive and has the virtue that it maps the parameter space of z, M_*, and specific SFR (sSFR) out to z = 5 using high-quality and uniform ancillary data from the COSMOS survey field. We thus can simultaneously constrain the functional dependencies on redshift, sSFR relative to the MS, and stellar mass. Our technique also does not suffer from the uncertainties introduced by variable excitation in the higher-J CO lines and the dissociation of CO in strong radiation fields. By contrast, the dust is much harder to destroy in strong radiation fields and is a 1% mass tracer as compared with CO, for which the abundance is ≲10⁻⁴.

A major uncertainty for both the dust and the CO line studies is, of course, their dependence on metallicity (Z). Both are probably more robust than is generally assumed. The CO line is heavily saturated in Galactic GMCs, which have typical τ_CO(1−0) ≳ 10. The ¹³CO emission line is typically ∼1/5 of the CO line flux in Galactic GMCs despite the much lower ¹³C/C abundance (∼1/60–1/90). Thus, the line luminosity must scale as the ∼1/3 power of the CO abundance and hence the metallicity (see Scoville & Solomon 1974, for a discussion of excitation by line photon trapping in optically thick lines).

With respect to the dust emission as a probe of gas, it is reassuring that the dust-to-gas-abundance ratio in low-redshift galaxies is approximately constant at ∼1% by mass from solar down to one-fifth solar metallicity (see Draine & Li 2007; Berta et al. 2016, Figure 16), although why this is the case is not understood.

Our finding that the gas-to-stellar-mass ratio and the accretion rates are both generally higher for lower-mass galaxies has implications for the gas-phase metallicities of galaxies. Assuming that the metallicity of freshly accreting gas is significantly lower than that of the internal gas in the galaxies, one would expect the gas-phase metallicity to increase in higher stellar mass galaxies. This is, of course, known to be true, and it is a major motivation for our focus on galaxies with relatively high M_*. It is also clear that a so-called "closed box" model for the evolution of metal content has little physical justification in light of the extremely large accretion, SF, and feedback rates.

N.S. wishes to thank his mentors (Phil Solomon and Peter Goldreich) and his colleagues who have made a career in astronomy and astrophysics so enjoyable. Each step in a new direction has been stimulating, with new insights and unsolved problems over the broad range of exploration. The support of the public through funding and their interest and appreciation of the expanding understanding are a constant encouragement.

We thank Zara Scoville for proofreading the manuscript and Alvio Renzini for several good suggestions. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada), NSC and ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ. The National Radio Astronomy Observatory is a facility of the National Science Foundation operated under cooperative agreement by Associated Universities, Inc. S.T. acknowledges support from the European Research Council (ERC) in the form of Consolidator grant, 648179, ConTExt. The Cosmic Dawn Center is funded by the Danish National Research Foundation under grant No. 140.

In this appendix we summarize the physical principles behind our use of the long-wavelength dust continuum as a tracer of gas mass. An empirical calibration of this technique is derived from a diverse sample of 128 galaxies both at low redshift and out to z = 4, all having global measurements of the long-wavelength dust continuum and the well-calibrated CO (1−0) line luminosities.

A.1. A Simple Physical Model for the Dust Radiative Equilibrium and Radiative Transfer

To illustrate the essential physics of the IR emission from dense molecular and dust cloud, we have computed the radiative equilibrium for a spherical dust cloud surrounding a central luminosity source, either a young star cluster or an AGN. In the calculation discussed below we specify the dust and gas in 15 shells logarithmically spaced in radius and with density decreasing as r⁻¹ out to an outer radius at 100 pc (as shown schematically in Figure 15, left panel). The central luminosity (presumably visible and UV photons) is taken as a point-source blackbody (T = 3 × 10⁴ K) with L = 10⁷ L_⊙. This luminosity source is surrounded by ionized gas out to ∼1 pc, at which point the optical and near-UV photons will be absorbed in a thin boundary dust shell and reradiated in the MID IR. To illustrate the effect of increasing cloud masses (and increasing dust optical depth), the total cloud masses were taken to be 10⁶, 10⁷, and 10⁸ M_⊙—thus spanning a factor 100 in optical depth. These parameters are chosen to be similar to those of massive clouds, which might exist in galaxies at z = 2–3, when there were generally much greater gas masses, or the massive gas concentrations in SB galaxy nuclei.

Zoom In Zoom Out Reset image size

For the dust opacity we adopt Bruce Draine's Milky Way dust absorption curve taken from his website file "Kext-albedo-WD-MW-3.1-60-DO3-all" (Figure 15, right panel). The scattering and absorption contributions in his table were summed at each wavelength to yield the κ values at 517 wavelengths from λ = 0.1 to 5000 μm. (At λ ≥ 250 μm we use a simple power law with spectral index β = 2, as expected from the Kramers–Kronig relations at long wavelengths.) To translate from dust to gas masses, we adopt a constant gas-to-dust-mass ratio of 105:1 (Draine 2011).

In our numerical integration of the radiative equilibrium and radiation transfer, the central input luminosity is conserved to within ∼20% by the emergent luminosity calculated at r = 100 pc. This is true even in the most optically thick model (e.g., having 10⁸ M_⊙ within r = 100 pc and radial τ_100μm ≃ 1 as shown in Figure 17).

A.2. Empirical Calibration of the RJ Luminosity-to-mass Conversion Factor

At long wavelengths on the RJ tail, the dust emission is almost always optically thin (see Figure 15, right panel), and the emission flux per unit mass of dust is linearly dependent on the dust temperature. Thus, the flux observed on the RJ tail provides a linear estimate (see Figure 16) of the dust mass and hence the gas mass, provided that the dust emissivity per unit mass and the dust-to-gas-mass ratio are constrained. Fortunately, both of these prerequisites are well established from observations of nearby galaxies (e.g., Draine et al. 2007; Galametz et al. 2011). Figure 17 shows the mass-weighted and luminosity-weighted dust temperatures of the emergent radiation for the three mass models and for the sample redshifts.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

On the optically thin RJ tail of the IR emission, the observed flux density is given by

$$\begin{eqnarray*}&&{S}_{\nu }\propto {\kappa }_{D}(\nu ){T}_{D}{\nu }^{2}\frac{{M}_{D}}{{d}_{L}^{2}},\end{eqnarray*}$$

where T _D is the temperature of the emitting dust grains, κ_D (ν) is the dust opacity per unit mass of dust, M _D is the total mass of dust, and d_L is the source luminosity distance. Thus, the mass of dust and gas can be estimated from observed specific luminosity L_ν on the RJ tail:

$$\begin{eqnarray}&&{M}_{\mathrm{dust}}\propto \displaystyle \frac{{L}_{\nu }}{\langle {T}_{D}{\rangle }_{M}{\kappa }_{D}(\nu )}\end{eqnarray} \tag{ A1 }$$

$$\begin{eqnarray}&&{M}_{\mathrm{gas}}\propto \displaystyle \frac{{L}_{\nu }}{{f}_{D}\langle {T}_{D}{\rangle }_{M}{\kappa }_{D}(\nu )}.\end{eqnarray} \tag{ A2 }$$

Here 〈T_D 〉_M is the mean mass-weighted dust temperature and f_D is the dust-to-gas-mass ratio (typically ∼1/100 for solar-metallicity gas).

The empirical calibration of the technique is based on both low- and high-redshift galaxy samples: (1) a sample of 30 local star-forming galaxies; (2) 12 low-z ultraluminous infrared galaxies (ULIRGs); and (3) 30 z ∼ 2 submillimeter galaxies (SMGs). These three samples with 128 galaxies are restricted to only those galaxies having good estimates of the total source-integrated, long-wavelength continuum and CO (1–0). (NB: Once these samples were selected for the calibration, no individual galaxies were selectively dropped owing to departures from the mean.)

We avoid using higher-J CO lines since only the 1–0 transition has been well calibrated using large samples of Galactic GMCs with viral mass estimates. The emissions in the J = 3 − 2 and higher CO lines originate only from high-excitation (both high-density and heated) cloud core regions. These higher CO lines will not be reliable as overall mass tracers of molecular gas contents. Unfortunately, ALMA does not currently have frequency coverage to enable observations of the J = 1 − 0 or 2 − 1 lines at z ≥ 2. There is also no physical justification for observers to claim a fixed ratio of the higher lines to J = 1 − 0 and to thereby infer overall gas contents from the higher-J CO emission. The 3 − 2, 4 − 3, and 5 − 4 lines originate from levels at E_u /k = 33, 55, and 82 K above the ground state, whereas most of the cloud mass is expected to be at ∼25 K. The higher-J lines originate largely from compact core regions, in contrast to the 1 − 0 line, which arises from the full cloud extent and thus samples the overall mass.

In calibrating the CO (1–0) masses, we have adopted ${\alpha }_{\mathrm{CO}(1-0)}=3\,\times \,{10}^{20}\,{\mathrm{cm}}^{-2}\,{\left({\rm{K}}\,\mathrm{km}\,{{\rm{s}}}^{-1}\right)}^{-1}$, which is derived from correlation of the CO line luminosities and virial masses for resolved Galactic GMCs. We believe that this is more correct than the value obtained from Galactic gamma-ray surveys (α = 2 × 10²⁰; see Bolatto et al. 2013), since the latter requires questionable assumptions: (1) that the cosmic rays that produce the ∼2 MeV gamma rays by interaction with the gas fully penetrate the GMCs, and (2) the cosmic-ray density with Galactic radius. (If one adopted the latter value of α_CO, the derived scaling for the dust-based gas masses would be reduced by a factor of 2/3.)

All galaxies in our RJ calibration sample were required to have global measurements of CO (1–0) and RJ dust continuum. The large range in apparent luminosities at 850 μm and in CO is due to the inclusion of high-redshift SMGs, many of which are strongly lensed. The samples are all processed in a common way: (1) all molecular gas masses were derived using the same CO (1−0) conversion factor X_CO = 3 × 10²⁰ N(H₂) cm⁻² (K km s⁻¹)⁻¹, and (2) the molecular gas masses all include a correction for the associated masses of He.

These diverse samples yield remarkably similar values for the dust-to-gas conversion factor α₈₅₀ = 5.6 × 10¹⁹ cgs per M_⊙ (see dashed line fit in the right panel of Figure 18). The Planck value for the Milky Way converted to the same X_CO is 6.2 × 10¹⁹ cgs M_⊙ ^–1 between the dust continuum flux and the molecular masses (including He). We have adopted a long-wavelength spectral index for the dust opacity of β_dust = 2 for all these data sets when converting the RJ flux at one wavelength to the reference λ = 850 μm.

Zoom In Zoom Out Reset image size

Figure 18 shows the ratio of specific luminosity at rest frame λ = 850 μm to that of the CO (1−0) line, and one clearly sees a similar ratio of RJ dust continuum to CO luminosity. Using a standard Galactic CO (1–0) conversion factor, we then obtain the relation by which we convert the RJ dust continuum to gas masses:

$$\begin{eqnarray}\begin{array}{rcl}{M}_{\mathrm{gas}} & = & 1.78\,{S}_{{\nu }_{\mathrm{obs}}}[\mathrm{mJy}]{\left(1+{\rm{z}}\right)}^{-4.8}\\ & & \times {\left(\displaystyle \frac{{\nu }_{850\mu {\rm{m}}}}{{\nu }_{\mathrm{obs}}}\right)}^{4.0}{\left({d}_{L}[\mathrm{Gpc}]\right)}^{2}\\ & & \times \,\left\{\displaystyle \frac{5.6\,\times \,{10}^{19}}{{\alpha }_{850}}\right\}\,\displaystyle \frac{{{\rm{\Gamma }}}_{0}}{{{\rm{\Gamma }}}_{\mathrm{RJ}}}{10}^{10}{M}_{\odot }.\end{array}\end{eqnarray} \tag{ A3 }$$

In Equation (A3), Γ is a correction for departures from strict ν² of the RJ continuum (Scoville et al. 2016), shown in the left panel of Figure 19. Here ${\alpha }_{850}=(5.6\pm 1.7)\,\times {10}^{19}\,\mathrm{erg}\,{{\rm{s}}}^{-1}\,{\mathrm{Hz}}^{-1}\,{M}_{\odot }^{-1}$ is the derived calibration constant between 850 μm luminosity and gas mass. We have adopted a dust opacity spectral index β = 2 (see Section A.1).

Zoom In Zoom Out Reset image size

In the present work, the conversion of the fluxes, measured with ALMA, to gas masses is done with Equation (A3). Using Equation (A3), the predicted fluxes for a fiducial gas mass of 10¹⁰ M_⊙ in the ALMA and Herschel-SPIRE bands are shown in Figure 19 as a function of redshift.

In the current work we restrict the observed galaxies to be relatively massive (M_* > 5 × 10¹⁰ M_⊙), since they should have close to solar metallicity and presumably close to solar dust-to-gas-abundance ratios. We note that for the first factor ∼5 down from solar metallicity the dust-to-gas ratio is constant for those galaxies with global measurements (Draine et al. 2007). Probing lower stellar mass galaxies, which presumably would have significantly subsolar metallicity, will require careful calibration as a function of metallicity or mass. We note that in Draine et al. (2007) there is little evidence of variation in the dust-to-gas-abundance ratio for the first factor of 4–5 down from solar metallicity. However, at lower metallicities the dust-to-gas-abundance ratio does clearly decrease (see Figure 17 in Draine et al. 2007 and Figure 16 in Berta et al. 2016).

A.3. A Caution for Higher Redshifts

In our analysis, we make use of a principle we refer to as the continuity of MS evolution—simply stated, the temporal evolution of the star-forming galaxy population may be followed by Lagrangian integration in the SFR/M_* plane. This follows from the fact that approximately 95% of star-forming galaxies at each epoch lie on the MS, with SFRs dispersed only a factor 2 above or below the MS (Rodighiero et al. 2011). A similar approach has been used by Noeske et al. (2007), Renzini (2009), and Leitner (2012) and references therein.

This continuity assumption ignores the galaxy buildup arising from major mergers of similar-mass galaxies, since they can depopulate the MS population in the mass range of interest. Major mergers may be responsible for some galaxies in the SB population above the MS (see Section 4). On the other hand, minor mergers may be considered simply as one element of the average accretion process considered in Section 3.

We are also neglecting the SF quenching processes in galaxies. This occurs mainly in the highest-mass galaxies (M_* > 2 × 10¹¹ M_⊙ at z > 1) and in dense environments at lower redshift (Peng et al. 2010; Darvish et al. 2016). At z > 1.2 the quenched red galaxies are a minor population (see Figure 12) and the quenching processes are of lesser importance.

The paths of evolution in the SFR versus M_* plane can be easily derived since the MS loci give dM_*/dt = SFR (M_*). One simply follows each galaxy in a Lagrangian fashion as it builds up its mass. In the Lagrangian integration, we move with the galaxies as they trace these paths, and the time derivatives of a mass component M are taken along the evolutionary track. Figure 10 shows these evolutionary tracks for a sample of galaxies. Using this continuity principle to evolve each individual galaxy over time, the evolution for MS galaxies across the SFR–M_* plane is as shown in Figure 10. Here we have assumed that 30% of the SFR is eventually put back into the gas by stellar mass loss. This is appropriate for the mass loss from a stellar population with a Chabrier initial mass function (IMF; see Leitner & Kravtsov 2011). In this figure, the curved horizontal lines are the MS at fiducial epochs or redshifts, while the downward curves are the evolutionary tracks for fiducial M_* from 1 to 10 × 10¹⁰ M_⊙. At higher redshift, the evolution tends toward increasing M_*, whereas at lower redshift the evolution is in both SFR and M_*. In future epochs, the evolution is more vertical as the galaxies exhaust their gas supplies. Thus, there are three phases in the evolution:

• 1.  
  the gas-accretion-dominated and stellar mass buildup phase at z > 2 corresponding to cosmic age less than 3.3 Gyr (see Section 3);
• 2.  
  the transition phase, where gas accretion approximately balances SF consumption and the evolution becomes diagonal; and
• 3.  
  the epoch of gas exhaustion at z ≲ 0.1 (age 12.5 Gyr), where the evolution is vertically downward in the SFR versus M_* plane.

These evolutionary phases are all obvious (and not a new development here). However, in Section 3 we make use of this continuity assumption to derive the accretion rates and hence substantiate the three phases as separated by their accretion rates relative to their SFRs. When these phases begin and end is a function of the galaxy stellar mass—the transitions in the relative accretion rates take place at higher redshift (i.e., earlier cosmic epoch) for the more massive galaxies. Figure 10 shows the evolution for three fiducial stellar masses.

At each epoch there exists a much smaller population (∼5% by number at z = 2) that has SFRs 2–100 times that of the MS at the same stellar mass. Do these SB galaxies quickly exhaust their supply of SF gas, thus evolving rapidly back to the MS? Or are their gas masses systematically larger so that their depletion times differ little from the MS galaxies? These SB galaxies must be either a short-duration, but common, evolutionary phase or of long duration but not undergone by the majority of the galaxy population. Their significance in the overall cosmic evolution of SF is certainly greater than 5%, since they have 2–50 times higher SFRs than the MS galaxies.

The major data sets used for our analysis are as follows:

• 1.  
  the latest COSMOS 2020 photometric redshift catalog (Weaver et al. 2022) for positions, stellar masses, and the unobscured SFRs;
• 2.  
  the far-infrared continuum fluxes from Herschel PACS and SPIRE for the IR-based obscured SFR rates; and
• 3.  
  ALMA Band 6 and 7 continuum imaging for estimating gas masses from the long-wavelength RJ fluxes.

C.1. Redshifts, Stellar Masses, and SFRs: Optical/UV and IR Star Formation Rates

C.2. IR Source Catalogs

C.3. ALMA Band 6 and 7 Continuum Data