EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 664 (August 2022) A&A, 664 (2022) A73 Full HTML
Open Access
Issue
A&A
Volume 664, August 2022
Article Number A73
Number of page(s) 39
Section Extragalactic astronomy
DOI https://doi.org/10.1051/0004-6361/202142554
Published online 05 August 2022

© D. Burgarella et al. 2022

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.

1. Introduction

Since the very first papers (e.g., Madau et al. 1996) on high redshift galaxies such as the Lyman break galaxies (LBGs), the issue of how much of their energy is emitted in the far-infrared (far-IR) has been an open question in the early universe. Today, a new question is coming to the forefront and we wonder what the dust cycle in high redshift galaxies is, that is how are large dust masses of dust formed at very high redshifts and efficiently destroyed or removed later, and what are the characteristics of these dust grains.

Thanks to deep observations with the Herschel space observatory (Herschel), the Atacama Large Millimeter/Submillimiter Array (ALMA) in the Southern Hemisphere, and the Northern Extended Millimeter Array (NOEMA) in the Northern Hemisphere, we have started to explore the dusty part of the spectral energy distributions (SEDs) of high redshift star-forming galaxies (Hiz-SFGs; Bouwens et al. 2016; Burgarella et al. 2020; Sugahara et al. 2021; Hashimoto et al. 2019; Koprowski et al. 2020; Faisst et al. 2020a,b, and others) to obtain an overall view of the energy budget of normal star-forming galaxies in the early universe.

Still, except for lensed objects, targeted studies are not very successful and most of the Hiz-SFGs are not individually detected in the rest-frame far-IR, that is in the observed submillimeter (submm) and millimeter ranges (e.g., Bouwens et al. 2016; Bethermin et al. 2020; Burgarella et al. 2020; Faisst et al. 2020b). Stacking remains the preferred tool to learn what the statistical dust properties of these galaxies are (e.g., Álvarez-Márquez et al. 2016, 2019; Carvajal et al. 2020).

Another interesting way to study these objects is through the far-IR fine-structure lines. These lines have been observed with the aim of understanding the interstellar medium properties and gas cooling in the neutral and ionized phases of local low-metallicity star-forming galaxies (e.g., Madden et al. 2013; Cormier et al. 2019; Fernández-Ontiveros et al. 2016). In the high redshift universe, several of these lines have also been identified (Harikane et al. 2014; De Breuck et al. 2019; Cunningham et al. 2020; Pavesi et al. 2019, for instance). We note that [OIII]88 μm is very strong and may be the most intense line at high redshift where we expect that most galaxies contain a large population of O stars and have a low metallicity (Arata et al. 2020). Attempts at detecting [OIII]88 μm have indeed been quite successful (e.g., Inoue et al. 2016; Hashimoto et al. 2019; Álvarez-Márquez et al. 2019). Since some of these lines (e.g., [OIII]88 μm and [CII]158 μm) are correlated with the star formation rate (SFR) of the galaxies, they can also be helpful when performing SED fitting by providing additional constraints on the recent SFR. Recently, the ALPINE-ALMA [CII] survey (Le Fèvre et al. 2020; Faisst et al. 2020a; Bethermin et al. 2020) initiated an effort to detect the [CII]158 μm line of 118 Hiz-SFGs in the redshift range 4.4 < z < 5.9. The [CII]158 μm line is one of the dominant gas coolants (Lagache et al. 2018). This line is detected in the ALPINE observations with an overall detection rate of 64% and a signal-to-noise ratio (S/N) threshold larger than 3.5σ. The ALPINE observations also allowed for the dust continuum of 23 of the 118 galaxies above 3.5σ to be detected (Bethermin et al. 2020).

The dust mass (Mdust) can be estimated from the observed flux density once we can assume a dust temperature (Tdust) and other dust-related physical parameters (e.g., Pozzi et al. 2021). However, even though it might seem obvious, estimating Mdust is much safer when using multiple data from the far-IR SEDs at different wavelengths. There is some degeneracy between Tdust and Mdust, and more than a few data points lying both on the Rayleigh-Jeans (RJ) side, above the wavelength of the peak of the dust emission and on the Wien side, and below the peak wavelength of the dust emission. This information is important to constrain the SED (Liang et al. 2019). Because the individual approach is still quite a difficult task when dealing with the dust emission of galaxies close to or in the epoch of reionization (EoR), we need to combine the information from individual galaxies – assuming they share characteristic dust emission – to improve our wavelength coverage.

In this paper, we adopt a statistical approach using a sample of ALMA-observed high-redshift star-forming galaxies (ALPINE, Le Fèvre et al. 2020) the methodology already described in Burgarella et al. (2020) and Nanni et al. (2020). The redshifts of the final sample cover a range 4.5 ≲ z ≲ 6.2 near or in the EoR that will be explored by the James Webb Space Telescope (JWST). The ALPINE sample is a unique source of data, which we combined to the one from Burgarella et al. (2020)1 to build an IR composite template corresponding to the characteristic far-IR dust emission of Hiz-SFGs in the early universe. Some original works (e.g., Ouchi et al. 1999) have tried to use local starburst galaxies to build such templates and therefore estimate upper limits of submm flux densities. Similar to in Shapley et al. (2003) where a high S/N spectrum of z ∼ 3 LBGs was derived, and even in a more similar way in Pearson et al. (2013) in which 40 H-ATLAS sources with previously measured redshifts in the range 0.5 < z < 4.2 were used to derive a suitable average template for high-redshift H-ATLAS sources, the observed data from our sample were used to build the template.

In addition to using the composite IR SED to derive the dust properties of the galaxies via SED fitting, these data are also unique to help understand the stellar populations of these Hiz-SFGs and to calibrate important diagnostic diagrams in the early universe. Once an IR composite template is safely estimated from the observed ALMA detection, the spectral information in the rest-frame ultraviolet (UV), optical, and near IR is combined and we can come back to the individual approach to fit each galaxy.

In this paper, we assume a Chabrier initial mass function (IMF, Chabrier 2003). We use WMAP7 cosmology (Komatsu et al. 2011).

2. The sample of studied galaxies

The ALPINE sample is representative of the overall Hiz-SFGs population in the redshift range 4.5 ≲ z ≲ 5.5. It is not dominated by IR-bright (i.e., in dust continuum) submm galaxies (SMGs). The ALPINE galaxies are mainly located on or near the main sequence (e.g., Rodighiero et al. 2011; Gruppioni et al. 2013; Donnari et al. 2019; Sherman et al. 2021) in the SFR versus stellar mass (Mstar) relation observed at these redshifts (Speagle et al. 2014; Pearson et al. 2018). It is mostly dominated by UV-selected galaxies (see Table 1 in Faisst et al. 2020a) with about 62% of the sample identified with the dropout technique, that is LBGs and 28% that are Lyman α emitters which were selected with narrow bands. We added, to the ALPINE galaxies, seven LBGs from Burgarella et al. (2020) and Nanni et al. (2020) to build our sample. The data of the latter seven LBGs were collected from several works (Bouwens et al. 2016; Capak et al. 2015; Faisst et al. 2017; Scoville et al. 2016; Willott et al. 2015; Hashimoto et al. 2018). This selection is certainly not complete and will very likely introduce biases on the results we obtained, but this is as good as we can currently do. The origins of the data are listed in Table 1.

Table 1.

Origins of the data used in this paper.

Using Band-7 of ALMA, the ALPINE project observed a sample of galaxies with spectroscopic redshifts in two well-observed fields: 105 galaxies in the Cosmic Evolution Survey field (COSMOS, Scoville et al. 2007) and the remaining 13 in the Extended Chandra Deep Field South (ECDFS, Giacconi et al. 2002), with extensive multiband Hubble space telescope (HST) data from CANDELS (Grogin et al. 2011; Koekemoer et al. 2011). In this survey, 64% of the galaxies have been detected in [CII] at 3.5σ above the noise, as well as 21% of the galaxies detected in the continuum (S/N threshold corresponding to a 95% purity). The sample is divided into two redshift ranges at 4.40 < z < 4.65 (median redshift ⟨z⟩ = 4.5, containing 67 galaxies) and 5.05 < z < 5.90 (median redshift ⟨z⟩ = 5.5, containing 51 galaxies), separated by a gap in the transmission of the atmosphere. We note that type I active galactic nuclei (AGN), identified from broad spectral lines, are excluded from the present sample.

As the galaxies included in the ALPINE sample belong to well-studied fields, they come with a rich ancillary data set (Faisst et al. 2020a). Due to the nature of the selection of these galaxies, they all have spectroscopic observations in rest-frame UV, performed with the Keck telescope and the European Very Large Telescope (VLT). A plethora of photometric observations are also available, from ground-based observatories in the UV to optical, from the HST in the UV, as well as from Spitzer above the Balmer break (all rest-frame features).

3. Methodology

The objective of the paper is to study the physical properties of the sample of Hiz-SFGs, and more specifically their stellar populations and dust grains. Because we do not have a wide wavelength coverage of their IR dust emission, we built an IR composite template from the subsample of ALMA-detected objects. To this aim, we applied the methodology described in Fig. 1 and it is detailed.

thumbnail Fig. 1.

Flow chart of the process followed to build the composite IR SED for the Hiz-SFGs sample. The various phases of this process are shown on the right side of the figures. They are also described in more detail in the text (see Sect. 3) In the first phase, we selected only ALMA-detected objects to build the first IR template (IRV1). Next, we added the data from the z ∼ 4.5 and ∼5.5 stacks from Bethermin et al. (2020) to build the second and final IR template (IRV2). Adding UV and optical data to IRV2, we fit the entire (detected and upper limits) sample.

3.1. Building the IR composite template

3.1.1. Phase 0: Selection of the galaxies used to build the IR composite template

From the ALPINE sample and the one compiled by Burgarella et al. (2020), we selected the objects for which the dust continuum was detected with ALMA. More details on the selection are provided in Table 1. The final sample of ALMA-detected objects contains 27 galaxies.

3.1.2. Phase 1: CIGALE initial individual fits of ALMA-detected galaxies to estimate the normalization at λ = 200 μm

The objective of this initial phase is to estimate a normalization factor for the SEDs of all the ALMA-detected galaxies. Galaxies with upper limits in the far-IR were not selected to build the IR composite SED. This normalization factor was derived from the estimated flux density at λ = 200 μm (rest-frame), which is approximately the maximum wavelength above which we have no observed data. The SEDs of all the ALMA-detected objects at λ = 200 μm were set to 1.0 to build the IR composite SED. All the other data were modified accordingly. The spectral models selected in this initial fit are described in detail in Appendix A. We selected a wide range of models, especially for the dust emission which is important for this work.

We fit all the ALMA-detected galaxies using several options with CIGALE (Burgarella et al. 2005; Boquien et al. 2009, see Appendix A for details on the priors for the models). CIGALE cannot model emission lines from photo-dissociation regions (PDRs). We used the [CII]158 μm fluxes to estimate the star formation rate (SFR[CII]158 μm), which we evaluated with the relation from Schaerer et al. (2020). These SFR[CII]158 μm values, along with the associated uncertainties derived from the line uncertainties, were added as properties to constrain the CIGALE SED fitting2. To perform this analysis, we selected three types of dust emission available in CIGALE: (i) a mid-IR power law and a general modified blackbody (Casey 2012): PL+G_MBB, (ii) a mid-IR power law and an optically thin modified blackbody (Casey 2012): PL+OT_MBB, and, (iii) the models from Draine et al. (2014): DL2014.

It is specifically worth checking whether the prior assumptions made on the emissivity index (βRJ) of the modified blackbody could bias the final result. We could not safely estimate βRJ for the individual objects but we could do so better with the IR composite template because more data points are used, thus we also fixed βRJ = 1.0, 1.5, and 2.0 to test the impact on the level of the luminosity at 200 μm (rest-frame) used for the normalization of the individual SEDs. Table 2 shows that the bias introduced by the initial prior on βRJ is much lower that the uncertainties coming from the values derived by an SED fitting. However, we would like to stress that a wide range of SED shapes were used to estimate the normalization, including several βRJ, Tdust, and DL2014 models.

Table 2.

Parameters of the dust emission.

3.1.3. Phase 2: Building of the IR composite template

In this second phase, we proceeded by actually building the observed composite IR template. For this, we made use of the normalized SEDs of each of the 27 ALMA-detected objects.

