Dust, Gas, and Metal Content in Star-forming Galaxies at z ∼ 3.3 Revealed with ALMA and Near-IR Spectroscopy

Tomoko L. Suzuki; Masato Onodera; Tadayuki Kodama; Emanuele Daddi; Masao Hayashi; Yusei Koyama; Rhythm Shimakawa; Ian Smail; David Sobral; Sandro Tacchella; Ichi Tanaka

doi:10.3847/1538-4357/abd4e7

1. Introduction

Molecular gas (H₂) is one of the fundamental physical quantities of galaxies because it is the fuel for star formation. It is well known that the gas surface density is correlated with the star formation rate (SFR) surface density (the Schmidt–Kennicutt relation; Schmidt 1959; Kennicutt 1998). The total gas mass is also connected with the total star-forming activity (e.g., Daddi et al. 2010; Genzel et al. 2010). The typical SFR of star-forming galaxies at a given stellar mass appears to monotonically increase with increasing redshifts (e.g., Whitaker et al. 2012; Sobral et al. 2014; Tomczak et al. 2016). More active star formation in galaxies at higher redshifts is considered to be supported by a larger amount of gas (e.g., Daddi et al. 2010; Genzel et al. 2010; Geach et al. 2011; Bothwell et al. 2013b; Tacconi et al. 2013; Birkin et al. 2021). Investigating the gas contents in galaxies at high redshifts is crucial for understanding the formation and evolution of galaxies in the universe (e.g., Walter et al. 2016; Riechers et al. 2019).

Observational studies in the past decade revealed gas properties not only for dusty starburst galaxies, such as submillimeter bright galaxies (SMGs), but also for ultraviolet (UV)/optical-selected star-forming galaxies on the stellar mass–SFR relation, the so-called "main sequence" of star-forming galaxies, at z ≳ 1 (e.g., Daddi et al. 2010; Genzel et al. 2010; Tacconi et al. 2010, 2013; Magdis et al. 2017). The increasing number of galaxies with individual measurements of the gas mass in a wide redshift range makes it possible to establish the scaling relations for gas mass fraction and gas depletion timescale (= M_gas/SFR) as a function of redshift, stellar mass, and SFR (e.g., Scoville et al. 2017; Tacconi et al. 2018; Freundlich et al. 2019; Liu et al. 2019). At z > 3, however, the number of UV/optical-selected star-forming galaxies with the individual measurements of gas content is still small (with CO emission lines: Magdis et al. 2017; Cassata et al. 2020, and with dust continuum: Schinnerer et al. 2016; Wiklind et al. 2019; Aravena et al. 2020). The evolution of the gas properties in UV/optical-selected galaxies at z > 3 is still unclear.

The atomic and/or molecular hydrogen gas mass is also said to be correlated with the gas-phase metallicity from both observations (e.g., Bothwell et al. 2013a, 2016b; Hunt et al. 2015; Seko et al. 2016a; Brown et al. 2018) and cosmological numerical simulations (e.g., Lagos et al. 2016; Torrey et al. 2019). It has been suggested that the gas mass is more fundamental than the SFR when the scatter of the mass–metallicity relation of star-forming galaxies is to be explained (e.g., Bothwell et al. 2013a; Zahid et al. 2014; Brown et al. 2018). Indeed, more gas-rich star-forming galaxies tend to be more metal poor and form stars more actively. At high redshifts, a direct comparison between gas mass and gas-phase metallicity is limited to a handful of galaxies at z ∼ 1–3 (Saintonge et al. 2013; Seko et al. 2016a; Shapley et al. 2020). Seko et al. (2016a) found a trend that the gas mass fraction decreases with increasing metallicities at a fixed stellar mass for star-forming galaxies at z ∼ 1.4. No such direct comparison between the two quantities is available at z > 3.

Galaxies evolve by interacting with the intergalactic medium (IGM). Gas accretes onto galaxies from the outside, chemical enrichment proceeds as stars form, and gas and metals are ejected from galaxies by outflow (e.g., Bouché et al. 2010; Davé et al. 2011; Lilly et al. 2013; Peng & Maiolino 2014; Tacchella et al. 2020). Gas mass fraction and gas-phase metallicity are often used to investigate the relative contributions between star formation, gas outflow, and inflow (Erb 2008; Cresci et al. 2010; Troncoso et al. 2014; Yabe et al. 2015; Seko et al. 2016b; Sanders et al. 2020; see also Burgarella et al. 2020; Nanni et al. 2020 using dust-to-stellar mass ratio). Most of these studies estimated gas mass fractions by converting the SFR surface density into gas surface density with the Schmidt–Kennicutt relation (Erb 2008; Cresci et al. 2010; Troncoso et al. 2014; Yabe et al. 2015; Sanders et al. 2020). Given that galaxies are more actively forming stars at higher redshifts, they are expected to interact more actively with the surrounding IGM via outflows and inflows (e.g., Yabe et al. 2015). At z > 3, it has been suggested that star-forming galaxies are no longer in equilibrium (e.g., Mannucci et al. 2010), where the gas consumption due to star formation and outflows is balanced with the gas acquisition by inflows (inflow = star formation + outflow), due to the intense gas inflows onto galaxies in the early universe (e.g., Padmanabhan & Loeb 2020). Estimating both the gas mass and gas-phase metallicity for star-forming galaxies at z > 3 allows tests of whether galaxies at z > 3 are out of equilibrium.

Several methods are commonly used to estimate gas masses. The first method is using CO emission line fluxes (e.g., Daddi et al. 2010; Genzel et al. 2010; Tacconi et al. 2010, 2013). This method has uncertainties on the CO-to-H₂ conversion factor, which changes depending on metallicity (Genzel et al. 2012), and on the CO excitation states when higher-J CO lines are used (e.g., Daddi et al. 2015; Riechers et al. 2020). Furthermore, observations of CO lines for galaxies at high redshifts are observationally expensive. The second method is converting a dust mass into a gas mass with an assumed gas-to-dust mass ratio (e.g., Santini et al. 2014; Béthermin et al. 2015). Because the gas-to-dust mass ratio depends on the metallicity (Leroy et al. 2011; Rémy-Ruyer et al. 2014), metallicity measurements are crucial for estimating the gas mass accurately. Gas masses can also be estimated with an empirically calibrated relation between single-band dust continuum flux at the Rayleigh–Jeans (R-J) tail and gas mass (e.g., Scoville et al. 2014, 2016; Groves et al. 2015). These scaling relations are calibrated with local galaxies or with local galaxies and SMGs up to z ∼ 2. In this method, the gas-to-dust mass ratio is included in the conversion factor, and therefore does not need to be considered. It remains unclear whether the empirical calibration methods are applicable to galaxies at z > 3 or how much scatter there is in the relationships. Given that dust continuum observations take much less time compared to the CO observations, using dust continuum as a tracer of gas has the advantage of increasing the number of galaxies at higher redshifts with individual gas estimates, but these will only be reliable when precise metallicities are available as well. Metallicity measurements based on rest-frame optical emission lines for dustier star-forming galaxies are thought to have larger uncertainties due to strong dust obscuration (e.g., Santini et al. 2010). Herrera-Camus et al. (2018) reported a discrepancy between metallicities derived with rest-frame optical emission lines and far-infrared (FIR) fine-structure lines for local (ultra) luminous infrared galaxies ((U)LIRGs).

In this paper, we present the results from submillimeter observations with the Atacama Large Millimeter/submillimeter Array (ALMA) of star-forming galaxies at z = 3–4. High-quality near-infrared (NIR) spectra obtained with Keck/MOSFIRE (McLean et al. 2010, 2012) are available for all of the targets, and their gas-phase metallicities have been measured (Onodera et al. 2016; Suzuki et al. 2017). By observing the dust continuum emission, we can estimate their dust masses and convert them into gas masses using the relation between the metallicity and gas-to-dust mass ratio. We then investigate the gas properties, namely, gas mass fractions and gas depletion timescales, of star-forming galaxies at z = 3–4. Comparing the gas contents with the gas-phase metallicities, we aim to understand how star-forming galaxies at this epoch interact with their surrounding IGM via gas inflows and outflows.

