The Dust Attenuation of Star-forming Galaxies at $z\sim3$ and Beyond: New Insights from ALMA Observations

TL;DR

This study investigates dust attenuation in high-redshift star-forming galaxies using ALMA data, revealing that dust correction methods valid at z~3 may not apply at higher redshifts, indicating evolution in dust properties.

Contribution

The paper provides new observational evidence of redshift evolution in the IRX-$eta_{ m UV}$ relation beyond z~3 using ALMA data, challenging the universality of local dust correction calibrations.

Findings

IRX-$eta_{ m UV}$ relation at z~3 aligns with local starburst galaxies

Large dispersion observed around the IRX-$eta_{ m UV}$ relation

Evidence for evolution of IRX-$eta_{ m UV}$ at z>3 towards lower IRX values

Abstract

We present results on the dust attenuation of galaxies at redshift $\sim3-6$ by studying the relationship between the UV spectral slope ($\beta_{\rm UV}$) and the infrared excess (IRX; $L_{\rm IR}$/$L_{\rm UV}$) using ALMA far-infrared continuum observations. Our study is based on a sample of 67 massive, star-forming galaxies with a median mass of $M_{\ast}\sim 10^{10.7}\,M_{\rm \odot}$ spanning a redshift range $z=2.6-3.7$ (median $z=3.2$) that were observed with ALMA at $\lambda_{rest}=300\,{\rm \mu m}$. Both the individual ALMA detections (41 sources) and stacks including all galaxies show the IRX-$\beta_{\rm UV}$ relationship at $z\sim3$ is mostly consistent with that of local starburst galaxies on average. However, we find evidence for a large dispersion around the mean relationship by up to $\pm0.5$ dex. Nevertheless, the locally calibrated dust correction factors based on the IRX-$\beta_{\rm UV}$ relation are on average applicable to main-sequence $z\sim3$ galaxies. This does not appear to be the case at even higher redshifts, however. Using public ALMA observations of $z\sim4-6$ galaxies we find evidence for a significant evolution in the IRX-$\beta_{\rm UV}$ and the IRX-$M_{\ast}$ relations beyond $z\sim3$ toward lower IRX values. We discuss several caveats that could affect these results, including the assumed dust temperature. ALMA observations of larger $z>3$ galaxy samples will be required to confirm this intriguing redshift evolution.

Tables

Table 1. Table 1: Results of the Stacking Analysis

${\beta }_{\mathrm{UV}}$	# of sources	log $M_{\ast }$a	log $L_{\mathrm{IR}}$b	logIRX
		[$M_{☉}$]	[$L_{☉}$]
All sources
	14	$10.7$	$11.8$	${0.97}_{-0.12}^{+0.12}$
$-1.5$ to $-1.0$	30	$10.6$	$11.9$	${1.11}_{-0.08}^{+0.07}$
$>-1.0$	23	$10.8$	$12.0$	${1.28}_{-0.09}^{+0.06}$
Non-detection sources
All	26	$10.7$	$11.1$	${0.37}_{-0.09}^{+0.11}$

Equations

$$\log (\mathrm{IRX})=0.72[\log ({\mathrm{M}}_{\ast })-10.35]+1.08$$

Full text

††pubyear: 2017††pagerange: The Dust Attenuation of Star-forming Galaxies at $z\sim 3$ and Beyond: New Insights from ALMA Observations–The Dust Attenuation of Star-forming Galaxies at $z\sim 3$ and Beyond: New Insights from ALMA Observations

1 Introduction

One of the keys to understand the formation and evolution of galaxies is to estimate the star-formation activity of galaxies across cosmic times. Recent optical/near-infrared (NIR) observations successfully measured the star formation rate density (SFRD) of galaxies out to $z\sim 10$ (e.g. Oesch et al., 2013, 2014; Ellis et al., 2013; Bouwens et al., 2015; McLeod et al., 2016), and showed that the SFRD steadily increases toward lower redshift, peaks around $z\sim 2$–3, and then declines by an order of magnitude over the last seven billion years (e.g. Hopkins & Beacom, 2006; Madau & Dickinson, 2014).

