EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 688 (August 2024) A&A, 688 (2024) A55 Full HTML
Open Access
Issue
A&A
Volume 688, August 2024
Article Number A55
Number of page(s) 14
Section Extragalactic astronomy
DOI https://doi.org/10.1051/0004-6361/202450081
Published online 01 August 2024

© The Authors 2024

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.

1. Introduction

One of the greatest achievements of the last 20 years in the field of galaxy evolution is the high precision measurement of the cosmic star formation rate density (SFRD) up to z ∼ 4 (Madau & Dickinson 2014). This accomplishment was made possible by combining rest-frame infrared (IR) and ultraviolet (UV) observations probing, respectively, the dust-attenuated and unattenuated star formation rate (SFR) of star-forming galaxies (SFGs). At z > 4, however, estimates of the cosmic SFRD become highly uncertain because current measurements rely primarily on rest-frame UV observations to measure the SFRs of galaxies, accounting for dust attenuation using the local relations between the UV spectral slope (βUV) and LIR/LUV (≡IRX; where LIR and LUV are the IR and UV luminosities, respectively) or between the stellar mass (M*) and IRX (Madau & Dickinson 2014). While such dust attenuation corrections are known to fail for a large fraction of massive galaxies at low redshifts (z ≲ 2.5; Wuyts et al. 2011), initial Lyman break studies argued that the low dust content of z > 4 galaxies made this method reliable in the early Universe (e.g., Bouwens et al. 2012). Observations from the Atacama Large Millimeter/submillimeter Array (ALMA) have since cast doubt on the ability of rest-frame UV data alone to provide an accurate view on the cosmic dust-attenuated SFRD, despite somewhat contradictory literature results. On the one hand, serendipitous detections in the ALMA large program ALPINE suggest that dust-attenuated star formation still dominates the cosmic SFRD up to z ∼ 6 and that initial rest-frame UV measurements have underestimated its contribution in the early Universe (Gruppioni et al. 2020; Khusanova et al. 2021). On the other hand, ALMA continuum follow-ups of UV-selected galaxies at z > 4 have revealed a redshift evolution of the dust attenuation, with galaxies at z > 4 having even lower IR luminosities than expected from their UV spectral slope (Fudamoto et al. 2020a), although opposite results have recently appeared (e.g., Bowler et al. 2024). These conflicting results provide a highly uncertain view on the dust attenuation and dust-attenuated star formation properties of galaxies at z > 4: the IRX–M* and IRX–βUV relations could be consistent with their low redshift incarnation or significantly steeper and flatter, respectively (Fudamoto et al. 2020a; Bowler et al. 2022, 2024); the SFR–M* relation is unknown within a factor ∼4 (Khusanova et al. 2021; Popesso et al. 2023; Goovaerts et al. 2024), although it is commonly used to select “normal” SFGs at these epochs; and the cosmic SFRD is uncertain by a factor ∼10 at z ∼ 4 (Gruppioni et al. 2020; Zavala et al. 2021; Algera et al. 2023; Traina et al. 2024).

Current confusion on the dust attenuation properties of high-redshift galaxies stems from the past use of biased samples of z > 4 SFGs (UV-, IR-, or zspec-selected), and from the fact that these past studies are also based on few objects, so that average measurements are dominated by the intrinsic scatter of the SFR–M* and IRX–M* relations. To constrain statistically these relations and the cosmic SFRD at z > 4, it is important to measure the mean IR and UV luminosities of a large sample of mass-selected SFGs at these redshifts. To this end, we are taking advantage of the deepest optical-to-near-IR data publicly available in the Cosmic Evolution Survey (COSMOS; i.e., COSMOS2020 catalog; Weaver et al. 2022) and we are capitalizing on the availability of a large number of ALMA archival data on this field, collected by the A3COSMOS project (Liu et al. 2019a; Adscheid et al. 2024). We selected SFGs at 4 < z < 5 from the COSMOS2020 catalog, a deep near-infrared izYJHKs-selected catalog (with 5σ limits of 27 and 24.7 magAB in the deepest i and shallowest Ks bands, respectively) that allows for detection of dusty galaxies that are absent from classical UV/optical-selected catalogs (Weaver et al. 2023). We applied an uv-domain stacking analysis to the archival ALMA data of this sample. Our stacking analysis is crucial as it provides measurements of the mean IR luminosities of these SFGs down to a stellar mass of ∼109.5 M⊙, even though most of these galaxies are not detected individually in the ALMA archival data.

The structure of this paper is as follows. In Sect. 2, we present the A3COSMOS database used in our analysis. In Sect. 3, we describe how we built our mass-complete sample of SFGs at 4 < z < 5. In Sect. 4, we describe how we measured the mean IR and UV properties of our galaxies using a uv-domain stacking analysis of the A3COSMOS database and the COSMOS2020 catalog, respectively. In Sects. 5.1–5.4, we present the main results of this study, that is, the IRX–βUV and IRX–M* relations, the SFR–M* relation, the cosmic dust-attenuated SFRD, and the “total” cosmic dust-attenuated SFRD at z ∼ 4.5, respectively. Finally, in Sect. 6, we summarize and conclude.

We assume a flat Λ cold dark matter cosmology with H0 = 70 km s−1 Mpc−1, ΩM = 0.30, and ΩΛ = 0.70. A Chabrier (2003) initial mass function (IMF) is assumed for all stellar masses and SFRs.

2. Data

We used all archival ALMA band 6 and band 7 observations publicly available in COSMOS (RA = 10h00m28.6s, Dec = +02°12′21.0″) as of September 1, 2022. We have limited our analysis to these bands, as they are close to the peak of the dust emission of z ∼ 4.5 SFGs and therefore provide good proxies of their infrared luminosities. In contrast, longer-wavelength bands probe the Rayleigh-Jeans emission of dust at these redshifts, which is a poor proxy for the infrared luminosity. Shorter-wavelength bands would in principle be better proxies of the infrared luminosity of SFGs at z ∼ 4.5, but the number of such observations in the ALMA archive is so small that they would have no impact on our analysis. This database of all archival ALMA band 6 and band 7 observations was assembled by the A3COSMOS project, which aims to process all ALMA projects targeting the COSMOS field in a homogeneous way, and to provide calibrated visibilities, cleaned images, and value-added source catalogs (targeted and serendipitously detected) via a single access portal (Liu et al. 2019a; Adscheid et al. 2024). We limited our analysis to projects from cycles ≥3 because measurement sets from earlier cycles have different definitions of visibility weights, which makes their use in our uv-domain stacking analysis problematic (Wang et al. 2022). Our final database contains 87 ALMA projects (53 in band 6 and 34 in band 7), and has a total of 2100 images (equivalently ALMA pointings; 946 in band 6 and 1154 in band 7).

3. Sample

To build our mass-selected sample of SFGs, we used the deepest data publicly available in COSMOS (i.e., COSMOS2020 catalog; Weaver et al. 2022), as this wide, deep survey is ideal for studying representative samples of massive, high-redshift galaxies. Following Weaver et al. (2023), we used the FARMER COSMOS2020 catalog, which is the areal union of the deep near-infrared UltraVISTA DR4 imaging and the Subaru Suprime-Cam intermediate bands PDR2 that covers 1.279 deg2 after removing contamination due to bright stars. Galaxies in the COSMOS2020 catalog were selected from the near-infrared izYJHKs co-added detection image, with 5σ limits of 27 and 24.7 magAB (2″ aperture) in the deepest i and shallowest Ks bands, respectively. As in Weaver et al. (2023), we selected SFGs from this FARMER COSMOS2020 catalog using a standard NUV − r − J criterion (Ilbert et al. 2013) and only kept SFGs with robust stellar mass estimates (mIRAC3.6 < 26 AB and χ reduced 2 < 10 $ \chi^{2}_{\mathrm{reduced}} < 10 $) and accurate photometric redshifts (Δz/(1 + z) < 0.1) at 4 < z < 5. As detailed in Weaver et al. (2023), these photometric redshift and stellar mass estimates were obtained by fitting with LePhare the FARMER photometry in 30 bands from UV to near-infrared. We adopted the photometric redshift and stellar mass estimates from LP_zPDF and LP_MASS_MED, defined as the median of their likelihood distributions as they are generally less susceptible to template fitting systematics than those taken at the minimum χ reduced 2 $ \chi^{2}_{\mathrm{reduced}} $. From this initial sample, we further excluded AGNs, identified through a X-ray Chandra detection within 0 . ″ $ \overset{\prime \prime }{.} $6 (Weaver et al. 2022) and/or a spectral energy distribution (SED) best fit with an AGN template ( χ reduced 2 (AGN)<0.5× χ reduced 2 (Gal.) $ \chi^2_{\rm reduced}({\rm AGN}) < 0.5\times\chi^2_{\rm reduced}({\rm Gal.}) $; Weaver et al. 2022). These AGN exclusions concern less than 0.1% of the SFGs in our z ∼ 4.5 parent sample, and has a negligible impact on our results. Finally, we kept only galaxies above the stellar mass completeness limit for SFGs in the COSMOS2020 catalog, that is, M* > 109.5 M⊙ at 4 < z < 5 (Weaver et al. 2023). This mass-complete sample, identical to that used in Weaver et al. (2023) to build the stellar mass function at z ∼ 4.5, contains 5810 galaxies. We note that, while this sample may still miss part of the population of optically-dark galaxies that are heavily dust attenuated and have recently been discovered through Spitzer-IRAC and ALMA observations (e.g., Wang et al. 2019; Xiao et al. 2023), their contribution to the dust-attenuated SFRD remains modest as discussed in Sect. 5.4.

Finally, we cross-matched this mass-complete sample with the A3COSMOS database, and excluded galaxies that were not covered by any ALMA observations (i.e., located at regions with a primary beam response of less than 0.2). This reduces our sample to 440 galaxies, with redshift, stellar mass, UV luminosity, and UV spectral slope distributions that are indistinguishable (Kolmogorov-Smirnov test) from those of their parent sample of 5810 galaxies (see Appendix A). Of these 440 galaxies, 64 can be considered as the primary target of a given ALMA observation, as they are located within 3″ of its phase center. In Appendix A, we found that these ALMA primary targets are biased toward bright, massive galaxies in a redshift range favorable to [CII] observations, but do not appear to be biased toward a particular UV spectral slope. These primary targets can thus potentially introduce complex and uncontrollable selection biases into our final ALMA-covered mass-complete galaxy sample. On the other hand, according to the analysis presented in Appendix A, excluding these ALMA primary targets provides us with a fairly representative sample of our parent sample of z ∼ 4.5 SFGs, albeit perhaps slightly biased against bright galaxies. In Sect. 5, we therefore investigate the impact on our results of including or excluding these ALMA primary targets in our stacks. The number of galaxies in each of our stellar mass bins, and the number of galaxies that can be considered as the ALMA primary target, are given in Table 1.

Table 1.

Measured dust-attenuated luminosities and derived properties of our six uv-domain stacks.

4. Methods

4.1. Measuring ⟨LIR⟩ by stacking ALMA data in the uv domain