– Phase 2.1: Checking the homogeneity in the galaxy sample → IRV1 template: By making use of the normalization factors and benefiting from the redshift coverage (4.3 < z < 6.2), we used the universe as a spectrograph and we combined the SEDs of 27 ALMA-detected objects to form an initial version (IRV1) of the observed composite SED (Figs. 1 and 2). We note that this process is different from the stacking method whose aim is to detect the average flux density at a given wavelength of a sample of detected and/or undetected objects.

thumbnail Fig. 2.

IR combined SED built from the ALMA-detected objects from the ALPINE sample, from Burgarella et al. (2020), and with the two stacks from Bethermin et al. (2020) and the associated uncertainty regions. We note that the four left-most points for the stacks are upper limits, as shown by the error bars reaching the bottom of the figure. The SEDs corresponding to optically thin modified black bodies are superimposed to the observed SEDs. Qualitatively, we can see that the IR-combined SED seems to be in agreement with greenish modified blackbody SEDs, that is Tdust ∼ 50–60 K. A qualitative analysis is performed later in this paper.

We began by verifying whether the subsample of ALMA-detected galaxies is homogeneous enough to build a single IR composite template. For instance, one single object or all the objects in a given redshift range might significantly impact the parameters of the dust emission and bias the IR composite template. Here, we only show the results of this test assuming PL+OT_MBB, but we checked that the results were the same regardless of the assumed dust emission option.

To reach this goal, we used the knife-jacking method, where we removed the SED of each galaxy (or a series of galaxies in a redshift range) sequentially from the list of objects, built the composite template, and fit it with CIGALE. Then, we put the removed objects back in the detected subsample and we reproduced this operation recursively as many times as there were objects (or series of galaxies in a redshift range) in the subsample. At the end, we could compare the parameters from the several composite templates built from Nobj–Nremoved-obj and we checked whether they were consistent.

From the analysis on individual objects (Fig. 3), we discovered that only when removing HZ10 at z = 5.659 did we get slightly different values for Tdust and βRJ. However, it does not significantly bias the parameters defining the shape of the composite template. We decided to keep it in the sample because it contributes to the representativity of the studied galaxies. The dust parameters are very stable with average values at Tdust ≈ 54 K and βRJ ≈ 0.9.

thumbnail Fig. 3.

Degeneracies in Tdust and βRJ. The knife jacking method allowed us to estimate whether the characteristics of the dust emission, i.e., Tdust and βRJ, vary when removing and adding back individual galaxies from the sample used to build the composite template. The ellipses show the uncertainties. The outlier point at βRJ ∼ 1.1 and Tdust ∼ 51 K was obtained when taking off HZ10 at z = 5.659 with five pieces of data in the far-IR+submm range. Only when removing it did the average locus move to slightly larger values of βRJ, but with larger uncertainties.

From the analysis on redshift ranges, we discovered that the redshift range 5.0 < z < 6.0 brings the objects with the widest wavelength range to the IR composite SED. Without the ALMA-detected galaxies in the range 5.0 < z < 6.0, there was a difference for Tdust and βRJ (Table 3). But again, the results are in agreement within the uncertainties, even though an evolution within the uncertainties cannot be ruled out.

Table 3.

Degeneracy of βRJ and Tdust.

This second phase provided an IR composite template IRV1 shown as green points in Fig. 2.

– Phase 2.2: Checking if the stacked data from Bethermin et al. (2020) are in agreement with the IRV1 → IRV2 template: In addition to the main ALPINE data, Bethermin et al. (2020) stacked data from two galaxy samples at z ∼ 4.5 and z ∼ 5.5 with Herschel (Pilbratt et al. 2010) data from the PEP (Lutz et al. 2011) and HerMES (Oliver et al. 2012) surveys and AzTEC/ASTE data from Aretxaga et al. (2011) at 1.1 mm. At 850 μm, they used the SCUBA2 data from Casey et al. (2013). This data set can be very useful because it extends the wavelength range to the mid-IR. However, because there are not enough ALPINE sources to obtain a sufficiently high S/N in the stacked Herschel data, the ALPINE selection was not used for the stacking. Bethermin et al. (2020) define two redshift bins (4 < z < 5 and 5 < z < 6) that match the redshift ranges probed by ALPINE. Before further continuing to fulfill our task of building the LBG composite IR template, we first need to assess whether the z ∼ 4.5 and/or z ∼ 5.5 stacks are consistent with the IRV1 composite IR template.

Our analysis allows us to conclude that both stacks are consistent (all χ ν 2 $ \chi^2_\nu $ from the SED fitting have ∼0.3 − 0.4) with the ALMA-detected galaxies’ composite IR template built by combining the full sample of the 27 ALMA-detected objects. The final IR composite template (IRV2, Fig. 2) was built from the ALMA-detected galaxies from Burgarella et al. (2020), the ALPINE sample, in addition to the two stacks from Bethermin et al. (2020).

In Table C.1, we provide the IR composite template based on the observed data and show the fits with the three dust emissions from 1 μm to 1000 μm fitted in Fig. 4. The modeled SEDs from 20 μm to 1 mm are listed in Table D.1.

thumbnail Fig. 4.

Comparison of the various fits of the IR composite SED built with data from Burgarella et al. (2020), from ALPINE, and from the z ∼ 4.5 and z ∼ 5.5 stacks from Bethermin et al. (2020). Top: fit with a modified general blackbody and power law in the mid-IR as in Casey (2012). Middle: fit with an optically thin modified blackbody and power law in the mid-IR as in Casey (2012). Bottom: Fit with a model from Draine et al. (2014). The fits are all statistically equivalently good with reduced χ2 ≈ 0.3–0.4.

3.1.4. Fitting the IR composite template

The main dust parameters derived from fitting the observed IRV2 IR composite template (Fig. 4) are listed in Table 4. The objects in the present sample are Hiz-SFGs with a low dust attenuation (e.g., Faisst et al. 2020a, and later in this paper). This assumption is valid for the ALPINE sample (Faisst et al. 2020a). Therefore, we assume that the best emission models for these objects should be optically thin. We checked this hypothesis by estimating the following (Eq. (2) in Jones et al. 2020):

τ ν = M dust A gal κ ν $$ \begin{aligned} \tau _\nu = \frac{M_{\rm dust}}{A_{\rm gal}}\kappa _{\nu } \end{aligned} $$

Table 4.

Main relevant physical parameters derived by fitting the IR template with the various assumptions of the dust emission.

where τν is the optical depth, Mdust is the dust mass, Agal is the area covered by the galaxy, and κν is the dust mass absorption coefficient:

κ ν = κ 0 ( ν / ν 0 ) β RJ $$ \begin{aligned} \kappa _{\nu }\ =\ \kappa _0\ (\nu \ / \nu _0)\ ^{\beta _{\rm RJ}} \end{aligned} $$

where ν0 is the frequency where the optical depth equals unity and βRJ is the spectral emissivity index from the Rayleigh-Jeans range. We tested the opacity both with the optically thin and with the general modified black bodies (Casey 2012). The median radius of the galaxies are re, [CII] = 2.1 ± 0.16 kpc (Fujimoto et al. 2020) to estimate Agal. Both provide values that are in the range [10−4–10−2], that is τν ≪ 1.0, which confirms that the optically thin assumption is valid here. In the rest of the paper, we only make use of DL2014 and PL+OT_MBB dust emissions.

After fitting the IR composite template, we derived values for Tdust and βRJ (Table 4). However, these dust parameters are known to be degenerate depending on the S/N and the wavelength sampling of the far-IR SED (e.g., Juvela et al. 2013; Tabatabaei et al. 2014). This degeneracy is a serious problem when using only one or a few data points in the far-IR. However, using the derived IR composite template helps to remove the degeneracy. We tested how well Tdust and βRJ can be estimated by performing fits with fixed Tdust and varying βRJ, then by fixing βRJ and keeping Tdust free. The results are compared to the parameters derived by keeping both parameters free in Fig. 5. Even though we can see the usual regular evolution of Tdust with βRJ (or vice versa), the quality of the fit with fixed parameters improves when getting closer to the Bayesian values derived with CIGALE when both parameters are free (Table 4).

thumbnail Fig. 5.

Results of the tests on the determination of Tdust and βRJ and their possible degeneracy are shown here. Dots show the evolution of Tdust when βRJ was fixed during the fits. Boxes show the evolution of βRJ when Tdust was fixed during the fits. Finally, the diamond symbol shows the result and uncertainties on both axes when Tdust and βRJ are both free. The color coding indicates the level of the goodness of fit.

The value of βRJ = 0.87  ±  0.28 for a power law plus optically thin modified blackbody found in this paper is low, but comparable to the minimum values found by Bendo et al. (2003), Galametz et al. (2012), for example, for nearby galaxies in the range 0.8 < βRJ < 2.5. Lower βRJ closer to 1.0 were estimated when fitting different types of galaxies and more specifically very low-metallicity galaxies such as SBS 0335-052 (Hunt et al. 2014), NGC 1705 (O’Halloran et al. 2010), and even low-metallicity regions in very nearby galaxies with excellent coverage to estimate βRJ similar to Messier 33, for instance (Tibbs et al. 2018; Tabatabaei et al. 2014). Tabatabaei et al. (2014) found an apparent decrease in both Tdust and βRJ with an increasing M33 radius. This corresponds to regions with low metallicities and they propose that the effect could help to find an origin in the different grain compositions and, possibly, different size distributions. Assuming such βRJ values for Hiz-SFGs would mean that the mean dust temperature could be higher than Tdust estimated with emissivities in the range 1.5 < βRJ < 2.0. This result (βRJ ∼ 1) needs to be further confirmed for the SEDs of similar low-Mstar, low-AFUV galaxies with a better RJ wavelength coverage.

The value of Tdust found from the IR composite template (54.1 ± 6.7 K) should be compared to other ALMA-based dust temperatures for objects at high redshift (z > 4.5). In a recent paper, Bakx et al. (2021) present the evolution of Tdust for “normal” (i.e., main sequence) galaxies from z ∼ 0 to z ∼ 8. A linear increase is suggested up to z = 6. At larger redshifts, a large dispersion in Tdust can be noticed even though the linear relation might hold, given the large uncertainties (Faisst et al. 2020b; Harikane et al. 2020; Sugahara et al. 2021; Laporte et al. 2019). Bakx et al. (2021) show that adding ALMA Band 9 significantly reduces the uncertainty on the dust temperature for a single object by further constraining the shape of the SED. At the mean redshift z = 4.94 ± 0.54 of our entire sample, the dust temperature that could be estimated assuming the linear relation from Fig. 4 in Bakx et al. (2021) would be Tdust ≈ 49 K, which is not significantly different from our estimation of 54.1 ± 6.7 K. However, Bakx et al. (2021) derived a value of βRJ = 1.61 ± 0.60 which is larger than the one found in the previous paragraph. Both are consistent if we account for the uncertainties. Lower values for βRJ mean higher dust temperatures, as shown in Fig. 5. Thus measuring βRJ with a good sampling of the RJ part of the spectrum is important to address the question about the decrease in βRJ in low-metallicity regions toward the outskirts of local galaxies and in some local low-metallicity galaxies. If so, the increase in dust temperature might be larger than when evaluated with βRJ ∼ 1.5–2.0.

3.2. Phase 3: Far-UV to far-IR fit of individual galaxies with the IR composite template

We now fit the complete sample of objects assuming the dust emission from the IRV2 composite SED. The rest-frame UV, optical, and near IR ranges are useful when constraining the stellar emission and, therefore, the properties of the stellar populations. However, attenuation also impacts the rest-frame UV spectral range, which is therefore useful to constrain the amount of dust attenuation. The assumptions used when fitting the entire galaxy sample with CIGALE are given in Table B.1. We focus on the properties of the stellar populations and the dust attenuation in these objects.

Even though the Bayes factor is conceptually more conclusive, its computation can still be very complex. We need to compute a quantity called the marginal likelihood or evidence (P(D/Mi)), where D are the photometric data and Mi are each model, which implies computing a very complicated and time-consuming integral because the likelihood of observation D under the model Mi must quantify over all possible parameters of that model. An analytic computation is almost never possible and direct evaluation by numerical quadrature is almost never feasible for models of real-world dimensionality and complexity. This is why a variety of approximations based on special properties of the models or their posterior distributions were developed (Kass & Raftery 1995).