However, a major caveat of current studies at $z>3$ in particular is that most measurements of the SFRD rely on the rest-frame UV light from young, massive stars to estimate the star formation rate of galaxies. This UV light is highly sensitive to obscuration by dust particles in the interstellar medium. The absorbed energy in the UV is re-emitted in infrared (IR) wavebands by the heated dust. As a result, a galaxy emits a considerable amount of energy in the IR, which together with the UV luminosity provides a good measurement of the total star-formation of a galaxy (e.g. Elbaz et al., 2011; Symeonidis et al., 2013; Ivison et al., 2016). Without a direct detection of dust emission, however, one requires a well calibrated relation between total energy output from star formation activity and observed energy in UV to recover the total star formation rate (SFR) of galaxies.

The useful quantities that parametrise the amount of obscuration are the UV spectral slope ($\beta_{\rm UV}$) and the infrared excess ({\rm IRX}\equiv$$L_{\rm IR}/$L_{\rm UV}$). $\beta_{\rm UV}$ is defined as a power law form of $f_{\lambda}\propto\lambda^{\beta_{\rm UV}}$ and thus describes galaxy colours in the UV waveband, which are modulated by dust attenuation. Thus one naturally expects a relation between the $\beta_{\rm UV}$ and the fraction of energy emitted by dust in the IR (i.e. ${\rm IRX}$). The IRX–$\beta_{\rm UV}$ relation has been well calibrated using local star bursting galaxies (Meurer et al., 1999; Overzier et al., 2011; Takeuchi et al., 2012; Casey et al., 2014). In particular, the local relation estimated in Meurer et al. (1999) is very widely used, and we will refer to it as M99 in the following.

The locally calibrated IRX–$\beta_{\rm UV}$ relation has been shown to be applicable to the general star-forming galaxy population at redshift $z=0$–2 (e.g. Reddy et al., 2006; Daddi et al., 2009; Pannella et al., 2009; Reddy et al., 2010; Overzier et al., 2011; Sklias et al., 2014). However, it is not clear if the relation is still valid at $z\gtrsim 3$, as the sample of individual IR detections is limited to the most luminous galaxies with extremely high SFR (e.g. Oteo et al., 2013; Riechers et al., 2013).

Despite the power of stacking analyses of high-redshift galaxies using single dish sub/mm telescopes (Heinis et al., 2013, 2014; Pannella et al., 2015; Álvarez-Márquez et al., 2016; Bourne et al., 2017), detailed studies require individual detections of IR emission from normal star-forming galaxies along the star-forming main sequence.

With currently accumulating observations with extremely sensitive sub/mm interferometers such as ALMA and NOEMA, the study of dust emission from high-redshift star forming galaxies is rapidly maturing. These observations hint at a potential evolution of the IRX–$\beta_{\rm UV}$ relation at $z\gtrsim 3$, with high-redshift star forming galaxies appearing to be ’infrared-dark’ (Ouchi et al., 2013; Ota et al., 2014; Capak et al., 2015; Schaerer et al., 2015; Bouwens et al., 2016; Dunlop et al., 2017). However, current studies are still based on small samples only, and the direct detection of IR emission from statistical samples of individual high-redshift galaxies is still missing.

In this paper, to improve our knowledge of the nature of dust emission from $z\gtrsim 3$ star-forming galaxies, we explore the IRX–$\beta_{\rm UV}$ relation at $z\sim 3.2$ by using our ALMA observations of UV selected, massive (M_{\ast}$$\gtrsim 10^{10.7}$$M_{\odot}), star-forming galaxies located on COSMOS field. We also include public ALMA observations of $z\sim 4$–6 galaxies, in order to investigate a possible redshift evolution of the dust attenuation between $z\sim 3$ and $z\sim 6$, within the first 2 Gyr of cosmic history.

This paper is organised as follows: In Section 2, we present our galaxy sample and the data used for the analysis. Section 3 describes the basic measurements performed, before our results are presented in Section 4 and we conclude in Section 5. Throughout this paper, we assume a cosmology with $(\Omega_{m},\Omega_{\Lambda},h_{0})=(0.3,0.7,0.7)$, and a Chabrier initial mass function (Chabrier, 2003) where applicable.

2 Observations

2.1 ALMA Observations of $z\sim 3.2$ Targets

The main data used for this study come from an ALMA Cycle-2 Program (2013.1.00151.S, PI E. Schinnerer), which targets the $\sim\,300\,\rm{\mu m}$ rest-frame continuum of a sample of main-sequence galaxies at $z\sim 3$. Details of the target selection and ALMA observations are presented in Schinnerer et al. (2016), and they are briefly summarised below.