Stacking the A3COSMOS dataset is a challenge, as this database is heterogeneous in terms of observed frequencies, angular resolution, and integration time. As demonstrated in Wang et al. (2022), this challenge can, however, be overcome. The problem of frequency heterogeneity is solved by prior rescaling of each individual dataset to a common rest-frame frequency using submillimeter SED templates, while the problem of angular resolution heterogeneity is solved by performing our stacking analysis in the uv domain. Finally, the problem of integration time heterogeneity is simply solved by keeping track of the weight of each galaxies in the final stacks.

To perform the uv-domain stacking analysis on the A3COSMOS dataset, we followed the same approach as Wang et al. (2022). Details regarding these steps can be found in their Sects. 3.1 and 3.2. First, we scaled the observed ALMA visibility amplitudes of each galaxy to its rest-frame luminosity at 150 μm (i.e., L 150 rest $ L_{\mathrm{150}}^{\mathrm{rest}} $) by using the redshift of the source, the SED template of Béthermin et al. (2017) and the CASA tasks gencal and applycal. This particular SED template was chosen to be consistent with the LIR-to- L 150 rest $ L_{\mathrm{150}}^{\mathrm{rest}} $ conversion used subsequently (see below) and because it provides a reasonable fit to the Herschel stacks of 4 < z < 5 SFGs performed in Béthermin et al. (2020). The reference wavelength of 150 μm was chosen to limit the impact of the assumed SED template on our results, as for galaxies located at z ∼ 4.5, it actually corresponds to the rest-frame wavelength probed by the ALMA band 7, which dominates in number the A3COSMOS database. For each galaxy, we then shifted the phase center of its visibilities to its coordinate using a CASA-based package, the STACKER (Lindroos et al. 2015). The visibilities of the galaxies in each of our stellar mass bins were then stacked together using the CASA task concat. The cleaned image was generated from the stacked measurement set with the CASA task tclean, using natural weighting and cleaning the image down to 3σ. Finally, we measured the stacked L 150 rest $ L_{\mathrm{150}}^{\mathrm{rest}} $ (i.e., ⟨ L 150 rest ⟩ $ \langle L_{\mathrm{150}}^{\mathrm{rest}}\rangle $) using the central pixel value, as at the angular resolution of our stacked measurement set (i.e., ∼ 0 . ″ 5 − 1 ″ $ {\sim}0{{\overset{\prime\prime}{.}}}5{-}1^{{\prime\prime}} $), our stacked population remains basically spatially unresolved. We verified, however, that using instead 2D Gaussian fits at the phase center of the stacked images gave consistent ⟨ L 150 rest ⟩ $ \langle L_{\mathrm{150}}^{\mathrm{rest}}\rangle $ measurements.

To assess the uncertainty on ⟨ L 150 rest ⟩ $ \langle L_{\mathrm{150}}^{\mathrm{rest}}\rangle $, we used a bootstrap analysis. For each stellar mass bin, we made 100 realizations of the stacking analysis, using for each realization a different sample, drawn from the original one, with the same number of sources but allowing for replacement (i.e., a galaxy can be picked several times). In what follows, we use the 50th percentile as the best estimate and the 16th and 84th percentiles of this bootstrap analysis as errors on ⟨ L 150 rest ⟩ $ \langle L_{\mathrm{150}}^{\mathrm{rest}}\rangle $.

These ⟨ L 150 rest ⟩ $ \langle L_{\mathrm{150}}^{\mathrm{rest}}\rangle $ measurements were converted into mean IR luminosities, ⟨LIR [8 − 1000 μm]⟩, following Khusanova et al. (2021),

⟨ L IR ⟩ L ⊙ = ( 3.80 ± 0.63 ) × 10 − 21 × ⟨ L 150 rest ⟩ erg s − 1 H z − 1 . $$ \begin{aligned} \frac{\langle L_{\rm IR}\rangle }{L_\odot } = (3.80\pm 0.63)\times 10^{-21}\times \frac{\langle L_{\rm 150}^\mathrm{rest}\rangle }{\mathrm{erg\,s}^{-1}\,Hz^{-1}}. \end{aligned} $$(1)

This conversion factor was calculated in Khusanova et al. (2021) by averaging the LIR-to- L 150 rest $ L_{\mathrm{150}}^{\mathrm{rest}} $ conversion factors of all SED templates providing a reasonable fit (i.e., χ reduced 2 <1.5 $ \chi^2_{\rm reduced} < 1.5 $) to the Herschel stacks of 4 < z < 5 SFGs performed in Béthermin et al. (2020). It corresponds approximately to that of a modified blackbody function with a dust temperature of 41 K (Khusanova et al. 2021). Uncertainties on the LIR-to- L 150 rest $ L_{\mathrm{150}}^{\mathrm{rest}} $ conversion factor are mitigated with respect to the assumed dust temperature, since at this rest-frame wavelength we are probing close to the dust SED peak. Errors in these LIR-to- L 150 rest $ L_{\mathrm{150}}^{\mathrm{rest}} $ conversion factors, as evaluated in Khusanova et al. (2021), are naturally propagated to all the measurements presented hereafter. We note that the SED template of Béthermin et al. (2017) used here to rescale each individual visibility amplitudes to a common rest-frame frequency is one of the templates providing the best χ reduced 2 $ \chi^2_{\rm reduced} $ from comparison to the Herschel stacks of Béthermin et al. (2020). Also, it yields a LIR-to- L 150 rest $ L_{\mathrm{150}}^{\mathrm{rest}} $ conversion factor that is within 5% to that of Eq. (1). Our rescaling to a common rest-frame frequency and our LIR-to- L 150 rest $ L_{\mathrm{150}}^{\mathrm{rest}} $ conversion are thus fully consistent.

Finally, these mean IR luminosities were converted into dust-attenuated SFRs following Madau & Dickinson (2014) for a Chabrier (2003) IMF,

⟨ SFR IR ⟩ M ⊙ yr − 1 = 1 × 10 − 10 × ⟨ L IR ⟩ L ⊙ . $$ \begin{aligned} \frac{\langle \mathrm{SFR}_{\rm IR}\rangle }{M_\odot \,\mathrm{yr^{-1}}} = 1\times 10^{-10}\times \frac{\langle L_{\rm IR}\rangle }{L_\odot }. \end{aligned} $$(2)

The results of this uv-domain stacking analysis (i.e., ⟨ L 150 rest ⟩ $ \langle L_{\mathrm{150}}^{\mathrm{rest}}\rangle $, ⟨LIR⟩, and ⟨SFRIR⟩) applied to our mass-complete sample of 4 < z < 5 SFGs divided into three stellar mass bins (i.e., 109.5 < M*/M⊙ < 1010, 1010 < M*/M⊙ < 1010.5, and 1010.5 < M*/M⊙ < 1011.5), are given in Table 1 and shown in Fig. 1. This particular choice of stellar mass bin was made to maximize the signal-to-noise ratio of the stacks while providing a good sampling of the stellar mass range probed by our sample. Figure 1 shows that indeed in each of our stellar mass bins we obtain robust detection in terms of signal-to-noise ratio (S/N). Naturally, errors on ⟨ L 150 rest ⟩ $ \langle L_{\mathrm{150}}^{\mathrm{rest}} \rangle $ from our bootstrap analysis are much larger, as they take into account not only the photometric noise in the stacks but also the intrinsic dispersion in the luminosities of the stacked population.

thumbnail Fig. 1.

Cutouts of the stacked images. These cutouts correspond to our 109.5 < M*/M⊙ < 1010 (left), 1010 < M*/M⊙ < 1010.5 (center), and 1010.5 < M*/M⊙ < 1011.5 (right) stellar mass bins, including (top) or excluding (bottom) the ALMA primary targets. The number of stacked galaxies is given in each cutout. The red dashed ellipse shows the ALMA synthesized beam. Color scales are similar for a given stellar mass bin. Each cutout has a size of 4″ × 4″.

We note that the dust-attenuated properties given in Table 1 are weighted values, as a galaxy with deeper ALMA coverage has a higher weight in our uv-domain stacking analysis. The weight of each galaxy is simply given by the number of visibilities to which it corresponds in our stacks (ALMA achieved greater sensitivity by increasing the number of observed visibilities, rather than by increasing the integration time per visibility). Therefore, in what follows, when we compare these dust-attenuated properties to other physical quantities (e.g., redshift, stellar mass, UV luminosity), these are also weighted using the number of visibilities of each galaxy.

4.2. Measuring ⟨LUV⟩ and ⟨βUV⟩ using the COSMOS2020 catalog

In our study, the UV luminosity refers to monochromatic luminosity at 1600 Å (i.e., LUV ≡ ν1600Lν1600), not corrected for dust attenuation. For each of our galaxies, LUV was taken from the COSMOS2020 catalog and corresponds to the UV luminosity of their best-fit SED found with LePhare assuming the Bruzual & Charlot (2003) model with a delayed star formation history (for more details see Weaver et al. 2022). The robustness of these LUV measurements is ensured by the fact that all our galaxies have at least one (but more for most; see below) photometric detection (i.e., S/N > 3) in the rest-frame UV (i.e., λrest < 3300 Å) and that they have accurate SED fits with LePhare (i.e., χ reduced 2 <10 $ \chi^2_{\rm reduced} < 10 $; see Sect. 3). From these LUV, we then measured ⟨LUV⟩, accounting for the weight of each galaxy in our uv-domain stacking analysis (see Sect. 4.1). We assumed that the uncertainties on ⟨LUV⟩ were dominated by the scatter within the population rather than by fitting uncertainties. Therefore, the uncertainties on ⟨LUV⟩ were estimated by bootstrapping over the distributions of LUV in each stellar mass bin. These mean UV luminosities were then converted into unattenuated SFRs following Madau & Dickinson (2014) for a Chabrier (2003) IMF,

⟨ SFR UV ⟩ M ⊙ yr − 1 = 1.5 × 10 − 10 × ⟨ L UV ⟩ L ⊙ . $$ \begin{aligned} \frac{\langle \mathrm{SFR}_{\rm UV}\rangle }{M_\odot \,\mathrm{yr^{-1}}} = 1.5\times 10^{-10}\times \frac{\langle L_{\rm UV}\rangle }{L_\odot }. \end{aligned} $$(3)

Finally, from the best-fit SED of each of our galaxies in COSMOS2020 catalog, we also inferred their UV spectral slope (i.e., βUV with fλ ∝ λβUV). To this end, we used their rest-frame luminosities at 1600 Å and 2300 Å (e.g., Calzetti et al. 1994). From these βUV, we then measured ⟨βUV⟩ and associated uncertainties using a bootstrap analysis and accounting for the weights of each galaxy in our uv-domain stacking analysis (see Sect. 4.1).

The mean UV luminosities ⟨LUV⟩, UV spectral slopes ⟨βUV⟩ and unattenuated SFRs ⟨SFRUV⟩ of the SFGs in our three stellar mass bins are given in Table 1. In addition, this table also provides their mean stellar mass and redshift, as inferred using the COSMOS2020 catalog, and accounting for the weights of each galaxy in our uv-domain stacking analysis (see Sect. 4.1).