This paper is organized as follows. In Section 2 we introduce our parent sample of star-forming galaxies at z = 3–4 and describe the observations conducted with ALMA. We also describe the reduction and analysis for the obtained data and stacking analysis. In Section 3 we present our estimates of the physical quantities, such as gas-phase metallicity, ionization parameter, and gas mass. In Section 4 we show our results on the dust and gas properties of the star-forming galaxies at z = 3–4 and discuss their metallicity dependences. We also compare our observational results with a gas regulator model to discuss how star-forming galaxies at this epoch interact with their surrounding IGM. We summarize this paper in Section 5.

Throughout this paper, we assume cosmological parameters of Ωおめが_m = 0.3, Ωおめが_Λらむだ = 0.7, and H₀ = 70 km s⁻¹ Mpc⁻¹. We use a Chabrier initial mass function (IMF; Chabrier 2003).

2. Observation and Reduction

2.1. Spectroscopically Confirmed Galaxies at z = 3–4

Our parent sample is constructed from the two different studies of star-forming galaxies at 3 < z < 4 in the COSMOS field using NIR spectroscopy with Keck/MOSFIRE. One study is based on a spectroscopic and photometric redshift selection (Onodera et al. 2016, Section 2.1.1), while the other study is based on narrow-band selection (Suzuki et al. 2017, Section 2.1.2). The parent sample for both studies consists of 53 galaxies with z_spec ∼ 3.0–3.8.

2.1.1. UV-selected Galaxies

Onodera et al. (2016, hereafter O16) selected targets for spectroscopic observation originally from the zCOSMOS-Deep redshift catalog (Lilly et al. 2007) and the 30-band COSMOS photometric redshift catalog (McCracken et al. 2012; Ilbert et al. 2013). O16 conducted H- and K-band spectroscopy and confirmed 43 galaxies at z_spec = 3.0–3.8 based on the rest-frame optical emission lines. The confirmed star-forming galaxies span a stellar mass range of $\mathrm{log}({M}_{* }/{M}_{\odot })\,\sim$ 8.5–11.0 and are distributed around the star-forming main sequence at z ∼ 3.3 (O16).

2.1.2. [O iii] Emission Line Galaxies

Suzuki et al. (2017, hereafter S17) selected targets for spectroscopic observation from a catalog of narrow-band (NB)-selected [O iii]λらむだ5007 emission line galaxies at z = 3.23, obtained by the High-Z Emission Line Survey (HiZELS; Best et al. 2013; Sobral et al. 2013; Khostovan et al. 2015). S17 conducted H- and K-band spectroscopy and confirmed 10 [O iii] emitters at z_spec = 3.23–3.28. The stellar mass range of the confirmed [O iii] emitters is $\mathrm{log}({M}_{* }/{M}_{\odot })\,\sim$ 9.1–10.2. The [O iii] emitters follow the star-forming main sequence at z ∼ 3.2 and the mass–metallicity relation established by O16 (S17).

2.2. ALMA Band-6 Observation

We selected galaxies with $\mathrm{log}({M}_{* }/{M}_{\odot })\geqslant 10$ and ≥ 3σしぐま detection of [O iii]λらむだ5007, Hβべーた, or [O ii]λらむだ3727 emission lines from the parent sample as targets for the ALMA observations. We excluded 2 galaxies classified as active galactic nuclei (AGNs) in O16. One has an X-ray counterpart detected with Chandra. The other shows strong [Ne iii]λらむだ3869 emission and a high [O iii]λらむだ4363/Hγがんま ratio, which is likely to be powered by the AGN (O16). As a result, the sample for the ALMA observations consists of 12 galaxies, 2 of which are [O iii] emitters from S17 (Table 1).

Table 1. Summary of the Targets of this Study and ALMA Band-6 Observations

ID ^a	R.A. ^a	Decl. ^a	z_spec	Exp. Time ^b	Central Freq. ^c	rms Level ^d	S_{1.3 mm} ^e	Reference
	(deg)	(deg)		(minutes)	(GHz)	(mJy beam⁻¹)	(mJy)
208681	149.90551	2.353990	3.267	5	232.1	0.05	1.02 ± 0.05	O16
214339	150.31607	2.372240	3.609	5	233.0	0.05	0.30 ± 0.05	''
406444	150.33032	2.072270	3.304	5	232.1	0.04	0.12 ± 0.04	''
3	149.97513	1.69375	3.230	41	235.6	0.01	0.11 ± 0.01	S17
434585	149.84702	2.373020	3.363	11	230.9	0.03	0.11 ± 0.03	O16
192129^f	150.30078	2.300540	3.495	49	237.3	0.01	0.05 ± 0.01	''
217753	149.89451	2.383700	3.254	5	232.1	0.05	<0.15	''
218783	149.92082	2.387060	3.297	5	232.1	0.04	<0.11	''
212298	150.34268	2.365390	3.108	5	246.2	0.04	<0.11	''
413391	149.78424	2.452890	3.365	11	230.9	0.03	<0.10	''
5	149.95568	1.68044	3.241	53	235.6	0.01	<0.04	S17
434618	149.89213	2.414710	3.285	88	233.5	0.01	<0.03	O16

Notes.

^aObject IDs and coordinates are extracted from the original papers.^bOn-source observing time for Band-6 observations.^cCentral frequency of the four spectral windows.^dMeasured in the tapered maps.^eMeasured in the tapered maps with imfit and corrected for the primary beam. 3σしぐま upper limits are shown for the ALMA nondetected sources.^fX-ray detected source (Section 2.2).

Download table as: ASCII Typeset image

Although the potential AGNs were excluded from our ALMA targets, we found that one of the ALMA targets, 192129, is detected in X-ray with Chandra (Elvis et al. 2009; Civano et al. 2012, 2016) and included in the X-ray selected AGN catalog of Kalfountzou et al. (2014) as a type-2 AGN. The optical–NIR spectral energy distribution (SED) of this source is not peculiar compared to the other galaxies (O16 and as shown in the best-fit SEDs in the Appendix) and its Hβべーた emission line is narrow, which is likely to be consistent with the classification by Kalfountzou et al. (2014). Although we expect that the optical–NIR emission is dominated by emission from the host galaxy, the physical quantities, such as stellar mass, SFR, gas-phase metallicity, and ionization parameter (Sections 2.4 and 3.1), may be affected by the emission from the AGN. On the other hand, the dust continuum observed at ALMA Band-6 (λらむだ_rest ∼ 280 μみゅーm) is expected to be dominated by cold dust emission from star-forming regions (T ∼ 20–40 K). We do not exclude this source in the following analyses but distinguish it from the other sources in each figure.

Our ALMA Cycle 6 program with Band 6 was conducted during 2018 December–2019 March (2018.1.00681.S, PI: T. Suzuki). The frequencies of the four spectral windows are slightly different for the targets depending on their spectroscopic redshifts (between 221.9 GHz and 254.4 GHz). We set the frequencies of the spectral windows so that we can cover the CO(9–8) line (νにゅー_rest = 1036.9 GHz) with one of the four spectral windows. The effective bandwidth of each spectral window is 1.875 GHz. The data were taken with the time division mode. The total on-source time is ∼5–90 minutes depending on the stellar mass, SFR, and gas-phase metallicity of the targets, as summarized in Table 1.

The brightest source at 1.3 mm, 208681, appears to be detected with the CO(9–8) line. Given that quasar host galaxies tend to have more extreme CO excitation states than normal star-forming galaxies (Carilli & Walter 2013), the CO(9–8) line detection may suggest that this source hosts an AGN. We will discuss the CO(9–8) line of this source in a forthcoming paper (T. L. Suzuki et al. 2021, in preparation). We note that the contribution of the CO(9–8) line to the dust continuum flux is negligible.

We used the Common Astronomy Software Application package (casa; McMullin et al. 2007) to calibrate the raw data. We ran the clean algorithm with natural weighting. When sources were detected with ≥5σしぐま level, we ran clean again by masking the sources. Because the synthesized beam sizes are slightly different for the targets, we created the uv-tapered maps for some of the sources to conduct flux measurement under similar beam sizes. The average beam size of the 12 ALMA maps is 1 farcs 52 × 1 farcs 32.