Generally, it is found that the conclusions drawn from the Bayes factor are more satisfying, but also more complex methodologies have not been qualitatively very different from those drawn from the simpler BIC method (Raftery 1998). However, we can elaborate on the use of the BIC versus Bayes Factor for our specific case. From Raftery (1998), we determined that in our case, the unit information prior should reasonably cover the range of observed data because we have a homogeneous sample of galaxies and we have a rough idea of the general range within which the data are likely to lie in advance. So we defined our priors to homogeneously cover the expected range of parameters. Volinsky (1997) showed, via simulations, that the performance of Bayes factors can be better than that of BIC if the prior used is spread out less than the unit information prior. In our case, the unit information prior is quite spread out (see Table 5) and BIC should be at least as good as the Bayes factor. Because of this spread out prior, BIC provides more conservative results. Also, if an effect is favored by BIC, this should be a solid result. On the other hand, if BIC does not favor an effect, it might still be possible that another, justifiable prior could change the conclusion. For us, that means that the most conclusive results would be that an SFH that includes a burst is solidly ruled out when compared to a delayed SFH with τ = 500 Myr SFH.

Table 5.

Comparison of the reference model (SFH delayed: t/τ2 exp(−t/τ) with τ = 500 Myr) to all the others using the BIC.

We decided to use the Bayesian information criterion (BIC) test to compare the various hypotheses made on the dust emission and on the star formation history in Fig. 6. To interpret the results from the BIC tests, we refer to Kass & Raftery (1995):

BIC = χ 2 + k ln ( n ) $$ \begin{aligned} \mathrm{BIC}\ = \ \chi ^{2}\ +\ k\ln (n) \end{aligned} $$

thumbnail Fig. 6.

Result of the ΔBIC test on our sample of Hiz-SFGs. A delayed SFH without a burst and τmain = 500 Myr is labeled as “tau500”. With an additional burst at the end of the SFH, it is labeled as “burst”. A constant SFH without a burst and τmain = 20 Gyrs is labeled as “tau20000”. And when several τmain could be selected in the SED fitting, it is labeled as “multitau” in the legend. Top: ΔBIC test that compares the influence of the DL2014 model and the PL+OT_MBB dust emissions in building the IR template. Center: ΔBIC test on the SFH assuming the PL+OT_MBB for the IR template. Bottom: ΔBIC test on the SFH assuming the DL2014 model for the IR template. The color band allows one to interpret the results of the evidence (ΔBIC) against the model with the higher BIC: red means “faint evidence”, orange means “positive evidence”, and green mean “strong evidence”. We do not see any strong evidence that DL2014 or PL+OT_MBB are better for fitting the data. An SFH that includes a burst is positively ruled out, while a delayed SFH with τ = 500 Myr is weekly favored.

where k is the number of model parameters in the test and n is the sample size, that is to say the number of pieces of photometric data in our case. When comparing several models, the one with the lowest BIC is preferred. To interpret the results, 0 < ΔBIC < 2 means that evidence for a difference between the two hypotheses is faint (2 < ΔBIC < 6 means positive, 6 < ΔBIC < 10 means strong, and ΔBIC > 10 means very strong).

We start by comparing the two dust emission models, DL2014 (Draine et al. 2014) and PL+OT_MBB (Casey 2012), available in CIGALE. Table 5 presents all the results of the BIC analysis on our sample of Hiz-SFGs, numerically. The histogram of ΔBIC values strongly peaks at 0 with a small tail extending to ΔBIC = −4 (Fig. 6), meaning that we do not see any strong difference between DL2014 and PL+OT_MBB when fitting our sample. However, for a minority of objects, PL+OT_MBB might be favored. We show the results using both DL2014 and PL+OT_MBB later in the paper. However, in general, no significant difference as to the quality of the fit is observed, confirming the BIC analysis.

We now focus on the comparison between the various SFHs. Figure 6 shows the results of the BIC analysis by comparing the various kinds of SFHs used in the CIGALE analysis. We note that we only used the initial observed data and not the IR composite template. The only significant difference is that a simple delayed SFH with τmain = 500 Myr without any final burst is preferred to a delayed SFH with a burst. A constant SFH (with τmain = 20 Gyrs) and multi-τmain are weakly disfavored when compared to a simple delayed SFH with τmain = 500 Myrs. These comparisons are compatible with both DL2014 and PL+OT_MBB. There is some dispersion in the results which suggests that the preferred modeling for the entire galaxy population might not be valid for each individual galaxy, as shown in Fig. 6.

In order to move a little bit further with our comparison of the models, we performed a mock analysis with CIGALE from the initial fit, again without using the IR composite template. We selected the PL+OT_MBB and DL2014 options for the dust emission in Fig. C.1. In brief, for each object, the CIGALE mock analysis consists in using the best fit model for each object to generate an input photometric catalog (called a mock catalog) to which we randomly added uncertainties drawn inside a Gaussian distribution with σ corresponding to the observed uncertainties (see Boquien et al. 2019). The SED fitting procedure was applied in the very same way as when fitting the observed SEDs. Then, we compared the derived output parameters to those of the known input (Fig. C.1). The results from this mock analysis (Table 6) show that both for the PL+OT_MBB and for the DL2014 dust emissions, we were able to recover the main physical parameters related to the SFH and to the dust emission. This means that we do not expect strong degeneracies (r2  =  0.80 shown in Fig. C.1, where r is Pearson’s correlation coefficient commonly used in linear regression). The mean value and standard dispersion of the main parameters derived from the individual fits of our galaxy sample are given in Table 7 and the individual values for the same parameters are listed in Tables E.1 and E.2 (available at the CDS) in their entirety.

Table 6.

Summary of the power of two for correlation coefficients, r2, from the mock analysis on the main physical parameters performed by CIGALE.

Table 7.

Mean value and standard deviations of the physical parameters derived from the fit analysis of the individual objects in our sample.

4. Analysis of the results

From the previous analysis, we derived a set of physical parameters for each of the galaxies in our sample. These parameters allowed us to define and build diagnostic diagrams that permitted us to characterize these galaxies and, more specifically, their SFH and their dust properties. We analyzed the locations of our sample in the IRX versus βUV diagram, the AFUV versus Mstar diagram, and in the specific dust mass versus specific star formation rate: sMdust versus sSFR = Mdust/Mstar versus SFR/Mstar (dust formation rate diagram or DFRD).

4.1. The SFR versus Mstar diagram

Figure 7 presents the SFR versus Mstar diagram3. The points corresponding to this work are found in the expected range when compared to Pearson et al. (2018) who also used CIGALE:

log 10 ( SFR ) = ( 1.00 ± 0.22 ) ( log 10 ( M star ) − 10.6 ) + ( 1.92 ± 0.21 ) $$ \begin{aligned} \log _{10}(\mathrm{SFR}) = (1.00 \pm 0.22)\ (\log _{10}(M_{\rm star}) - 10.6) + (1.92 \pm 0.21)\end{aligned} $$

thumbnail Fig. 7.

In the log10 (SFR) vs. log10 (Mstar) diagram, the main part of the sample is found within the limits for the fits of the main sequence at z = 5.0 by Speagle et al. (2014; green shading), by Pearson et al. (2018) at z ∼ 5.2 (purple shading), and by Faisst et al. (2020a), who found that the galaxies are in agreement with Speagle et al. (2014). Is is important to note, however, that some of our objects are at redshifts larger than the 4.5–5.6 ALPINE sample (see color code of the markers). For Speagle et al. (2014), we used the “mixed” (preferred fit)” function (as defined by Speagle et al. 2014). The results from the two fits with DL2014 and PL+OT_MBB are presented. The two types of dust emission do not significantly modify the location of the points in the diagram. Detections are shown as dots and boxes and upper limits are shown as downward- and upward-pointing triangles. The uncertainty range for upper limits extends to the bottom of the plot. In addition to having mainly upper limits, at log10 (Mstar) < 9.5, the sample is very likely incomplete which means that it is difficult to estimate a trend from these data over the entire mass range. We also added the objects and stacks from Khusanova et al. (2021) with open markers. Finally, the selection of TNG50 galaxies used in Schulz et al. (2020) is also provided (red-shaded area).

and compared to Speagle et al. (2014) with their “mixed” (preferred fit)”:

log 10 ( SFR ) = [ ( 0.73 ± 0.02 ) − ( 0.027 ± 0.006 ) × t [ Gyr ] ] log 10 ( M star ) − [ ( 5.42 ± 0.22 ) + ( 0.42 ± 0.07 ) × t [ Gyr ] ] $$ \begin{aligned}\log _{10}(\mathrm{SFR})&= [(0.73 \pm 0.02) - (0.027 \pm 0.006) \times t[\mathrm{Gyr}] ] \log _{10}(M_{\rm star}) \nonumber \\&\quad - [(5.42 \pm 0.22) + (0.42 \pm 0.07)\times t[\mathrm{Gyr}]] \end{aligned} $$

function evaluated at z = 5.0 (i.e., tuniverse = 1.186 Gyr), which was converted to a Chabrier IMF by subtracting 0.03 to log10(SFR) and to log10(Mstar). The first result is that regardless of the dust emission used (DL2014 or PL+OT_MBB), we found about the same main sequence. Our data are in good agreement with Faisst et al. (2020a) and Khusanova et al. (2021) and also they generally follow the relations derived by Speagle et al. (2014) and Pearson et al. (2018). Most of the detections are found at a relatively large stellar mass (log10(Mstar)∼10). We confirm the evolution of the main sequence to z = 4.5 − 5.5 with our sample of galaxies (mainly objects not detected in dust continuum), which extend to log10(Mstar)∼8.5.

4.2. The AFUV versus Mstar diagram

The relation between the dust attenuation (AFUV or its proxy, IRX = Ldust/LFUV) and the stellar mass (Mstar) is another way to estimate the dust attenuation in galaxies without far-IR data. This relation between the stellar mass and dust attenuation has been the focus of numerous studies (as early as Xu et al. 2007; Buat et al. 2009, and references therein). However, even if this relation could be very useful in addition to the IRX versus βFUV one (fλ ∝ λβFUV), the link between the far-UV (FUV) dust attenuation and the stellar mass (AFUV versus Mstar) is not well established at all redshifts.

Bouwens et al. (2016) define what a consensus relation could be: log10(IRX) = log10(Mstar/M⊙) – 9.17 assuming the dust temperature evolves with the redshift. This relation is linear in the plane log10(IRX) versus log10 (Mstar/M⊙), and in the range 9.0 < log10 (Mstar/M⊙) < 11.0. In a recent paper, Carvajal et al. (2020) used the stacking method for 1582 UV LBGs with photometric redshifts in the range z ∼ 2–8 to reach down to log10 (Mstar/M⊙) = 6.0. However, the constraints from this stacking are only upper limits, which are less useful than detections (see their Fig. 16). Fudamoto et al. (2020) made use of the ALPINE data and show that the log10(IRX) versus Mstar relation derived from their observations is inconsistent with the previously determined relations at z ≤ 4. They found a fast decrease in IRX at z ∼ 4 in massive galaxies which suggests an evolution of the average amount of dust attenuation in star forming galaxies. Bernhard et al. (2014) assume an evolving normalization of the log10(IRX) versus log10(Mstar) relation in the low redshift universe at z < 1. We need an opposite evolution (lower in IRX) in the high redshift universe at z > 4 as in Bogdanoska & Burgarella (2020).

Figure 8 shows that the two Mstar stacks from Fudamoto et al. (2020) at z ∼ 4.5 are marginally (accounting for the uncertainties) in agreement with Bogdanoska & Burgarella (2020). They are in better agreement at z ∼ 5.5. However, the low stellar mass range is crucial, especially in the early universe where galaxies are expected to have lower stellar masses because they were still building their stars at that time. From our analysis, we confirm the effect found by Fudamoto et al. (2020): our galaxies have lower AFUV or IRX values at a fixed mass, compared to previously studied IRX–Mstar relations at z < 4, with quite a large scatter.

thumbnail Fig. 8.