To test the robustness of our UV luminosity and UV spectral slope measurements and especially to ensure that they are not biased because they are based on SED fits, we turned to the observed photometry of our galaxies, after applying to the COSMOS2020 magnitudes the appropriate zero-point offsets and foreground extinction corrections (see Laigle et al. 2016; Weaver et al. 2022). To measure the UV spectral slope (i.e., β UV phot . $ \beta_{\mathrm{UV}}^{\mathrm{phot.}} $), we required for each galaxy at least two detections (i.e., S/N > 3) in the rest-frame UV (i.e., λrest < 3300 Å) separated by 500 Å. For our two lower stellar mass bins, this requirement is met by 100% of our galaxies, while in our highest stellar mass bin it is met by only 80% of our galaxies, probably due to higher dust attenuation. We then calculated the UV luminosity (i.e., L UV phot . $ L_{\mathrm{UV}}^{\mathrm{phot.}} $) of these galaxies using the photometry closest to the rest-frame 1600 Å and applying a k-correction based on the previously measured UV spectral slope. For all galaxies for which this analysis was possible (i.e., ∼95%), we find very good agreements between LUV and L UV phot . $ L_{\mathrm{UV}}^{\mathrm{phot.}} $. The median and 16th and 84th percentiles of the L UV phot . $ L_{\mathrm{UV}}^{\mathrm{phot.}} $-to-LUV ratio is 1 . 03 − 0.12 + 0.11 $ 1.03_{-0.12}^{+0.11} $, with no significant dependence on stellar mass. Our fiducial UV luminosities LUV can therefore be considered accurate to within ∼10% and only slightly overestimated by ∼3%. As these systematics and uncertainties are much smaller than those introduced by the dispersion in LUV within the galaxy population of each of our stellar mass bins, we decided, for the sake of completeness, to use LUV instead of L UV phot . $ L_{\mathrm{UV}}^{\mathrm{phot.}} $ in the rest of our analysis. Similarly, we find good, although poorer, agreements between β UV phot . $ \beta_{\mathrm{UV}}^{\mathrm{phot.}} $ and βUV, with a median β UV phot . $ \beta_{\mathrm{UV}}^{\mathrm{phot.}} $-to-βUV ratio and 16th and 84th percentiles of 0 . 74 − 0.51 + 0.35 $ 0.74_{-0.51}^{+0.35} $. Despite these slight differences, for the sake of completeness and consistency, we decided to use βUV instead of β UV phot . $ \beta_{\mathrm{UV}}^{\mathrm{phot.}} $ in the rest of our analysis. However, in Sect. 5.1, we illustrate the impact of using β UV phot . $ \beta_{\mathrm{UV}}^{\mathrm{phot.}} $ instead of βUV on the observed IRX–βUV relation.

5. Results

Using our uv-domain stacking analysis, we are able to measure for the first time the mean dust attenuation properties (i.e., ⟨LIR⟩ and ⟨SFRIR⟩) of a mass-complete sample of SFGs at 4 < z < 5 divided in three stellar mass bins, that is, 109.5 < M*/M⊙ < 1010, 1010 < M*/M⊙ < 1010.5, and 1010.5 < M*/M⊙ < 1011.5 (see Table 1). Combining this information with the emission of these galaxies in the rest-frame UV (i.e., ⟨LUV⟩, ⟨βUV⟩, and ⟨SFRUV⟩), we now study the IRX–βUV, IRX–M*, and SFR–M* relations at z ∼ 4.5 and infer the total, dust-attenuated, and unattenuated cosmic SFRD at this epoch.

5.1. The IRX–βUV and IRX–M* relations

The IR excess (IRX≡LIR/LUV) and its relation with βUV and M* is often used to account for the dust-attenuated fraction of the SFR seen in the UV. Indeed, in the absence of IR photometry, the wealth of optical-to-near-IR data available for high-redshift galaxies allows us in most cases to accurately measure LUV, βUV, and/or M*, which, combined with the (locally) calibrated IRX–βUV and IRX–M* relations, can be used to estimate LIR (and consequently SFRIR).

The IRX–βUV relation is the most fundamental of these two relations. Indeed, assuming a universal, dust unattenuated UV spectral slope β0 for SFGs, the frequency dependence of the attenuation curve directly determines the position of these galaxies in the IRX–βUV plane (see, e.g., Meurer et al. 1999; Fudamoto et al. 2020a),

log 10 IRX ≡ log 10 ( L IR L UV ) = 1.5 × ( 10 0.4 × d A 1600 d β UV × ( β UV − β 0 ) − 1 ) , $$ \begin{aligned} \mathrm{log_{10}\,IRX} \equiv \mathrm{log_{10}\,}\Biggl (\frac{L_{\rm IR}}{L_{\rm UV}}\Biggr ) = 1.5 \times \biggl ( 10^{0.4\times \frac{\mathrm{d}A_{1600}}{\mathrm{d}\beta _{\rm UV}} \times (\beta _{\rm UV} - \beta _0)} -1 \biggr ), \end{aligned} $$(4)

where dA1600/dβUV depends on the frequency-dependence of the attenuation curve (and thus mostly on the dust grain size distribution and composition), and where the factor 1.5 is the bolometric correction of LUV and can be obtained by comparing Eq. (3) and Eq. (2), as LIR is approaching a bolometric luminosity (i.e., for dust LIR [8 − 1000 μm] ∼ Lbol [0 − ∞]).

The distribution in the IRX–βUV plane of the mean ⟨IRX⟩ and ⟨βUV⟩ measurements for our mass-complete sample of SFGs at 4 < z < 5 divided in three stellar mass bins is shown in the left panel of Fig. 2. This distribution is overall very consistent with the local reference for starbursts from Meurer et al. (1999, i.e., with β0 = −2.22 and dA1600/dβUV = 1.99; see also Calzetti et al. 2000) and its updated version from Overzier et al. (2011, i.e., with β0 = −1.96 and dA1600/dβUV = 1.96). On the contrary, our measurements disfavor a flatter IRX–βUV relation, such as that observed in the Small Magellanic Cloud (SMC; here parametrized with β0 = −2.22 and dA1600/dβUV = 1.10, following Fudamoto et al. 2020a). Considering ⟨ L UV phot . ⟩ $ \langle L_{\mathrm{UV}}^{\mathrm{phot.}}\rangle $ and ⟨ β UV phot . ⟩ $ \langle\beta_{\mathrm{UV}}^{\mathrm{phot.}}\rangle $ instead of ⟨LUV⟩ and ⟨βUV⟩, and excluding instead of including the ALMA primary targets in our stacks, does not qualitatively change our results. These various measurements are indeed consistent with each other, except for our lowest stellar mass bin where the photometric estimates slightly shift toward redder UV spectral slopes. Due to these uncertainties, we are unfortunately unable to impose strict constraints on the value of β0.

thumbnail Fig. 2.

Dust attenution properties of SFGs at z ∼ 4.5. Left: IRX–βUV relation derived at 4 < z < 5 by stacking all SFGs with a A3COSMOS coverage (dark-green hexagons) or by stacking all SFGs except the ALMA primary targets (lime-green circles). Hexagons and circles with faded colors correspond to our ⟨ L UV phot . ⟩ $ \langle L_{\mathrm{UV}}^{\mathrm{phot.}}\rangle $ and ⟨ β UV phot . ⟩ $ \langle\beta_{\mathrm{UV}}^{\mathrm{phot.}}\rangle $ estimates, i.e., inferred from a COSMOS2020 photometry as opposed to their best-fit SED (see Sect. 4.2). For each stellar mass bin, the fraction of galaxies for which we were able to make such photometric measurement (see text for details) is given in the bottom part of the panel. In parenthesis, we provide the same number but after excluding the ALMA primary targets. The gray shaded squares present measurements for UV-selected SFGs at z ∼ 4 − 8 from Bowler et al. (2024). The black solid, dash-dotted, dashed, and dotted lines show the IRX–βUV relation observed in local starbursts by Meurer et al. (1999) and by Calzetti et al. (2000), an updated version for local galaxies (Overzier et al. 2011), and an SMC-like dust attenuation relation, respectively. The orange dash-dot-dotted line shows the relation observed in a UV-selected sample of SFGs at z ∼ 4.5 (Fudamoto et al. 2020a). Right: IRX–M* relation at 4 < z < 5. Symbols are the same as in the left panel. Dark-green and lime-green solid lines are linear fits to these data points, while the associated shaded regions show the 1σ uncertainties on these fits (i.e., the 16th to 84th ranges). The pink dotted line shows the redshift-independent IRX–M* relation derived by Pannella et al. (2015) for z ≲ 4 SFGs. The black dash-dotted line shows the IRX–M* relation at z ∼ 4 inferred by Koprowski et al. (2018). The blue dashed line shows the IRX–M* relation at z ∼ 3 and z ∼ 4.5 found in Fudamoto et al. (2020a,b).

Overall, our results suggest that the grain size distribution and composition of the dust in z ∼ 4.5 SFGs are very similar to those of local starbursts, extending to a higher redshift the results for z ∼ 3 SFGs obtained by Fudamoto et al. (2020b). The inconsistency with the rather flat IRX–βUV relation found at z ∼ 4.5 in Fudamoto et al. (2020a) could lie on the fact that this latter study is based on a UV-selected sample of SFGs that missed a fraction of highly dust-attenuated galaxies. However, more recently, by re-analysing the ALPINE measurements used in Fudamoto et al. (2020a) and combining them with z ∼ 7 measurements from REBEL, Bowler et al. (2024) find that using a consistent methodology, UV-selected SFGs at z ∼ 4 − 8 appear to follow the same IRX–βUV as local starbursts, in perfect agreement with our findings.

The IRX–M* is less fundamental than the IRX–βUV relation, as it involves two quantities that are only indirectly related, the stellar mass being a proxy for the integrated past star formation and hence dust production. Therefore, unlike the IRX–βUV relation that mainly depends on the grain size distribution and composition of the dust, the IRX–M* relation also depends on its mass per unit of SFR and its geometry relative to the newly-formed stars. Strikingly, the IRX–M* relation is found to evolve very weakly from z ∼ 0 to z ∼ 3 (Pannella et al. 2015), despite the significant evolution in this redshift range in the star formation efficiency (SFR/Mgas; e.g., Liu et al. 2019b; Wang et al. 2022), metallicity (e.g., Bellstedt et al. 2021), and size of SFGs (e.g., van der Wel et al. 2014) at a given stellar mass (see also discussion in Shapley et al. 2022).