We used imfit to fit a 2D Gaussian to each target. The central position was fixed at the centroid determined in the Ks-band image from UltraVISTA.¹¹ Our detection criterion was that the peak flux obtained by imfit has >3σしぐま significance. As a result, 6 out of 12 galaxies satisfy this criterion, as summarized in Table 1. We ran imfit several times with different parameter settings to check the robustness of the obtained fluxes. We confirmed that the fitting results for the 6 detected sources are not affected by the parameter settings. In the following sections, we use peak fluxes measured with imfit as the total fluxes of the detected sources. The obtained peak fluxes broadly agree with the aperture fluxes (r = 1 farcs 5) measured at the position of the Ks-band centroid within the uncertainties, which suggests that our targets are not spatially resolved in the ALMA maps. As for the nondetected sources, we assigned 3σしぐま upper limit fluxes. The measured fluxes and limits are summarized in Table 1. The listed fluxes were corrected for the primary beam. Because the ALMA targets are located at the center of each ALMA map, the primary beam correction is lower than 1% for all of the targets.

Figure 1 shows the ALMA maps of our target galaxies together with the Ks-band centroids. Some of the ALMA-detected sources, such as 406444 and 434585, show a clear spatial offset (up to ∼0 farcs 5) between the dust continuum emission and the Ks-band centroid. Their relatively large positional offsets are probably due to the lower signal-to-noise ratios of their dust continuum emission. Indeed, we found a trend in the ALMA-detected sources that the positional offset between the dust emission peak and the Ks-band centroid becomes larger with decreasing signal-to-noise ratio of the dust continuum emission.

**Figure 1.** ALMA Band-6 images before tapering of 12 targets (image size: 5'' × 5''). A black circle shows the beam size of each image. Black contours correspond to 4σしぐま, 8σしぐま, 12σしぐま, and 16σしぐま. A plus represents the centroid determined in the Ks-band image. Some of the ALMA-detected sources, such as 406444 and 434585, show a larger spatial offset between the dust continuum emission and the Ks-band centroid than the other detected sources. This is probably due to the lower signal-to-noise ratios of their dust continuum emission. The stacked image of the five individually ALMA nondetected sources with $\mathrm{log}({M}_{* }/{M}_{\odot })\,=$ 10.0–10.4 is also shown. The central position is shown with a plus. The stacked emission is detected at 5σしぐま.
Download figure:
Standard image High-resolution image

2.3. Stacking Analysis

We stacked the Band-6 images of the ALMA nondetected sources to investigate their average dust continuum flux. As a result of SED fitting in Section 2.4, one source, 434618, turns out to be only $\mathrm{log}({M}_{* }/{M}_{\odot })={9.39}_{-0.01}^{+0.12}$ . This is probably due to using a different SED fitting code and/or different photometric catalog with deeper NIR and Spitzer data from the previous estimate (O16). Because the stellar mass of 434618 is ≳0.6 dex lower than those of the other nondetected sources, we excluded this source so that we were able to investigate the average dust and gas properties of star-forming galaxies with similar stellar masses of $\mathrm{log}({M}_{* }/{M}_{\odot })\,=$ 10.0–10.4.

We cut out the tapered 20'' × 20'' ALMA images centered on the Ks-band position. Then, we stacked the cutout images by weighting with the rms levels in the tapered maps (Table 1). The stacked image is shown in Figure 1. We measured the total flux of the stacked image with imfit as done in Section 2.2. The obtained flux is 0.06 ± 0.01 mJy, which satisfies our detection criterion.

2.4. SED Fitting

We conducted SED fitting including the dust continuum flux or limit at 1.3 mm obtained with ALMA. We used the SED fitting code magphys (da Cunha et al. 2008, 2015; Battisti et al. 2019), which can fit the SEDs from the optical to radio wavelengths consistently. magphys combines the emission by stellar populations with the attenuation and emission by dust in galaxies based on the energy balance technique. We used the updated version of magphys for galaxies at high redshifts (da Cunha et al. 2015).

magphys adopts stellar population synthesis models of Bruzual & Charlot (2003) with the Chabrier (2003) IMF. The metallicity range is set to be from 0.2 to 2 times solar. Star formation history is parameterized as a continuous delayed exponential function, i.e., the SFR rises at the earlier epoch and then declines exponentially with a certain timescale between 0.075 and 1.5 Gyr⁻¹. The age is randomly drawn between 0.1 and 10 Gyr. magphys also includes star bursts of random duration and amplitude to account for the stochasticity on star formation history. The current SFR is determined by averaging SFH over the last 100 Myr. As for the dust attenuation, magphys uses the two-component model of Charlot & Fall (2000). A number of tests of the application of magphys to dust-obscured galaxies, including simulated galaxies from EAGLE, are discussed in Dudzevičiūtė et al. (2020).

magphys takes into account four main dust components, namely, a polycyclic aromatic hydrocarbons, hot dust at mid-infrared, warm dust, and cold dust in the thermal equilibrium. The warm and cold dust components in the thermal equilibrium emit as modified blackbodies with an emissivity index βべーた of 1.5 for the warm components and 2 for the cold components (da Cunha et al. 2015). magphys assumes a dust mass absorption coefficient at 850 μみゅーm of κかっぱ_abs = 0.77 g⁻¹ cm².

We combined the flux densities at 1.3 mm from ALMA with the optical–NIR broadband photometries (u, B, V, r, i_p, z_pp, Y, J, H, Ks, 3.6, 4.5, 5.8, and 8.0 μみゅーm) from the COSMOS2015 catalog (Laigle et al. 2016). Because magphys does not include emission lines from the ionized gas, we subtracted the emission line fluxes measured with the NIR spectra from the H-band ([O ii]) and Ks-band ([O iii] doublet, Hβべーた) fluxes. We took into account the upper limits for the optical–NIR photometries by giving 0 to the flux column and a 3σしぐま value to the photometric error column according to da Cunha et al. (2015). As for the 1.3 mm flux of the ALMA nondetected sources, we gave a 1.5σしぐま ± 1σしぐま value as done in Dudzevičiūtė et al. (2020). Using a 1.5σしぐま ± 1σしぐま value provides a better weighting of the submillimeter constraint in the best-fit model returned by magphys than using a 3σしぐま upper limit. The derived physical parameters, such as stellar masses and SFRs, do not significantly change depending on the adopted flux and error values (Dudzevičiūtė et al. 2020).

We also conducted SED fitting for the stacked sample with the obtained 1.3 mm flux in Section 2.3. When taking an average of the optical–NIR photometries, we used the same weights as used in the ALMA image stacking (Section 2.3).

The best-fit SEDs of the individual galaxies and the stacked sample are shown in the Appendix. We use the median values of the probability distribution function (PDF) for stellar mass, SFR, and dust mass in the following analyses. These physical quantities are summarized in Table 2. The uncertainties correspond to the 16–84th percentile values of the PDF. As for the dust masses of the ALMA nondetected sources, we use the 97.5th percentile values of the PDF as the upper limits.

Table 2. Summary of the Physical Quantities of the Star-forming Galaxies at z ∼ 3.3 and the Stacked Sample

ID	$\mathrm{log}({M}_{* }/{M}_{\odot })$ ^a	$\mathrm{log}(\mathrm{SFR})$ ^a	$12+\mathrm{log}({\rm{O}}/{\rm{H}})$ ^b	log(q)	δでるた_GDR	$\mathrm{log}({M}_{\mathrm{dust}}/{M}_{\odot })$ ^a ^, ^c	$\mathrm{log}({M}_{\mathrm{gas}}/{M}_{\odot })$
		(M_⊙ yr⁻¹)		(cm s⁻¹)