The main purpose of the ALMA program was to observe massive, star-forming galaxies at $z\sim 3$ on the main sequence. In this paper, we only include sources with reliable photometric redshifts in the range 2.6–3.7 in the latest COSMOS15 catalogue (Laigle et al., 2016) and we exclude sources that potentially host an AGN (as discussed in Schinnerer et al., 2016). Our final $z\sim 3$ sample thus consists of 67 galaxies, which have a median redshift of 3.18, and a median stellar mass of ${\rm log(M_{\ast}[M_{\odot}])\sim 10.7}$, with the minimum/maximum stellar mass being ${\rm log(M_{\ast}[M_{\odot}])=9.2/11.5}$. While most sources only have photometric redshifts, 7 galaxies of our sample have high-quality spectroscopic redshift measurements from VIMOS Ultra Deep Survey (VUDS; Le Fèvre et al., 2015).

ALMA observations of these targets were performed to detect the continuum emission at $240\,\rm{GHz}$, which corresponds to rest-frame $\lambda_{\rm rf}\sim\,300\,\rm{\mu m}$ at the redshift of our sample. Targets were observed with typically 38 antennas between December, 25th and 30th, 2015 for an average on-source time of 2 min. The achieved beam-size and rms of the ALMA images are $1.8$′′$\times$1.1*′′(1.7′′$\times$1.1′′*) and $66(77)\,\rm{\mu Jy/beam}$ for 19(48) fields using the NATURAL weighting scheme.

A total of 41 of our 67 targets were detected at more than $>3\,\sigma$ at the counterpart positions of the latest COSMOS Ultra-Vista NIR images (McCracken et al., 2012). Four out of these sources have reliable spectroscopic redshifts from VUDS. For the 26 targets with ALMA non-detections, we compute $3\sigma$ upper limits on the ALMA fluxes that are used throughout the paper.

Stellar masses for our $z\sim 3$ sample are calculated based on multi-wavelength spectral energy distribution (SED) fitting using the MAGPHYS code (da Cunha et al., 2008). In addition to the COSMOS optical+NIR photometry (Laigle et al., 2016), we also exploit all available longer wavelength data in these fits. This includes the Spitzer/MIPS $24\,{\rm\mu m}$ fluxes from Le Floc’h et al. (2009), the Herschel/PACS $100\,\mu m$ fluxes from Lutz et al. (2011), and our ALMA 1.3-mm fluxes (see Schinnerer et al., 2016, for details on the SED fitting).

2.2 $z\geq 4$ Galaxy Samples

While our main analysis is based on the $z\sim 3$ sample from Schinnerer et al. (2016), we also investigate a potential redshift evolution of the IRX–$\beta_{\rm UV}$ relation in section 4.4. To do this, we exploit the published ALMA fluxes of $z\sim 4$–5 galaxies from Scoville et al. (2016), and the IRX–$\beta_{\rm UV}$ data points of $z\sim 5$–6 galaxies from Capak et al. (2015).

ALMA band-6 fluxes of the $z\sim 4-5$ sample are taken from Table 8 of Scoville et al. (2016). This sample consists of 19 galaxies with median $\log(M_{\ast}[M_{\odot}])$ of 10.45 and median redshift of 4.3. For this sample, the observed ALMA fluxes probe rest-frame wavelengths $\lambda_{\rm rf}\sim 240\,{\rm\mu m}$. The UV luminosities and continuum slopes are derived in the same manner as for our $z\sim 3$ sample based on the photometry of the COSMOS15 catalogue. Stellar masses for the $z\sim 4-5$ sources are directly taken from that catalogue.

For further constraints at even higher redshift, we exploit the data set by Capak et al. (2015). That sample consists of 9 galaxies (4 detections and 5 non-detections) with median mass $\log(M_{\ast}[M_{\odot}])$ of 9.86 and median redshift of 5.55. We utilise the $L_{\rm IR}$, $L_{\rm UV}$, and $\beta_{\rm UV}$ estimates of the targets from Table 1 of Barisic et al. (2017), which updates the original measurements by using additional HST photometry (and also provides measurements for individual clumps in two sources). We excluded the X-ray detected quasar from their sample.

3 Analysis

The analysis of the IRX–$\beta_{\rm UV}$ relation requires us to measure the galaxies’ luminosities at UV and IR wavelengths, as well as the UV continuum slope $\beta_{\rm UV}$, which we describe in the following section.

3.1 $L_{\rm IR}$, $L_{\rm UV}$, and $\beta_{\rm UV}$