The distribution in the IRX–M* plane of the ⟨IRX⟩ and ⟨M*⟩ measurements for our mass-complete sample of SFGs at 4 < z < 5 is shown on the right panel of Fig. 2. Overall, this distribution (including or excluding the ALMA primary targets in our stacking analysis) is consistent with a steepening of the IRX–M* relation, our measurements being only marginally consistent with the z = 1 − 4 relation of Pannella et al. (2015) at low mass. These findings reinforces those of Fudamoto et al. (2020b), who also find a Meurer-like IRX–βUV relation but a steeper IRX–M* relation in a mass-complete sample of SFGs at z ∼ 3. These results are also very consistent with the z ∼ 4 − 8 measurements of Bowler et al. (2024), obtained by combining the UV-selected sample of SFGs from the ALPINE and REBELS surveys. Therefore, although the dust in z ∼ 4.5 SFGs appears to have the same composition as in local starbursts, its mass and geometry result in lower attenuation in ≲1010 M⊙ SFGs at z ∼ 4.5 than at z ≲ 3 (e.g., −0.5 dex in IRX at M* ∼ 109.75 M⊙). We note that at these low stellar masses, Shapley et al. (2023) do not find significantly lower dust attenuation as measured by the JWST using the Balmer decrement of SFGs, perhaps suggesting a difference between continuum and nebular reddening.

A linear fit to our data points yields the following IRX–M*,

log 10 IRX = 1 . 1 − 0.4 + 0.4 × log 10 ( M ∗ 10 10 M ⊙ ) + 0 . 67 − 0.14 + 0.13 , $$ \begin{aligned} \mathrm{log_{10}\,IRX} = 1.1^{+0.4}_{-0.4}\times \mathrm{log_{10}}\,\Biggl (\frac{M_*}{10^{10}M_{\odot }}\Biggr )+0.67^{+0.13}_{-0.14}\,, \end{aligned} $$(5)

when including the ALMA primary targets in our stacks, or,

log 10 IRX = 1 . 8 − 0.7 + 1.0 × log 10 ( M ∗ 10 10 M ⊙ ) + 0 . 45 − 0.42 + 0.30 , $$ \begin{aligned} \mathrm{log_{10}\,IRX} = 1.8^{+1.0}_{-0.7}\times \mathrm{log_{10}}\,\Biggl (\frac{M_*}{10^{10}\,M_{\odot }}\Biggr )+0.45^{+0.30}_{-0.42}, \end{aligned} $$(6)

when excluding the ALMA primary targets. These two fits are very consistent overall, but the latter gives a less constrained and slightly steeper IRX–M* relation, due to the fact that by excluding the ALMA primary targets, our lowest stellar mass bin becomes an upper limit. In the rest of our analysis, we consider the whole set of IRX–M* relations found in these two fits (shaded regions in Fig. 2), and use as the fiducial IRX–M* relation at z ∼ 4.5 the linear combination in the log-space of these two relations,

log 10 IRX = 1.45 × log 10 ( M ∗ 10 10 M ⊙ ) + 0.56 . $$ \begin{aligned} \mathrm{log_{10}\,IRX} = 1.45\times \mathrm{log_{10}}\,\Biggl (\frac{M_*}{10^{10}\,M_{\odot }}\Biggr )+0.56. \end{aligned} $$(7)

5.2. The SFR–M* relation

In the SFR–M* plane, SFGs lie on a tight locus, known as the main sequence (MS) of SFGs and characterized by an almost linear increase in SFR with stellar mass (i.e., SFR ∝ M*; e.g., Noeske et al. 2007; Elbaz et al. 2007; Schreiber et al. 2015; Leslie et al. 2020; Daddi et al. 2022; Popesso et al. 2023; Goovaerts et al. 2024; Koprowski et al. 2024, and references therein). The small scatter of the MS suggests that secular evolution is the dominant mode of growth in SFGs, where gas inflow, outflow, and star formation are in equilibrium (e.g., Bouché et al. 2010; Davé et al. 2012; Lilly et al. 2013; Peng & Maiolino 2014; Rathaus & Sternberg 2016; Magnelli et al. 2020). Over the last decade, the MS has been a cornerstone of galaxy evolution studies, enabling the selection of “normal” SFGs at a given (z, M*) and the study of galaxy properties as a function of their position relative to the MS. While well established up to z ∼ 3, the MS (i.e., its normalization, slope, and dispersion) is not well constrained at higher redshifts, even though it is commonly used to select normal SFGs at these epochs. Indeed, current mass-complete samples of SFGs, including those from JWST, lack systematic detections in the IR and thus rely on highly uncertain dust-attenuated SFRs, while even state-of-the-art stacking analysis with Herschel are limited to the most massive end of the SFG population (i.e., ≳1011 M⊙; e.g., Schreiber et al. 2015; Koprowski et al. 2024). These latter stacking estimates are also hampered by clustering signal, which has to be removed beforehand due to the coarse angular resolution of Herschel-SPIRE (> 14″). Our ALMA stacks are not affected by such clustering bias thanks to their high angular resolution.

The distribution in the SFR–M* plane of the ⟨SFR⟩ and ⟨M*⟩ measurements for our mass-complete sample of SFGs at 4 < z < 5 is shown in Fig. 3. The measurements including the ALMA primary targets lie almost perfectly on the MS of Schreiber et al. (2015, i.e., MSAll,  Fit ≡  MSS15). This agreement is somewhat surprising as at these redshifts, Schreiber et al. (2015) are only able to probe the most massive SFGs (i.e., ≳1011 M⊙) with their Herschel stacks. In contrast to the previous section, the exclusion of the ALMA primary targets in our stacks has a significant impact on the inferred MS relation. Indeed, although these two measurements formally agree within their uncertainties, those obtained by excluding the ALMA primary targets give systematically lower SFRs; this tendency increasing at higher stellar mass where the weight of the ALMA primary targets in our stacks increases (see numbers at the bottom of Fig. 3). By fitting these latter measurements with the MS of Schreiber et al. (2015) leaving its normalization as a free parameter, we find MSOffcenter,  Fit ≡ MSS15 − 0.3 dex. Naturally, this difference stems from the fact that a fraction of the principal investigators of ALMA has targeted high-redshift starbursts (e.g., submillimeter-selected galaxies) or at least galaxies located on the upper envelop of the MS dispersion. Consequently, the SFR–M* relation obtained here by including the ALMA primary targets is probably slightly biased toward high values, while the SFR–M* relation obtained by excluding the ALMA primary targets is probably biased toward low values. In the rest of our analysis, we take into account this uncertainty on the normalization of the MS (i.e., ±0.15 dex) and use as our fiducial MS the linear combination in the log-space of these two relations, that is, MSfid. ≡ MSS15 − 0.15 dex.

thumbnail Fig. 3.

SFR–M* relation at 4 < z < 5. Dark-gray hexagons show the dust-attenuated SFRs (i.e., SFRIR) inferred by stacking all SFGs with a A3COSMOS coverage, while dark-green hexagons show their total SFRs, i.e., adding the contribution of their unattenuated SFRs measured in the UV (i.e., SFRUV; symbols are slightly offset along the stellar mass axis for clarity). Light-gray and lime-green circles display the same quantities but for our stacking analysis which excludes the ALMA primary targets (symbols are slightly offset along the stellar mass axis for clarity). The dark-green solid and red dotted lines represent the MS and its ±0.3 dex dispersion from Schreiber et al. (2015), whose normalization happens to be perfectly consistent with our dark-green hexagons, i.e., MSAll,  Fit ≡ MSS15. The lime-green line is obtained by renormalizing the MS of Schreiber et al. (2015) to fit our lime-green circles, i.e., MSOffcenter,  Fit ≡ MSS15 − 0.3 dex. For each stellar mass bin, we give the fraction of ALMA primary targets and in parenthesis their weights in our stacks. Finally, pink squares and diamonds show the total and dust-attenuated SFRs found in Khusanova et al. (2021) by stacking in three stellar mass bins, ALMA observations of all (detected and undetected) z ∼ 4.5 SFGs from the ALMA ALPINE survey.

5.3. The cosmic SFRD

Combining the results of the previous sections with the stellar mass function (SMF) of SFGs at z ∼ 4.5, we can calculate the dust-attenuated and unattenuated cosmic SFRD at this epoch and the evolution of these quantities with stellar mass. In the top left panel of Fig. 4, we display the SMF of SFGs used in our analysis and inferred in Weaver et al. (2023). The choice of this particular SMF was relatively straightforward as it is arguably the best estimate available to date and, more importantly, it was measured from the same parent sample as that used here. As shown in Fig. 4 and also noted in Weaver et al. (2023), at the massive end, their data points deviate, however, from a canonical Schechter function (Schechter 1976) and thus from their best fit. After an in-depth analysis of this population of massive SFGs, Weaver et al. (2023) conclude that their luminosity, redshift, mass, and number density are not compatible with AGN-dominated galaxies, but rather with the population of optically-dark galaxies recently discovered thanks to Spitzer-IRAC and ALMA observations (e.g., Wang et al. 2019; Xiao et al. 2023). They argue that the Ks ∼ 24.7 AB depth of the UltraVISTA DR4 used to build the izYJHKs detection image on which the COSMOS2020 catalog was selected could indeed be sufficient to reach out into this population missed by previous optical or even H-band selections. To account analytically for this massive population, we fit simultaneously all the SMF data points (3.5 < z < 4.5 and 4.5 < z < 5.5) of Weaver et al. (2023) with the combination of a Schechter function dominating at low stellar masses and a power-law dominating at high stellar masses. The result of this fit and its associated uncertainties is represented by the solid line and the gray shaded region in the top-left panel of Fig. 4. In the rest of the analysis, we use this fit as our fiducial SMF, but also discuss the impact of using instead the Schechter fit of Weaver et al. (2023) on our calculations.

thumbnail Fig. 4.

Evolution with stellar mass of the cosmic SFRD at 4 < z < 5. Top left: fiducial SMF used in our calculations is represented by the gray line and shaded region. This is the combination of a Schechter and power-law function fit to the SMF of Weaver et al. (2023, gray and black circles). A Schechter function fit to these data points is shown by the gray dotted line. Bottom left: cosmic SFRD (thick black line) and dust-attenuated SFRD (red line) above a given stellar mass, as inferred using the fiducial SMF. The dotted gray and red lines present the same quantities but using a simple Schechter SMF in our calculations. Top right: fraction of the M* > 108 M⊙ SFRD (gray histogram) and dust-attenuated SFRD (red histogram) that is attributed to a given stellar mass bin. Bottom right: fraction of the cosmic SFRD above a given stellar mass that is dust-attenuated (red region) and unattenuated (blue region). The dotted gray line present the same quantity but using a simple Schechter SMF in our calculations. In all panels, the vertical line represents the 0.03 × M⋆ lower integration limits commonly used to infer the “total” cosmic SFRD (e.g., Madau & Dickinson 2014). All quantities in the bottom-left, top-right, and bottom-right panels were calculated using our fiducial SMF, IRX–M*, MS relations. The propagation of the 1-σ uncertainties on the SMF and the IRX–M* relation are represented by shaded regions. The uncertainties on the MS would move the black and red lines of the bottom-left panel up and down by 0.15 dex, but would have no effect on the other quantities displayed in this figure.