208681	10.82_−0.03 ^+0.00	2.10_−0.01 ^+0.08	8.59_−0.05 ^+0.04	7.64 ± 0.03	126₋₁₄ ⁺¹³	8.89_−0.08 ^+0.11	11.29_−0.09 ^+0.11
214339	10.48_−0.05 ^+0.04	1.76 ± 0.13	8.30_−0.12 ^+0.09	7.82_−0.06 ^+0.05	264₋₇₀ ⁺⁵⁴	8.22_−0.16 ^+0.17	10.94_−0.20 ^+0.19
406444	10.87_−0.01 ^+0.08	2.31_−0.11 ^+0.08	8.39_−0.09 ^+0.08	7.77_−0.09 ^+0.08	211₋₄₉ ⁺⁴⁰	7.64_−0.13 ^+0.12	10.26_−0.17 ^+0.15
3	10.43_−0.08 ^+0.06	1.73_−0.15 ^+0.18	8.40_−0.06 ^+0.05	7.69 ± 0.06	205₋₂₉ ⁺²⁶	7.71_−0.13 ^+0.15	10.32_−0.14 ^+0.16
434585	10.13_−0.12 ^+0.04	1.73_−0.01 ^+0.17	8.47_−0.14 ^+0.10	7.69 ± 0.25	172₋₆₄ ⁺⁴⁶	7.67_−0.22 ^+0.20	10.21_−0.28 ^+0.23
192129	10.45_−0.00 ^+0.12	1.35_−0.08 ^+0.00	8.41_−0.08 ^+0.07	7.83_−0.04 ^+0.03	199₋₄₁ ⁺³³	7.40_−0.28 ^+0.18	10.00_−0.29 ^+0.19
217753	10.39_−0.06 ^+0.05	1.67_−0.18 ^+0.12	8.57_−0.05 ^+0.04	7.55 ± 0.06	131₋₁₇ ⁺¹⁵	<8.00	<10.42
218783	10.12_−0.07 ^+0.08	1.70_−0.17 ^+0.16	8.42_−0.07 ^+0.06	7.67 ± 0.03	197₋₃₅ ⁺³⁰	<7.84	<10.43
212298	10.38_−0.08 ^+0.14	1.74_−0.20 ^+0.26	8.39_−0.08 ^+0.07	7.76_−0.04 ^+0.03	213₋₄₃ ⁺³⁵	<7.84	<10.47
413391	10.08_−0.01 ^+0.00	1.88 ± 0.00	8.33_−0.10 ^+0.08	7.69_−0.07 ^+0.06	246₋₅₆ ⁺⁴⁷	<7.73	<10.42
5	10.14_−0.05 ^+0.09	1.37 ± 0.15	8.37 ± 0.05	7.59 ± 0.07	220₋₂₇ ⁺²⁵	<7.33	<9.97
434618	9.39_−0.01 ^+0.12	1.18_−0.09 ^+0.00	8.26_−0.09 ^+0.08	7.78_−0.04 ^+0.03	284₋₅₈ ⁺⁴⁹	<7.28	<10.04
stack^d	10.18_−0.06 ^+0.05	1.52_−0.14 ^+0.15	8.38 ± 0.05	7.62 ± 0.06	216 ± 27	7.33_−0.15 ^+0.17	9.97_−0.16 ^+0.18

Notes.

^aMedian value of the PDF obtained from magphys. Error bars correspond to the 16th–68th percentiles.^bEstimated using an empirical calibration method by Curti et al. (2017).^cThe 97.5th percentile values from magphys are given as upper limits.^dStacking result for the five ALMA nondetected sources with log(M_*/M_⊙) = 10.0–10.4 (Section 2.3).

Download table as: ASCII Typeset image

In order to evaluate whether the upper limits on the dust masses are reasonable, we estimated dust mass upper limits with a different method. We calculated a ratio between the dust mass and the luminosity density at 997.6 GHz in the rest-frame, M_dust/L_{997.6 GHz}, for each detected source. We then converted the 3σしぐま upper limits of 1.3 mm fluxes into the dust mass upper limits with the median M_dust/L_{997.6 GHz} ratio. The difference of the rest-frame frequencies among the sources was corrected for assuming L_νにゅー ∝ λらむだ^−3.7. The estimated dust mass upper limits are similar to the 97.5th percentile values of the PDF from magphys.

One of the uncertainties on the dust mass is the assumed dust mass absorption coefficient, κかっぱ_abs. It is reported that dust masses obtained with magphys are lower by a factor of two than those estimated based on the Draine & Li (2007) models, which assume the smaller dust mass absorption coefficient of κかっぱ_abs = 0.38 g⁻¹ cm² (Hunt et al. 2019).

3. Analysis

3.1. Gas-phase Metallicity and Ionization Parameter

We recalculated the gas-phase metallicities with the following four relations, which are locally calibrated in Curti et al. (2017):

$\begin{eqnarray}&&\mathrm{log}\ {R}_{2}=0.418-0.961x-3.505{x}^{2}-1.949{x}^{3},\end{eqnarray} \tag{ 1 }$

$\begin{eqnarray}&&\mathrm{log}\ {R}_{3}=-0.277-3.549x-3.593{x}^{2}-0.981{x}^{3},\end{eqnarray} \tag{ 2 }$

$\begin{eqnarray}&&\mathrm{log}\ {{\rm{O}}}_{32}=-0.691-2.944x-1.308{x}^{2},\end{eqnarray} \tag{ 3 }$

$\begin{eqnarray}&&\mathrm{log}\ {R}_{23}=0.527-1.569x-1.652{x}^{2}-0.421{x}^{3},\end{eqnarray} \tag{ 4 }$

where R₂ = [O ii]/Hβべーた, R₃ = [O iii]λらむだ5007/Hβべーた, O₃₂ = [O iii]λらむだ5007/[O ii], R₂₃ = ([O iii]λらむだλらむだ4959,5007 + [O ii])/Hβべーた, and x is $12+\mathrm{log}({\rm{O}}/{\rm{H}})$ normalized to the solar value. The emission line fluxes of the sources are available in O16 and S17.

We also estimated the ionization parameter, log(q), for the galaxies observed with ALMA as done in O16. The ionization parameter is described as the ratio of the number of the ionizing photons and the hydrogen atoms to be ionized. The definition of q is as follows:

$\begin{eqnarray}&&q=\displaystyle \frac{{Q}_{{{\rm{H}}}^{0}}}{4\pi {R}_{s}^{2}{n}_{{n}_{{\rm{H}}}}},\end{eqnarray} \tag{ 5 }$

where ${Q}_{{{\rm{H}}}^{0}}$ is the flux of the ionizing photons produced by the existing stars above the Lyman limit, R_s is the Strömgren radius, and n_H is the local density of hydrogen atoms (Kewley & Dopita 2002).

We use the following relation by Kobulnicky & Kewley (2004) to estimate the ionization parameter from the [O iii]λらむだλらむだ4959,5007/[O ii] ratio (O₃₂) and gas-phase metallicity;

$\begin{eqnarray}&&\begin{array}{l}\mathrm{log}(q)=\{32.81-1.153{y}^{2}+[12+\mathrm{log}({\rm{O}}/{\rm{H}})]\\ \times \ (-3.396-0.025y+0.1444{y}^{2})\}\\ \times \ \{4.63-0.3119y-0.163{y}^{2}+[12+\mathrm{log}({\rm{O}}/{\rm{H}})]\\ \times \ (-0.48+0.0271y+0.02037{y}^{2})\}{}^{-1},\end{array}\end{eqnarray} \tag{ 6 }$

where $y=\mathrm{log}\ {{\rm{O}}}_{32}$ .

3.2. M_–SFR and M_–Gas-phase Metallicity Diagram

Figure 2 shows the star-forming main sequence and the mass–metallicity relation for the star-forming galaxies at z ∼ 3.3. In the left panel, we also show star-forming galaxies and SMGs at z = 3–4 from the literature, which are introduced in Section 3.4. In the right panel, we show stacking results at z ∼ 3.3 from the MOSDEF survey (Sanders et al. 2020). We use the line ratios given in Sanders et al. (2020) and the same metallicity calibration method as shown in Section 3.1. Our targets distribute around the star-forming main sequence and the mass–metallicity relation, and thus, are not biased in terms of the star-forming activity and gas-phase metallicity. The stacked sample is also close to the star-forming main sequence and the mass–metallicity relation, indicating that the stacked sample has a typical SFR and gas-phase metallicity for its stellar mass.

3.3. Gas Mass

We converted the dust masses from magphys into gas masses with the relation between the gas-phase metallicity and gas-to-dust mass ratio. We used the relation shown in Magdis et al. (2012) as follows:

$\begin{eqnarray}&&\begin{array}{l}\mathrm{log}({\delta }_{\mathrm{gdr}})=(10.54\pm 1.0)\\ \quad -\ (0.99\pm 0.12)\times (12+\mathrm{log}({\rm{O}}/{\rm{H}})),\end{array}\end{eqnarray} \tag{ 7 }$