Infrared luminosities ($L_{\rm IR}$) are estimated by fitting a single component modified black body function to the ALMA $240\,{\rm GHz}$ fluxes of our galaxies and integrating over the wavelength range between 8 and $1000\,{\rm\mu m}$. We assume a dust emissivity index of $\beta_{\rm{d}}=1.5$ (Dunne & Eales, 2001).

Dust temperatures ($T_{d}$) are estimated based on a stacking analysis of Spitzer, Herschel, and ALMA images of our full $z\sim 3$ galaxy sample (see Magnelli et al., 2014, for a similar analysis). From the full IR SED of these stacked data, we obtain a best-fit dust temperature of $T_{d}=41.5\pm 2.5\,{\rm K}$ by fitting a standard modified black body spectrum (Casey, 2012). The derived dust temperature is consistent with recently obtained values at $z\sim 4$ galaxies (e.g. Schreiber et al., 2016). Hereafter, we will use $T_{d}=40\,{\rm K}$ as our fiducial value throughout the paper, but we will comment on the potential impact of assuming a lower dust temperature of $T_{d}=35$ K.

Uncertainties on our ${L_{\rm IR}}$ measurements are estimated using the ALMA flux measurement errors. For the 26 non-detected targets, we estimated upper limits on their $L_{\rm IR}$  using the $3\sigma$ ALMA flux limits.

The $L_{\rm IR}$ of our detected targets range between $5.9\times 10^{11}$ and $5.8\times 10^{12}\,L_{\odot}$, thus straddling the range between luminous infrared galaxies (LIRGs, $10^{11}<L_{\rm IR}<10^{12}\,L_{\odot}$) and ultra luminous infrared galaxies (ULIRGs, $10^{12}<L_{\rm IR}<10^{13}\,L_{\odot}$).

While the increased CMB temperature can affect the luminosity measurements of higher redshift galaxies (e.g. da Cunha et al., 2013; Zhang et al., 2016), we estimate that the expected effect on our $z\sim 3$ sample is still less than 5% and thus negligible compared to our typical flux measurement uncertainties.

UV spectral slopes ($\beta_{\rm{UV}}$) are calculated by fitting a power law $f_{\rm{\lambda}}\propto\lambda^{\beta_{\rm{UV}}}$ to the broad- and narrow-band photometry of the galaxies over the rest-frame wavelength range $1500-2500\,$Å. For this purpose, we utilise the photometry measured within 3*′′* apertures listed in the COSMOS15 catalogue (Laigle et al., 2016). Typically more than 10 data points are available covering this waveband for each of our targets. The monochromatic UV luminosities are then calculated at rest-frame $1600\,$Å using the above power law fits to the UV SEDs and defining $L_{\rm{UV}}=\nu_{1600}\,L_{\nu_{1600}}$. The range of obtained $L_{\rm UV}$ covers $1.2\times 10^{10}<\nu_{1600}\,L_{\nu_{1600}}<1.5\times 10^{11}\,L_{\odot}$.

3.2 Stacking Analysis

While 41 out of our sample of 67 galaxies are detected with ALMA, the remaining targets still contain important information on the average IRX relation at $z\sim 3$. We therefore perform a stacking analysis of all target fields. In particular, we are interested in the average IRX of our sample as a function of UV continuum slope $\beta_{\rm UV}$. We thus split our full galaxy sample in three bins: $\beta_{\rm UV}$ <-1.5,~{}$$\beta_{\rm UV}$$\,=-1.5 to $-1.0$, and $\beta_{\rm UV}$ $>-1.0$, which contain 14, 30, and 23 galaxies, respectively.

Thanks to the very uniform synthesised beam full-width-at-half-maximum (FWHM) of our ALMA images, stacks can be created by simply taking the median of each pixel in the ALMA maps. During the stacking procedure, the ALMA images are weighted by the inverse of the galaxy’s $L_{\rm UV}$, in order to compute the median IRX of the sample going into the stack. All three ALMA stacks reveal a significant detection, as shown in Fig. 1.

We also compute a stacked image of all 26 ALMA non-detections using the same procedure as above. The stack reveals a tentative detection with a significance of $\sim 4.3\,\sigma$ (lower right panel of Fig. 1). We note that this tentative detection is the most significant positive/negative signal within the central $23^{\prime\prime}\times 23^{\prime\prime}$ of the stacked image, corresponding to the FWHM of the primary beam of each observation at 240 GHz.