First, we calculated the cosmic SFRD above a given stellar mass (i.e., SFRD( > M*)) and the contribution of each stellar mass bin to SFRD(> 108 M⊙) by multiplying and then integrating our fiducial SMF with our fiducial MS (bottom-left and top-right panels of Fig. 4, respectively). Then, using our fiducial IRX–M* relation, we calculated the cosmic dust-attenuated SFRD above a given stellar mass (i.e., SFRDIR(> M*)) and the contribution of each stellar mass bin to SFRDIR(> 108 M⊙; bottom-left and top-right panels of Fig. 4, respectively). Finally, we calculated the fraction of the cosmic SFRD above a given stellar mass that is dust attenuated (fattenuated; bottom-right panel of Fig. 4). In these calculations, the complex propagations of the uncertainties on the SMF and the IRX–M* relation were inferred using Monte Carlo realizations and are represented by shaded regions. In contrast, using MSAll,  Fit or MSOffcenter,  Fit instead of our fiducial MS simply move the black and red lines of the bottom-left panel of Fig. 4 (i.e., SFRD( > M*) and SFRDIR(> M*)) up or down by 0.15 dex, respectively, but would have no effect on the other quantities displayed in this figure. Finally, we note that our calculations does not account for the contribution of starbursts (i.e., SFGs located > 0.6 dex above the MS) to the cosmic SFRD, since by construction we have assumed that at a given mass, all SFGs are MS galaxies. However, this rare population of galaxies is known to have a modest ∼10% contribution to the SFRD (Rodighiero et al. 2011; Schreiber et al. 2015, i.e., at most +0.05 dex to our estimates) and part of this contribution is actually taken into account by our MSAll,  Fit relation that is supposedly biased towards this population (see Sect. 5.2).

Three main conclusions can be drawn from Fig. 4. Firstly, the combination of the MS and SMF slopes implies that the cosmic SFRD converges at low stellar masses: the number of SFGs at low stellar mass does not increase sufficiently to counterbalance the decrease in their SFRs. This is highlighted by the fact that, based on our extrapolations, SFGs at 108 < M*/M⊙ < 109 already contribute less than 15% of the SFRD(> 108 M⊙). The steep slope of the IRX–M* relation makes this result even more pronounced for SFRDIR(> M*), with SFGs at 108 < M*/M⊙ < 109 contributing less than 5% of the SFRDIR(> 108 M⊙). This implies that the study of SFGs much less massive than ≲109 M⊙ is almost irrelevant for our understanding of the cosmic SFRD at z ∼ 4.5. Secondly, the cosmic SFRD is mainly dominated by SFGs with a stellar mass of 109.5 − 10.5 M⊙: galaxies that are more massive and therefore more star-forming, are too rare to contribute significantly to the SFRD; while if less massive galaxies are numerous, they do not form enough stars to contribute significantly to the SFRD. This characteristic stellar mass of the SFG population is consistent with that observed at lower redshifts (∼1010.3 M⊙; e.g., Karim et al. 2011; Leslie et al. 2020). Thirdly, the fraction of the cosmic SFRD that is attenuated by dust remains significant even at such an early cosmic epoch, with fattenuated ≡ SFRDIR(> M*)/SFRD(> M*) converging to a value of ∼60 ± 10% for M* = 108 M⊙.

Using a Schechter SMF instead of our fiducial SMF would not qualitatively change any of these results (see dotted lines in Fig. 4). Indeed, while a massive population is responsible for a large deviation of our fiducial SMF from a Schechter function, its number density remains too low to have a significant impact on the global star formation activity in the Universe. In fact, excluding this population from our calculation shifts the characteristic mass of the SFG population to lower values by only ∼0.3 dex and decrease the cosmic SFRD by ∼0.1 dex, this latter value being in perfect agreement with the finding of, for example, Wang et al. (2019) for optically-dark galaxies. As these variations are negligible compared to those introduced by uncertainties on the MS and IRX–M* relations, in what follows we restrict our calculation to our fiducial SMF (i.e., accounting for this massive population of SFGs).

Finally, as detailed in Appendix B, we also tested the robustness of our calculations against the most recent JWST-derived SMF from Weibel et al. (2024). Because the SMF of Weibel et al. (2024) differs from that of Weaver et al. (2023) only at M* ≲ 109 M⊙ where it has a slightly steeper slope (−1.79 versus −1.56), our results remain qualitatively unchanged using this alternative SMF: the convergence of the cosmic SFRD occurs but at lower stellar mass (107 M⊙ vs. 108 M⊙); at the convergence, the fraction that is dust attenuated is lower but still significant (43% vs. 60%); the contribution of very low stellar mass galaxies increases but remains sub-dominant (35% vs. 15% for M* < 109 M⊙); while the “total” cosmic SFRD (i.e. at 108.9 M⊙; see Sect. 5.4) and the fraction of it that is attenuated by dust remains mostly unchanged from those deduced using the COSMOS2020-derived SMF of Weaver et al. (2023). Because of these very good agreements over the stellar mass range for which we are able to measure the dust-attenuated properties of z ∼ 4.5 SFGs, and because using the COSMOS2020-derived SMF ensure consistency with our stack samples, we have decided to retain the COSMOS2020-derived SMF from Weaver et al. (2023) as our fiducial SMF. We defer to future work the possibility of extending the constraints on the dust-attenuation of SFGs to lower stellar masses using JWST-based catalogs.

5.4. The “total” cosmic SFRD

In the literature, the IR and UV luminosity functions are commonly integrated down to 0.03 × L⋆ to infer the “total” cosmic SFRD (e.g., Madau & Dickinson 2014), where L⋆ is the characteristic luminosity of the Schechter function. The slope of our fiducial MS being equal to one (Schreiber et al. 2015), integrating the SMF down to 0.03 × M⋆ (i.e., 108.9 M⊙; vertical lines in Fig. 4) should provide consistent results. In Fig. 5, we compare the total cosmic SFRD and dust-attenuated SFRD at z ∼ 4.5 calculated here by integrating the SMF down to this limit with the literature results. The symbols correspond to the calculation using our fiducial SMF, MS and IRX–M* relations, while the range on the cosmic SFRD corresponds to the ±0.15 dex uncertainties on the MS relation (i.e., excluding or not the ALMA primary targets) and the range on the cosmic dust-attenuated SFRD takes into account the uncertainties on both the MS relation (i.e., excluding or not the ALMA primary targets) and the IRX–M* relations. Hereafter these total cosmic SFRD and total cosmic dust-attenuated SFRD are referred as SFRDUV + IR and SFRDIR, respectively. We find SFRDUV + IR  = 0.012 − 0.028 M⊙ yr−1 Mpc−3 and SFRDIR = 0.007 − 0.021 M⊙ yr−1 Mpc−3.

thumbnail Fig. 5.

Redshift evolution of the “total” cosmic SFRD. The cosmic SFRD (UV+IR) and dust-attenuated SFRD (IR) derived in this work are shown by black and orange circles, respectively (offset slightly along the redshift axis for clarity). To be consistent with Madau & Dickinson (2014), these measurements were obtained by integrating the SMF down to 0.03 × M⋆ (i.e., 108.9 M⊙). Orange open symbols are a collection of dust-attenuated SFRDs (i.e., SFRDIR) from the literature (Sanders et al. 2003; Takeuchi et al. 2003; Schiminovich et al. 2005; Dahlen et al. 2007; Reddy & Steidel 2009; Magnelli et al. 2011, 2013; Cucciati et al. 2012; Bouwens et al. 2012; Gruppioni et al. 2013; Finkelstein et al. 2015; Bouwens et al. 2015; Liu et al. 2018; Khusanova et al. 2020, 2021; Zavala et al. 2021; Algera et al. 2023; Traina et al. 2024). Light blue open symbols are a collection of unattenuated SFRDs (i.e., SFRDUV,  uncorr) from the literature, while dark blue open symbols correspond to the same measurements after correction for dust attenuation (i.e., SFRDUV,  corr). The dark-gray line is the redshift evolution of the cosmic SFRD inferred in Madau & Dickinson (2014). Orange dotted and dashed lines show the redshift evolution of the dust-attenuated SFRD inferred in Koprowski et al. (2017) and Maniyar et al. (2018), respectively.

Our estimate of SFRDUV + IR is consistent with Madau & Dickinson (2014), in particular the upper value of our range, which corresponds to MSAll,  Fit ≡ MSS15. Our estimate of SFRDIR is also in agreement with previous measurements based on Herschel (Liu et al. 2018), SCUBA (Maniyar et al. 2018) or ALMA (e.g., Zavala et al. 2021; Traina et al. 2024).

In contrast, our estimate of SFRDIR disagrees with the recent ALMA measurements by Gruppioni et al. (2020) and Khusanova et al. (2021), both of which found large SFRDIR values, about two times and four times larger than the SFRDUV + IR and SFRDIR inferred in Madau & Dickinson (2014), respectively. This disagreement (although our estimates are formally consistent within the uncertainties) could naturally come from our sides, owing to the fact that our method is based on a near-infrared-selected sample (i.e., COSMOS2020) that may miss a population of heavily dust-attenuated SFGs only recoverable using an ALMA-selected sample such as that used in, e.g., Gruppioni et al. (2020). The existence of a population of optically-dark galaxies (e.g., Wang et al. 2019; Xiao et al. 2023) supports this view. However, the impact of this population, while undeniable, remains modest with the most recent estimates ranging from SFRDIR ∼ 0.002 M⊙ yr−1 Mpc−3 to ∼0.012 M⊙ yr−1 Mpc−3 (e.g. Xiao et al. 2023; Barrufet et al. 2023; Williams et al. 2024). In addition, their impact on our measurements would be further mitigated because at the Ks ∼ 24.7 AB depth of the UltraVISTA DR4 used to build the izYJHKs detection image, the COSMOS2020 catalog should detect part of this population (Weaver et al. 2023). Therefore, while the inclusion of optically-dark galaxies could lead to a increase of our SFRDUV + IR and SFRDIR measurements by at most ∼0.012 M⊙ yr−1 Mpc−3 (+0.25 dex compared to our fiducial SFRDIR), such increase is well within the uncertainties associated with our measurements and not enough to explain the disagreements (∼ + 0.7 dex compared to our fiducial SFRDIR value) with Gruppioni et al. (2020) and Khusanova et al. (2021). Future ALMA stacks of JWST-selected SFGs at z ∼ 4.5 will enable us to assess the exact contribution of this population of galaxies to the cosmic SFRD.