which is based on the relation of Leroy et al. (2011) and uses the metallicity estimated with the [N ii]/Hαあるふぁ ratio of Pettini & Pagel (2004). Note that the dust mass estimation in Leroy et al. (2011) and Magdis et al. (2012) is based on the Draine & Li (2007) models. The scatter of this relation is 0.15 dex (Magdis et al. 2012).

We need to convert the gas-phase metallicity in Table 2 into that based on the Pettini & Pagel (2004) calibration. We estimated [N ii]/Hαあるふぁ ratios using the relation between $12+\mathrm{log}({\rm{O}}/{\rm{H}})$ and [N ii]/Hαあるふぁ of Curti et al. (2017), and then converted the estimated [N ii]/Hαあるふぁ ratios into the gas-phase metallicities using the Pettini & Pagel (2004) calibration. The empirical relation between $12+\mathrm{log}({\rm{O}}/{\rm{H}})$ and [N ii]/Hαあるふぁ has a scatter of 0.1 dex in the metallicity direction (Curti et al. 2017). This scatter causes ∼0.19 dex uncertainty on average on the estimated [N ii]/Hαあるふぁ ratios for our sample. Given that the Pettini & Pagel (2004) calibration has a scatter of 0.18 dex, the converted gas-phase metallicities have a typical uncertainty of 0.26 dex.

Then, we estimated gas masses as follows:

$\begin{eqnarray}&&{M}_{\mathrm{gas}}={M}_{\mathrm{dust}}\times {\delta }_{\mathrm{gdr}},\end{eqnarray} \tag{ 8 }$

where M_gas includes both molecular and atomic hydrogen. We multiplied our dust masses by a factor of two when converting them into gas masses with Equation (8) to correct for the systematic difference of the dust mass estimation between magphys and Draine & Li (2007) models (Section 2.4). The estimated gas masses and limits of the individual sources and the stacked sample are summarized in Table 2.

Given the ∼0.30 dex uncertainty coming from the assumed κかっぱ_abs (Hunt et al. 2019), ∼0.26 dex uncertainty on the converted gas-phase metallicities and the ∼0.15 dex scatter of the relation between the gas-phase metallicity and gas-to-dust mass ratio (Magdis et al. 2012), the systematic uncertainty of our gas mass estimate is roughly 0.42 dex.

3.4. Comparison Sample from the Literature

We next introduce samples from the literature to which we compare our data in Section 4. Because different works use different approaches to estimate dust and/or gas masses, these comparisons must be interpreted with care. We refer to the papers cited below for more details about the samples selection, observations, and methods used to estimate dust and/or gas masses.

1.
Magdis et al. (2017) investigated the dust and gas masses of two massive Lyman break galaxies (LBGs) at z ∼ 3. The dust masses were estimated with the Draine & Li (2007) models. They used several independent methods to estimate the gas masses, namely, CO(3–2) line, dust mass from the IR SED, and the empirical relation of Groves et al. (2015).
2.
Wiklind et al. (2019) targeted star-forming galaxies at z ∼ 3. They used the empirical relation of Scoville et al. (2016) to estimate molecular gas masses. We show 11 galaxies with individual molecular gas estimates.
3.
ASPECS: We extract two galaxies at z ∼ 3.6 from the ASPECS 1.2 mm continuum source catalog (Aravena et al. 2020). The dust masses are estimated from SED fitting with magphys. We use the gas masses estimated with the dust mass and the fixed gas-to-dust mass ratio of 200 in their catalog.
4.
Cassata et al. (2020) observed CO emission lines for massive LBGs at z ∼ 3–4. They used the CO(5–4) emission lines to estimate molecular gas masses.
5.
AS2UDS is an ALMA survey targeting 700 SMGs (Dudzevičiūtė et al. 2020). Here we show the AS2UDS galaxies at z = 3–4. Dudzevičiūtė et al. (2020) estimated the dust masses of the AS2UDS galaxies from the SED fitting with magphys. They converted the dust mass into the gas mass assuming a fixed gas-to-dust mass ratio of 100.
6.
Tan et al. (2014) investigated the dust masses of three SMGs at z = 4.05 with multiband photometry in the IR regime. They used the Draine & Li (2007) dust models to estimate the dust masses.

We also introduce the following galaxy samples at lower redshifts, which have individual measurements of gas mass and gas-phase metallicity.

1.
Saintonge et al. (2013) investigated the dust and molecular gas masses of 17 gravitationally lensed star-forming galaxies at z ∼ 1.4–3. They used the Draine & Li (2007) models to estimate dust masses. Molecular gas masses were estimated from the CO(3–2) lines. We show 4 galaxies at z ∼ 2–3 in Section 4.1 and 4.4.
2.
Seko et al. (2016a) investigated the dust and molecular gas masses of star-forming galaxies at z ∼ 1.4. They converted the dust continuum flux into a dust mass assuming modified blackbody with fixed T_dust = 30 K and βべーた = 1.5. They used the CO(5–4) line to estimate the molecular gas masses.
3.
xCOLD-GASS is a CO(1–0) line survey for local Sloan Digital Sky Survey (SDSS) galaxies. We used the public catalog in Saintonge et al. (2017) and combined it with the catalog of the xGASS project (Catinella et al. 2018) to estimate the total gas masses. The stellar masses of the xCOLD-GASS galaxies are in the range of $\mathrm{log}({M}_{* }/{M}_{\odot })\,=$ 9.1–11.2.
4.
ALLSMOG is a CO(2–1) line survey for local SDSS galaxies (Bothwell et al. 2014). Most of the ALLSMOG galaxies have the measurements of the atomic hydrogen gas by different studies (see Cicone et al. (2017) for more details). The stellar mass range of the ALLSMOG galaxies is $\mathrm{log}({M}_{* }/{M}_{\odot })\,=$ 9.3–10.0.

4. Results and Discussion

4.1. Dust Mass and Its Metallicity Dependence

The dust masses of the star-forming galaxies at z ∼ 3.3 are estimated to be $\mathrm{log}({M}_{\mathrm{dust}}/{M}_{\odot })\,\sim$ 7.4–8.9 (Table 2). The dust mass of the stacked sample is $\mathrm{log}({M}_{\mathrm{dust}}/{M}_{\odot })={7.33}_{-0.15}^{+0.17}$ . The left panel of Figure 3 shows the comparison of dust masses between the star-forming galaxies at z ∼ 3.3 and the galaxies at z ∼ 1.4–4 in the literature (Section 3.4). As mentioned in Section 3.4, Tan et al. (2014), Magdis et al. (2017), and Saintonge et al. (2013) estimated dust masses with the Draine & Li (2007) models. To correct for the systematic difference of the dust mass estimate between magphys and Draine & Li (2007) models, the dust masses of the galaxies in these studies are divided by a factor of two in the left panel of Figure 3.

The right panel of Figure 3 shows the relation between the gas-phase metallicity and dust-to-stellar mass ratio for the star-forming galaxies at z ∼ 3.3. Given that dust is produced from metals, we would expect galaxies with higher metallicities to have higher dust masses at a given stellar mass. We find no statistically significant trend between the dust-to-stellar mass ratio and gas-phase metallicity in our sample.

The brightest source at 1.3 mm in our sample, 208681 (Table 1 and Figure 3), has a dust mass of $\mathrm{log}({M}_{\mathrm{dust}}/{M}_{\odot })={8.89}_{-0.08}^{+0.11}$ , which is comparable to those of SMGs at z ∼ 3–4 (Tan et al. 2014; Dudzevičiūtė et al. 2020). This source can be classified as an SMG in terms of its dust content. The other five galaxies and the stacked sample have ∼1 dex lower dust masses than the SMGs at z ∼ 3–4 with similar stellar masses. This brightest source, 208681, is also the most metal-rich galaxy with $12+\mathrm{log}({\rm{O}}/{\rm{H}})={8.59}_{-0.05}^{+0.04}$ in our sample and appears to be distributed apart from the other galaxies in the right panel of Figure 3. This may suggest that SMGs have a different relation between the gas-phase metallicity and dust-to-stellar mass ratio from UV/optical-selected star-forming galaxies.