We then measure total fluxes for each of the stacks using the peak fluxes and applying a correction to total flux derived from the mean of the ALMA detected sources. Uncertainties on all stacked quantities, including the IRX, are estimated using bootstrap resampling. The results of this stacking analysis are summarised in Table 1, where we also tabulate the resulting median IRX of the different stacks that are used in our analysis later on.

4 Results and Discussion

In the following, we discuss our results on the IRX–$\beta_{\rm UV}$ and the IRX–$M_{\ast}$ relations at $z\sim 3.2$, as well as their possible redshift evolution to $z\sim 4-6$.

4.1 IRX–$\beta_{\rm UV}$ Relation at $z\sim 3.2$

The IRX–$\beta_{\rm UV}$ diagram of our $z\sim 3.2$ galaxy sample is shown in Fig. 2 together with the expected relations of two different local dust extinction curves, the M99 relation and SMC-type extinction (e.g. Prevot et al., 1984). These two attenuation laws have both been discussed in the past to represent different sub-samples of $z\sim 2-3$ galaxies (e.g. Reddy et al., 2012; Bouwens et al., 2016).

The figure shows that our ALMA detected targets typically lie above the M99 line. Similarly, the $3\sigma$ upper limits of individually un-detected sources often lie above the M99 relation. However, when we combine these samples through stacks in bins of $\beta_{\rm UV}$, we obtain three measurements of the IRX–$\beta_{\rm UV}$ relation that lie at the low end of the individually detected IRX values, showing that the typical $z\sim 3$ galaxy is mostly consistent with the M99 relation, in particular at $\beta>-1.5$.

Since the sample with $\beta<-1.5$ only contains 14 sources, its stacked IRX value is quite uncertain. It lies significantly above the M99 relation. However, in the redder UV slope bins, where we have better statistics, both median stacked IRX values are in excellent agreement with the M99 relation.

Interestingly, the median stacks show that an SMC-like attenuation is clearly not applicable for the general galaxy sample at $z\sim 3.2$. In the $\beta_{\rm UV}$ range probed by our stacks, the median IRX values lie $\gtrsim 0.7\,{\rm dex}$ above an SMC-like attenuation curve. However, considering the stack of ALMA non-detected sources only, we find a very low IRX value that lies within $<1\sigma$ of the SMC curve.

This provides some evidence for a large dispersion in IRX values at a given $\beta_{\rm UV}$ from one galaxy to another at $z\sim 3$, with certain subsamples of galaxies that show significantly lower IRX values at a given $\beta_{\rm UV}$ than the average population. This confirms previous analyses of $z\sim 2-3$ galaxies which find that the local M99 relation is typically applicable to the general galaxy population, but that there are a few galaxies (that are typically younger) that are more consistent with an SMC-type attenuation (Siana et al., 2009; Reddy et al., 2012). Within our sample, we do not find a clear separation in physical parameters, however, between our ALMA detections and non-detections (e.g. $z_{\rm phot}$, age, $M_{\ast}$).

Note that our conclusions would not change significantly under the assumption of a lower dust temperature, $T_{d}=35\,{\rm K}$ that has been used in some previous analyses of high-redshift galaxies (e.g., Bouwens et al., 2016). Under this assumption, the values of $L_{\rm IR}$ as well as the IRX decrease by $\sim 0.25\,{\rm dex}$ for our sample (see arrow in Fig 2).

4.2 Comparison to Previous $z\sim 3$ Studies

In Fig. 3, we compare our results with several previous studies of the IRX–$\beta_{\rm UV}$ relation at $z\sim 3$ for samples covering a similar mass range (M_{\ast}$$\gtrsim 10^{10}$$M_{\odot}). In particular, the comparison samples include stacking analyses of Lyman Break Galaxies (LBGs) and/or IR-selected galaxies (Coppin et al., 2015; Bourne et al., 2017; Álvarez-Márquez et al., 2016; Reddy et al., 2017), as well as individual ALMA detections of $z\gtrsim 1.5$–2.7 galaxies from Bouwens et al. (2016). Despite similar mass ranges probed, these previous analyses found IRX values that are sometimes different by $\gtrsim 1\,{\rm dex}$ at a fixed $\beta_{\rm UV}$.

Interestingly, our data points and stacks are in good agreement with the full range of these previous analyses. For instance, the stack of our non-detections shows an IRX value similar to individual galaxies found by Bouwens et al. (2016). On the other hand, our stacked measurements and our individual detections are in good agreement with the values found by Coppin et al. (2015), who utilised a stacking analysis of SCUBA-2 and Herschel observations of $\sim 4200$ LBGs. Similar agreement is found with the $z\sim 3$ high-mass stacks from Álvarez-Márquez et al. (2016).