AFUV − Mstar diagram at z ∼ 4.5 (left) and z ∼ 5.5 (right). The gray areas correspond to the expected relation at z = 4.5 and 5.5 from Bogdanoska & Burgarella (2020). This relation was formed by a broken line, which is flat at log10Mstar ≤ 9.8 and rises at log10 (Mstar) > 9.8. The ALPINE data are very dispersed at z ∼ 4.5, while this flatness is supported by the data at z ∼ 5.5. The conversion from IRX to AFUV is from Burgarella et al. (2005): AFUV = −0.028[log10 (IRX)]3 + 0.392 [log10 (IRX)]2 + 1.094 [log10 (IRX)] + 0.5.

In order to constrain the amount of dust attenuation of these low-mass galaxies, it is possible to adopt a global approach and compare what the hypothesis on the AFUV versus log10 (Mstar) implies for the redshift evolution of the average dust attenuation in the universe as presented in Cucciati et al. (2012), Burgarella et al. (2013), Madau & Dickinson (2014), for example. Using this approach, Bogdanoska & Burgarella (2020) conclude that it is not possible to extend a linear relation to the lowest mass range without strongly underpredicting the average dust attenuation in the universe at all redshifts. We need a flattening of the relation at low stellar mass (log10(Mstar/M⊙) < 9.0). This means that the apparent dust attenuation of the low-mass galaxies is significantly higher than 0. This very interesting and unexplained point is also suggested by other observational and theoretical works (e.g., Salim et al. 2016; Takeuchi et al. 2010; Cousin et al. 2019; Ma et al. 2016). They propose modeling the AFUV versus Mstar relation with the following broken law:

A FUV = { a × 1.1 at log 10 M star < 9.8 a × ( log 10 M star − 8.7 ) at log 10 M star ≥ 9.8 $$ \begin{aligned} A_{\rm FUV} = {\left\{ \begin{array}{ll} {a \times 1.1}&{\mathrm{at}\,\log _{10}\,M_{\rm star} < 9.8}\\ {a \times (\log _{10} M_{\rm star} - 8.7)}&{\mathrm{at} \log _{10}\,M_{\rm star} \ge 9.8}\\ \end{array}\right.} \end{aligned} $$(1)

where a = (z + γ)×αβ − (z + γ) with z being the redshift and the constants α = 1.84  ±  0.11, β = 1.84  ±  0.12, and γ = 0.14 ± 0.04. Although there is some dispersion and the completeness is very likely small at log10Mstar < 10.0, Fig. 8 suggests that the relation from Bogdanoska & Burgarella (2020) at z = 5.0 is in broad agreement with the data. However, the large number of upper limits at low mass could suggest that the level of the flat relation might be overestimated with respect to the present data. Alternate ways to measure AFUV for low-mass galaxies should be explored, maybe via the Balmer decrements with JWST or with radio data using the radio-to-IR ratio qIR that allows one to constrain the far-IR emission of star-forming galaxies from radio data (e.g., Helou et al. 1985; Delvecchio et al. 2021).

4.3. The IRX – βFUV diagram

To estimate the UV slope βFUV, we used the definition given by Calzetti et al. (1994) in the range 125–260 nm within ten selected windows (see their Table 2) designed to remove all absorption features and the 217.5 nm dust bump (βCalzetti − 1994) from the fitting procedure. The IRX – βFUV is a classical tool to estimate the dust attenuation of galaxies. Figure 9 shows the location of the galaxies in the IRX versus βFUV diagram. We also show the classical positions assuming a Calzetti law and a small Magellanic cloud (SMC) law. Our galaxies appear to be systematically shifted to the left of the diagram, that is to bluer values of βFUV compared to the locus estimated by Overzier et al. (2010) for galaxies at z = 0. This might be related to the evolution of the stellar population at z ∼ 4.5 − 5.0 as they were younger and bluer (and/or an intrinsic evolution of other parameters similar to the IMF). Age effects were already found at low redshift (e.g., Kong et al. 2004; Boquien et al. 2009).

thumbnail Fig. 9.

IRX – βCFUV diagram. IRX values were estimated with CIGALE and βFUV were fitted on the data directly. They are not model-dependant. The black continuous line corresponds to the original Calzetti law, under the assumption that the underlying dust curve follows the Calzetti et al. (2000) attenuation. The dashed line to the predicted law assumes the SMC extinction law (e.g., Gordon et al. 2003). Both are from McLure et al. (2018). The color of the symbols are related to the axis to the right of the figure, the age of the stellar population. Both results with dust emission, DL2014 and PL+OT_MBB, are presented with different symbols with ALMA upper limits and ALMA detections. In the figures, we also plotted the laws from Schulz et al. (2020) at z = 0.0 (blue continuous line), 4.5 (green dotted line), and 5.5 (red dashed line). A comparison with the IRX – βFUV plot in Fudamoto et al. (2020) with the same sample suggests that our (especially ALMA-detected) galaxies extend less to very blue βFUV. This is mainly true for the subsample of galaxies not detected in continuum with ALMA. This is probably due to two effects: first, in Fudamoto et al. (2020), 3σ upper limits are plotted. Another possible effect could be because we used a unique IR composite template in the SED fitting, which reduces the uncertainties on Ldust.

The IllustrisTNG Project (TNG hereafter: Nelson et al. 2018; Pillepich et al. 2018) is a suite of cosmological simulations for the formation of galaxies. TNG50 does not model dust directly, but Schulz et al. (2020) assume that the diffuse dust content of a galaxy is traced by the gas-phase metal distribution assigned to this galaxy. The dust density distribution is derived from the TNG50 gas density distribution with assumptions on the dust-to-metal ratio, the gas metallicity, the gas temperature, and the instantaneous SFR. Finally, SKIRT (Baes et al. 2011) is used to model the emission of the galaxies. Their fiducial model is a multicomponent dust mix, which models a composition of graphite, silicate, and polycyclic aromatic hydrocarbon (PAH) grains, with various grain size bins for each grain type which reproduce the properties of Milky Way (MW), large Magellanic cloud (LMC), and SMC type dust. In SKIRT, the dust of the molecular birth clouds mentioned is treated separately from the diffuse ISM dust.

Schulz et al. (2020) used the output of the TNG50 simulation and suggest a redshift-dependent systematic shift toward lower βFUV with increasing redshift modeled by adding a component βz(z) = 0.142z  −  0.081 to the βFUV value at z = 0. The trend was calibrated up to z = 4. If we extrapolate their trend out to the redshift of our galaxies, the corresponding locus in Fig. 9 would appear too blue when compared to most of our galaxies. In order to extend the redshift range, we computed βz = β − βOverzieret al.(2010) because Overzier et al. (2010) is the reference at z = 0 used by Schulz et al. (2020).

An evolution in redshift of this effect is observationally confirmed. One of the main parameters that can be at the origin of the shift is likely to be the evolution of the stellar population, as mentioned above. So, even though using the redshift as the independent variable might be a simple solution, a more realistic and physical approach should be related to the stellar populations themselves. Schulz et al. (2020) present, in their Fig. 7 (middle), the variations of the intrinsic UV slopes β0 of the galaxy stellar population (i.e., before the attenuation by dust) against their mass weighted mean stellar population ages. This figure shows that β0 correlates with the stellar population age: the younger the stellar population, the lower β0 is. From Fig. 7 (middle) from Schulz et al. (2020), we measured the mean and uncertainties over each axis of the intrinsic UV–slopes β0 and of the mass-weighted mean stellar population ages (age_Mstar) for each redshift range. CIGALE provides mass-weighted mean stellar population ages for each of the galaxies and for our two best models defined above (i.e., delayed SFH with τ = 500 Myr plus DL2014 and PL+OT_MBB for the IR emission). These values are plotted in Fig. 10. The average in the four bins of age_Mstar from Schulz et al. (2020) and our data points seem to follow a linear relation. Although the points from Schulz et al. (2020) are not fully compatible with ours, there is a general trend. Statistical tests with the LINMIX library, which uses a hierarchical Bayesian approach for the linear regression with an error in both X and Y (Kelly 2007), suggest that the trend is highly significant (see correlation coefficient in Fig. 10), given the number of points used. The equations in Fig. 10 allow one to quantify this systematic shift along βFUV at large redshift.

thumbnail Fig. 10.

Linear relation between the shift of the IRX-βFUV relation at z = 0 (namely that of Overzier et al. 2010) vs. the mass-weighted age in years derived from the analysis of the galaxies studied here. Red symbols correspond to ALMA-detected objects, while blue objects correspond to upper limits from ALMA. We note, however, that these upper limits apply to ALMA, but not to the parameters presented in this figure: ageMstar and βz.

4.4. The dust formation rate diagram (DFRD)

Burgarella et al. (2020) and Nanni et al. (2020) built a diagram to follow the evolution of the dust mass called the dust formation rate diagram (DFRD) which shows the specific dust mass (Mdust/Mstar = sMdust) versus the specific star formation rate (SFR/Mstar = sSFR). The specific dust mass was already identified by Calura et al. (2017) as a quantity that represents a true measure of how much dust per unit stellar mass survives the various destruction processes in galaxies. However, this is also a quantity that allows one to quantify the various dust formation processes. The interpretation of this DFRD is that galaxies would follow a dust cycle where they start to build their dust grains at high sSFR leading to a fast rise in sMdust. After this phase, the galaxies would reach a maximum in sMdust before losing their dust grains and going down in the DFRD. At the end of this dust cycle, Burgarella et al. (2020) and Nanni et al. (2020) found that about 80% of the mass fraction of the total baryons had been removed by the outflow and the rest had mainly been destroyed by supernovae (SNe). However, we still have difficulty understanding the evolutionary status of galaxies with very high sMdust > 0.01. Calura et al. (2017) propose that spiral galaxies are characterized by a nearly constant sMdust as a function of the stellar mass and cosmic time, whereas proto-spheroids present an early steep increase of the sMdust, which stops at a maximal value and decreases in the latest stages.

We built this same diagram from the large galaxy sample studied in this paper and for our two main models (delayed SFH with τ = 500 Myr plus DL2014 and PL+OT_MBB for the IR emission). When using the module DL2014 in CIGALE, the dust mass is provided with the models. We note that the models of DL2014 assume an optically thin dust. The CIGALE module casey2012_OT, that is PL+OT_MBB, does not give Mdust, and we need to compute it.

The IR SED

S ν ∝ ν β RJ B ν ( T dust ) $$ \begin{aligned} S_\nu \propto \nu ^{\beta _{\rm RJ}} B_{\nu }(T_{\rm dust})\end{aligned} $$

was computed with a modified blackbody

B ν ( T dust ) = 2 h c 2 ν 3 e h ν k T dust − 1 $$ \begin{aligned} B_{\nu }(T_{\rm dust}) = \frac{2h}{c^2}\frac{\nu ^3}{e\ ^{\frac{h\nu }{kT_{\rm dust}}}-1}\end{aligned} $$

and the dust mass was derived with the following formula:

M dust = L ν 4 π κ ν B ν ( T dust ) $$ \begin{aligned} M_{\rm dust}\ =\ \frac{L_{\nu }}{4\pi \kappa _{\nu } B_{\nu }(T_{\rm dust})}\end{aligned} $$

where ν = c/λ200 μm, h is the Planck constant, and c is the speed of light. We note that the emission of the CMB was neglected given the dust temperatures T ≳ 40 K found in this paper (da Cunha et al. 2013).

In the models for DL2014, the dust composition (Weingartner & Draine 2001) corresponds to 30% of graphite and 70% of silicate (MgFeSiO4). We first adopted the constant given in Draine (2003) at λ = 200 μm: κ0 = 0.637 m2 kg−1 because we also assumed λ0 = 200 μm, κν = κ0 (λ0 /λ) β  = κ0. When comparing Mdust(DL2014) to Mdust(PL+OT_MBB), we found the latter to be systematically smaller by a factor of 0.37 ± 0.01. This was already noticed in several works (Magdis et al. 2013; Santini et al. 2014; Bianchi et al. 2019; Magnelli et al. 2012). The origin of this shift might be found in the single-temperature models that are unable to account for the wide range in the temperature of dust grains that are exposed to different intensities of the interstellar radiation field (see Liang et al. 2019, for a very pedagogical description).

Figure 11 presents the DFRD with data points color-coded with the parameters agemain. Only individual results obtained using DL2014 are shown. However, the trends are also presented for other dust emissions (PL+OT_MBB with various assumptions on κ0, as explained below). The lines are based on a fit with only the ALMA-detected objects. We see an age sequence from right to left, that is with decreasing sSFR.

thumbnail Fig. 11.

Comparison of observed DFRDs with models. Top: DFRDs color-coded with log10agemain. The individual symbols were computed assuming a PL+OT_MBB emission. We clearly observe an age sequence from right to left. The lines show the trend assuming different dust emissions and different values for κ0. The top one assumes DL2014 dust emission. The second one from the top was derived when using an optically thin modified blackbody with κ0 = 0.637 m2 kg−1 (PL+OT_MBB). We note that the factor 0.37 needed to match DL2014 to PL+OT_MBB emission was not applied. The third one from the top assumes optically thin modified blackbodies with κ0 corresponding to SNe dust mass absorption coefficients. It corresponds to the situation before the SNe reverse shock destruction, and the bottom line represents after the reverse shock destruction. Both from Hirashita et al. (2017). Bottom: same models as in as (a). The same objects plotted in (a) are color-coded in Mstar. A sample of DSFGs (Table 8), also color-coded in Mstar, was added to the plot with the code as follows: the upper triangle shows the maximum value, and the lower triangle shows the minimum value. These high redshift DSFGs are found on the same sequence as the other objects. However, these DSFGs have stellar masses larger than the underlying galaxy population that we study. It is important to note that the physical parameters (especially Mstar because the IR emission is dominant) of these high redshift DSFGs have very large uncertainties and their position location in the diagram can almost cover the entire plot. Better estimates of these parameters coming from JWST would help.

We note that κ0 = 0.637 m2 kg−1 corresponds to the dust opacity in the models of DL2014. In the early universe, most if not all of the dust could only be produced by SNe (e.g., Burgarella et al. 2020; Nanni et al. 2020). Hirashita et al. (2017) propose values for κ0 for dust condensed in SNe before reverse shock destruction, κ158 μm = 0.557 m2 kg−1, and for dust ejected from SNe after reverse shock destruction, κ158 μm = 0.894 m2 kg−1. The dust composition and grain size distribution assumed in Hirashita et al. (2017) are from Nozawa et al. (2003). After correction to get these dust mass absorption coefficients at λ = 200 μm using κν  =  κ0 (λ0 /λ) βRJ with βRJ from the SED fitting (Table 4), we obtained κ0 = 0.45 and 0.72 m2 kg−1. This means that, when adopting the SNe value for κ0, the dust masses would be roughly of the same order as PL+OT_MBB with κ0 = 0.637 m2 kg−1. This is certainly a crucial point: given the age of these galaxies, it is very likely that the dust grains have been produced by SNe. So we could wonder how to reconcile this assumption with the fact that the present models are far from being able to reproduce diagnostic diagrams such as Fig. 11. This point is also debated from a laboratory standpoint. For instance, Fanciullo et al. (2020), Ysard et al. (2019) suggest that current dust masses are overestimated by up to a factor of 10–20 or 2–5, depending on the assumptions on grain structure (porous or compact, respectively). Laboratory measurements of dust analogs show that FIR opacities, that is to say mass absorption coefficients, are usually higher than the values used in models and that they depend on several parameters, including temperature, composition, shape, and morphology, for instance. This would mean that dust mass estimates may be overestimated. The properties (e.g., dust composition and grain size) are still unknown for galaxies at a very large redshift. This is important when deriving the dust mass (e.g., Ysard et al. 2019; Hirashita et al. 2017; Inoue et al. 2020).

One interesting point, which is beyond of the scope of this work, is related to the nature of objects that are usually considered to belong to another class: the high redshift dusty star-forming galaxies (DSFG) such as ADFS-27 (Riechers et al. 2021) or HFLS3 (Riechers et al. 2013). The very massive and very dusty galaxies have been found in the high redshift universe thanks to wide-field far-IR and submm observations with, for example, the South Pole Telescope or the European Space Agency’s Herschel Space Observatory. High redshift DSFGs are at the high mass end, which suggests that they are not similar to the early phases of the SFGs analyzed here. In an attempt to check their location in the DFRD, we added the DFSGs from Table 8 to Fig. 12. All of them fall on the same sequence identified for our studied sample; however, they do not show the same stellar mass as our sample. We can wonder whether the stellar mass of these high redshift DSFGs are badly estimated or whether they follow a different path. It should be noticed that CIGALE with similar assumptions on the SFH and on the dust emission has been used for some of these objects, such as ADFS-27.

thumbnail Fig. 12.

DFRD for models with dust masses computed by CIGALE with DL2014 dust models. Only models with a delayed star formation history are plotted. A delayed star formation history only was consistently used in the SED fitting with CIGALE to estimate sMdust and sSFR. The models shown do not have any grain growth in the ISM. We also show the other models in light gray to show which parameter space the models cover. On this figure, we superimposed the three lines that show the respective linear fits as in Fig. 11. The models shown are from Burgarella et al. (2020) and Nanni et al. (2020). This figure shows that most of the objects (ALMA detections and ALMA upper limits) are approximately located at the same positions as the models. Only HZ9 is found at larger sSFR and sMdust. If confirmed at this location in the DFRD, these objects would be very difficult to explain with these models.

Table 8.

Physical parameters for the sample of high redshift dusty star-forming galaxies.

In Fig. 12, we compare the entire sample with the models built in Burgarella et al. (2020) and in Nanni et al. (2020). The models are fully explained in Nanni et al. (2020). In brief, we performed the calculations for the metal evolution using the One-zone Model for the Evolution of GAlaxies (OMEGA) code (Côté et al. 2017). We assumed the metal yields for type II SNe from Kobayashi et al. (2006) computed up to 40 M⊙, and from the FRUITY database for low-mass stars with M > 1.3M⊙ evolving through the thermally pulsing AGB phase and developing stellar winds (Cristallo et al. 2011, 2015; Piersanti et al. 2013). The yields for population III stars are from Heger & Woosley (2010) and are limited to the mass range 10 < M/M⊙ < 30 in OMEGA. Dust removal from the galaxy through galactic outflow follows a rate proportional to the SFR through the “mass-loading factor”: ML × SFR. This assumes that the galactic outflow is generated by the feedback of stars on the gas in the ISM (e.g., Murray et al. 2005). Two kinds of IMFs were tested: top-heavy IMFs and a Chabrier IMF. We assume that a fraction of silicates (olivine and pyroxene), iron, and carbon grains ejected in the ISM are condensed. The best models selected by Burgarella et al. (2020) generally agree with the Alpine sample. Only one object (HZ9 at sSFR ∼ 5 – 6 × 10−8 yr−1) lies at high sSFR. This is consistent with a fast evolution of these young objects, regardless of their sMdust. With the present data, we cannot rule out that the trend could keep rising to the right of the diagram. Galaxies have to start building dust grains early in their lifetime, at log10 sSFR ≥ 10−8 − 10−7 yr−1, to the right of the DFRD. It therefore seems important to measure dust and stellar masses for galaxies at high sSFR and young ages to locate the locus of this rising sMdust and thus bring some observational data to constrain the models.

5. Conclusions

This paper studies a sample of star-forming galaxies observed with ALMA in the redshift range 4.5 ≲ z ≲ 6.2. Some of them, detected in the continuum with ALMA, were used to build a composite IR SED covering a large wavelength range from about 15 μm to 218 μm. Using CIGALE, we modeled the IR composite template and derived a set of physical properties. Building on the assumption that this IR composite template is valid for an ensemble of galaxies selected in the same way, such as ALPINE’s SFGs and LBGs, we use it to fit the SEDs of the ALPINE sample plus the sample of LBGs from Burgarella et al. (2020) and Nanni et al. (2020) from the far-UV to the submm ranges.

The following results were found:

  • We built a unique z > 4 IR SED. It is compatible with the ALPINE sample completed with galaxies from Burgarella et al. (2020) plus the z ∼ 4.5 and z ∼ 5.5 stacks from Bethermin et al. (2020). We provide a table to the readers with the raw composite IR SED and the modeled SEDs.

  • Except for an SFH based on a delayed with a final burst history, all other SFHs (constant, delayed with τmain = 500 Myr, and delayed with several τmain) are in agreement with the data. For simplicity, we selected the delayed SFH with τmain = 500 Myr.

  • We checked the position of the sample in the SFR versus Mstar diagram. The objects in the sample follow the same trend as those previously derived in papers at the same redshift.

  • When comparing the position of the sample with the evolution in redshift of the IRX versus Mstar relation found by Bogdanoska & Burgarella (2020) at z ∼ 4.5 and at z ∼ 5.5, we found a reasonable agreement. This agreement especially holds at z ∼ 5.5. However, the absolute level of this flat relation must be reevaluated using detections of low-mass galaxies.

  • The sample studied in this paper was placed in the IRX versus βFUV diagram. An evolution in age is observed with younger galaxies having bluer βFUV.

  • Moreover, our sample was found to be shifted to bluer βFUV with respect to a z = 0 reference. We modeled this evolution with the mass-weighted age and provide an equation that allows one to predict the position of the sequence as a function of the redshift.

  • We plotted our galaxies in the DFRD (Mdust/Mstar versus SFR/Mstar, DFRD). The objects form a sequence that can be described by evolution in various time-dependent parameters, and more specifically the age of the stellar populations and the stellar mass. This suggests that the sequence is due to the evolution of the galaxies over cosmic time.

  • The models built by Burgarella et al. (2020), Nanni et al. (2020) are generally in agreement with the data, although the sample is much larger than in these previous papers.

  • High redshift, high mass DSFGs over-plotted in the DFRD are found on the same sequence followed by our sample. However, the stellar mass estimated for these high redshift DSFGs are larger than the one found at the same locus for our sample. Uncertainties on the stellar mass of these objects are notoriously high because the stellar emission is faint, but the disagreement seems large. Better estimates of their stellar masses should be available if they are observed with JWST.

1

The combined ALPINE plus Burgarella et al. (2020) and Nanni et al. (2020) samples studied in this paper are referred to as our sample hereafter.

2

By comparing the results of the fits with and without [CII]158 μm, we get the following: Mstar([w/CII]) / Mstar([w/CII]) = 0.90 ± 0.27, SFR([w/CII]) / SFR([w/CII]) = 1.22 ± 0.32, Ldust([w/CII]) / Ldust([w/CII]) = 1.22 ± 0.27, and LFUV([w/CII]) / LFUV([w/CII]) = 1.02 ± 0.02 with no significant difference for the selected dust emission type.

3

We only present the results using the reference SFH derived in the present section, i.e., delayed SFR(t) = t/ τ main 2 $ t/\tau_{\rm main}^2 $ exp(t/τ) with τmain = 500 Myr, but with the following two options for dust emission: Draine et al. (2014) and a power law in a mid-IR and optically thin blackbody.

Acknowledgments

This program receives funding from the CNRS national program Cosmology and Galaxies. D.R. acknowledges support from the National Science Foundation under grant No. AST-1910107. D.R. also acknowledges support from the Alexander von Humboldt Foundation through a Humboldt Research Fellowship for Experienced Researchers. M.T. acknowledges the support from grant PRIN MIUR 2017 20173ML3WW 001. A.N. acknowledges support from the Narodowe Centrum Nauki (UMO-2018/30/E/ST9/00082 and UMO-2020/38/E/ST9/00077). G.C.J. acknowledges ERC Advanced Grant 695671 “QUENCH” and support by the Science and Technology Facilities Council (STFC). Y.F. further acknowledges support from NAOJ ALMA Scientific Research Grant number 2020-16B. M.R. acknowledges support from the Narodowe Centrum Nauki (UMO-2020/38/E/ST9/00077). M.B. gratefully acknowledges support by the ANID BASAL project FB210003 and from the FONDECYT regular grant 1211000. E.I. acknowledges partial support from FONDECYT through grant N° 1171710.

References

  1. Álvarez-Márquez, J., Burgarella, D., Heinis, S., et al. 2016, A&A, 587, A122 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  2. Álvarez-Márquez, J., Burgarella, D., Buat, V., Ilbert, O., & Pérez-González, P. G. 2019, A&A, 630, A153 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  3. Arata, S., Yajima, H., Nagamine, K., Abe, M., & Khochfar, S. 2020, MNRAS, 498, 5541 [NASA ADS] [CrossRef] [Google Scholar]
  4. Aretxaga, I., Wilson, G. W., Aguilar, E., et al. 2011, MNRAS, 415, 3831 [Google Scholar]
  5. Baes, M., Verstappen, J., De Looze, I., et al. 2011, ApJS, 196, 22 [Google Scholar]
  6. Bakx, T. J. L. C., Sommovigo, L., Carniani, S., et al. 2021, MNRAS, 508, L58 [NASA ADS] [CrossRef] [Google Scholar]
  7. Bellstedt, S., Robotham, A. S. G., Driver, S. P., et al. 2021, MNRAS, 503, 3309 [NASA ADS] [CrossRef] [Google Scholar]
  8. Bendo, G. J., Joseph, R. D., Wells, M., et al. 2003, AJ, 125, 2361 [NASA ADS] [CrossRef] [Google Scholar]
  9. Bernhard, E., Béthermin, M., Sargent, M., et al. 2014, MNRAS, 442, 509 [NASA ADS] [CrossRef] [Google Scholar]
  10. Bethermin, M., Fudamoto, Y., Ginolfi, M., et al. 2020, A&A, 643, A2 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  11. Bianchi, S., Casasola, V., Baes, M., et al. 2019, A&A, 631, A102 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  12. Bogdanoska, J., & Burgarella, D. 2020, MNRAS, 496, 5341 [NASA ADS] [CrossRef] [Google Scholar]
  13. Boquien, M., Calzetti, D., Kennicutt, R., et al. 2009, ApJ, 706, 553 [NASA ADS] [CrossRef] [Google Scholar]
  14. Boquien, M., Burgarella, D., Roehlly, Y., et al. 2019, A&A, 622, A103 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  15. Bouwens, R. J., Aravena, M., Decarli, R., et al. 2016, ApJ, 833, 72 [NASA ADS] [CrossRef] [Google Scholar]
  16. Bruzual, G., & Charlot, S. 2003, MNRAS, 344, 1000 [NASA ADS] [CrossRef] [Google Scholar]
  17. Buat, V., Takeuchi, T. T., Burgarella, D., Giovannoli, E., & Murata, K. L. 2009, A&A, 507, 693 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  18. Burgarella, D., Buat, V., & Iglesias-Páramo, J. 2005, MNRAS, 360, 1413 [NASA ADS] [CrossRef] [Google Scholar]
  19. Burgarella, D., Buat, V., Gruppioni, C., et al. 2013, A&A, 554, A70 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  20. Burgarella, D., Nanni, A., Hirashita, H., et al. 2020, A&A, 637, A32 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  21. Calura, F., Pozzi, F., Cresci, G., et al. 2017, MNRAS, 465, 54 [Google Scholar]
  22. Calzetti, D., Kinney, A. L., & Storchi-Bergmann, T. 1994, ApJ, 429, 582 [Google Scholar]
  23. Calzetti, D., Armus, L., Bohlin, R. C., et al. 2000, ApJ, 533, 682 [NASA ADS] [CrossRef] [Google Scholar]
  24. Capak, P. L., Carilli, C., Jones, G., et al. 2015, Nature, 522, 455 [NASA ADS] [CrossRef] [Google Scholar]
  25. Carvajal, R., Bauer, F. E., Bouwens, R. J., et al. 2020, A&A, 633, A160 [EDP Sciences] [Google Scholar]
  26. Casey, C. M. 2012, MNRAS, 425, 3094 [Google Scholar]
  27. Casey, C. M., Chen, C.-C., Cowie, L. L., et al. 2013, MNRAS, 436, 1919 [NASA ADS] [CrossRef] [Google Scholar]
  28. Castellano, M., Sommariva, V., Fontana, A., et al. 2014, A&A, 566, A19 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  29. Chabrier, G. 2003, PASP, 115, 763 [Google Scholar]
  30. Ciesla, L., Boquien, M., Boselli, A., et al. 2014, A&A, 565, A128 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  31. Cooray, A., Calanog, J., Wardlow, J. L., et al. 2014, ApJ, 790, 40 [NASA ADS] [CrossRef] [Google Scholar]
  32. Cormier, D., Abel, N. P., Hony, S., et al. 2019, A&A, 626, A23 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  33. Cortese, L., Boselli, A., Franzetti, P., et al. 2008, MNRAS, 386, 1157 [NASA ADS] [CrossRef] [Google Scholar]
  34. Côté, B., O’Shea, W., Ritter, C., Herwig, F., & Venn, K. A. 2017, ApJ, 835, 128 [CrossRef] [Google Scholar]
  35. Cousin, M., Buat, V., Lagache, G., & Bethermin, M. 2019, A&A, 627, A132 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  36. Cristallo, S., Piersanti, L., Straniero, O., et al. 2011, ApJS, 197, 17 [NASA ADS] [CrossRef] [Google Scholar]
  37. Cristallo, S., Straniero, O., Piersanti, L., & Gobrecht, D. 2015, ApJS, 219, 40 [Google Scholar]
  38. Cucciati, O., Tresse, L., Ilbert, O., et al. 2012, A&A, 539, A31 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  39. Cunningham, D. J. M., Chapman, S. C., Aravena, M., et al. 2020, MNRAS, 494, 4090 [NASA ADS] [CrossRef] [Google Scholar]
  40. da Cunha, E., Groves, B., Walter, F., et al. 2013, ApJ, 766, 13 [Google Scholar]
  41. De Breuck, C., Weiß, A., Béthermin, M., et al. 2019, A&A, 631, A167 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  42. Delvecchio, I., Daddi, E., Sargent, M. T., et al. 2021, A&A, 647, A123 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  43. Donnari, M., Pillepich, A., Nelson, D., et al. 2019, MNRAS, 485, 4817 [Google Scholar]
  44. Douglas, L. S., Bremer, M. N., Lehnert, M. D., Stanway, E. R., & Milvang-Jensen, B. 2010, MNRAS, 409, 1155 [NASA ADS] [CrossRef] [Google Scholar]
  45. Draine, B. T. 2003, ARA&A, 41, 241 [NASA ADS] [CrossRef] [Google Scholar]
  46. Draine, B. T., Aniano, G., Krause, O., et al. 2014, ApJ, 780, 172 [Google Scholar]
  47. Faisst, A. L., Capak, P. L., Yan, L., et al. 2017, ApJ, 847, 21 [Google Scholar]
  48. Faisst, A. L., Schaerer, D., Lemaux, B. C., et al. 2020a, ApJS, 247, 61 [NASA ADS] [CrossRef] [Google Scholar]
  49. Faisst, A. L., Fudamoto, Y., Oesch, P. A., et al. 2020b, MNRAS, 498, 4192 [NASA ADS] [CrossRef] [Google Scholar]
  50. Fanciullo, L., Kemper, F., Scicluna, P., Dharmawardena, T. E., & Srinivasan, S. 2020, MNRAS, 499, 4666 [NASA ADS] [CrossRef] [Google Scholar]
  51. Fernández-Ontiveros, J. A., Spinoglio, L., Pereira-Santaella, M., et al. 2016, ApJS, 226, 19 [CrossRef] [Google Scholar]
  52. Fudamoto, Y., Oesch, P. A., Faisst, A., et al. 2020, A&A, 643, A4 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  53. Fujimoto, S., Silverman, J. D., Bethermin, M., et al. 2020, ApJ, 900, 1 [Google Scholar]
  54. Galametz, M., Kennicutt, R. C., Albrecht, M., et al. 2012, MNRAS, 425, 763 [NASA ADS] [CrossRef] [Google Scholar]
  55. Giacconi, R., Zirm, A., Wang, J., et al. 2002, VizieR Online Data Catalog: II/13 [Google Scholar]
  56. Gordon, K. D., Clayton, G. C., Misselt, K. A., et al. 2003, ApJ, 594, 279 [NASA ADS] [CrossRef] [Google Scholar]
  57. Grogin, N. A., Kocevski, D. D., Faber, S. M., et al. 2011, ApJS, 197, 35 [NASA ADS] [CrossRef] [Google Scholar]
  58. Gruppioni, C., Pozzi, F., Rodighiero, G., et al. 2013, MNRAS, 432, 23 [Google Scholar]
  59. Harikane, Y., Ouchi, M., Yuma, S., et al. 2014, ApJ, 794, 129 [NASA ADS] [CrossRef] [Google Scholar]
  60. Harikane, Y., Ouchi, M., Inoue, A. K., et al. 2020, ApJ, 896, 93 [Google Scholar]
  61. Hashimoto, T., Laporte, N., Mawatari, K., et al. 2018, Nature, 557, 392 [NASA ADS] [CrossRef] [Google Scholar]
  62. Hashimoto, T., Inoue, A. K., Mawatari, K., et al. 2019, PASJ, 71, 71 [Google Scholar]
  63. Heger, A., & Woosley, S. E. 2010, ApJ, 724, 341 [Google Scholar]
  64. Helou, G., Soifer, B. T., & Rowan-Robinson, M. 1985, ApJ, 298, L7 [Google Scholar]
  65. Hirashita, H., Burgarella, D., & Bouwens, R. J. 2017, MNRAS, 472, 4587 [NASA ADS] [CrossRef] [Google Scholar]
  66. Hunt, L. K., Testi, L., Casasola, V., et al. 2014, A&A, 561, A49 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  67. Inoue, A. K., Tamura, Y., Matsuo, H., et al. 2016, Science, 352, 1559 [NASA ADS] [CrossRef] [Google Scholar]
  68. Inoue, A. K., Hashimoto, T., Chihara, H., & Koike, C. 2020, MNRAS, 495, 1577 [NASA ADS] [CrossRef] [Google Scholar]
  69. Jones, G. C., Maiolino, R., Caselli, P., & Carniani, S. 2020, MNRAS, 498, 4109 [NASA ADS] [CrossRef] [Google Scholar]
  70. Juvela, M., Montillaud, J., Ysard, N., & Lunttila, T. 2013, A&A, 556, A63 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  71. Kass, R. E., & Raftery, A. E. 1995, J. Am. Stat. Assoc., 90, 773 [Google Scholar]
  72. Kelly, B. C. 2007, ApJ, 665, 1489 [Google Scholar]
  73. Khusanova, Y., Bethermin, M., Le Fèvre, O., et al. 2021, A&A, 649, A152 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  74. Kobayashi, C., Umeda, H., Nomoto, K., Tominaga, N., & Ohkubo, T. 2006, ApJ, 653, 1145 [NASA ADS] [CrossRef] [Google Scholar]
  75. Koekemoer, A. M., Faber, S. M., Ferguson, H. C., et al. 2011, ApJS, 197, 36 [NASA ADS] [CrossRef] [Google Scholar]
  76. Komatsu, E., Smith, K. M., Dunkley, J., et al. 2011, ApJS, 192, 18 [Google Scholar]
  77. Kong, X., Charlot, S., Brinchmann, J., & Fall, S. M. 2004, MNRAS, 349, 769 [NASA ADS] [CrossRef] [Google Scholar]
  78. Koprowski, M. P., Coppin, K. E. K., Geach, J. E., et al. 2020, MNRAS, 492, 4927 [NASA ADS] [CrossRef] [Google Scholar]
  79. Lagache, G., Cousin, M., & Chatzikos, M. 2018, A&A, 609, A130 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  80. Laporte, N., Katz, H., Ellis, R. S., et al. 2019, MNRAS, 487, L81 [Google Scholar]
  81. Le Fèvre, O., Béthermin, M., Faisst, A., et al. 2020, A&A, 643, A1 [Google Scholar]
  82. Leitherer, C., Calzetti, D., & Martins, L. P. 2002, ApJ, 574, 114 [NASA ADS] [CrossRef] [Google Scholar]
  83. Liang, L., Feldmann, R., Kereš, D., et al. 2019, MNRAS, 489, 1397 [NASA ADS] [CrossRef] [Google Scholar]
  84. Lutz, D., Poglitsch, A., Altieri, B., et al. 2011, A&A, 532, A90 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  85. Ma, X., Hopkins, P. F., Faucher-Giguère, C.-A., et al. 2016, MNRAS, 456, 2140 [NASA ADS] [CrossRef] [Google Scholar]
  86. Madau, P., & Dickinson, M. 2014, ARA&A, 52, 415 [Google Scholar]
  87. Madau, P., Ferguson, H. C., Dickinson, M. E., et al. 1996, MNRAS, 283, 1388 [Google Scholar]
  88. Madden, S. C., Rémy-Ruyer, A., Galametz, M., et al. 2013, PASP, 125, 600 [NASA ADS] [CrossRef] [Google Scholar]
  89. Magdis, G. E., Rigopoulou, D., Helou, G., et al. 2013, A&A, 558, A136 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  90. Magnelli, B., Lutz, D., Santini, P., et al. 2012, A&A, 539, A155 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  91. Marrone, D., Spilker, J. S., Hayward, C. C., et al. 2018, Nature, 553, 51 [NASA ADS] [CrossRef] [Google Scholar]
  92. McLure, R. J., Dunlop, J. S., Cullen, F., et al. 2018, MNRAS, 476, 3991 [Google Scholar]
  93. Murray, N., Quataert, E., & Thompson, T. A. 2005, ApJ, 618, 569 [NASA ADS] [CrossRef] [Google Scholar]
  94. Nanni, A., Burgarella, D., Theulé, P., Côté, B., & Hirashita, H. 2020, A&A, 641, A168 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  95. Nelson, D., Pillepich, A., Springel, V., et al. 2018, MNRAS, 475, 624 [Google Scholar]
  96. Nozawa, T., Kozasa, T., Umeda, H., Maeda, K., & Nomoto, K. 2003, ApJ, 598, 785 [NASA ADS] [CrossRef] [Google Scholar]
  97. O’Halloran, B., Galametz, M., Madden, S. C., et al. 2010, A&A, 518, L58 [CrossRef] [EDP Sciences] [Google Scholar]
  98. Oliver, S. J., Bock, J., Altieri, B., et al. 2012, MNRAS, 424, 1614 [NASA ADS] [CrossRef] [Google Scholar]
  99. Ouchi, M., Yamada, T., Kawai, H., & Ohta, K. 1999, ApJ, 517, L19 [NASA ADS] [CrossRef] [Google Scholar]
  100. Overzier, R. A., Heckman, T. M., Schiminovich, D., et al. 2010, ApJ, 710, 979 [NASA ADS] [CrossRef] [Google Scholar]
  101. Pavesi, R., Riechers, D. A., Capak, P. L., et al. 2016, ApJ, 832, 151 [Google Scholar]
  102. Pavesi, R., Riechers, D. A., Faisst, A. L., Stacey, G. J., & Capak, P. L. 2019, ApJ, 882, 168 [Google Scholar]
  103. Pearson, E. A., Eales, S., Dunne, L., et al. 2013, MNRAS, 435, 2753 [NASA ADS] [CrossRef] [Google Scholar]
  104. Pearson, W. J., Wang, L., Hurley, P. D., et al. 2018, A&A, 615, A146 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  105. Piersanti, L., Cristallo, S., & Straniero, O. 2013, ApJ, 774, 98 [NASA ADS] [CrossRef] [Google Scholar]
  106. Pilbratt, G. L., Riedinger, J. R., Passvogel, T., et al. 2010, A&A, 518, L1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  107. Pillepich, A., Nelson, D., Hernquist, L., et al. 2018, MNRAS, 475, 648 [Google Scholar]
  108. Pozzi, F., Calura, F., Fudamoto, Y., et al. 2021, A&A, 653, A84 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  109. Raftery, A. E. 1998, Bayes Factors and BIC: Comment on Weakliem [Google Scholar]
  110. Riechers, D. A., Bradford, C. M., Clements, D. L., et al. 2013, Nature, 496, 329 [Google Scholar]
  111. Riechers, D. A., Hodge, J. A., Pavesi, R., et al. 2020, ApJ, 895, 81 [Google Scholar]
  112. Riechers, D. A., Nayyeri, H., Burgarella, D., et al. 2021, ApJ, 907, 62 [NASA ADS] [CrossRef] [Google Scholar]
  113. Rodighiero, G., Daddi, E., Baronchelli, I., et al. 2011, ApJ, 739, L40 [Google Scholar]
  114. Salim, S., Lee, J. C., Janowiecki, S., et al. 2016, ApJS, 227, 2 [NASA ADS] [CrossRef] [Google Scholar]
  115. Santini, P., Maiolino, R., Magnelli, B., et al. 2014, A&A, 562, A30 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  116. Schaerer, D., Ginolfi, M., Béthermin, M., et al. 2020, A&A, 643, A3 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  117. Schulz, S., Popping, G., Pillepich, A., et al. 2020, MNRAS, 497, 4773 [NASA ADS] [CrossRef] [Google Scholar]
  118. Scoville, N., Aussel, H., Brusa, M., et al. 2007, ApJS, 172, 1 [Google Scholar]
  119. Scoville, N., Sheth, K., Aussel, H., et al. 2016, ApJ, 820, 83 [Google Scholar]
  120. Shapley, A. E., Steidel, C. C., Pettini, M., & Adelberger, K. L. 2003, ApJ, 588, 65 [Google Scholar]
  121. Sherman, S., Jogee, S., Florez, J., et al. 2021, MNRAS, 505, 947 [CrossRef] [Google Scholar]
  122. Strandet, M. L., Weiss, A., De Breuck, C., et al. 2017, ApJ, 842, L15 [Google Scholar]
  123. Speagle, J. S., Steinhardt, C. L., Capak, P. L., & Silverman, J. D. 2014, ApJS, 214, 15 [Google Scholar]
  124. Sugahara, Y., Inoue, A. K., Hashimoto, T., et al. 2021, ApJ, 923, 5 [NASA ADS] [CrossRef] [Google Scholar]
  125. Tabatabaei, F. S., Braine, J., Xilouris, E. M., et al. 2014, A&A, 561, A95 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  126. Takeuchi, T. T., Buat, V., Heinis, S., et al. 2010, A&A, 514, A4 [CrossRef] [EDP Sciences] [Google Scholar]
  127. Tibbs, C. T., Israel, F. P., Laureijs, R. J., et al. 2018, MNRAS, 477, 4968 [CrossRef] [Google Scholar]
  128. Volinsky, C. T. 1997, Ph.D. Thesis, University of Washington [Google Scholar]
  129. Weingartner, J. C., & Draine, B. T. 2001, ApJ, 548, 296 [Google Scholar]
  130. Willott, C. J., Carilli, C. L., Wagg, J., & Wang, R. 2015, ApJ, 807, 180 [NASA ADS] [CrossRef] [Google Scholar]
  131. Witt, A. N., & Gordon, K. D. 2000, ApJ, 528, 799 [Google Scholar]
  132. Xu, C. K., Shupe, D., Buat, V., et al. 2007, ApJS, 173, 432 [NASA ADS] [CrossRef] [Google Scholar]
  133. Ysard, N., Koehler, M., Jimenez-Serra, I., Jones, A. P., & Verstraete, L. 2019, A&A, 631, A88 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  134. Yuan, F.-T., Burgarella, D., Corre, D., et al. 2019, A&A, 631, A123 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]

Appendix A: CIGALE parameters for the initial fits

A.1. Star formation histories (SFHs)

  • *

    A delayed SFH: SFR(t) = t/ τ main 2 $ t/\tau_{\rm main}^2 $ exp(−t/τmain) with τmain in the range 25 to 10000 Myrs and main ages (t or agemain) in the range 2 - 1200 Myrs is the first option. This type of SFH allows for SFRs to increase or decrease, depending on the age of the stellar population.

  • *

    The second option is a delayed SFH and a final young burst with the burst age (ageburst) are in the range 2 to 50 Myrs, but main ages are in the range 100 to 1200 Myrs, and the fraction of burst (fburst) is as follows: 0.0, 0.001, 0.01, 0.10, 0.50 and τburst = 20000Myrs, that is to say it is similar to a constant SFH given the burst ages.

  • *

    Fixed parameters: because it was not possible to estimate τmain with confidence, we also sequentially fixed and tried τmain = 25, 250, and 2500 Myrs both for the delayed and burst runs before comparing the results to the above runs.

A.2. Dust characteristics

  • *

    The models from Draine et al. (2014): because no data are available in the rest-frame near-IR and mid-IR, we could not constrain the mass fraction of PAH (qPAH). However, it is generally accepted (e.g., Douglas et al. 2010; Castellano et al. 2014; Yuan et al. 2019; Bellstedt et al. 2021) that the metallicity of these galaxies are likely subsolar. Ciesla et al. (2014) show a relation between qPAH and the metallicity. Without any other constraint, we arbitrarily chose a low value for the mass fraction of PAH: qPAH = 0.47. The impact of this parameter on the dust mass or luminosity is very low. We used the full range of allowed parameters for the minimum of the distribution of the starlight intensity relative to the local interstellar radiation field, Umin, for α, the power law slope dU/dM = Uα, as well as for γ, the fraction of the dust heated by starlight above the lower cutoff Umin. For the same reason as for τmain, in the SFH parameters, we also sequentially tried fixed values for α = 1.0, 2.0, and 3.0.

  • *

    Modified blackbody plus a power law similar to Casey (2012): Our initial priors are the dust temperature in the range 30K ≤ Tdust ≤ 85K and the emissivity in the RJ part of the SED 0.5 ≤ βRJ ≤ 2.0. We kept the MIR power slope at the default αMIR = 2.0 because no data are available to constrain αMIR. These values are consistent with Faisst et al. (2020b), that is 40 < Tdust[K] < 60 (called TSED, i.e., SED dust temperature in their paper), with a median at 48 K and emissivity indices (βd in their paper) between 1.6 and 2.4 for all galaxies, and a median of 2.0. As before, because it could not be safely estimated with the data available, we tried several fixed values for βRJ= 1.0, 1.5, and 2.0.

  • *

    Other parameters: In addition to these parameters, we selected a Chabrier IMF and a wide range of dust attenuation. We selected CIGALE’s dust attenuation law from the Calzetti and Leitherer module. It is based on the Calzetti et al. (2000) starburst attenuation curve, extended with the Leitherer et al. (2002) curve between the Lyman break and 120 nm. A nebular contribution was also added (see Boquien et al. 2019) where Dλ is the Drude profile, and the last term renormalizes the curve so that E(B-V) remains equal to the input E(B-V) when δ is not 0:

    k λ = k λ starburst × ( λ / λ V ) δ × D λ ) E ( B − V ) δ = 0 E ( B − V ) δ . $$ \begin{aligned}k_\lambda =k_\lambda ^{starburst}\times (\lambda /\lambda _V)^\delta \times D_\lambda )\frac{E(B-V)_{\delta =0}}{E(B-V)_\delta }.\end{aligned} $$