The star-forming galaxies at z ∼ 3–4 except for the SMGs show a positive correlation between the stellar mass and dust mass. The dust-to-stellar mass ratio takes a value of between 1 × 10⁻³ and 5 × 10⁻³ (median value is 2 × 10⁻³) in the stellar mass range of $\mathrm{log}({M}_{* }/{M}_{\odot })\,\sim$ 10.1–11.4.

When they are compared with star-forming galaxies at z ∼ 1.4 (Seko et al. 2016a) and z ∼ 2–3 (Saintonge et al. 2013), star-forming galaxies at lower redshifts have similar dust-to-stellar mass ratios (1 × 10⁻³–5 × 10⁻³) as our galaxies at z ∼ 3.3 with similar stellar masses. Although a fair comparison is difficult due to different sample selection of the three studies, the evolution of the dust-to-stellar mass ratios seems to be mild since z ∼ 1.4–3.3 as shown in Béthermin et al. (2015).

4.2. Gas Properties of Star-forming Galaxies at z = 3–4

Figure 4 shows the gas mass fraction, f_gas = M_gas/(M_gas + M_*), and gas depletion timescale, t_dep = M_gas/SFR, as a function of stellar mass for the star-forming galaxies at z ∼ 3.3. Our estimated gas mass fractions are 0.20–0.75, and the gas depletion timescales are 0.09–1.55 Gyr. The typical uncertainties of f_gas and t_dep are ±0.09 and ±0.22 dex, respectively. Given that the gas masses have the systematic uncertainty of ∼0.42 dex (Section 3.3), f_gas and t_dep have an additional error of ∼0.19 and 0.42 dex (1σしぐま), respectively. As for the stacked sample, the gas mass fraction and gas depletion timescale are estimated to be ${f}_{\mathrm{gas}}={0.38}_{-0.09}^{+0.10}$ and ${t}_{\mathrm{dep}}={0.28}_{-0.14}^{+0.15}$ Gyr, respectively.

We show star-forming galaxies and SMGs at z ∼ 3–4 from the literature (Section 3.4) in Figure 4. The solid line in each panel represents the scaling relation for galaxies on the star-forming main sequence at z ∼ 3.3 from Tacconi et al. (2018). The dashed lines correspond to the case when galaxies are at±0.3 dex from the star-forming main sequence. Our sample including the stacking result reaches down to f_gas ∼ 0.2–0.3, which is lower by a factor of ≳2 than the scaling relation. We also find a scatter of ≳1 dex for the gas depletion timescale at a fixed stellar mass. Such a large scatter of the gas properties is also seen in the samples of Wiklind et al. (2019) and Dudzevičiūtė et al. (2020).

It has been reported that the gas mass fraction and depletion timescale of the main-sequence galaxies gradually change depending on the deviation from the star-forming main sequence (Δでるた_MS; e.g., Saintonge et al. 2012; Sargent et al. 2014; Tacconi et al. 2018). We checked how the offset of the gas mass fraction or depletion timescale from the scaling relation for galaxies on the main sequence changes depending on Δでるた_MS. We find a trend consistent with Tacconi et al. (2018) when we combine our sample with the samples of the literature (Magdis et al. 2017; Wiklind et al. 2019; Aravena et al. 2020; Cassata et al. 2020). However, at a fixed Δでるた_MS, our sample shows a large scatter of the gas mass fraction and depletion timescale. The observed scatter of the gas mass fraction and depletion timescale in our sample cannot be explained by Δでるた_MS alone. These results suggest that the fundamental gas properties of galaxies are highly diverse even when they have similar stellar masses and SFRs (Elbaz et al. 2018). Given that the scaling relations are possibly biased toward dusty and gas-rich galaxies especially at higher redshifts, the scaling relations may not be representative of the majority of the galaxy populations at z ∼ 3–4.

Given that we use the metallicities to derive the gas properties, the observed trends as a function of stellar mass in Figure 4 may be partly caused by the mass–metallicity relation. However, the distribution of our sample in Figure 4 does not change significantly when a constant gas-to-dust mass ratio is assumed. This means that our results are not affected by the fact that we use the gas-phase metallicities to derive the gas properties.

4.3. Gas Mass Fraction versus Physical Conditions of the Ionized Gas

We investigate the relation between the gas mass fraction and the physical conditions of the ionized gas, namely, gas-phase metallicity and ionization parameter (Section 3.1), for our galaxy sample at z ∼ 3.3. Figure 5 shows the comparison between the gas mass fraction and gas-phase metallicity. We note that the gas mass fraction and gas-phase metallicity (and also ionization parameter) are not independent, as mentioned in the previous section. In Figure 5 we show the dependence of the two quantities on each other when the dust-to-stellar mass ratio is fixed with dashed lines.

**Figure 5.** Gas mass fraction vs. gas-phase metallicity for our galaxy sample at z ∼ 3.3. The horizontal bar in the bottom right corner shows an additional ±1σしぐま error on the gas mass fraction coming from the systematic uncertainty on M_gas. The dashed lines show the dependence of the two quantities on each other when the dust-to-stellar mass ratio is fixed at 1 × 10⁻³ and 3 × 10⁻³. We find no statistically significant correlation between the gas mass fraction and gas-phase metallicity in our sample.
Download figure:
Standard image High-resolution image

Given that the abundance of oxygen with respect to hydrogen changes depending on the amount of the hydrogen gas in galaxies, the gas-phase metallicity, 12+log(O/H), would be expected to decrease as the gas mass fraction increases (e.g., Bothwell et al. 2013a, 2016a; Zahid et al. 2014; Seko et al. 2016a). However, we find no statistically significant correlation between the gas mass fraction and gas-phase metallicity for the star-forming galaxies at z ∼ 3.3, which is the same as the result obtained from the comparison between the gas-phase metallicity and dust-to-stellar mass ratio (Figure 3).

We find no clear correlation between the gas mass fraction and ionization parameter either. According to Kashino & Inoue (2019), the gas mass fraction and ionization parameter are related indirectly via three parameters, namely, the specific SFR (sSFR), the gas-phase metallicity, and the electron density. When the gas-phase metallicity increases or the sSFR decreases, both the gas mass fraction and ionization parameter decreases. When the electron density increases, the gas mass fraction increases but the ionization parameter decreases (Kashino & Inoue 2019). Because the gas mass fraction and ionization parameter depend on the three parameters in a different way, the correlation of the gas mass fraction with the ionization parameter is not straightforward. We would need to fix some of the parameters to investigate the trend between the two quantities.

The lack of a clear correlation between gas mass fractions and the physical conditions of the ionized gas may reflect stochastic star formation histories for star-forming galaxies at high redshifts. Star formation in galaxies at higher redshifts is suggested to be burstier than in local galaxies (Guo et al. 2016; Faucher-Giguère 2018; Tacchella et al. 2020). When the star-forming activity in galaxies changes on a short timescale, it becomes more difficult to identify a global trend between the physical quantities.

Note that our sample size may be too small to find any correlation. We need a larger sample of galaxies covering a wider range of stellar mass to confirm whether the gas mass fraction correlates with gas-phase metallicity and ionization parameter.

4.4. Comparison of Galaxies at z = 0–3.3 in the f_gas versus 12+log(O/H) Diagram

Figure 6 shows the star-forming galaxies from z = 0 to 3.3 (Section 3.4) in the gas mass fraction versus metallicity diagram.

In the literature (Saintonge et al. 2013, 2017; Seko et al. 2016a; Cicone et al. 2017), the gas-phase metallicities are estimated with the [N ii]/Hαあるふぁ ratios. In order to compare with the gas-phase metallicities of our sample, which are estimated based on [O iii], Hβべーた, and [O ii] lines (Curti et al. 2017), we convert the given [N ii]/Hαあるふぁ ratios in the previous studies into the gas-phase metallicities using the empirical relation between $12+\mathrm{log}({\rm{O}}/{\rm{H}})$ and [N ii]/Hαあるふぁ of Curti et al. (2017).