We also investigate whether an IRX–$\beta$ relation with bluer intrinsic UV slope could explain our data. Such a trend to bluer intrinsic slopes at higher redshift is expected due to, e.g., lower metallicities and younger stellar population ages (e.g. Dayal & Ferrara, 2012; Schaerer et al., 2013; Alavi et al., 2014; Castellano et al., 2014; de Barros et al., 2014; Sklias et al., 2014; Cullen et al., 2017), and evidence for such bluer intrinsic slopes is continuously emerging from several studies. For instance, Smit et al. (2016) derive a similar best fit $\beta_{\rm UV,int}=-2.5$ by comparing UV-based SFRs with those inferred from ${\rm H_{\alpha}}$ fluxes at $z\sim 4$. Reddy et al. (2017) provide theoretical derivations of IRX–$\beta$ relationships based on different assumptions of the intrinsic SED and dust curves. A relation based on a young SED with $\beta_{\rm UV,int}=-2.62$, reddened by the Reddy et al. (2015) dust curve indeed provides an adequate representation of our stacked points, including the bluest bin, in better agreement than the M99 relation. It will thus be interesting to further investigate the overall shape of the IRX–$\beta_{\rm UV}$ relation and the $\beta_{\rm UV,\,int}$ based on larger ALMA surveys of high redshift galaxies spanning a wide range in UV continuum slopes in the future.

4.3 The IRX-$M_{\ast}$ Relation

Several studies have shown that more massive star-forming galaxies show larger dust attenuation than their lower mass counterparts, which results in a mass dependence of the IRX values (e.g. Reddy et al., 2010; Buat et al., 2012; Heinis et al., 2013; Coppin et al., 2015; Pannella et al., 2015; Bouwens et al., 2016; Bourne et al., 2017; Dunlop et al., 2017).

In a recent analysis, Álvarez-Márquez et al. (2016) discuss the IRX–$M_{\ast}$ relation based on a stacking analysis of $z\sim 3$ LBGs using Herschel and AzTEC data. Although their mean IRX–$\beta_{\rm UV}$ relation exhibits a lower IRX than M99 at fixed $\beta_{\rm UV}$ as shown in our Fig. 3, their stacking analysis of a mass-matched sub-sample with M_{\ast}$$\sim 10^{10.5}–$10^{11}$ results in a $\sim 0.5$ dex higher IRX at fixed $\beta_{\rm UV}$ than their mean relation, agreeing with the M99 relation and our results (see Fig. 3), and it clearly demonstrates the importance of comparing mass-matched galaxy samples.

In Fig. 4 right we analyse the IRX–$M_{\ast}$ relation of our sample. Although the dynamic range in mass is small, there is generally a correlation with a slope that is consistent with previous studies. However, our targets show lower IRX values than previous studies at fixed $M_{\ast}$. We note that this could potentially be due to different $\beta_{\rm UV}$ distributions, in particular the lack of redder objects in our sample (i.e., with $\beta_{\rm UV}$ $>-0.5$).

We quantify the offset by fitting for the normalisation but keeping the slope fixed to the one found in Heinis et al. (2013). This results in a best-fit IRX–$M_{\ast}$ relation following:

[TABLE]

The offset in normalisation relative to Heinis et al. (2013) is thus $0.24\,{\rm dex}$. A similar offset is found relative to the IRX–$M_{\ast}$ relation from Álvarez-Márquez et al. (2016).

Interestingly, the sub-sample of our ALMA non-detections shows an even larger offset to lower IRX values compared to previous studies when we consider the location of their median stacked IRX. This suggests that there is also a significant scatter around the mean IRX–$M_{\ast}$ relation at $z\gtrsim 3$. In Reddy et al. (2010), the authors discussed the IRX–$M_{\ast}$ relation of $z\sim 2$ UV selected star forming galaxies, and estimated a scatter in the relation of $\sim 0.46\,{\rm dex}$ about a linear fit. Consistent with this, our result suggests that individual galaxies indeed can span a $\sim 0.5\,{\rm dex}$ range in IRX values at $z\sim 3$ around the IRX–$\beta_{\rm UV}$ relation. Such a large dispersion will result in large uncertainties on the inferred dust corrections of UV-based SFRs of individual galaxies, and it will be crucial to search for correlations with other physical parameters in an attempt to decrease these dispersions in the future with larger samples. For instance, Nordon et al. (2013) show that the scatter in the IRX is correlated with the position of a galaxy with respect to the main-sequence of star-formation.