Instead, we argue that the large SFRDIR values found in Gruppioni et al. (2020) and Khusanova et al. (2021) could be due to low number statistic and extrapolation of the IRX–M* relation, respectively. On the one hand, using the same approach as Gruppioni et al. (2020) but using a much larger number of ALMA observations drawn from the A3COSMOS database, Traina et al. (2024) measure a SFRDIR that is more consistent with our measurements. This suggests that a larger statistic and the combination of Herschel and ALMA in Traina et al. (2024) have lowered the inferred IR luminosity function and, subsequently, SFRDIR at z ∼ 4.5. On the other hand, for ≥108.35 M⊙ where they have direct constraints on the LIR–M* relation, Khusanova et al. (2021) inferred a SFRDIR of 0.023 M⊙ yr−1 Mpc−3, very consistent with our measurements. This suggests that the disagreement with Khusanova et al. (2021) stems from their extrapolation to lower stellar mass, as their total cosmic SFRDIR is for SFGs with > 106 M⊙. The reason for these very low integration limits and for the increase of SFRDIR from 108.35 M⊙ to 106 M⊙ when our measurements have already converged at these masses (see Fig. 4) is threefold: Khusanova et al. (2021) used a MS relation with a sub-unity slope of ∼0.85, which implies that the 0.03 × L⋆ integration limit is converted to a very low-mass limit of 106 M⊙; Next, Khusanova et al. (2021) used a LIR–M* relation with a sub-unity slope of ∼0.85, which implies an unrealistic mass-independent value for IRX; Finally, all these elements combined imply that as stellar mass decreases, the SFR of their SFGs decreases more slowly than in our case (flatter MS relation) and more importantly the dust-attenuated fraction of these SFRs remains significant. In brief, an flat IRX–M* relation combined with a shallow MS relation result in unrealistic extrapolations to lower stellar masses and an increase of SFRDIR(> M*) from 0.023 M⊙ yr−1 Mpc−3 for M* = 108.35 M⊙ to 0.085 M⊙ yr−1 Mpc−3 for M* = 106 M⊙.

Finally, in Fig. 6, we present the cosmic dust attenuation in magnitudes, AUV ≡ 2.5 × log10(1 + SFRDIR/SFRDUV) (see, e.g., Madau & Dickinson 2014) and compare it to various measurements from the literature. Our estimate that about ∼68% of the total cosmic SFRD is attenuated by dust at M* ≥ 108.9 M⊙ (and thus AUV = 1.14) is in line with previous results from the literature, which range from ∼35% (Khusanova et al. 2020), ∼45% (Zavala et al. 2021), to ∼70% (Finkelstein et al. 2015). The fraction of the total cosmic SFRD that is attenuated by dust remains thus significant at this early cosmic epoch, even though this fraction has clearly decreased from z ∼ 1 − 2, where it peaks, to z ∼ 4.5.

thumbnail Fig. 6.

Cosmic dust attenuation in magnitudes (AUV ≡ 2.5 × log10(1 + SFRDIR/SFRDUV)) as a function of redshift. The cosmic dust attenuation inferred in this work is shown by an orange circle. Blue and orange symbols are a collection of cosmic dust attenuation measurements from the literature.

6. Summary

Using the wealth of multiwavelength observations publicly available over the COSMOS field, we measured for the first time the mean UV and IR emission properties of a mass-complete sample of SFGs at 4 < z < 5 divided in three stellar mass bins, that is, 109.5 < M*/M⊙ < 1010, 1010 < M*/M⊙ < 1010.5, and 1010.5 < M*/M⊙ < 1011.5. We built this mass-complete sample of SFGs and measured their rest-frame UV properties (i.e., ⟨LUV⟩, ⟨βUV⟩, and ⟨SFRUV⟩) using the deepest optical-to-near-IR data publicly available in COSMOS (i.e., the COSMOS2020 catalog; Weaver et al. 2022). Their mean IR properties (i.e., ⟨LIR⟩ and ⟨SFRIR⟩) were then measured by stacking in the uv domain all archival ALMA band 6 and 7 observations publicly available for these galaxies (A3COSMOS; Liu et al. 2019a; Adscheid et al. 2024). With this unique approach, we find the following:

  1. The relation between the IR excess (IRX ≡ LIR/LUV) and the UV spectral slope (βUV) in our mass-complete sample of SFGs at 4 < z < 5 is overall very consistent with the local IRX–βUV relation from Meurer et al. (1999) and Calzetti et al. (2000), and disfavors a flatter IRX–βUV relation, such as that observed in the SMC. This result, consistent with the most recent z ∼ 4 − 8 measurements from Bowler et al. (2024), suggests that the grain size distribution and composition of the dust in SFGs with M* > 109.5 M⊙ and at z ∼ 4.5 are very similar to those of local starbursts.

  2. Our measurements favor a slight steepening of the IRX–M* relation at z ∼ 4.5, when compared to the redshift-independent IRX–M* relation observed at z ∼ 1 − 3 (e.g., Pannella et al. 2015). Thus, although the dust in z ∼ 4.5 SFGs appears to have the same composition than in local starbursts, its mass and geometry result in lower attenuation in ≲1010 M⊙ SFGs at z ∼ 4.5 than at z ≲ 3 (e.g., −0.5 dex in IRX at M* ∼ 109.75 M⊙).

  3. In the SFR–M* plane, our galaxies lie almost perfectly on the MS of Schreiber et al. (2015), while they lie 0.3 dex below this relation when excluding from our stacks the ALMA primary targets (i.e., sources within 3″ from the ALMA phase center). ALMA primary targets are probably slightly biased toward galaxies located on the upper envelope of the MS dispersion. We set our fiducial MS to MSS15 − 0.15 dex.

  4. The combination of the MS and SMF slopes implies that the SFRD( > M*) converges at M* ≲ 109 M⊙, as the number of SFGs at lower stellar mass does not increase sufficiently to counterbalance the decrease in their SFRs. For example, SFGs with 108 < M*/M⊙ < 109 contribute already less than 15% of the SFRD(> 108 M⊙), and less than 5% of the SFRDIR(> 108 M⊙).

  5. The “total” cosmic SFRD inferred here at z ∼ 4.5 for SFGs with M* > 0.03 × M⋆ (M⋆ being the characteristic stellar mass of SFGs at this epoch) is consistent with Madau & Dickinson (2014) and is dominated by SFGs with a stellar mass of 109.5 − 10.5 M⊙. At z ∼ 4.5, SFRDUV + IR = 0.012 − 0.028 M⊙ yr−1 Mpc−3 and SFRDIR  = 0.007 − 0.021 M⊙ yr−1 Mpc−3. The population of optical-dark galaxies potentially missed by our study could add at most 0.012 M⊙ yr−1 Mpc−3 to these values (Xiao et al. 2023; Barrufet et al. 2023; Williams et al. 2024). Future ALMA stacks of JWST-selected SFGs at z ∼ 4.5 will assess the exact contribution of this population of galaxies to the cosmic SFRD.

  6. The fraction of the cosmic SFRD that is attenuated by dust remains significant even at such an early cosmic epoch, with fattenuated ≡ SFRDIR(> M*)/SFRD(> M*) having a value of 68 ± 10% for M* ≥ 108.9 M⊙ (i.e., 0.03 × M⋆) and converging to a value of 60 ± 10% for M* ≥ 108 M⊙.

Acknowledgments

We would like to thank the referee for their comments that have helped to improve our paper. BM, DE, and CGG acknowledges support from CNES. S.A. gratefully acknowledges the Collaborative Research Center 1601 (SFB 1601 sub-project C2) funded by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) – 500700252. E.F.-J.A. acknowledge support from UNAM-PAPIIT project IA102023, and from CONAHCyT Ciencia de Frontera project ID: CF-2023-I-506. MF acknowledges support from NSF grant AST-2009577 and NASA JWST GO program 1727. CG acknowledges the support from grant PRIN MIUR 2017-20173ML3WW_001 “Opening the ALMA window on the cosmic evolution of gas, stars, and supermassive black hole”, and co-funding by the European Union – NextGenerationEU within PRIN 2022 project n.20229YBSAN “Globular clusters in cosmological simulations and in lensed fields: from their birth to the present epoch”. BM acknowledges the following open source software used in the analysis: Astropy (Astropy Collaboration 2022), photutils (Bradley et al. 2022), and NumPy (Harris et al. 2020).