The gas mass fraction of the galaxies at z ∼ 0 and 3.3 is the total (molecular+atomic) gas mass fraction to compare with a gas regulator model, in which the atomic and molecular hydrogen are indistinguishable, in the following sections. The gas mass fraction of the galaxies in Seko et al. (2016a) and Saintonge et al. (2013) is the molecular gas mass fraction. We expect that the comparison in Figure 6 is not significantly affected by the fact that we do not include the atomic gas for the two samples. The fraction of the molecular gas is suggested to increase with increasing redshifts because of the higher surface density of galaxies at higher redshifts (e.g., Popping et al. 2015). Popping et al. (2015) suggest the fraction of the molecular gas in the total gas is ∼0.6–0.8 at z ∼ 1.5–3.0 based on their simulations.

Focusing on the local galaxies in Figure 6, the gas-phase metallicity gradually decreases with increasing gas mass fraction, as shown in Bothwell et al. (2013a) and Hunt et al. (2015). Such a gradual decrease of the gas-phase metallicity with increasing gas mass fraction indicates that we need to cover a wide range of gas mass fractions to identify the correlation between the two quantities.

Whereas the star-forming galaxies at z ∼ 1.4 from Seko et al. (2016a) appear to be located at the gas-rich end of the distribution of the local star-forming galaxies, the galaxies at z ≳ 2 from this study and Saintonge et al. (2013) show an offset toward the lower gas-phase metallicities (∼0.2 dex) with respect to the distribution of the local galaxies. This result suggests that star-forming galaxies at z ≳ 2 are less chemically enriched than those at z = 0 and even at z ∼ 1.4 with similar gas mass fractions.

The molecular gas mass is estimated with the CO lines in the literature (Section 3.4). Although the systematic difference caused by using different methods to estimate gas mass could change the relative distribution of the galaxies in the horizontal direction, it cannot explain the offset of the galaxies at z ∼ 3.3 toward the low gas-phase metallicity with respect to the local galaxies. As for the gas-phase metallicity, Curti et al. (2017) showed that the offset of gas-phase metallicities calibrated with different line ratios is 0.04 dex on average. The systematic uncertainty caused by using different line ratios to calibrate the metallicity is also unlikely to affect our results. We have a caveat on our gas-phase metallicity measurement for the ALMA-detected, dusty star-forming galaxies in our sample. Herrera-Camus et al. (2018) reported that gas-phase metallicities calibrated with rest-frame optical emission lines tend to be lower than those calibrated with FIR fine-structure lines by up to a factor of two for local (U)LIRGs. The FIR lines are likely to trace the ionized gas in the dense and dusty star-forming regions, which are no longer traced by the optical emission lines, and such dense and dusty star-forming regions would be more metal enriched (e.g., Santini et al. 2010). When this is also the case for our ALMA-detected galaxies at z ∼ 3.3, the offset toward the low metallicity with respect to the local galaxies in Figure 6 could be explained by the underestimated gas-phase metallicities for our sample. However, when we compare the dust extinction values, A_V, between our sample and local (U)LIRGs in Rupke et al. (2008), the median A_V of our ALMA-detected galaxies (∼0.6 mag) is much smaller than that of local (U)LIRGs (∼3.6 mag). This implies that the ALMA-detected galaxies in our sample are not as dusty as the local (U)LIRGs, and thus, that the metallicities calibrated with the optical emission lines can be regarded as representative values for our sample at z ∼ 3.3.

4.4.1. Comparison with a Gas Regulator Model

"Equilibrium," "bathtub", or "gas regulator" models are used to track the evolution of the fundamental physical quantities of galaxies, such as gas mass, SFR, and metallicity, by considering gas inflows, outflows, star formation, and metal production in galaxies (e.g., Finlator & Davé 2008; Bouché et al. 2010; Davé et al. 2012; Dayal et al. 2013; Lilly et al. 2013; Peng & Maiolino 2014; Tacchella et al. 2020). We compare the observational data at z = 0–3.3 with the gas regulator model by Peng & Maiolino (2014).

Peng & Maiolino (2014) derived the analytic formula to track the evolution of the physical quantities, such as gas mass, SFR, metallicity, and stellar mass. The input parameters of this model are gas inflow rate (Φふぁい), star formation efficiency (εいぷしろん = SFR/M_gas), mass-loading factor (λらむだ = outflow rate/SFR), and return mass fraction (R). The gas accretion onto the galaxy is assumed to scale with the growth rate of the dark matter halo. The dark matter halo growth rate is derived from the cosmological hydrodynamic simulations (Faucher-Giguère et al. 2011). The outflow rate is assumed to be proportional to SFR. The return mass fraction takes on values ∼0.2 to ∼0.5 depending on the IMF. This model assumes that these input parameters are constant with time or change with longer timescales than the equilibrium timescale. The equilibrium timescale is the timescale to reach the equilibrium state, where gas acquisition by inflows balance with gas consumption by star formation and outflows. The equilibrium timescale is expressed as follows:

$\begin{eqnarray}&&{\tau }_{\mathrm{eq}}=\displaystyle \frac{1}{\varepsilon (1-R+\lambda )}.\end{eqnarray} \tag{ 9 }$

The time evolution of the gas mass fraction and gas-phase metallicity is described as follows:

$\begin{eqnarray}&&{f}_{\mathrm{gas}}(t)=\displaystyle \frac{1}{1+\varepsilon (1-R)\left(\tfrac{t}{1-{e}^{-\tfrac{t}{{\tau }_{\mathrm{eq}}}}}-{\tau }_{\mathrm{eq}}\right)},\end{eqnarray} \tag{ 10 }$

$\begin{eqnarray}&&{Z}_{\mathrm{gas}}(t)=[{Z}_{0}+y{\tau }_{\mathrm{eq}}\varepsilon (1-{e}^{-\tfrac{t}{{\tau }_{\mathrm{eq}}}})][1-{e}^{-\tfrac{t}{{\tau }_{\mathrm{eq}}(1-{e}^{-t/{\tau }_{\mathrm{eq}}})}}],\end{eqnarray} \tag{ 11 }$

where Z₀ is the metallicity of the infalling gas, and y is the average yield per stellar generation.

In Figure 6 we show the model tracks obtained from the gas regulator model of Peng & Maiolino (2014). We assume R = 0.4 (for the Chabrier IMF; Madau & Dickinson 2014) and y = 1.5Z_⊙ (e.g., Yabe et al. 2015). The gas depletion timescale (1/εいぷしろん) is set to be t_dep = 0.8 Gyr. Note that the normalization of model tracks in Figure 6 does not depend on the absolute value of t_dep. We assume five different mass-loading factors between λらむだ = 0.5 and 2.5 (Figure 6).

We find that the distribution of the star-forming galaxies at z ∼ 3.3 and those from Saintonge et al. (2013) can be broadly explained by the model tracks with the high mass-loading factor of λらむだ ∼ 2.0–2.5 rather than lower values such as λらむだ ∼ 0.5 or 1. We need stronger outflow with larger λらむだ to achieve a lower gas-phase metallicity for star-forming galaxies at z ≳ 2 than the local ones with similar gas mass fractions. This result may suggest a redshift evolution of the mass-loading factor λらむだ from z = 0 to 3.3, as we discuss below.

Yabe et al. (2015) showed the increasing outflow rate normalized by SFR with increasing redshifts up to z ∼ 2 by comparing the observational data (stellar mass, gas mass fraction, and gas-phase metallicity) with a simple chemical evolution model (see also Troncoso et al. 2014). Some theoretical studies based on analytic models or numerical simulations showed the redshift evolution of the mass-loading factor (e.g., Barai et al. 2015; Mitra et al. 2015; Muratov et al. 2015; Hayward & Hopkins 2017). Observationally, Sugahara et al. (2017) showed a trend that the star-forming galaxies at higher redshift (up to z = 2) have larger mass-loading factor at a fixed circular velocity. Our results obtained from the comparison between the observational data, and the model tracks support the idea that star-forming galaxies at higher redshifts have larger mass-loading factors, and thus, more massive outflow.

4.4.2. Equilibrium Timescale

We estimate the equilibrium timescales (Equation (9)) for the star-forming galaxies at z ∼ 3.3 with the gas depletion timescale obtained from the observation and the mass-loading factor inferred from the comparison with the model tracks in Figure 6. The equilibrium timescales of the ALMA-detected sources are estimated to be 0.03–2.21 Gyr (average value: 0.52 Gyr) assuming R = 0.4 for the Chabrier IMF (Madau & Dickinson 2014).