4.4 The Evolution of IRX–$\beta_{\rm UV}$ Beyond $z=3$

Our finding that the IRX–$\beta_{\rm UV}$ relation at $z\sim 3.2$ is consistent with the local M99 relation is in contrast to some recent results on high redshift galaxies.

For instance, Capak et al. (2015) find that their sample of $z\sim 5.5$ galaxies generally lies $\sim 1\,{\rm dex}$ below the M99 relation (see also Barisic et al., 2017). Similarly, Bouwens et al. (2016) find that the $z>4$ IRX–$\beta_{\rm UV}$ relation lies below the M99 relation based on a stacking analysis of Lyman Break Galaxies (LBGs) in the ASPECS data over the Hubble Ultra-Deep Field.

To investigate this evolution further, we compare our $z\sim 3.2$ galaxies with the $z\sim 4-5$ sample from Scoville et al. (2016) and with the $z\sim 5.5$ galaxies from Capak et al. (2015); Barisic et al. (2017) in Fig. 4. While the $z>4$ sample is currently very small, we see tentative evidence for an evolutionary trend at $z>3$. Between $z\sim 3$–6, the median IRX gradually drops by $\gtrsim 1\,{\rm dex}$, and a large fraction of the $z\geq 4$ sample lies close to or even below the SMC like extinction curve.

This evolution to lower IRX values at a given $\beta_{\rm UV}$ is very intriguing, and could potentially hint at an evolution of the dust properties within the first 2 Gyr of cosmic history. This is expected theoretically due to the younger stellar population ages of galaxies in the early universe for which the AGB star population has not yet had enough time to produce a significant amount of dust. At very high redshift, the main dust production mechanism is thus thought to be supernovae explosions, which can be expected to lead to a different dust composition and attenuation curve than at later times (Todini & Ferrara, 2001; Gallerani et al., 2010). Alternative solutions as to why higher redshift galaxies show lower dust emission could involve a higher fraction of gas to be locked up in a cold, molecular component compared to lower redshift (Ferrara et al., 2016). Furthermore, it will be important to understand the relative morphological distributions and the mixing of the stars and the dust in multi-component high redshift (e.g. Casey et al., 2014; Koprowski et al., 2016).

However, before making strong claims about a possible redshift evolution, we have to control for the differences in stellar masses of the different samples. To do this, we explore the IRX–$M_{\ast}$ relation of the three high-redshift samples in the right panel of Fig. 4. As is evident from the plot, there is indeed a difference in mass among the different samples. However, the higher redshift samples generally show a larger offset from the $z\sim 3$ stacking analysis of Heinis et al. (2014) at lower masses. The median offsets from our best-fit IRX–$M_{\ast}$ relation are $\sim 0.31\,{\rm dex}$ and $\sim 1.2\,{\rm dex}$ for the samples at $z\sim 4.3$ and $z\sim 5.5$ (including $3\,\sigma$ upper limits), respectively.

Note that this rapid redshift evolution is in general agreement with the best-fit IRX–$M_{\ast}$ relation presented in Bouwens et al. (2016), where the authors studied star forming galaxies with ${\rm log(}M_{\ast})<9.5$ at $z\sim 4$–10, and found lower IRX values by $\sim 0.5\,{\rm dex}$ at fixed $M_{\ast}$ compared to a consensus relation estimated from previous analyses at $z\sim 2-3$ (Reddy et al., 2010; Whitaker et al., 2014; Álvarez-Márquez et al., 2016). However, this was based on the assumption of $T_{d}=35\,{\rm K}$. When assuming that the dust temperature increases with redshift, their results are in better agreement with the previous relations (Fig 4).

It is thus clear that it is critical to constrain the dust temperatures at $z>4$. Changes in $T_{d}$ can have a significant impact on the $L_{\rm IR}$ and IRX values estimated from single-band observations, and could affect all samples included in Fig.4. An evolution in the $T_{\rm d}$ as a function of redshift is has been observed by many studies using stacking analyses (e.g. Magdis et al., 2012; Magnelli et al., 2014; Béthermin et al., 2015; Ivison et al., 2016), typically finding an evolution following $T_{d}\propto(1+z)^{0.32}$. Assuming $T_{d}=40\,{\rm K}$ at $z=3.2$, this evolution then predicts $T_{d}\sim\,43\,{\rm K}$ and $T_{d}\sim\,46\,{\rm K}$ at $z\sim 4.3$ and $z\sim 5.5$, respectively.