References

  1. Adscheid, S., Magnelli, B., Liu, D., et al. 2024, A&A, 685, A1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  2. Algera, H. S. B., Inami, H., Oesch, P. A., et al. 2023, MNRAS, 518, 6142 [Google Scholar]
  3. Astropy Collaboration (Price-Whelan, A. M., et al.) 2022, ApJ, 935, 167 [NASA ADS] [CrossRef] [Google Scholar]
  4. Bagley, M. B., Finkelstein, S. L., Koekemoer, A. M., et al. 2023, ApJ, 946, L12 [NASA ADS] [CrossRef] [Google Scholar]
  5. Barrufet, L., Oesch, P. A., Weibel, A., et al. 2023, MNRAS, 522, 449 [NASA ADS] [CrossRef] [Google Scholar]
  6. Bellstedt, S., Robotham, A. S. G., Driver, S. P., et al. 2021, MNRAS, 503, 3309 [NASA ADS] [CrossRef] [Google Scholar]
  7. Béthermin, M., Wu, H.-Y., Lagache, G., et al. 2017, A&A, 607, A89 [Google Scholar]
  8. Béthermin, M., Fudamoto, Y., Ginolfi, M., et al. 2020, A&A, 643, A2 [Google Scholar]
  9. Bouché, N., Dekel, A., Genzel, R., et al. 2010, ApJ, 718, 1001 [Google Scholar]
  10. Bouwens, R. J., Illingworth, G. D., Oesch, P. A., et al. 2012, ApJ, 754, 83 [Google Scholar]
  11. Bouwens, R. J., Illingworth, G. D., Oesch, P. A., et al. 2015, ApJ, 803, 34 [Google Scholar]
  12. Bowler, R. A. A., Cullen, F., McLure, R. J., Dunlop, J. S., & Avison, A. 2022, MNRAS, 510, 5088 [NASA ADS] [CrossRef] [Google Scholar]
  13. Bowler, R. A. A., Inami, H., Sommovigo, L., et al. 2024, MNRAS, 527, 5808 [Google Scholar]
  14. Bradley, L., Sipőcz, B., Robitaille, T., et al. 2022, https://zenodo.org/records/6825092 [Google Scholar]
  15. Bruzual, G., & Charlot, S. 2003, MNRAS, 344, 1000 [NASA ADS] [CrossRef] [Google Scholar]
  16. Calzetti, D., Kinney, A. L., & Storchi-Bergmann, T. 1994, ApJ, 429, 582 [Google Scholar]
  17. Calzetti, D., Armus, L., Bohlin, R. C., et al. 2000, ApJ, 533, 682 [NASA ADS] [CrossRef] [Google Scholar]
  18. Chabrier, G. 2003, PASP, 115, 763 [Google Scholar]
  19. Cucciati, O., Tresse, L., Ilbert, O., et al. 2012, A&A, 539, A31 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  20. Daddi, E., Delvecchio, I., Dimauro, P., et al. 2022, A&A, 661, L7 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  21. Dahlen, T., Mobasher, B., Dickinson, M., et al. 2007, ApJ, 654, 172 [Google Scholar]
  22. Davé, R., Finlator, K., & Oppenheimer, B. D. 2012, MNRAS, 421, 98 [Google Scholar]
  23. Eisenstein, D. J., Willott, C., Alberts, S., et al. 2023, ApJS, submitted [arXiv:2306.02465] [Google Scholar]
  24. Elbaz, D., Daddi, E., Le Borgne, D., et al. 2007, A&A, 468, 33 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  25. Finkelstein, S. L., Ryan, R. E., Jr, Papovich, C., et al. 2015, ApJ, 810, 71 [NASA ADS] [CrossRef] [Google Scholar]
  26. Fudamoto, Y., Oesch, P. A., Faisst, A., et al. 2020a, A&A, 643, A4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  27. Fudamoto, Y., Oesch, P. A., Magnelli, B., et al. 2020b, MNRAS, 491, 4724 [Google Scholar]
  28. Goovaerts, I., Pello, R., Burgarella, D., et al. 2024, A&A, 683, A184 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  29. Gruppioni, C., Pozzi, F., Rodighiero, G., et al. 2013, MNRAS, 432, 23 [Google Scholar]
  30. Gruppioni, C., Béthermin, M., Loiacono, F., et al. 2020, A&A, 643, A8 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  31. Harris, C. R., Millman, K. J., Walt, S. J. V. D., et al. 2020, Nature, 585, 357 [NASA ADS] [CrossRef] [Google Scholar]
  32. Ilbert, O., McCracken, H. J., Le Fèvre, O., et al. 2013, A&A, 556, A55 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  33. Karim, A., Schinnerer, E., Martínez-Sansigre, A., et al. 2011, ApJ, 730, 61 [Google Scholar]
  34. Khusanova, Y., Le Fèvre, O., Cassata, P., et al. 2020, A&A, 634, A97 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  35. Khusanova, Y., Bethermin, M., Le Fèvre, O., et al. 2021, A&A, 649, A152 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  36. Koprowski, M. P., Dunlop, J. S., Michałowski, M. J., et al. 2017, MNRAS, 471, 4155 [Google Scholar]
  37. Koprowski, M. P., Coppin, K. E. K., Geach, J. E., et al. 2018, MNRAS, 479, 4355 [NASA ADS] [CrossRef] [Google Scholar]
  38. Koprowski, M. P., Wijesekera, J. V., Dunlop, J. S., et al. 2024, A&A, submitted [arXiv:2403.06575] [Google Scholar]
  39. Laigle, C., McCracken, H. J., Ilbert, O., et al. 2016, ApJS, 224, 24 [Google Scholar]
  40. Leslie, S. K., Schinnerer, E., Liu, D., et al. 2020, ApJ, 899, 58 [Google Scholar]
  41. Lilly, S. J., Carollo, C. M., Pipino, A., Renzini, A., & Peng, Y. 2013, ApJ, 772, 119 [NASA ADS] [CrossRef] [Google Scholar]
  42. Lindroos, L., Knudsen, K. K., Vlemmings, W., Conway, J., & Martí-Vidal, I. 2015, MNRAS, 446, 3502 [NASA ADS] [CrossRef] [Google Scholar]
  43. Liu, D., Daddi, E., Dickinson, M., et al. 2018, ApJ, 853, 172 [Google Scholar]
  44. Liu, D., Lang, P., Magnelli, B., et al. 2019a, ApJS, 244, 40 [Google Scholar]
  45. Liu, D., Schinnerer, E., Groves, B., et al. 2019b, ApJ, 887, 235 [Google Scholar]
  46. Madau, P., & Dickinson, M. 2014, ARA&A, 52, 415 [Google Scholar]
  47. Magnelli, B., Elbaz, D., Chary, R. R., et al. 2011, A&A, 528, A35 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  48. Magnelli, B., Popesso, P., Berta, S., et al. 2013, A&A, 553, A132 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  49. Magnelli, B., Boogaard, L., Decarli, R., et al. 2020, ApJ, 892, 66 [NASA ADS] [CrossRef] [Google Scholar]
  50. Maniyar, A. S., Béthermin, M., & Lagache, G. 2018, A&A, 614, A39 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  51. Meurer, G. R., Heckman, T. M., & Calzetti, D. 1999, ApJ, 521, 64 [Google Scholar]
  52. Noeske, K. G., Weiner, B. J., Faber, S. M., et al. 2007, ApJ, 660, L43 [CrossRef] [Google Scholar]
  53. Overzier, R. A., Heckman, T. M., Wang, J., et al. 2011, ApJ, 726, L7 [NASA ADS] [CrossRef] [Google Scholar]
  54. Pannella, M., Elbaz, D., Daddi, E., et al. 2015, ApJ, 807, 141 [Google Scholar]
  55. Peng, Y.-J., & Maiolino, R. 2014, MNRAS, 443, 3643 [NASA ADS] [CrossRef] [Google Scholar]
  56. Popesso, P., Concas, A., Cresci, G., et al. 2023, MNRAS, 519, 1526 [Google Scholar]
  57. Rathaus, B., & Sternberg, A. 2016, MNRAS, 458, 3168 [NASA ADS] [CrossRef] [Google Scholar]
  58. Reddy, N. A., & Steidel, C. C. 2009, ApJ, 692, 778 [Google Scholar]
  59. Rodighiero, G., Daddi, E., Baronchelli, I., et al. 2011, ApJ, 739, L40 [Google Scholar]
  60. Sanders, D. B., Mazzarella, J. M., Kim, D. C., Surace, J. A., & Soifer, B. T. 2003, AJ, 126, 1607 [Google Scholar]
  61. Schechter, P. 1976, ApJ, 203, 297 [Google Scholar]
  62. Schiminovich, D., Ilbert, O., Arnouts, S., et al. 2005, ApJ, 619, L47 [Google Scholar]
  63. Schreiber, C., Pannella, M., Elbaz, D., et al. 2015, A&A, 575, A74 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  64. Shapley, A. E., Sanders, R. L., Salim, S., et al. 2022, ApJ, 926, 145 [NASA ADS] [CrossRef] [Google Scholar]
  65. Shapley, A. E., Sanders, R. L., Reddy, N. A., Topping, M. W., & Brammer, G. B. 2023, ApJ, 954, 157 [NASA ADS] [CrossRef] [Google Scholar]
  66. Takeuchi, T. T., Yoshikawa, K., & Ishii, T. T. 2003, ApJ, 587, L89 [Google Scholar]
  67. Traina, A., Gruppioni, C., Delvecchio, I., et al. 2024, A&A, 681, A118 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  68. van der Wel, A., Franx, M., van Dokkum, P. G., et al. 2014, ApJ, 788, 28 [Google Scholar]
  69. Wang, T., Schreiber, C., Elbaz, D., et al. 2019, Nature, 572, 211 [Google Scholar]
  70. Wang, T.-M., Magnelli, B., Schinnerer, E., et al. 2022, A&A, 660, A142 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  71. Weaver, J. R., Kauffmann, O. B., Ilbert, O., et al. 2022, ApJS, 258, 11 [NASA ADS] [CrossRef] [Google Scholar]
  72. Weaver, J. R., Davidzon, I., Toft, S., et al. 2023, A&A, 677, A184 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  73. Weibel, A., Oesch, P. A., Barrufet, L., et al. 2024, MNRAS, submitted [arXiv:2403.08872] [Google Scholar]
  74. Williams, C. C., Alberts, S., Ji, Z., et al. 2024, ApJ, 968, 34 [NASA ADS] [CrossRef] [Google Scholar]
  75. Wuyts, S., Förster Schreiber, N. M., van der Wel, A., et al. 2011, ApJ, 742, 96 [NASA ADS] [CrossRef] [Google Scholar]
  76. Xiao, M. Y., Elbaz, D., Gómez-Guijarro, C., et al. 2023, A&A, 672, A18 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  77. Zavala, J. A., Casey, C. M., Manning, S. M., et al. 2021, ApJ, 909, 165 [CrossRef] [Google Scholar]

Appendix A: Representativeness of the stacked samples

As our stacked samples include only a small fraction of the sources in our parent sample of z ∼ 4.5 SFGs (those with an ALMA coverage), their representativeness may not be optimal. To check this, we compared the distributions of key physical properties (redshift, stellar mass, UV luminosity, and UV spectral slope) in our parent and stacked samples, i.e., all galaxies with an ALMA coverage (“All”) and galaxies with an ALMA coverage but excluding the ALMA primary targets (“Off”; Fig. A.1). To highlight the potential bias of the ALMA primary target subsample, we also examined its distributions, although this subsample has not been stacked separately in our analysis.

thumbnail Fig. A.1.

Representativeness of the stacked samples. Normalized distribution of some key physical properties, i.e., redshift (top left), stellar mass (top right), UV luminosity (bottom left), and UV spectral slope (bottom left) in our parent sample (green; “Par.”), our stack samples, including (red; “All”) and excluding the ALMA primary targets (“Off”), and the sub-sample of only the ALMA primary targets (yellow; “Targ.”). In each panel, we provide the p-value from the two-sample Kolmogorov-Smirnov (KS) analysis performed between these different samples and our parent sample.

For all these physical parameters, our stacked samples (“All” and “Off”) have very similar distributions to our parent sample. A two-sample Kolmogorov-Smirnov (KS) analysis indicates that our “All” and “Off” samples are both consistent with being randomly drawn from the same distribution as our parent sample (p-value  >  0.05). Moreover, the similarity between our parent and the “Off” samples is consistently better than between our parent and the “All” samples. This is naturally explained by the fact that the “All” sample contains the small subsample of ALMA primary targets (64 out of 440), which can, with high confidence, be ruled out as being randomly drawn from our parent sample (p-value  <  0.05). As expected, these ALMA primary targets are biased toward bright, massive galaxies, in a redshift range favorable for [CII] observations. However, this subsample does not appear to be biased toward galaxies with a particular UV spectral slope. This last finding probably explains why there is no significant difference between the IRX–βUV and IRX–M* relations inferred by including or excluding these ALMA primary targets from our stacks.

In short, this analysis reinforces our choice to perform our stacking analysis by including and excluding these ALMA primary targets, as these latter are slightly biased towards bright, massive galaxies. Excluding these ALMA primary targets gives us a fairly representative sample of our parent sample of z ∼ 4.5 SFGs.

Appendix B: An alternative fiducial SMF

Here we test the robustness of our calculations against the most recent SMF of Weibel et al. (2024), derived by combining the public JWST/NIRCam imaging programs from CEERS (Bagley et al. 2023), PRIMER (PI: J. Dunlop), and JADES (Eisenstein et al. 2023), covering a total area of ∼500 arcmin2. Using this alternative SMF, we repeated the analysis performed in Sect. 5.3, that is, combining it with our IRX–M* and MS relations to calculate the dust-attenuated and unattenuated cosmic SFRD at z ∼ 4.5 and the evolution of these quantities with stellar mass. The result of this analysis is shown in Fig. B.1.

thumbnail Fig. B.1.

Evolution with stellar mass of the cosmic SFRD at 4 < z < 5. This evolution was inferred using the JWST-derived SMF of Weibel et al. (2024, red and orange circles) instead of the COSMOS2020-derived SMF of Weaver et al. (2023). The symbols are the same as in Fig. 4. However, to allow a direct comparison with the results displayed in Fig. 4, in the top left panel we show with black and gray circles the SMF of Weaver et al. (2023), in the bottom left panel we show with black and red circles the cosmic SFRD and dust-attenuated SFRD at 0.03 × M⋆ (i.e., “total” cosmic SFRDs) measured using the SMF of Weaver et al. (2023), and in the bottom right panel we show with the black line the fraction of the SFRD that is attenuated measured using the SMF of Weaver et al. (2023). Also, compared to Fig. 4, the x-axis has been extended down to 107 M⊙ and the top right panel is now the fraction of the 107 M⊙ SFRD.