According to Peng & Maiolino (2014), when the equilibrium timescale is much shorter than the Hubble time, galaxies are expected to be in the equilibrium state, where gas acquisition by inflows and gas consumption by star formation and outflows are balanced. On the other hand, when the equilibrium timescale is comparable to the Hubble time, galaxies are considered to have a much larger gas reservoir and to be out of equilibrium.

The average equilibrium timescale of the detected sources is roughly one order of magnitude shorter than the Hubble time at z = 3.3 (2.82 Gyr). However, not all of the galaxies necessarily start forming stars at the beginning of the universe. Given that the age of galaxies must be younger than the Hubble time, the equilibrium timescale should probably be compared with the age of the galaxies rather than the Hubble time.

We here use the ratio of M_*/SFR, which can be regarded as the minimum age of a galaxy. The star-forming galaxies at z ∼ 3.3 have M_*/SFR = 0.25–1.25 Gyr (average: 0.57 Gyr), which is closer to the equilibrium timescales than the Hubble time. Especially the galaxies with relatively larger gas mass fractions, f_gas ∼ 0.6–0.8, in our sample tend to have the equilibrium timescales comparable to the minimum ages. This result may suggest that normal star-forming galaxies at z ∼ 3 with relatively large gas mass fractions have not yet reached the equilibrium state, as suggested in Mannucci et al. (2010).

In the future, it will be of interest to study how our results are affected when the assumption is relaxed that galaxies are in equilibrium and when bursty star formation histories are considered (Tacchella et al. 2020). Direct measurements of gas outflow and inflow rates would be also important to further investigate whether the star-forming galaxies at z ∼ 3.3 are out of equilibrium. The spatially resolved emission line maps for the individual galaxies will enable us to search for outflow signatures and estimate the mass outflow rates (e.g., Genzel et al. 2011; Davies et al. 2019). Furthermore, a simulation study suggests a correlation between metallicity gradients and gas accretion rates (Collacchioni et al. 2020). We would be able to investigate the gas inflow rates by obtaining metallicity gradients from the spatially resolved emission line maps.

5. Summary

We conducted ALMA Band-6 observations of star-forming galaxies at z ∼ 3.3, whose metallicity measurements are based on rest-frame optical spectroscopy. Thus we were able to directly compare the metallicities with the dust and inferred gas properties from our ALMA observations for star-forming galaxies at z ∼ 3.3. We detected the dust continuum emission individually from six out of 12 galaxies. We stacked the ALMA maps of the five ALMA nondetected sources with $\mathrm{log}({M}_{* }/{M}_{\odot })\,=$ 10.0–10.4 and obtained a ∼5σしぐま detection of this sample.

We estimated dust masses from SED fitting with magphys including the 1.3 mm fluxes from ALMA. We converted the dust mass into the gas mass with a relation between the gas-phase metallicity and gas-to-dust mass ratio. With the estimates of dust mass, gas mass, and the physical conditions of the ionized gas, we conclude the following:

1.
The median value of the dust-to-stellar mass ratios is M_dust/M_* ∼ 3.0 ± 2.0 × 10⁻³. The dust-to-stellar mass ratio of the stacked sample is ∼1.4 ± 0.5 × 10⁻³. We find no clear trend between the dust-to-stellar mass ratio and gas-phase metallicity.
2.
The estimated gas mass fractions and gas depletion timescales are f_gas = 0.20–0.75 and t_dep = 0.09–1.55 Gyr, respectively. The stacked sample shows ${f}_{\mathrm{gas}}={0.38}_{-0.09}^{+0.10}$ and ${t}_{\mathrm{dep}}={0.28}_{-0.14}^{+0.15}$ Gyr. The gas mass fractions and gas depletion timescales of the galaxies at z ∼ 3.3 show a wider spread at a fixed stellar mass than the scaling relations of galaxies on the main sequence at z ∼ 3.3. Given that most of our galaxies at z ∼ 3.3 distribute around the star-forming main sequence with ±0.3 dex, the large scatter of the gas mass fraction and depletion timescale may suggest a significant diversity of these fundamental properties within the so-called main sequence.
3.
We find no clear correlation between the gas mass fraction and the physical conditions of the ionized gas, namely, gas-phase metallicity and ionization parameter, at z ∼ 3.3. We may require a larger sample of galaxies covering a wider range of the physical quantities to confirm whether gas mass fractions correlate with the ionized gas conditions.
4.
Comparing star-forming galaxies at different redshifts in the gas mass fraction versus metallicity diagram, we find that the star-forming galaxies at z ≳ 2 show an offset toward lower metallicities compared to the distribution of local star-forming galaxies, in the sense that star-forming galaxies at z ≳ 2 appear to be more metal poor than the local galaxies with similar gas mass fractions.
5.
We find that the distribution of star-forming galaxies at z ∼ 3.3 in the gas mass fraction versus gas-phase metallicity diagram can be broadly explained by models assuming higher mass-loading factors in outflows of λらむだ ∼ 2.0–2.5 from the gas regulator model of Peng & Maiolino (2014). This result supports the idea that star-forming galaxies at higher redshifts have powerful outflows with higher mass-loading factors.
6.
Comparing the equilibrium timescales (Peng & Maiolino 2014) and the minimum ages of the galaxies (M_*/SFR), we find that the equilibrium timescale of the relatively gas-rich galaxies (f_gas ∼ 0.7) is comparable to their minimum ages, suggesting that they may be out of equilibrium.

It remains unclear whether star-forming galaxies at high redshifts follow the same relation between gas-phase metallicity and gas-to-dust mass ratio as local galaxies (Saintonge et al. 2013; Seko et al. 2016a, but see also Magdis et al. 2012; Shapley et al. 2020). Observations of independent gas tracers, such as CO or [C i] emission lines, are required to investigate the relation between the gas-phase metallicity and gas-to-dust mass ratio at z > 3.

Another caveat is whether the metallicities derived from the rest-frame optical emission lines are applicable to dusty star-forming galaxies at high redshifts. Metallicity measurements with FIR fine-structure lines are required to investigate this further.

High-resolution integral-field-unit observation with the James Webb Space Telescope will enable us to investigate the metallicity gradients within the individual galaxies and to search for the outflow signatures within them. The spatially resolved emission line maps would be useful to investigate the effects of gas inflows and outflows more directly.

We thank the anonymous referee for a careful reading and comments that improved the clarity of this paper. T.L.S. would like to thank Ken-ichi Tadaki and Nao Fukagawa for useful comments. I.R.S. acknowledges support from STFC (ST/T000244/1). This work was supported by NAOJ ALMA Scientific Research grant No. 2018-08A. This paper makes use of the following ALMA data: ADS/JAO.ALMA#2018.1.00681.S. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada), MOST 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. Some of the data presented herein were obtained at the W. M. Keck Observatory, which is operated as a scientific partnership among the California Institute of Technology, the University of California and the National Aeronautics and Space Administration. The Observatory was made possibility by the generous financial support of the W. M. Keck Foundation. The authors wish to recognize and acknowledge the very significant cultural role and reverence that the summit of Maunakea has always had within the indigenous Hawaiian community. We are most fortune to have the opportunity to conduct observations from this mountain. Data analyses were in part carried out on the open use data analysis computer system at the Astronomy Data Center, ADC, of the National Astronomical Observatory of Japan (NAOJ).

Facilities: ALMA - Atacama Large Millimeter Array, Keck:I (MOSFIRE). -

Software: astropy (Astropy Collaboration et al. 2013; Price-Whelan et al. 2018), casa (McMullin et al. 2007), topcat (Taylor et al. 2005).

Appendix: magphys Best-fit SEDs

Figure 7 shows the best-fit SEDs from magphys for 12 galaxies observed with ALMA. We also show the best-fit SED for the stacked sample, including the 5 ALMA nondetected sources with $\mathrm{log}({M}_{* }/{M}_{\odot })\,=$ 10.0–10.4 (Section 2.3).

Dust, Gas, and Metal Content in Star-forming Galaxies at z ∼ 3.3 Revealed with ALMA and Near-IR Spectroscopy

Article metrics

Permissions

Author e-mails

Author affiliations

ORCID iDs

Dates

Abstract

1. Introduction

2. Observation and Reduction