Relative to our fiducial value of $T_{\rm d}=40\,{\rm K}$, these assumptions would increase the IRX values of the observed samples by $\sim 0.1\,{\rm dex}$, and $\sim 0.2\,{\rm dex}$ at $z\sim 4.3$ and $z\sim 5.5$, respectively. Clearly, a $T_{\rm d}$ evolution could thus partially compensate for the observed IRX deficit at $z>4$. However, a much more extreme evolution of the dust temperature would be required to fully account for the median decrement in IRX values at $z>4$. Nevertheless, it is clear that better constraints on the dust temperatures of high redshift galaxies will be crucial for future studies.

5 Conclusions

In this paper, we have presented an investigation of the IRX-$\beta_{\rm UV}$ relation of massive (M_{\ast}$$\sim 10^{10.7}$$M_{\odot}) star-forming galaxies at $z\sim 3.2$, through FIR continuum observations with ALMA. We also explored a possible redshift evolution of the IRX–$\beta_{\rm UV}$ and the IRX–$M_{\ast}$ relations by utilising public ALMA observations of star forming galaxies between $z\sim 3$ and $z\sim 6$. Our findings can be summarised as follows:

• •

The IRX–$\beta_{\rm UV}$ relation of our sample is on average mostly consistent with that of local starburst galaxies derived by Meurer et al. (1999) and we can exclude an SMC-type extinction curve for the general star-forming $z\sim 3$ galaxy (see Fig. 2, 3). However, as the stack of our ALMA non-detections shows, certain sub-samples of galaxies exist that are in better agreement with the SMC relation.

• •

Our individually detected sources and our stacks span a large range in IRX values at fixed $\beta_{\rm UV}$, consistent with previous analyses at $z\sim 3$. This indicates the presence of a very large dispersion around the average IRX–$\beta_{\rm UV}$ relation of up to $\pm 0.5$dex. The same is observed for the IRX–$M_{\ast}$ relation (see Fig. 3,4).

• •

With public ALMA observations of galaxies at $z\sim 4$–6, we explored the redshift evolution of the IRX–$\beta_{\rm UV}$ and the IRX–$M_{\ast}$ relations, which we find to be significant. Even though the current samples are still small, we find a strong, gradual decrease to lower IRX values at fixed $\beta_{\rm UV}$ or $M_{\ast}$ with increasing redshift, with the $z\sim 6$ sample generally lying an order of magnitude below our $z\sim 3$ relations on average (see Fig. 4).

The fact that our $z\sim 3$ sample is consistent with the local M99 dust law, yet higher redshift samples appear to fall off this relation, is particularly interesting, as it could point to significant changes in the ISM and dust properties of high-redshift galaxies. It will be important to confirm the suggested evolution between $z\sim 3$ and $z\sim 6$ with future datasets. In particular, our study demonstrates that individual IR detections of the dust emission from a statistical high-redshift sample are crucial due to considerable scatter around the average IRX–$\beta_{\rm UV}$ that can not be fully captured by stacking analyses. Our analysis thus motivates a larger investment of ALMA time to study the dust emission of a statistical sample of normal high-redshift galaxies.

Acknowledgements

The authors thank A. Ferrara, R. Bouwens, and N. Reddy for very helpful discussions related to this work. This paper makes use of the following ALMA data: ADS/JAO.ALMA#2013.1.00151.S, ADS/JAO.ALMA#2013.1.00034.S, and ADS/JAO.ALMA#2012.1.00523.S. ALMA is a partnership of ESO (representing its member states), NSF (USA) and NINS (Japan), together with NRC (Canada) and NSC and ASIAA (Taiwan) and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO and NAOJ. Based on data products from observations made with ESO Telescopes at the La Silla Paranal Observatory under ESO programme ID 179.A-2005 and on data products produced by TERAPIX and the Cambridge Astronomy Survey Unit on behalf of the UltraVISTA consortium. B.G. gratefully acknowledges the support of the Australian Research Council as the recipient of a Future Fellowship (FT140101202). MTS acknowledges support from a Royal Society Leverhulme Trust Senior Research Fellowship (LT150041). VS acknowledges support from the European Union’s Seventh Frame-work program under grant agreement 337595 (ERC Starting Grant, ’CoSMass’). OLF acknowledges support from the European Research Council Advanced Grant ERC-2010-AdG-268107-EARLY.