In the common stellar mass range, there is very good agreement in the SMFs of Weaver et al. (2023) and Weibel et al. (2024). In particular, as in Weaver et al. (2023), the JWST-derived SMF deviates at the massive end from the canonical Schechter function and requires a combination of a Schechter function and a power-law function to be fit. Moreover, such fit yields to a characteristic mass, M⋆, of 1010.4 M⊙, which is in perfect agreement with that deduced from the Weaver et al. (2023) data points. In fact, the only significant difference between the SMF of Weibel et al. (2024) and that of Weaver et al. (2023) appears at stellar masses not probed by the COSMOS2020 catalog (i.e., at ≲109 M⊙), where the JWST-derived SMF has a slightly steeper faint-end slope, −1.79 versus −1.56.

The very good agreement between the SMF of Weibel et al. (2024) and that of Weaver et al. (2023), especially down to 0.03 × M⋆ where the “total” cosmic SFRD is commonly defined, implies that our results remain qualitatively unchanged using this alternative SMF: the convergence of the cosmic SFRD occurs but at lower stellar mass (107 M⊙ vs. 108 M⊙); at the convergence, the fraction that is dust attenuated is lower but still significant (43% vs. 60%); the contribution of very low stellar mass galaxies increases but remains sub-dominant (35% vs. 15% for M* < 109 M⊙); while the “total” cosmic SFRD (i.e. at M* = 108.9 M⊙) and the fraction of it that is attenuated by dust remains mostly unchanged from those deduced using the COSMOS2020-derived SMF of Weaver et al. (2023).

Because this alternative SMF does not alter our results over the range of stellar mass for which we are able to measure the dust-attenuated properties of z ∼ 4.5 SFGs (i.e., ≳109 M⊙), we have decided to retain as our fiducial SMF the COSMOS2020-derived SMF from Weaver et al. (2023). This ensures also consistency between the SMF and our stack samples. We defer to future work the possibility of extending the constraints on the dust-attenuation of SFGs to lower stellar masses using JWST-based catalogs.

All Tables

Table 1.

Measured dust-attenuated luminosities and derived properties of our six uv-domain stacks.

In the text

All Figures

thumbnail Fig. 1.

Cutouts of the stacked images. These cutouts correspond to our 109.5 < M*/M⊙ < 1010 (left), 1010 < M*/M⊙ < 1010.5 (center), and 1010.5 < M*/M⊙ < 1011.5 (right) stellar mass bins, including (top) or excluding (bottom) the ALMA primary targets. The number of stacked galaxies is given in each cutout. The red dashed ellipse shows the ALMA synthesized beam. Color scales are similar for a given stellar mass bin. Each cutout has a size of 4″ × 4″.
In the text
thumbnail Fig. 2.

Dust attenution properties of SFGs at z ∼ 4.5. Left: IRX–βUV relation derived at 4 < z < 5 by stacking all SFGs with a A3COSMOS coverage (dark-green hexagons) or by stacking all SFGs except the ALMA primary targets (lime-green circles). Hexagons and circles with faded colors correspond to our ⟨ L UV phot . ⟩ $ \langle L_{\mathrm{UV}}^{\mathrm{phot.}}\rangle $ and ⟨ β UV phot . ⟩ $ \langle\beta_{\mathrm{UV}}^{\mathrm{phot.}}\rangle $ estimates, i.e., inferred from a COSMOS2020 photometry as opposed to their best-fit SED (see Sect. 4.2). For each stellar mass bin, the fraction of galaxies for which we were able to make such photometric measurement (see text for details) is given in the bottom part of the panel. In parenthesis, we provide the same number but after excluding the ALMA primary targets. The gray shaded squares present measurements for UV-selected SFGs at z ∼ 4 − 8 from Bowler et al. (2024). The black solid, dash-dotted, dashed, and dotted lines show the IRX–βUV relation observed in local starbursts by Meurer et al. (1999) and by Calzetti et al. (2000), an updated version for local galaxies (Overzier et al. 2011), and an SMC-like dust attenuation relation, respectively. The orange dash-dot-dotted line shows the relation observed in a UV-selected sample of SFGs at z ∼ 4.5 (Fudamoto et al. 2020a). Right: IRX–M* relation at 4 < z < 5. Symbols are the same as in the left panel. Dark-green and lime-green solid lines are linear fits to these data points, while the associated shaded regions show the 1σ uncertainties on these fits (i.e., the 16th to 84th ranges). The pink dotted line shows the redshift-independent IRX–M* relation derived by Pannella et al. (2015) for z ≲ 4 SFGs. The black dash-dotted line shows the IRX–M* relation at z ∼ 4 inferred by Koprowski et al. (2018). The blue dashed line shows the IRX–M* relation at z ∼ 3 and z ∼ 4.5 found in Fudamoto et al. (2020a,b).
In the text
thumbnail Fig. 3.

SFR–M* relation at 4 < z < 5. Dark-gray hexagons show the dust-attenuated SFRs (i.e., SFRIR) inferred by stacking all SFGs with a A3COSMOS coverage, while dark-green hexagons show their total SFRs, i.e., adding the contribution of their unattenuated SFRs measured in the UV (i.e., SFRUV; symbols are slightly offset along the stellar mass axis for clarity). Light-gray and lime-green circles display the same quantities but for our stacking analysis which excludes the ALMA primary targets (symbols are slightly offset along the stellar mass axis for clarity). The dark-green solid and red dotted lines represent the MS and its ±0.3 dex dispersion from Schreiber et al. (2015), whose normalization happens to be perfectly consistent with our dark-green hexagons, i.e., MSAll,  Fit ≡ MSS15. The lime-green line is obtained by renormalizing the MS of Schreiber et al. (2015) to fit our lime-green circles, i.e., MSOffcenter,  Fit ≡ MSS15 − 0.3 dex. For each stellar mass bin, we give the fraction of ALMA primary targets and in parenthesis their weights in our stacks. Finally, pink squares and diamonds show the total and dust-attenuated SFRs found in Khusanova et al. (2021) by stacking in three stellar mass bins, ALMA observations of all (detected and undetected) z ∼ 4.5 SFGs from the ALMA ALPINE survey.
In the text
thumbnail Fig. 4.

Evolution with stellar mass of the cosmic SFRD at 4 < z < 5. Top left: fiducial SMF used in our calculations is represented by the gray line and shaded region. This is the combination of a Schechter and power-law function fit to the SMF of Weaver et al. (2023, gray and black circles). A Schechter function fit to these data points is shown by the gray dotted line. Bottom left: cosmic SFRD (thick black line) and dust-attenuated SFRD (red line) above a given stellar mass, as inferred using the fiducial SMF. The dotted gray and red lines present the same quantities but using a simple Schechter SMF in our calculations. Top right: fraction of the M* > 108 M⊙ SFRD (gray histogram) and dust-attenuated SFRD (red histogram) that is attributed to a given stellar mass bin. Bottom right: fraction of the cosmic SFRD above a given stellar mass that is dust-attenuated (red region) and unattenuated (blue region). The dotted gray line present the same quantity but using a simple Schechter SMF in our calculations. In all panels, the vertical line represents the 0.03 × M⋆ lower integration limits commonly used to infer the “total” cosmic SFRD (e.g., Madau & Dickinson 2014). All quantities in the bottom-left, top-right, and bottom-right panels were calculated using our fiducial SMF, IRX–M*, MS relations. The propagation of the 1-σ uncertainties on the SMF and the IRX–M* relation are represented by shaded regions. The uncertainties on the MS would move the black and red lines of the bottom-left panel up and down by 0.15 dex, but would have no effect on the other quantities displayed in this figure.
In the text
thumbnail Fig. 5.

Redshift evolution of the “total” cosmic SFRD. The cosmic SFRD (UV+IR) and dust-attenuated SFRD (IR) derived in this work are shown by black and orange circles, respectively (offset slightly along the redshift axis for clarity). To be consistent with Madau & Dickinson (2014), these measurements were obtained by integrating the SMF down to 0.03 × M⋆ (i.e., 108.9 M⊙). Orange open symbols are a collection of dust-attenuated SFRDs (i.e., SFRDIR) from the literature (Sanders et al. 2003; Takeuchi et al. 2003; Schiminovich et al. 2005; Dahlen et al. 2007; Reddy & Steidel 2009; Magnelli et al. 2011, 2013; Cucciati et al. 2012; Bouwens et al. 2012; Gruppioni et al. 2013; Finkelstein et al. 2015; Bouwens et al. 2015; Liu et al. 2018; Khusanova et al. 2020, 2021; Zavala et al. 2021; Algera et al. 2023; Traina et al. 2024). Light blue open symbols are a collection of unattenuated SFRDs (i.e., SFRDUV,  uncorr) from the literature, while dark blue open symbols correspond to the same measurements after correction for dust attenuation (i.e., SFRDUV,  corr). The dark-gray line is the redshift evolution of the cosmic SFRD inferred in Madau & Dickinson (2014). Orange dotted and dashed lines show the redshift evolution of the dust-attenuated SFRD inferred in Koprowski et al. (2017) and Maniyar et al. (2018), respectively.
In the text
thumbnail Fig. 6.

Cosmic dust attenuation in magnitudes (AUV ≡ 2.5 × log10(1 + SFRDIR/SFRDUV)) as a function of redshift. The cosmic dust attenuation inferred in this work is shown by an orange circle. Blue and orange symbols are a collection of cosmic dust attenuation measurements from the literature.
In the text
thumbnail Fig. A.1.

Representativeness of the stacked samples. Normalized distribution of some key physical properties, i.e., redshift (top left), stellar mass (top right), UV luminosity (bottom left), and UV spectral slope (bottom left) in our parent sample (green; “Par.”), our stack samples, including (red; “All”) and excluding the ALMA primary targets (“Off”), and the sub-sample of only the ALMA primary targets (yellow; “Targ.”). In each panel, we provide the p-value from the two-sample Kolmogorov-Smirnov (KS) analysis performed between these different samples and our parent sample.
In the text
thumbnail Fig. B.1.

Evolution with stellar mass of the cosmic SFRD at 4 < z < 5. This evolution was inferred using the JWST-derived SMF of Weibel et al. (2024, red and orange circles) instead of the COSMOS2020-derived SMF of Weaver et al. (2023). The symbols are the same as in Fig. 4. However, to allow a direct comparison with the results displayed in Fig. 4, in the top left panel we show with black and gray circles the SMF of Weaver et al. (2023), in the bottom left panel we show with black and red circles the cosmic SFRD and dust-attenuated SFRD at 0.03 × M⋆ (i.e., “total” cosmic SFRDs) measured using the SMF of Weaver et al. (2023), and in the bottom right panel we show with the black line the fraction of the SFRD that is attenuated measured using the SMF of Weaver et al. (2023). Also, compared to Fig. 4, the x-axis has been extended down to 107 M⊙ and the top right panel is now the fraction of the 107 M⊙ SFRD.
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.