In Boquien et al. (submitted), a variation of the dust attenuation law is studied. However, changes in the shape of the dust attenuation law do not impact the IRX versus AFUV relation for star forming galaxies because it is almost completely independent of the extinction mechanisms (i.e., dust and star geometry, attenuation law (e.g., Witt & Gordon 2000; Cortese et al. 2008, and references therein).

Finally, from these initial fits, we examined and discarded the fits for which the reduced χ ν 2 ≥5.0 $ \chi_\nu^2 \geq 5.0 $ for each of our detected objects. These bad fits only represent a few of them with respect to all of the attempts. For each object in the sample, we computed the mean normalization factor from the CIGALE Bayesian luminosity at 200 μm: L200μm. After applying this normalization, all the observed SEDs started to share the same flux density at λ = 200 μm: fν(200 μm) = 1.0 (see Fig.1) and we could then proceed to the next phase to try and combine all the observed SEDs into a single one.

Appendix B: CIGALE parameters for the final fit

Table B.1.

CIGALE modules and input parameters used for all the fits. BC03 means Bruzual & Charlot (2003), and the Chabrier IMF refers to Chabrier (2003).

Appendix C: Results for the mock analysis performed with CIGALE

thumbnail Fig. C.1.

Results from the mock analysis performed with CIGALE where we compared input parameters (i.e., parameters from the best fit for each object) to the same parameters estimated by CIGALE with the very same priors used for the true analysis (see Boquien et al. 2019). We note that some parameters were estimated from another one via a linear relation, e.g., Mdust from L200μm, which explains why the figures are identical.

Table C.1.

IR composite built from the sample of ALMA-detected objects. S/N = -1.00 and negative errors mean that the data in the given band is an upper limit.

Appendix D: The observed IR composite template

Table D.1.

IR template for each of the models, normalized to 200 μm. Only the first rows are given here.

Appendix E: Physical parameters for each galaxy of the studied sample

Table E.1.

Physical parameters. Only the first rows are shown. For each object in the sample, the table contains the following information for PL+OT_MBB: id, redshift, βFUV, sSFR, sM_dust, M_dust, A_FUV, IRX, age_main, L_dust, L_FUV, SFR, and M_star.

Table E.2.

Physical parameters. Only the first rows are shown. For each object in the sample, the table contains the following information for PL+OT_MBB: id, redshift, βFUV, sSFR, sM_dust, M_dust, A_FUV, IRX, age_main, L_dust, L_FUV, SFR, and M_star.

All Tables

Table 1.

Origins of the data used in this paper.

In the text
Table 2.

Parameters of the dust emission.

In the text
Table 3.

Degeneracy of βRJ and Tdust.

In the text
Table 4.

Main relevant physical parameters derived by fitting the IR template with the various assumptions of the dust emission.

In the text
Table 5.

Comparison of the reference model (SFH delayed: t/τ2 exp(−t/τ) with τ = 500 Myr) to all the others using the BIC.

In the text
Table 6.

Summary of the power of two for correlation coefficients, r2, from the mock analysis on the main physical parameters performed by CIGALE.

In the text
Table 7.

Mean value and standard deviations of the physical parameters derived from the fit analysis of the individual objects in our sample.

In the text
Table 8.

Physical parameters for the sample of high redshift dusty star-forming galaxies.

In the text
Table B.1.

CIGALE modules and input parameters used for all the fits. BC03 means Bruzual & Charlot (2003), and the Chabrier IMF refers to Chabrier (2003).

In the text
Table C.1.

IR composite built from the sample of ALMA-detected objects. S/N = -1.00 and negative errors mean that the data in the given band is an upper limit.

In the text
Table D.1.

IR template for each of the models, normalized to 200 μm. Only the first rows are given here.

In the text
Table E.1.

Physical parameters. Only the first rows are shown. For each object in the sample, the table contains the following information for PL+OT_MBB: id, redshift, βFUV, sSFR, sM_dust, M_dust, A_FUV, IRX, age_main, L_dust, L_FUV, SFR, and M_star.

In the text
Table E.2.

Physical parameters. Only the first rows are shown. For each object in the sample, the table contains the following information for PL+OT_MBB: id, redshift, βFUV, sSFR, sM_dust, M_dust, A_FUV, IRX, age_main, L_dust, L_FUV, SFR, and M_star.

In the text

All Figures

thumbnail Fig. 1.

Flow chart of the process followed to build the composite IR SED for the Hiz-SFGs sample. The various phases of this process are shown on the right side of the figures. They are also described in more detail in the text (see Sect. 3) In the first phase, we selected only ALMA-detected objects to build the first IR template (IRV1). Next, we added the data from the z ∼ 4.5 and ∼5.5 stacks from Bethermin et al. (2020) to build the second and final IR template (IRV2). Adding UV and optical data to IRV2, we fit the entire (detected and upper limits) sample.
In the text
thumbnail Fig. 2.

IR combined SED built from the ALMA-detected objects from the ALPINE sample, from Burgarella et al. (2020), and with the two stacks from Bethermin et al. (2020) and the associated uncertainty regions. We note that the four left-most points for the stacks are upper limits, as shown by the error bars reaching the bottom of the figure. The SEDs corresponding to optically thin modified black bodies are superimposed to the observed SEDs. Qualitatively, we can see that the IR-combined SED seems to be in agreement with greenish modified blackbody SEDs, that is Tdust ∼ 50–60 K. A qualitative analysis is performed later in this paper.
In the text
thumbnail Fig. 3.

Degeneracies in Tdust and βRJ. The knife jacking method allowed us to estimate whether the characteristics of the dust emission, i.e., Tdust and βRJ, vary when removing and adding back individual galaxies from the sample used to build the composite template. The ellipses show the uncertainties. The outlier point at βRJ ∼ 1.1 and Tdust ∼ 51 K was obtained when taking off HZ10 at z = 5.659 with five pieces of data in the far-IR+submm range. Only when removing it did the average locus move to slightly larger values of βRJ, but with larger uncertainties.
In the text
thumbnail Fig. 4.

Comparison of the various fits of the IR composite SED built with data from Burgarella et al. (2020), from ALPINE, and from the z ∼ 4.5 and z ∼ 5.5 stacks from Bethermin et al. (2020). Top: fit with a modified general blackbody and power law in the mid-IR as in Casey (2012). Middle: fit with an optically thin modified blackbody and power law in the mid-IR as in Casey (2012). Bottom: Fit with a model from Draine et al. (2014). The fits are all statistically equivalently good with reduced χ2 ≈ 0.3–0.4.
In the text
thumbnail Fig. 5.

Results of the tests on the determination of Tdust and βRJ and their possible degeneracy are shown here. Dots show the evolution of Tdust when βRJ was fixed during the fits. Boxes show the evolution of βRJ when Tdust was fixed during the fits. Finally, the diamond symbol shows the result and uncertainties on both axes when Tdust and βRJ are both free. The color coding indicates the level of the goodness of fit.
In the text
thumbnail Fig. 6.

Result of the ΔBIC test on our sample of Hiz-SFGs. A delayed SFH without a burst and τmain = 500 Myr is labeled as “tau500”. With an additional burst at the end of the SFH, it is labeled as “burst”. A constant SFH without a burst and τmain = 20 Gyrs is labeled as “tau20000”. And when several τmain could be selected in the SED fitting, it is labeled as “multitau” in the legend. Top: ΔBIC test that compares the influence of the DL2014 model and the PL+OT_MBB dust emissions in building the IR template. Center: ΔBIC test on the SFH assuming the PL+OT_MBB for the IR template. Bottom: ΔBIC test on the SFH assuming the DL2014 model for the IR template. The color band allows one to interpret the results of the evidence (ΔBIC) against the model with the higher BIC: red means “faint evidence”, orange means “positive evidence”, and green mean “strong evidence”. We do not see any strong evidence that DL2014 or PL+OT_MBB are better for fitting the data. An SFH that includes a burst is positively ruled out, while a delayed SFH with τ = 500 Myr is weekly favored.
In the text
thumbnail Fig. 7.

In the log10 (SFR) vs. log10 (Mstar) diagram, the main part of the sample is found within the limits for the fits of the main sequence at z = 5.0 by Speagle et al. (2014; green shading), by Pearson et al. (2018) at z ∼ 5.2 (purple shading), and by Faisst et al. (2020a), who found that the galaxies are in agreement with Speagle et al. (2014). Is is important to note, however, that some of our objects are at redshifts larger than the 4.5–5.6 ALPINE sample (see color code of the markers). For Speagle et al. (2014), we used the “mixed” (preferred fit)” function (as defined by Speagle et al. 2014). The results from the two fits with DL2014 and PL+OT_MBB are presented. The two types of dust emission do not significantly modify the location of the points in the diagram. Detections are shown as dots and boxes and upper limits are shown as downward- and upward-pointing triangles. The uncertainty range for upper limits extends to the bottom of the plot. In addition to having mainly upper limits, at log10 (Mstar) < 9.5, the sample is very likely incomplete which means that it is difficult to estimate a trend from these data over the entire mass range. We also added the objects and stacks from Khusanova et al. (2021) with open markers. Finally, the selection of TNG50 galaxies used in Schulz et al. (2020) is also provided (red-shaded area).
In the text
thumbnail Fig. 8.

AFUV − Mstar diagram at z ∼ 4.5 (left) and z ∼ 5.5 (right). The gray areas correspond to the expected relation at z = 4.5 and 5.5 from Bogdanoska & Burgarella (2020). This relation was formed by a broken line, which is flat at log10Mstar ≤ 9.8 and rises at log10 (Mstar) > 9.8. The ALPINE data are very dispersed at z ∼ 4.5, while this flatness is supported by the data at z ∼ 5.5. The conversion from IRX to AFUV is from Burgarella et al. (2005): AFUV = −0.028[log10 (IRX)]3 + 0.392 [log10 (IRX)]2 + 1.094 [log10 (IRX)] + 0.5.
In the text
thumbnail Fig. 9.

IRX – βCFUV diagram. IRX values were estimated with CIGALE and βFUV were fitted on the data directly. They are not model-dependant. The black continuous line corresponds to the original Calzetti law, under the assumption that the underlying dust curve follows the Calzetti et al. (2000) attenuation. The dashed line to the predicted law assumes the SMC extinction law (e.g., Gordon et al. 2003). Both are from McLure et al. (2018). The color of the symbols are related to the axis to the right of the figure, the age of the stellar population. Both results with dust emission, DL2014 and PL+OT_MBB, are presented with different symbols with ALMA upper limits and ALMA detections. In the figures, we also plotted the laws from Schulz et al. (2020) at z = 0.0 (blue continuous line), 4.5 (green dotted line), and 5.5 (red dashed line). A comparison with the IRX – βFUV plot in Fudamoto et al. (2020) with the same sample suggests that our (especially ALMA-detected) galaxies extend less to very blue βFUV. This is mainly true for the subsample of galaxies not detected in continuum with ALMA. This is probably due to two effects: first, in Fudamoto et al. (2020), 3σ upper limits are plotted. Another possible effect could be because we used a unique IR composite template in the SED fitting, which reduces the uncertainties on Ldust.
In the text
thumbnail Fig. 10.

Linear relation between the shift of the IRX-βFUV relation at z = 0 (namely that of Overzier et al. 2010) vs. the mass-weighted age in years derived from the analysis of the galaxies studied here. Red symbols correspond to ALMA-detected objects, while blue objects correspond to upper limits from ALMA. We note, however, that these upper limits apply to ALMA, but not to the parameters presented in this figure: ageMstar and βz.
In the text
thumbnail Fig. 11.

Comparison of observed DFRDs with models. Top: DFRDs color-coded with log10agemain. The individual symbols were computed assuming a PL+OT_MBB emission. We clearly observe an age sequence from right to left. The lines show the trend assuming different dust emissions and different values for κ0. The top one assumes DL2014 dust emission. The second one from the top was derived when using an optically thin modified blackbody with κ0 = 0.637 m2 kg−1 (PL+OT_MBB). We note that the factor 0.37 needed to match DL2014 to PL+OT_MBB emission was not applied. The third one from the top assumes optically thin modified blackbodies with κ0 corresponding to SNe dust mass absorption coefficients. It corresponds to the situation before the SNe reverse shock destruction, and the bottom line represents after the reverse shock destruction. Both from Hirashita et al. (2017). Bottom: same models as in as (a). The same objects plotted in (a) are color-coded in Mstar. A sample of DSFGs (Table 8), also color-coded in Mstar, was added to the plot with the code as follows: the upper triangle shows the maximum value, and the lower triangle shows the minimum value. These high redshift DSFGs are found on the same sequence as the other objects. However, these DSFGs have stellar masses larger than the underlying galaxy population that we study. It is important to note that the physical parameters (especially Mstar because the IR emission is dominant) of these high redshift DSFGs have very large uncertainties and their position location in the diagram can almost cover the entire plot. Better estimates of these parameters coming from JWST would help.
In the text
thumbnail Fig. 12.

DFRD for models with dust masses computed by CIGALE with DL2014 dust models. Only models with a delayed star formation history are plotted. A delayed star formation history only was consistently used in the SED fitting with CIGALE to estimate sMdust and sSFR. The models shown do not have any grain growth in the ISM. We also show the other models in light gray to show which parameter space the models cover. On this figure, we superimposed the three lines that show the respective linear fits as in Fig. 11. The models shown are from Burgarella et al. (2020) and Nanni et al. (2020). This figure shows that most of the objects (ALMA detections and ALMA upper limits) are approximately located at the same positions as the models. Only HZ9 is found at larger sSFR and sMdust. If confirmed at this location in the DFRD, these objects would be very difficult to explain with these models.
In the text
thumbnail Fig. C.1.

Results from the mock analysis performed with CIGALE where we compared input parameters (i.e., parameters from the best fit for each object) to the same parameters estimated by CIGALE with the very same priors used for the true analysis (see Boquien et al. 2019). We note that some parameters were estimated from another one via a linear relation, e.g., Mdust from L200μm, which explains why the figures are identical.
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.