UV-CONTINUUM SLOPES AT z  ∼  4–7 FROM THE HUDF09+ERS+CANDELS OBSERVATIONS: DISCOVERY OF A WELL-DEFINED UV COLOR–MAGNITUDE RELATIONSHIP FOR z ⩾ 4 STAR-FORMING GALAXIES^*

The first set of observations we utilize is from the HUDF09 program (GO11563: PI: Illingworth) and involves ultra-deep WFC3/IR observations over three fields that already have ultra-deep ACS observations. These fields include the HUDF (Beckwith et al. 2006) and two HUDF05 flanking fields (Oesch et al. 2007). These observations are primarily of use in establishing the distribution of UV-continuum slopes β for very faint z  ∼  4–7 galaxies.

At present, the full two years of WFC3/IR observations (all 192 orbits from the GO11563 program) have been obtained over the three HUDF09 fields. The WFC3/IR observations over the HUDF include some 111 orbits of observations—while the HUDF09-1 and HUDF09-2 fields include 33 orbits and 48 orbits of observations, respectively. Our reductions of the WFC3/IR observations over these fields are already described in Bouwens et al. (2011b) and are conducted using standard procedures. Particularly important to this process was ensuring that the geometric distortion and velocity aberration were treated properly so that the registration between bands and with the ACS observations was very accurate (≲001). Accurate registration is absolutely essential not only for maximizing the accuracy of our color measurements, but also for minimizing the effect that misregistration might have on scatter in these same color measurements. The depths of our WFC3/IR observations reach to ≳29 AB mag at 5σ and are presented in detail in Table 1. The FWHM of the point spread function (PSF) in the WFC3/IR observations is ∼016.

Our reductions of the ACS observations over the HUDF09-2 area include all available data—including the 82 orbits that were taken during the execution of the original HUDF05 program and an additional 111 orbits that were taken in parallel with the HUDF09 WFC3/IR observations over the HUDF. The total number of orbits per filter over the HUDF09-2 are 10 orbits F435W, 32 orbits F606W, 46 orbits F775W, 16 orbits F814W, and 89 orbits F850LP. The total integration time is ≳50% of that available for the HUDF, allowing us to reach within ≲0.4 mag of the HUDF. This provides us with sufficiently high signal-to-noise (S/N) levels required to keep contamination in z ≳ 6 selections to a minimum. A more detailed description of this data set is provided by Bouwens et al. (2011b).

Wide-area WFC3/IR observations allow us to establish the UV-continuum slope β distribution for the rarer, higher luminosity sources. These observations are available over the CDF-South GOODS field from the Early Release Science (Windhorst et al. 2011) and CANDELS (Grogin et al. 2011; Koekemoer et al. 2011) programs.

The WFC3/IR observations from the ERS program are distributed over the northern part of the CDF South GOODS field (Windhorst et al. 2011; Figure 1). These observations consist of 60 orbits of F098M, F125W, and F160W observations distributed over 10 distinct 21 × 23 pointings, with two orbits F098M, two orbits F125W, and two orbits F160W per pointing. The WFC3/IR observations extend over ∼43 arcmin² in total, although only ∼39 arcmin² of those observations overlap with the ACS GOODS data and can be used here.

WFC3/IR observations from the CANDELS program are distributed over the southern two-thirds of the CDF-South GOODS field (106 arcmin²), with the deepest section planned for the central ∼66 arcmin² portion. Three orbits of F105W, four orbits of F125W, and four orbits of the F160W observations are available over the central sections of the CDF-South (representing all 10 SN search epochs to 2012 February 18) while only one orbit of the F105W, F125W, and F160W observations are available over the southern portion. The layout of the CANDELS and ERS observations with the CDF-South GOODS field is shown in Figure 1.

Our reductions of the WFC3/IR observations over the CDF-South are performed in exactly the same manner as the HUDF09 observations (see Bouwens et al. (2011b) for a detailed description). As in our HUDF09 reductions, a 006 pixel size was used. For the ACS observations over the CDF-South, we used the Bouwens et al. (2007) reductions for the F435W, F606W, F775W, and F850LP bands. These reductions are comparable to the GOODS v2.0 reduction (Giavalisco et al. 2004) but take advantage of the substantial supernova (SN) follow-up observations over the CDF-South (Riess et al. 2007) which add ∼0.1–0.3 mag of depth to the z₈₅₀ band.

We also reduced the new ACS I₈₁₄-band observations over the CDF-South GOODS fields from the CANDELS+ERS programs to take advantage of the superb depth available with these data at ∼8000 Å (typically 13 orbits over the CANDELS-DEEP region or 0.8 mag deeper than the GOODS i₇₇₅-band exposures). The F814W observations are valuable for controlling for contamination in our z  ∼  7 selections and for minimizing the uncertainties in our β determinations at z  ∼  4. We included all F814W observations from the ERS+CANDELS programs over the CDF-South. After using public codes (e.g., Anderson & Bedin 2010) to correct the raw frames for charge transfer efficiency defects and row-by-row banding artifacts, we performed the alignment, cosmic-ray rejection, and drizzling with the ACS GTO apsis pipeline (Blakeslee et al. 2003).

The procedure for source detection and photometry is discussed in many of our previous papers (e.g., Bouwens et al. 2007); we summarize the key steps here. We generate catalogs for our samples using the SExtractor (Bertin & Arnouts 1996) software run in double image mode, with the detection image taken to be the square root of χ² image (Szalay et al. 1999; similar to a co-added image of the WFC3/IR observations) and the measurement image to be a PSF-matched image from our ACS+WFC3/IR image sets (see below).

We select each of our high-redshift samples (i.e., z  ∼  4, z  ∼  5, z  ∼  6, and z  ∼  7) from separately generated source catalogs. This allows us to take advantage of the fact that z  ∼  4–6 samples can be selected entirely (or almost entirely) based on the ACS data (involving generally smaller apertures and therefore reducing the uncertainties in our color measurements) while z  ∼  7 selections require PSF matching to the WFC3/IR data.⁷ For our z  ∼  4–6 catalogs, the first step is to match PSFs by smoothing the shorter wavelength images to the resolution of the ACS z₈₅₀-band image. The square root of χ² images we use for the SExtractor detection images are constructed from the V₆₀₆i₇₇₅z₈₅₀, i₇₇₅z₈₅₀, and z₈₅₀ images for our z  ∼  4, z  ∼  5, and z  ∼  6 selections, respectively. For our z  ∼  7 catalogs, all the observations are PSF matched to the WFC3/IR H₁₆₀ band. The square root of χ² images are constructed from the available WFC3/IR imaging observations, e.g., Y₁₀₅J₁₂₅H₁₆₀ for our HUDF09 fields, Y₀₉₈J₁₂₅H₁₆₀ images for the ERS observations, and J₁₂₅H₁₆₀ for the CANDELS observations.

Colors are measured within small scalable apertures using a Kron (1980) factor of 1.2. These small-aperture color measurements are then corrected to total magnitudes in two steps. First, a correction is made for the light in a larger scalable aperture, with Kron factor equal to 2.5. This correction is made from the square root of χ² image to optimize the S/N. Second, a correction is made for light outside this large scalable aperture—using the encircled energies measured for point sources (typically a ∼0.1–0.15 mag correction).

Photometry was done using the latest WFC3/IR zero-point calibrations (2012 March 6) which differ by ∼0.01–0.02 mag from the earlier zero points. In addition, a correction to the photometry (i.e., E(B − V) = 0.009) was performed to account for foreground dust extinction from our own galaxy. This correction was based on the Schlegel et al. (1998) dust map.

To select star-forming galaxies at high redshift, we will use the well-established Lyman-break galaxy (LBG) technique. This technique takes advantage of the unique spectral characteristics of high-redshift star-forming galaxies, which show a very blue spectrum overall but a sharp cutoff blueward of Lyα. It has been shown to be very robust through extensive spectroscopic follow-up (Steidel et al. 1996; Steidel et al. 2003; Bunker et al. 2003; Dow-Hygelund et al. 2007; Popesso et al. 2009; Vanzella et al. 2009; Stark et al. 2010).

We will base our high-redshift samples on selection criteria from previous work on z ⩾ 4 galaxies. For our z  ∼  4 B₄₃₅ and z  ∼  5 V₆₀₆ dropout samples, we will apply almost the same criteria as Bouwens et al. (2007), namely,

$$\begin{eqnarray*} (B_{435}-V_{606} > 1.1) &\wedge& (B_{435}-V_{606} > (V_{606}-z_{850})+1.1) \\ &\wedge& (V_{606}-z_{850}<1.6) \end{eqnarray*}$$

for our B-dropout sample and

$$\begin{eqnarray*} &&[(V_{606}-i_{775} > 0.9(i_{775}-z_{850})+1.5) \vee (V_{606}-i_{775} > 2))] \wedge \\ &&\qquad\qquad\qquad\quad(V_{606}-i_{775}>1.2) \wedge (i_{775}-z_{850}<0.8) \end{eqnarray*}$$

for our V₆₀₆-dropout selection (but note that we use a i₇₇₅ − z₈₅₀ < 0.8 color selection instead of the i₇₇₅ − z₈₅₀ < 1.2 selection used by Bouwens et al. 2007 to minimize contamination). For our z  ∼  6 i₇₇₅-dropout selection, we expanded the criteria used in Bouwens et al. (2007) to take advantage of the deep near-IR observations from WFC3/IR to set limits on the color redward of the break. Our z  ∼  6 i₇₇₅-dropout criterion is

$$\begin{eqnarray*} &&(i_{775}-z_{850}>1.3)\wedge (z_{850}-J_{125}<0.9). \end{eqnarray*}$$

Finally, our z  ∼  7 z₈₅₀-dropout criterion is

$$\begin{eqnarray*} &&(z_{850} - Y_{105} > 0.7)\wedge (Y_{105}-J_{125} < 0.8) \wedge \\ &&\qquad(z_{850}-Y_{105}>1.4(Y_{105}-J_{125})+0.42), \end{eqnarray*}$$

or

$$\begin{eqnarray*} &&\qquad\qquad\qquad\quad(z_{850}-J_{125}>0.9)\wedge \\ &&\quad(z_{850}-J_{125}>0.8+1.1(J_{125}-H_{160})) \wedge \\ &&\qquad(z_{850}-J_{125}>0.4+1.1(Y_{098}-J_{125})), \end{eqnarray*}$$

depending upon whether our search field is from the HUDF09/CANDELS data sets or the ERS data set, respectively. The above selection criteria closely match those used in our previous study of the UV-continuum slope β at z  ∼  7 (Bouwens et al. 2010a). These criteria are slightly more inclusive than those considered by Bouwens et al. (2011b), but this is to allow us to maximize the size of our samples and to extend our selection to the highest redshift possible without suffering significant contamination from Lyα emission or IGM absorption. In cases of a non-detection in the dropout band, we set the flux in the dropout band to be equal to the 1σ upper limit. To take advantage of the deep I₈₁₄-band observations over the CDF-South to keep contamination in our z  ∼  7 samples to a minimum, we required all z  ∼  7 z₈₅₀-dropouts in our selection to have I₈₁₄ − J₁₂₅ colors greater than 2.0 or to be undetected in the I₈₁₄ band (at <1.5σ). No attempt is made to select or measure UV-continuum slopes for star-forming galaxies at z  ∼  8 due to the fact that the J₁₂₅ − H₁₆₀ color (the only available color providing high S/N information on the rest-frame UV SED at z  ∼  8) is affected by Lyman-series absorption and Lyα emission at z > 8.1 (but, see Taniguchi et al. 2010).

Figure 2 provides a convenient illustration of the approximate range in UV-continuum slopes β and redshifts selected by the above two-color criteria.

Zoom In Zoom Out Reset image size

We utilize several additional selection criteria in defining our final samples. Sources are required to have SExtractor stellarity parameters less than 0.8 (i.e., they show evidence of being extended) to ensure our samples are free of contamination from active galactic nuclei (AGNs) or low-mass stars. Given that >97% of bright sources in Lyman-break selections show clear evidence of having extended profiles (i.e., are not pointlike), this criterion has almost no effect on the overall composition of our samples.

We also require that sources be undetected in all bands blueward of the dropout bands. Sources are rejected if they are detected at >2σ in a single band, >1.5σ in two bands, or have a χ²_opt > 3 in the optical bands. We take χ²_opt = Σ_iSGN(f_i)(f_i/σ_i)² where f_i is the flux in band i in our smaller scalable apertures, σ_i is the uncertainty in this flux, and SGN(f_i) is equal to 1 if f_i > 0 and −1 if f_i < 0 (Bouwens et al. 2011b). As Bouwens et al. (2011b) illustrate, a χ²_opt criterion can be particularly useful for minimizing contamination in high redshift samples (see also Oesch et al. 2012a, 2012b).

Table 2 summarizes the properties of the z  ∼  4, z  ∼  5, z  ∼  6, and z  ∼  7 samples we derive from the HUDF09+ERS+CANDELS observations. In total, 308, 137, 70, and 57 z  ∼  4, z  ∼  5, z  ∼  6, and z  ∼  7 galaxies are found in the HUDF09 fields and 1514, 277, 101, and 44 z  ∼  4, z  ∼  5, z  ∼  6, and z  ∼  7 galaxies are found in the ERS+CANDELS fields. The approximate redshift distributions for our z  ∼  4, z  ∼  5, z  ∼  6, and z  ∼  7 selections are shown more explicitly in Figure 3. These redshift distributions are as calculated in Bouwens et al. (2007) and Bouwens et al. (2011b).

Zoom In Zoom Out Reset image size

Establishing the UV-continuum slope β for a source (where β is defined such that f_λ∝λ^β) requires that we have at least two flux measurements of a source in the UV continuum not significantly affected by Lyα or IGM absorption (unless both the redshift and Lyα flux are already known). In practice, this means that we should not make use of the flux information in the band immediately redward of the break since it is frequently contaminated by Lyα emission or absorption from the IGM. Assuming we have wavelength coverage in the B₄₃₅V₆₀₆i₇₇₅z₈₅₀Y₁₀₅J₁₂₅H₁₆₀ bands, we can obtain good estimates of the UV-continuum slope β for high-redshift sources from z  ∼  4 to at least z  ∼  7.5.

There are a few different approaches we could adopt in using the available flux measurements in the UV continuum to estimate these slopes. One approach is to derive the UV-continuum slopes from a fit to all available flux measurements in the rest-frame UV. An alternate approach is to derive a slope from the flux measurements at both ends of a fixed rest-frame wavelength range (e.g., as utilized by Bouwens et al. 2009).

Each approach has its advantages. The first approach uses the full flux information available on each source and also takes advantage of a much more extended wavelength baseline. This results in much smaller uncertainties in our β determinations, particularly at z  ∼  4 and z  ∼  5, where fluxes in four separate bandpasses are available. The simulations we perform in Appendix B.3 suggest factors of ∼1.5 improvement in the uncertainties on β at z  ∼  4 and z  ∼  5. On the other hand, the second approach has the advantage that all UV-continuum slope β determinations are made using a similar wavelength baseline at all redshifts. As a result, even if the spectral energy distribution (SED) of galaxies in the rest-frame UV is not a perfect power law (e.g., from the well-known bump in dust extinction law at 2175 Å; Stecher 1965), we would expect β determinations made at one redshift to show no systematic offset relative to those made at another, allowing for more robust measurements of evolution across cosmic time.

After some testing, we adopted the approach that utilizes all the available flux information to determine the UV-continuum slopes β. Figure 4 provides an example of such a determination for a z  ∼  4 galaxy in our HUDF sample and an example of such a determination for a z  ∼  6 galaxy. The passbands we will consider in deriving the slopes include the i₇₇₅z₈₅₀Y₁₀₅J₁₂₅ bands for z  ∼  4 galaxies, the z₈₅₀Y₁₀₅J₁₂₅H₁₆₀ bands for z  ∼  5 galaxies, the Y₁₀₅J₁₂₅H₁₆₀ bands for z  ∼  6 galaxies, and the J₁₂₅H₁₆₀ bands for z  ∼  7 galaxies. We include the I₈₁₄ band in these fits for z  ∼  4 galaxies over the ERS/CANDELS fields. We replace the Y₁₀₅ band with the Y₀₉₈ band in these fits for galaxies over the ERS field. In selecting these passbands, we explicitly excluded passbands which could be contaminated by Lyα emission, the Lyman-continuum break, or flux redward of the Balmer break. These choices are motivated by our expected redshift distributions for these samples (Figure 3). Table 3 includes a list of all the bands we use to perform these fits. The mean rest-frame wavelengths for our derived UV-continuum slopes β at z  ∼  4, z  ∼  5, z  ∼  6, and z  ∼  7 are 2041 Å, 1997 Å, 1846 Å, and 1731 Å, respectively.

Zoom In Zoom Out Reset image size

In determining the UV-continuum slopes β from the flux information just discussed, we find the UV-continuum slope β that minimizes the value of χ²:

$$\begin{equation} \chi ^2 = \Sigma _{i} \left(\frac{f_i - f_0 \lambda ^{\beta }}{\sigma _{i}} \right)^2, \end{equation} \tag{ 1 }$$

where f_i and σ_i are the observed fluxes and uncertainties, respectively, β is the best-fit UV-continuum slope, and f₀ is the best-fit normalization factor. The fluxes f_i ± σ_i in the above fits are from photometry of the full ACS+WFC3/IR data set PSF matched to the WFC3/IR H₁₆₀ band, using typical Kron-style apertures.

This photometry is therefore distinct from that used to select our z  ∼  4, z  ∼  5, and z  ∼  6 samples, since sources in our z  ∼  4, z  ∼  5, and z  ∼  6 Lyman-break samples were selected based on catalogs where the PSF matching was done to the ACS z₈₅₀ band. This approach offers two clear advantages. (1) By selecting sources based on catalogs where the PSF matching is done to the z₈₅₀ band, therefore making our photometric apertures smaller, we utilize much higher S/N photometry for the selection process, improving the robustness of our z  ∼  4, 5, and 6 selections. (2) By using different photometry (and smaller apertures) to select sources than we use to estimate β (typically involving larger apertures), we ensure that our β measurements are not subject to exactly the same noise as affects source selection (because the photometric apertures are different and will have somewhat different noise characteristics). This makes our z  ∼  4–6 β measurements less susceptible (by ∼50%) to the biases that arise from a coupling of errors in our β measurements with similar errors in our photometric selection ("photometric error coupling bias"; Appendix B.1.2). For our z  ∼  7 Lyman-break sample, however, the photometry used to estimate β is the same as what we use to select the sources. See Section 3.1 for a description of the photometry.

In Appendix A, we show that the present approach produces similar results to that based on UV colors spanning the wavelength range 1600 Å to 2200 Å (but with much smaller errors). For typical sources, the derived β's agree to within δβ ∼ 0.1 which is about as well as we can determine β using UV colors (where possible systematics in β estimates are on the order of ∼0.1). The present UV-continuum slope β estimates should therefore be directly comparable with others in the literature (Bouwens et al. 2006, 2009, 2010a; Meurer et al. 1999; Ouchi et al. 2004a; Stanway et al. 2005; Finkelstein et al. 2010) where the wavelength baseline 1600 Å to 2200 Å was typical for β determinations.

In deriving β's for sources in our samples from the available photometry, we tested for significant deviations from a pure power-law shape. This is important for ensuring that the model we use to characterize the UV-continuum SEDs of sources is meaningful. For the average source, we found that the observed photometry and the best-fit power-law SED differed by less than ∼0.02 mag. Such deviations are not large enough to have a sizeable impact on our β determinations, shifting the observed β determinations by ≲ 0.1. Such changes are well within the uncertainties we quote on our β determinations (Section 3.5).

To establish the actual distribution of UV-continuum slopes β from the observations, we must account for the effect that object selection and measurement have on the observed distribution. These effects can cause the observed distribution to look very different from what it is in reality. Examples of such effects include (1) our LBG selection criteria preferentially including those galaxies in our samples with the bluest UV-continuum slopes and (2) photometric noise increasing the spread in the measured UV-continuum slopes β.

To correct for these selection and measurement biases, we follow a very similar procedure to what we used in our first major study of the UV-continuum slope β distribution at high redshift (Bouwens et al. 2009). We create mock catalogs of galaxies, generate artificial images for each source in our catalogs, add these sources to the real observations, and then reselect the sources and measure their properties in the same way as with the real observations. We then determine the approximate biases that we would expect in the mean UV-continuum slope and 1σ scatter due to source selection or photometric scatter. We compute these biases as a function of magnitude for each of our data sets. To ensure that we are able to determine the relevant biases to high precision, we repeat these experiments for >10⁵ artificial galaxies.

Appendix B provides a detailed discussion of the simulations we use to estimate the likely biases in the UV-continuum slope β distribution and establish the needed corrections.

Figure 5 shows the results of these simulations for z  ∼  7 selections. The effect of biases related to the LBG selection itself (Figure 21) and due to a possible coupling between the selection and measurement processes ("photometric error coupling bias") is shown explicitly. It is immediately clear from this figure that we can successfully recover the input values of the UV-continuum slopes β and do not expect especially large biases in the derived slopes. We include similar panels for our z  ∼  4, z  ∼  5, and z  ∼  6 selections in Figure 22 of Appendix B.

Zoom In Zoom Out Reset image size

An important check to perform in assessing the overall quality of our UV-continuum slope determinations here is to compare the results from our deep data with the results from our wide data. Such a check is provided in Figure 25 of Appendix B, and it seems clear that the results from the HUDF09 and ERS+CANDELS data sets are in good agreement (Δβ ≲ 0.1–0.2) over the luminosity range where they overlap and have good statistics. This suggests that the UV-continuum slope β distributions we recover are accurate and free of any sizeable biases.

We can utilize the UV-continuum slopes β measured for individual sources (Section 3.3) to establish the distribution of UV-continuum slopes β as a function of luminosity for each of our z  ∼  4–7 LBG samples.

To establish this distribution while properly accounting for the relevant selection and measurement biases, we analyze the UV-continuum slope distribution for each data set and Lyman-break sample separately. We determine the mean UV-continuum slope and 1σ scatter. We correct the mean and scatter for the relevant biases (Appendix B and Section 3.4). Typical corrections in the mean UV-continuum slope β are ≲0.1 in general; the corrections do, however, reach values of Δβ ∼ 0.1 near the selection limit of our shallow and deep probes.

When determining the mean and scatter for each distribution, we use the biweight mean C_BI

$$\begin{eqnarray*} &&C_{{\rm BI}} = M + \frac{\sum _{|u_i|<1} (x_i - M)\big(1-u_i ^2\big)^2}{\sum _{|u_i|<1} \big(1-u_i ^2\big)^2} \end{eqnarray*}$$

and biweight scale S_BI

$$\begin{eqnarray*} &&S_{{\rm BI}} = n^{1/2} \frac{\big[ \sum _{|u_i|<1} (x_i - M)^2 \big(1-u_i^2\big)^4\big]^{1/2}}{|\sum _{|u_i|<1} \big(1-u_i ^2\big)\big(1- 5u_i ^2\big)|}, \end{eqnarray*}$$

where u_i = (x_i − M)/(c MAD), MAD = median(|x_i − M|), M is the median of the x_i's, x_i is the data, n is the number of sources summed over, and c is the "tuning constant" (Beers et al. 1990). c is taken to be equal to 6 in computing the biweight mean and 9 in computing the biweight scale. Use of robust statistics like the biweight mean and scale is valuable, given that a small fraction of the sources in our samples may be contaminated by light from nearby sources or may lie outside the target redshift range. Biweight mean determinations of β are similar to median determinations of β (median difference Δβ ∼0.01 with ∼0.04 scatter), but somewhat bluer (Δβ ∼ 0.1) than mean determinations (due to the β distribution having a somewhat extended tail toward red β's and the biweight mean de-weighting the tail).

For the absolute magnitude of each galaxy in the UV, we use the mean absolute magnitude of that source in all the HST bands that contribute to its UV-continuum slope determinations (Table 3). By using the geometric mean for this luminosity, we aim to avoid giving too much weight to the bluer or redder bands in defining this luminosity and hence artificially creating a correlation between the derived UV-continuum slope β and UV luminosity.

Figure 6 shows the UV-continuum slope β distribution versus the magnitude for our z  ∼  4, z  ∼  5, z  ∼  6, and z  ∼  7 samples. Scatter in the β distribution is relatively modest for our z  ∼  4 and z  ∼  5 samples, but is larger for our z  ∼  6 and z  ∼  7 samples. This is the result of the fact that our z  ∼  4–5 samples use a larger number of passbands and larger wavelength baseline to measure the UV-continuum slope β.

Zoom In Zoom Out Reset image size

In Table 4, we present our corrected determinations of the biweight mean UV-continuum slope and 1σ scatter as a function of UV luminosity, for our z  ∼  4, z  ∼  5, z  ∼  6, and z  ∼  7 samples. We also include an estimate of the uncertainties in the biweight mean UV-continuum slopes that result from the measurement uncertainties, intrinsic scatter in the β distribution, and small number statistics.

For the systematic error on our UV-continuum slope measurements, we conservatively adopt values of 0.10–0.28. The dominant component of this error comes from uncertainties in our photometry. Allowing for small errors in the photometric zero points, aperture corrections, and PSF matching, we estimate a maximum error of 0.05 mag in our flux measurements. This translates into possible systematic errors of 0.10, 0.10, 0.14, and 0.28 in the β measurements at z  ∼  4, z  ∼  5, z  ∼  6, and z  ∼  7.

We also allow for small systematic errors in our UV-continuum slope β measurements as a result of differences in the wavelength baseline we use to derive β (Table 3). These differences would have an effect on the β's we derive from the available photometry, if the UV SEDs differed substantially from a pure power-law form (i.e., f^β_λ). The tests we perform in Appendix A suggest that the relevant systematic errors are small Δβ ≲ 0.13. Similar β's are found from fits using the full UV continuum (Section 3.3 and Figure 4) as are found using a smaller wavelength range.

We observe a clear trend in the UV-continuum slope β as a function of UV luminosity in all four LBG samples, such that the UV-continuum slope β becomes progressively bluer at fainter luminosities (Figures 6 and 25). Such trends in β had already been identified by Bouwens et al. (2009) in their analyses of z  ∼  2.5 and z  ∼  4 galaxy samples (see also Meurer et al. 1999; Labbé et al. 2007; Figure 8 of Overzier et al. 2008) and by Bouwens et al. (2010a) in their analyses of z  ∼  5–7 galaxies (see also Wilkins et al. 2011). The observed luminosity dependence of β is thought to be due to a change in the dust content and perhaps the age of galaxies as a function of luminosity (Labbé et al. 2007; Bouwens et al. 2009). It is likely to be a manifestation of the well-established mass–metallicity relationship seen at a wide variety of redshifts (e.g., Tremonti et al. 2004; Erb et al. 2006a; Maiolino et al. 2008).

The UV-continuum slope β shows an approximately linear relationship on the magnitudes M_UV of the sources (Figure 6), and therefore it makes sense for us to model this relationship using a first-order polynomial and determine the best-fit slopes dβ/dM_UV and intercepts $\beta _{M_{\rm UV}=-19.5}$.

In fitting a line to our UV-continuum slopes, we use a finer binning (i.e., 0.5 mag) than shown in Figure 6 to minimize the impact of the binning scheme on the best-fit slopes and intercept. The results of our fits to the biweight means are shown as blue lines in Figure 6, and it is clear that the mean β's are well fit by the lines.⁸ The fits exhibit a similar dependence on the UV magnitude M_UV at each redshift. The intercept to the lines is also similar at all redshifts but appears to evolve monotonically with cosmic time. The best-fit parameters—slope and intercept—we derive for these lines are presented in Table 5 and Figure 7. These determinations are in excellent agreement with previous work, as we discuss in Section 4.

Zoom In Zoom Out Reset image size

Table 5. Best-fit Slopes and Intercepts to the UV-continuum Slope β to UV Luminosity Relationship

Dropout	Mean
Sample	Redshift	$\beta _{M_{\rm UV}=-19.5}$	dβ/dMUV
B435	3.8	−2.00 ± 0.02 ± 0.10	−0.11 ± 0.01a
V606	5.0	−2.08 ± 0.03 ± 0.10	−0.16 ± 0.03a
i775	5.9	−2.20 ± 0.05 ± 0.14	−0.15 ± 0.04a
z850	7.0	−2.27 ± 0.07 ± 0.28	−0.21 ± 0.07
U300b	2.5	−1.70 ± 0.07 ± 0.15	−0.20 ± 0.04

Notes. Section 3.6, see also Figure 7. ^aWhile the dβ/dM_UV slopes for our z  ∼  4, 5, and 6 samples are quite similar overall, these slopes show even better agreement if we consider the correlations over the same luminosity range, i.e., excluding β determinations faintward of −17 AB mag. In this case, we find a dβ/dM_UV slope of −0.13 ± 0.02 for our z  ∼  4 sample. ^bFrom Bouwens et al. (2009).

Download table as:  ASCIITypeset image

Before the WFC3/IR camera on HST became operational, Bouwens et al. (2009) made use of the ACS+NICMOS observations to quantify the distribution of UV-continuum slopes β for star-forming galaxies over range in redshift (z  ∼  2.5–6) and luminosities. Bouwens et al. (2009) derive β directly from the UV colors (e.g., as in Appendix A). How do the present UV-continuum slope determinations compare with those from Bouwens et al. (2009)? Both old and new results are shown in Figure 8. Comparing the UV-continuum slope measurements made on the identical sources in the two data sets (old ACS+NICMOS data versus the new WFC3/IR observations), we find reasonable agreement, with mean offsets β_WFC3/IR − β_{ACS + NICMOS} of only −0.10, 0.04, 0.09 in the derived β's at z  ∼  4, z  ∼  5, and z  ∼  6, respectively. The biweight mean UV-continuum slopes β found here are also in good agreement with the slopes derived by Bouwens et al. (2009). At z  ∼  4 and z  ∼  5, the present biweight mean β's are just ∼0.26 bluer and ∼0.16 redder than the values found in Bouwens et al. (2009); the present biweight mean β's at z  ∼  6 show no average shift at all relative to the values found by Bouwens et al. (2009). Δβ ∼ 0.1 of the differences result from our use of the more robust biweight means to express the central β's, so the agreement is quite good overall.

Bouwens et al. (2010a) took advantage of the first-year WFC3/IR observations over the HUDF and the ERS WFC3/IR observations to quantify the UV-continuum slope β distribution to higher redshifts (z  ∼  7) and lower luminosities. How do the present determinations of the UV-continuum slope β compare with Bouwens et al. (2010a)? At high luminosities, we find excellent agreement, with both studies preferring mean UV-continuum slopes β of −2 (see Figure 8). At lower luminosities, however, we now find somewhat redder values of the UV-continuum slope β, i.e., Δβ ∼ 0.2–0.4 than found by Bouwens et al. (2010a). The observed differences in the mean UV-continuum slope β seem to have resulted from both the small number of sources in previous samples and uncertainties in the photometry of faint sources in the early HUDF09 data. We remark that the current WFC3/IR observations of the HUDF from the HUDF09 WFC3/IR program are approximately twice as deep as what Bouwens et al. (2010a) used.

In the present analysis, we find that the UV-continuum slope β shows a correlation with the redshift of galaxies, with higher redshift galaxies being bluer. This evolution is evident both in Figure 6 (compare the solid lines in the z  ∼  5, z  ∼  6, z  ∼  7 panels with the dotted lines) and in Figure 7.

How does the correlation we find compare with other studies? A brief summary of the evidence for this redshift dependence is provided in Figure 8, and there is a clear trend from bluer UV-continuum slopes β at z  ∼  6–7 to redder slopes at z  ∼  2–4. Essentially all studies of the UV-continuum slope over the range z  ∼  3–7 (Lehnert & Bremer 2003; Stanway et al. 2005; Bouwens et al. 2006, 2009, 2010a; Wilkins et al. 2011; Finkelstein et al. 2012) find evidence for this evolution of β with redshift. The only exception to this is the recent study of Dunlop et al. (2012) who find no dependence (but Dunlop et al. 2012 do note that a comparison of their results with those at z  ∼  3 does argue for evolution); we discuss the Dunlop et al. (2012) results in Section 4.5.

As we discuss in Section 3.6, we find that the UV-continuum slope β shows a clear correlation with the rest-frame UV luminosity of galaxies, with lower luminosity galaxies being bluer at all redshifts. A good illustration of this correlation is provided in Figure 6. The best-fit relationship (solid blue line) shows almost exactly the same dependence for each high-redshift sample. The uniformity of the slope and modest variation as a function of redshift is also clearly illustrated in Figure 7. The uniformity extends even to our z  ∼  4 samples. While β in these samples appear to show a somewhat weaker dependence on luminosity than our other samples, a slightly steeper dependence is found, i.e., dβ/dM_UV = −0.13 ± 0.02, if we use the same luminosity baseline as our z  ∼  5 and z  ∼  6 samples (excluding sources faintward of −17 AB mag).

In general, the correlation we find agrees very well with most previous studies. Figure 9 shows different β determinations as a function of luminosity at z  ∼  4–5 (Ouchi et al. 2004a; Overzier et al. 2008; Bouwens et al. 2009; Lee et al. 2011; Wilkins et al. 2011; Dunlop et al. 2012).⁹ Despite some dispersion in the precise values of the β measurements and some small variance in the best-fit slopes dβ/dM_UV, nearly all published determinations of the UV-continuum slope β at z  ∼  4–5 find bluer values for β at lower luminosities, with the same luminosity dependence. Two recent analyses that found minimal or no correlation of β with luminosity are those of Dunlop et al. (2012) at z  ∼  5–7 and Finkelstein et al. (2012) at z  ∼  4–5. We discuss the Dunlop et al. (2012) results in Section 4.5 and the Finkelstein et al. (2012) results in Section 4.7.

Zoom In Zoom Out Reset image size

Similar correlations with luminosity are found in the UV-continuum slope results at z  ∼  2–3 (Labbé et al. 2007; Bouwens et al. 2009; Sawicki 2012) and in the UV–optical colors (Papovich et al. 2001; Labbé et al. 2007; González et al. 2012). Again, two analyses did not find a correlation of β with luminosity, and those are the Adelberger & Steidel (2000) and Reddy et al. (2008) analyses at z  ∼  2–3. In these cases, not only was the luminosity baseline too short to provide much leverage for quantifying this correlation, but also the luminosity range probed was around L* which is where the dispersion in dust properties relative to UV luminosity is at a maximum (e.g., Figure 13 of Reddy et al. 2010). However, when one adds the Adelberger & Steidel (2000) samples to the faint (∼26–27 mag) z  ∼  2.5 samples observed by Bouwens et al. (2009), a strong correlation is present (Figure 3 of Bouwens et al. 2009; see also Figure 5 of Sawicki 2012). Taken together these results indicate that there is also a clear trend with luminosity at z  ∼  2–3.

Dunlop et al. (2012) use the WFC3/IR observations over the HUDF09 and ERS observations to quantify the UV-continuum slope distribution at z  ∼  5–7. They select sources using a photometric redshift procedure and then measure their UV-continuum slopes β from the UV colors. In both respects, their procedure differs from that followed here (Sections 3.2 and 3.3). Overall, the individual UV-continuum slope β measurements of Dunlop et al. (2012) are in reasonable agreement with the present results (see Figure 8). The mean β they derive for bright z  ∼  5 galaxies is Δβ ∼ 0.3 bluer and the mean β they derive for faint z  ∼  6 galaxies is Δβ ∼ 0.4 redder. Given the quoted uncertainties on the measurements, the differences are not particularly significant.

The main differences arise when we look at the trends in the UV-continuum slope β with redshift and luminosity. Dunlop et al. (2012) find a mean UV-continuum slope β of galaxies equal to ∼ − 2.1, with no significant dependence on luminosity or redshift. This is in contrast to the strong correlation we find of β with both luminosity (Figures 6 and 7) and redshift (Figure 7). Their results are also inconsistent with the trends reported in the literature (Sections 4.3– 4.4; Figures 7 and 9).

We have attempted to understand the source of this difference both by a qualitative assessment of some issues that we know are important for obtaining reliable results and by a quantitative assessment using simulations (see Appendix D and Section 4.6). First, we remark that Dunlop et al. (2012) make no attempt to correct their mean β's for the fact that galaxies with bluer UV-continuum slopes β are easier to select than galaxies with redder UV-continuum slopes β. This effectively biases their UV-continuum slopes β to bluer values. It is a fairly straightforward process to correct for this issue (e.g., we describe such a correction in Appendix B.1.1 and Figure 5 where we show its effect).¹⁰ Second, Dunlop et al. (2012) use an overlapping set of information both to select sources and to measure the UV-continuum slopes β. This biases their results (see Appendix D), though the magnitude of this bias is mitigated by their consideration of only the brightest sources.

Third, one further limitation of the Dunlop et al. (2012) analysis is their exclusion of the lowest luminosity sources. This significantly reduces the leverage they have to quantify trends in the mean UV-continuum slope as a function of luminosity. Dunlop et al. (2012) exclude faint sources because their simulations suggested that β could not be measured in an unbiased way. However, as we show in Appendices B and D, we are able to obtain reliable measurements (see also Section 4.6). Not only can we recover the mean UV-continuum slope for our samples to very faint magnitudes with excellent accuracy, but we can recover these slopes for samples selected using a photometric redshift procedure, if the information used for source selection is clearly separated from that used to measure β (see Figure 10).

Zoom In Zoom Out Reset image size

As mentioned above, these three issues led us to consider a more quantitative evaluation of the Dunlop et al. (2012) procedure so that we could better understand what was happening with their measurement of β. This is discussed in more detail in the next section.

As discussed above, a key question that has arisen in recent papers is whether it is possible to determine the UV-continuum slope β distribution to very low luminosities with small biases. Dunlop et al. (2012), in particular, have suggested that it is not possible and have supported this suggestion with a series of simulations where they add noise to model galaxies and reselect these galaxies with their photometric redshift code. Dunlop et al. (2012) argue that noise in the photometry combined with a preference for selecting sources with blue colors would result in highly biased estimates for the mean UV-continuum slope β. Dunlop et al. (2012) show that such a bias toward bluer UV-continuum slopes at faint magnitudes is present in their mock data sets. Dunlop et al. (2012) argue that similar biases are likely present in Lyman-break selections, without further substantiating this claim.

We agree that selection biases can affect the distribution of UV-continuum slopes β. However, the size of these selection biases is extremely dependent upon how one selects the galaxies and measures their UV-continuum slopes β. As we show in Figure 10, our Lyman-break selections yield much smaller selection biases overall. The reason we expect biases in our selections to be small is that we select galaxies using a different part of the rest-frame UV SED (the bluest two passbands redward of the break) than we use to measure UV-continuum slopes (the second bluest passband and redder). As a result, photometric scatter in the bands we use to measure the UV-continuum slope β is largely independent of similar scatter in the bands we use for selection. Therefore, we would not expect our β measurements to be significantly biased as a result of the object selection process.¹¹

Achieving similarly small biases to faint magnitudes with a photometric redshift technique is also possible. However, one must again be careful to use different information to select sources from what one uses to measure the UV-continuum slope β. To illustrate this, we consider the situation at z  ∼  7 in the current HST ACS + WFC3/IR data set. We have run simulations where we attempt to measure the mean UV-continuum slope β for a set of z  ∼  7 galaxies selected using a photometric redshift procedure. Both the simulations and results are discussed in Appendix D. We consider (1) the case where five HST bands are used to select sources and determine redshifts (this excludes those bands used to measure β), (2) the case where six bands are used (and so now including one of the two bands used to measure β), and (3) the case where all seven HST bands are used for selection and redshift determination (and so both bands used to measure β are also included to measure redshifts and select the sources). This latter approach is basically what Dunlop et al. (2012) do.

The results are shown in Figure 26 of Appendix D. While β measurements show substantial biases when all seven bands are used for the photometric redshift estimates (similar to the procedure of Dunlop et al. 2012) β measurements made using five or six bands show much smaller biases. This demonstrates that the UV-continuum slope β can be measured with very small biases to faint magnitudes.

In an independent analysis, Finkelstein et al. (2012) also use the recent WFC3/IR observations over the HUDF and CDF-South GOODS field to quantify the UV-continuum slope distribution for star-forming galaxies at z  ∼  4–8. Finkelstein et al. (2012) split sources into five different redshift samples using a photometric redshift procedure and then estimate β by finding the model SED which best fits the photometry of a source and deriving β from this model. Finkelstein et al. (2012) find that β shows a clear dependence on redshift, but report only a limited dependence on the UV luminosity.

β versus luminosity trends. While the redshift dependence Finkelstein et al. (2012) find for β is in excellent agreement with what we find (compare the solid black squares and large blue circles in the lower panel of Figure 8), the luminosity dependence Finkelstein et al. (2012) observe would appear to be considerably weaker. After all, Finkelstein et al. (2012) report no significant correlation of β with UV luminosity in their baseline analyses of their five redshift samples—seemingly very different than the clear correlation of β with luminosity we report.

Despite these apparent differences, the overall results from our two studies are actually in fairly good agreement. For example, with regard to the z  ∼  6 and z  ∼  7 samples, the best-fit dβ/dM_UV values Finkelstein et al. (2012) find, i.e., −0.10 ± 0.07 and −0.20 ± 0.11, respectively, are strikingly similar to the values we find, i.e., −0.15 ± 0.04 and −0.21  ±  0.07, respectively.

For the z  ∼  4 and z  ∼  5 samples, Finkelstein et al. (2012) do not report a significant correlation between β and UV luminosity in what they identify as their baseline analysis (finding dβ/dM_UV values of 0.01 ± 0.03 and 0.00 ± 0.06, respectively). However, in this analysis, Finkelstein et al. (2012) only make use of a 2 mag baseline in luminosity (due to their binning scheme), with their low luminosity anchor point largely coming from the relatively shallow ∼1.6 orbit CANDELS data (with only a small contribution from the HUDF data). The CANDELS observations are clearly poorly suited to determine the trend in β to very low luminosities, given the very low S/N's and potentially large biases expected for the faintest sources in the CANDELS fields. However, the situation changes if we take advantage of the additional leverage in luminosity provided by the β measurements they provide for faint sources in the ultra-deep HUDF observations.

We can check this by extracting the HUDF measurements from Figure 5 of their paper. By comparing their median β measurements for brighter galaxies with their median β measurements for fainter galaxies in the HUDF, we find evidence for a significant correlation with luminosity. We find dβ/dM_UV trends of −0.06 ± 0.02 and −0.13 ± 0.04, respectively. In the final version of their paper, Finkelstein et al. (2012) also note a similar correlation with luminosity making use of the faintest HUDF sources, finding dβ/dM_UV trends of −0.07 ± 0.01 and −0.09 ± 0.03, respectively. While not in exact agreement with the trends we derive based on our own β measurements, i.e., −0.11 ± 0.01 and −0.16 ± 0.03, respectively, the agreement is much better. Use of the faint sources in the HUDF is important to take full advantage of the available leverage in luminosity to quantify the β versus M_UV trend.

In Figure 11, we show the β versus M_UV trend that Finkelstein et al. (2012) find in their baseline z  ∼  4 and z  ∼  5 analyses (magenta lines) and the trend we find from their measurements making exclusive use of sources from the HUDF to constrain β to fainter magnitudes (dashed red lines). In addition, we show the median β's Finkelstein et al. (2012) find for two fainter z  ∼  4 and z  ∼  5 subsamples within the HUDF (large solid red circles: we can extract β measurements for individual sources within the HUDF from their Figure 5) and the bootstrap uncertainties on these medians. In both samples there is a clear trend in the median β's toward bluer values at the lowest luminosities. It is striking how well the median β's Finkelstein et al. (2012) derive from the HUDF agree with our own β measurements, particularly in the luminosity interval [−20 mag, −18 mag]. While it is true that the median β's Finkelstein et al. (2012) derive for their entire ERS+CANDELS+HUDF09 sample are redder in general than what we find at these luminosities, these median β's receive their largest weight from the shallower ERS+CANDELS samples and therefore may be subject to substantial selection, measurement, or contamination biases (see the discussion at the end of this section).

Zoom In Zoom Out Reset image size

Expected trend in luminosity. Finkelstein et al. (2012) defend the weak correlation of β they report versus luminosity (particularly as derived in their baseline analysis), arguing that β should show a stronger correlation with stellar mass than UV luminosity. We do not dispute this assertion; however, it would be most surprising if a correlation of β with stellar mass did not also appear as a correlation with UV luminosity. Given the correlation found between SFR and stellar mass in high-redshift galaxies (e.g., Stark et al. 2009; González et al. 2011; McLure et al. 2011; Lee et al. 2012; Reddy et al. 2012b), we would expect the UV luminosity to be broadly correlated with stellar mass. A similar conclusion can be drawn from high-redshift angular correlation function results. Higher luminosity galaxies are consistently found to be more clustered than lower luminosity galaxies (e.g., Ouchi et al. 2004b; Lee et al. 2006). Such would not be the case if UV luminosity was not correlated with mass (in this case halo mass). Independent of these considerations, we remark that β also shows a clear correlation with UV luminosity in various cosmological hydrodynamical simulations (e.g., Finlator et al. 2011; Dayal & Ferrara 2012), and the predicted trends (e.g., dβ/dM_UV ∼ −0.10 is expected in the Finlator et al. 2011 simulations) are comparable to what we find (Figure 7).

Possible biases in the Finkelstein et al. measurements. Finkelstein et al. (2012) have suggested that the β–M_UV trends we find may be stronger than what they find due to the fact that we measure the UV luminosity at a different rest-frame wavelength than they do, and β versus M_UV trends may depend on this rest-frame wavelength. For our z  ∼  4 samples, for example, we measure the rest-frame UV luminosity at 2041 Å (see Table 3) while Finkelstein et al. (2012) measure it at 1500 Å. We are in full agreement that dβ/dM_UV will depend on the rest-frame wavelength. However, the results from Figure 3 in Labbé et al. (2007) suggest that one would find an even stronger dβ/dM_UV trend at bluer wavelengths than one would find at redder wavelengths, which is different from what Finkelstein et al. (2012) find. It is therefore not clear if this explains the differences.

Instead, one might be concerned that the β versus dM_UV trend Finkelstein et al. (2012) find may be biased as a result of the rest-frame wavelength Finkelstein et al. (2012) use to measure the M_UV luminosity.¹² By determining luminosity at the blue end (at 1500 Å) of the wavelength baseline they use to derive β, Finkelstein et al. (2012) effectively introduce a coupling between the errors that affect both their β measurements and their determinations of the UV luminosity M_UV. This could be problematic since any errors in the flux measurements of sources would cause sources to be either fainter and redder or brighter and bluer, causing the dβ/dM_UV trend derived by Finkelstein et al. (2012) to be biased toward too high of values. Repeating the determination of dβ/dM_UV at z  ∼  4–5 based on our own flux measurements but basing M_UV on the flux measurement at the blue end of the wavelength baseline to derive β and using only the wide-area CANDELS+ERS sources, we estimate that this could bias the derived dβ/dM_UV trend too high by Δ(dβ/dM_UV) ∼ 0.05. This bias would be analogous to the photometric error coupling bias we discuss in Appendix D. However, instead of the coupling being between source selection and the β measurements, it would be between the measurement of β and the measurement of the UV luminosity M_UV.¹³ Contamination in the shallower CANDELS+ERS samples (from lower redshift galaxies) could also be an issue for Finkelstein et al. (2012) in deriving the trend in β to lower luminosities M_UV.

Given the much smaller flux uncertainties for sources in the HUDF (and smaller contamination rates), we would expect the β's Finkelstein et al. (2012) measure there to show significantly smaller biases than β's measured for sources in ERS+CANDELS fields at the same luminosities. Encouragingly enough, the median β's Finkelstein et al. (2012) derive from the HUDF are in good agreement with our own, particularly in the luminosity interval [−20 mag, −18 mag] (compare the blue squares and large red circles in both the upper and lower panels of Figure 11).

Finally, we remark Finkelstein et al. (2012) use the same photometry both to select sources and measure β. Given the argumentation in the previous section and Appendix D (see also Dunlop et al. 2012), we might expect the bias in the derived β's to be non-zero. However, in practice, given the large number of passbands used to select sources and measure their redshifts, the bias is likely to be quite small except at z  ∼  7 (similar to our study: but see Section 4.8). Finkelstein et al. (2012) are aware of this issue and explicitly discuss it in their paper.

Galaxies with the most extreme UV properties are expected to lie at very high redshift and have low luminosities—given the early cosmic times in which they are observed and likely low masses. It has therefore been of considerable interest to establish the UV-continuum slopes β for the faintest observable z  ∼  7–8 galaxies. Early observations of such galaxies in the HUDF gave tantalizingly steep values of the UV-continuum slope β, i.e., β ∼ −3.

In the present study, we find a mean UV-continuum slope β of −2.7 ± 0.2 for faint z  ∼  7 galaxies. This is slightly redder than the mean UV-continuum slope β (−3.0 ± 0.24) found in our earlier study of β for z  ∼  7 galaxies in the HUDF. It is also slightly redder than that (−3.0 ± 0.5) found by Finkelstein et al. (2010) using the same data.¹⁴ Dunlop et al. (2012) have argued that it is not possible to estimate the UV-continuum slope β at such low luminosities, but as we discuss in Section 4.6, such measurements are possible if care is taken to minimize biases by ensuring the information used for selection is independent of that used to measure β.

How robust are our measurements of the mean UV-continuum slope β for faint z  ∼  7 galaxies? While our simulations suggest that the biases are not large, it is useful to check this result by obtaining an independent estimate of the mean β. For this estimate, we use two completely independent data sets to select sources and to measure β. Source selection is done using the first-year WFC3/IR observations over the HUDF (18 orbits Y₁₀₅, 16 orbits J₁₂₅, 28 orbits H₁₆₀) from the HUDF09 program while the measurement of β is done using the second-year WFC3/IR observations over the HUDF (18-orbit J₁₂₅, 25-orbit H₁₆₀ data). Since the second-year data were not used to select the sample, we can make an unbiased measurement of the UV-continuum slope β for this sample using the new observations.¹⁵ For comparison with the Bouwens et al. (2010a) study on β, we use the same z  ∼  7 sample. The biweight mean β we derive for our lowest luminosity (M_{UV, AB} ∼ −18.5) subsample is −2.8 ± 0.2. While this is slightly redder than the β = −3.0 ± 0.2 we find from the first-year observations (Bouwens et al. 2010a), this completely independent and unbiased estimate does suggest that the UV-continuum slopes β for lower luminosity galaxies at z  ∼  7 are very blue.

As one final check on the mean β for faint z  ∼  7 galaxies, we considered one variation on the previous test. We divided the J₁₂₅-band data for each HUDF09 field (HUDF09, HUDF09-1, HUDF09-2) into two disjoint subsets and produced separate J₁₂₅-band reductions from each. The Y₁₀₅-band data and first third of the J₁₂₅-band data were used for the selection of z  ∼  7 sources, and the H₁₆₀-band data and final two-thirds of the J₁₂₅-band data were used to measure the UV-continuum slopes β. As in the previous test, this is to ensure that the information used for source selection is completely independent of that used for the β measurements, and therefore the photometric error coupling bias (discussed in Section 4.6 and Appendix B.1.2) must be zero. The biweight mean β we derive for the faintest z  ∼  7 sample (M_{UV, AB} ∼ −18.3) based on the three HUDF09 fields (HUDF09, HUDF09-1, and HUDF09-2) is −2.7 ± 0.2, again consistent with our other estimates.

This new determination of the mean UV-continuum slope β for very low luminosity z  ∼  7 galaxies is very blue (i.e., β ∼ −2.7  ±  0.2). However, this does not appear to be especially anomalous. In fact, it appears to be consistent with what one might expect extrapolating the z  ∼  4–6 UV-continuum slope β relationship to z  ∼  7 (see, e.g., dashed line in the z  ∼  7 panel to Figure 6). The observed correlations of β with both redshift and luminosity are such that one would expect faint z  ∼  7 sources to be very blue. This blue β is also not inconsistent with what one can achieve with standard stellar population modeling (where the UV-continuum slope β can become as steep as ∼−2.7; Schaerer 2003; Bouwens et al. 2010a; Robertson et al. 2010). The observed UV-continuum slopes β at z  ∼  7 therefore seem to provide no particularly compelling evidence for exotic stellar populations, i.e., very low metallicity stellar populations or a high escape fraction (see also Finkelstein et al. 2012). Bouwens et al. (2010a) briefly speculated as to what such blue β's might imply, if future observations confirmed that β was really as blue as −3 with small uncertainties.

In the previous section, we presented evidence that the UV-continuum slopes β of star-forming galaxies at high redshift were distributed along a well-defined sequence in UV luminosity. The intrinsic scatter in β along the sequence is small (σ_β ∼ 0.34), and the dependence of β is such that galaxies become bluer toward lower luminosities. Such a sequence is particularly prominent in our z  ∼  4 sample, but all of our higher redshift selections (z  ∼  4, z  ∼  5, z  ∼  6, z  ∼  7) show strong evidence for such a sequence as well (Figure 6).

The fact that we observe the same dependence of β on UV luminosity in each of four Lyman-break selections suggests the trend we are recovering from the observations is real. Such a color–magnitude relationship was already evident in many studies of galaxies over the redshift range z  ∼  2–5, though the clearest evidence was presented by Bouwens et al. (2009) for z  ∼  2.5 and z  ∼  4 samples and Labbé et al. (2007) for z  ∼  1.05, z  ∼  1.8, and z  ∼  2.7 samples (see Section 4.4). The present work confirms these trends and extends them to z  ∼  5, z  ∼  6, and z  ∼  7.

Consistent with these trends, the UV-continuum slopes β we measure for the most luminous (and presumably most massive) galaxies have similarly red values for β of −1.6 to −1.9 in all four redshift samples examined here (Table 4; similar values were also found by Lee et al. 2011 and Willott et al. 2012) while the lowest luminosity galaxies probed have relatively similar blue values for β of −2.2 to −2.7 for all four samples.

The existence of well-defined sequences for both star-forming and evolved galaxies at lower redshift is now very well established. In the local universe (z  ∼  0.1), for example, Salim et al. (2007) find that galaxies fall along a well-defined sequence in SFR versus stellar mass. Similar star-forming sequences were found by Noeske et al. (2007) and Martin et al. (2007) at somewhat higher redshifts, from z  ∼  0.2 to z  ∼  1 (0.3 dex scatter in the SFRs). find evidence for such a sequence in star-forming galaxies at z  ∼  2, and Elbaz et al. (2011) show that such a sequence exists for even more luminous systems from recent Herschel observations. While the present color–magnitude sequence we observe is not an SFR versus stellar mass sequence, the existence of such a sequence suggests that galaxies build up in a relatively well-defined way versus cosmic time.

Theoretically, we would expect such a sequence due to the build up in metals and dust anticipated to occur as galaxies grow in luminosity and mass. One useful illustration of this can be seen in some recent work by Davé et al. (2006) and Finlator et al. (2011) who use smooth particle hydrodynamics to model the evolution of galaxies to z  ∼  6. In Figure 6 of Davé et al. (2006), for example, we see a clear mass dependence in the metallicity of galaxies, with ∼0.3 dex change in metallicity per ∼1 dex change in mass.

Figure 7 of Finlator et al. (2011) shows the expected trends in UV-continuum slopes β as a function of luminosity, including the effects of starlight, dust, and emission lines. Finlator et al. (2011) predict mean UV-continuum slopes of ∼−2.02 for luminous (M_{UV, AB} ∼ −20.5) z  ∼  7 galaxies and ∼−2.28 for lower luminosity (M_{UV, AB} ∼ −18.5) z  ∼  7 galaxies—equivalent to an approximate slope to the β–M_UV relationship of just ∼−0.13. A fit to the β versus M_UV relationship for all the z  ∼  4, z  ∼  5, z  ∼  6, and z  ∼  7 sources in the Finlator et al. (2011) simulations yields dβ/dM_UV slopes of ∼−0.10, ∼−0.08, ∼−0.13, and ∼−0.09, respectively (K. Finlator 2011, private communication). The mean slope to this relationship dβ/dM_UV of ∼ − 0.10 is in excellent agreement with what is observed (Figure 7 and Table 4). To better illustrate this, we include a comparison of the UV-continuum slopes β observed with that predicted from the Finlator et al. (2011) simulations (Figure 12).

Zoom In Zoom Out Reset image size

Rest-frame optical studies of z  ∼  4–6 galaxies with Spitzer IRAC have provided evidence for a similar sequence at z ⩾ 4 for star-forming galaxies (Stark et al. 2009; Labbé et al. 2010b; González et al. 2011). Typical star-forming galaxies at z  ∼  5 have rest-frame UV–optical colors of ∼0.6 mag, with a scatter of 0.5 dex (González et al. 2011, 2012). These UV–optical colors appear to show a slight dependence on luminosity, in the sense that brighter galaxies are redder and fainter galaxies are bluer (González et al. 2012). So far there is no evidence for evolution in the UV–optical colors over the redshift range 4 < z < 6 (Stark et al. 2009; Labbé et al. 2010b; González et al. 2011). Both trends in the UV–optical colors parallel those found in the UV-continuum slopes. Together these findings suggest that the evolution of galaxies at high redshift may be self-similar (see also González et al. 2012).

In the previous section, we briefly discussed the well-defined sequence in β's we observed versus luminosity in our z  ∼  4–7 samples as another instance of a "star-forming" sequence for galaxies. How shall we interpret the changes we observe in the mean UV-continuum slope β of galaxies on this sequence versus their UV luminosity?

Given the approximate correlation of galaxy mass with UV luminosity (e.g., Ouchi et al. 2004b; Lee et al. 2006, 2009; Stark et al. 2009) and gradual build up of galaxies in mass, we would expect the mean properties of galaxies to change gradually as a function of their UV luminosity. While we can imagine many properties of galaxies driving changes in the UV-continuum slope β as a function of UV luminosity, i.e., dust, age, metallicity, AGN content, changes in the dust content of galaxies would likely have the largest effect. Figure 13 provides a simple illustration of this. 0.3 dex changes in the dust content have a much larger effect on the UV-continuum slope β than similarly sized changes in the age, metallicity, or the stellar IMF.

Zoom In Zoom Out Reset image size

Moreover, given the mass–metallicity relationship observed at z  ∼  0–4 (e.g., Tremonti et al. 2004; Erb et al. 2006a; Maiolino et al. 2008), we would expect galaxy metallicity—and hence dust content—to show a strong correlation with mass (and luminosity). Higher luminosity galaxies would largely be redder because of their greater dust content while lower luminosity galaxies would be bluer due to a scarcity of dust. Such a correlation of dust content with mass has been explicitly shown (e.g., Figure 18 from Reddy et al. 2010 and Figure 5 from Pannella et al. 2009). We would, of course, also expect changes in the metallicity and age of star-forming galaxies to contribute to the observed trends in β; however, their effect on the UV-continuum slope would likely be much smaller in general (see also Bouwens et al. 2009, Section 4.4); and Labbé et al. 2007).

Based upon the above tests and discussion, we will assume that dust extinction is the dominant variable in setting up these observed trends. We therefore use our results on the UV-continuum slope β distribution to estimate a mean dust extinction for high-redshift galaxies. We will make use of well-known infrared excess (IRX)–β relationships known to work well at z  ∼  0 (e.g., Meurer et al. 1999; Burgarella et al. 2005; Overzier et al. 2011) and z  ∼  2 (e.g., Reddy & Steidel 2004; Reddy et al. 2006b, 2010, 2012a; Daddi et al. 2007). The canonical z = 0 IRX–β relationship (Meurer et al. 1999) is

$$\begin{equation} A_{1600} = 4.43 + 1.99\beta, \end{equation} \tag{ 2 }$$

where A₁₆₀₀ is the dust extinction at 1600 Å. The Meurer et al. (1999) approach is functionally equivalent to correcting for dust extinction based upon the Calzetti et al. (2000) dust curve.

Of course, use of the Meurer et al. (1999) IRX–β relation at z ⩾ 3 has not been without controversy, and there has been suggestions that the dust extinction in high-redshift galaxies may be either higher or lower than that implied by the Meurer et al. (1999) relationship. Certainly we might expect some change given that it is likely that AGB stars—thought to be the principal sites for the formation of dust—will not be present in the universe until the universe is at least 1 Gyr old, and therefore the dust that existed in the first 1 Gyr of the universe must have formed in another way, e.g., in the winds of SNe (e.g., Maiolino 2006; Maiolino et al. 2008). The attenuation curve for dust from AGB stars may be very different from dust of other origin (e.g., from SNe).

Dust obscuration would be higher in the high-redshift universe if the attention curve were flatter than Calzetti et al. (2000) while the obscuration would be lower if the attenuation curve were steeper than Calzetti et al. (2000), i.e., much more like that from the Small Magellanic Cloud. A flatter (steeper) attenuation curve implies a higher (lower) dust extinction for a given UV slope. Arguments for its being flatter come from efforts to derive the dust properties of QSOs (Gallerani et al. 2010) while arguments for its being steeper follow from studies of very young galaxies at z  ∼  2–3. Both Reddy et al. (2006b) and Siana et al. (2008, 2009) find that dust corrections implied by the UV-continuum slopes β of young galaxies are much too large for the Meurer et al. (1999) IRX–β relation to apply. Chary & Pope (2011) argue for low dust extinction in high-redshift galaxies based upon extragalactic background light stacking results. It also been argued that the Carilli et al. (2008) stacking results of the radio emission in z  ∼  3 galaxies from COSMOS also suggest lower values for the dust extinction, but Reddy et al. (2012a) dispute this, arguing that the Carilli et al. (2008) results are consistent with previous results (which support the Meurer et al. 1999 IRX–β relationship in the mean).

The above arguments aside, there is circumstantial evidence that dust obscuration in high-redshift galaxies is likely at least as large as implied by the Meurer et al. (1999) relationship (see also discussion in Section 6.3). Perhaps the strongest piece of evidence for substantial dust extinction is provided by the large number of high mass (∼1–3× 10¹⁰ M_☉) galaxies found at z  ∼  5–6 in the GOODS fields (Eyles et al. 2005; Yan et al. 2005, 2006).¹⁶ Building up these stellar masses by z  ∼  6 requires SFRs of ∼20–30 M_☉ yr⁻¹, even assuming constant SFRs for 1 Gyr. Without dust extinction, the progenitors to these massive galaxies would need to be very luminous indeed, i.e., −22 AB mag (Yan et al. 2006). However, such galaxies are not observed at z ≳ 7 in the requisite numbers (e.g., Bouwens et al. 2011b), arguing that the progenitors must be moderately dust obscured or have existed in a smaller form (i.e., having subsequently merged).

One other possible uncertainty in the dust corrections regards the possible impact of scatter in the IRX–β relationship. Smit et al. (2012) found that scatter in this relationship could potentially have a modest effect, i.e., ∼0.09 dex, on the estimated dust extinction, but it depended in detail on cross correlations between dust extinction, β, UV luminosity, and the SFR. Since the effect of scatter in the IRX–β relation on the mean dust extinction is not at all clear (and plausibly consistent with no net change), we will not consider a correction at this time.

Utilizing Equation (2) and the observed distribution of UV-continuum slopes (Table 4), we can estimate the mean extinction corrections (L_IR/L_UV + 1) as a function of UV luminosity. We have plotted the results in Figure 14 for the four different redshift intervals considered here. As in other work (Bouwens et al. 2009; Reddy & Steidel 2009; Sawicki 2012), we find that the typical dust extinction in galaxies increases systematically as a function of the UV luminosity. Note that in calculating these mean extinction factors we integrate over the full UV-continuum slope distribution. This distribution is approximated as a normal distribution with the biweight means given in Table 4 and 1σ scatter of 0.34 (the median scatter presented in Table 4). We take A₁₆₀₀ = 0 when A₁₆₀₀ < 0 in Equation (2) above.

Zoom In Zoom Out Reset image size

To determine the actual extinction corrections that are appropriate for real samples, we must weight these extinction corrections according to the UV LFs determined at z  ∼  4, z  ∼  5, z  ∼  6, and z  ∼  7 integrated to specific limiting luminosities. We will make use of the z  ∼  4, z  ∼  5, z  ∼  6, and z  ∼  7 LFs of Bouwens et al. (2007) and Bouwens et al. (2011b). The results are presented in Table 6 and Figure 15. We also include the dust extinction we would derive using an alternate IRX–β relationship A₁₆₀₀ = 4.01 + 1.81β derived by Overzier et al. (2011). Overzier et al. (2011) derive this relationship based upon a small sample of z  ∼  0 galaxies with similar properties to z  ∼  2–3 LBGs.

Zoom In Zoom Out Reset image size

The extinction corrections given in Table 6 and Figure 15 are an update to the estimates we previously provided in Bouwens et al. (2009) based upon the best UV-continuum slope β estimates available at that time. Relative to the Bouwens et al. (2009) extinction estimates, the most significant change is in the extinction estimates we find at z  ∼  4 which are now ∼1.5–2× lower (see also Castellano et al. 2012). This is the direct result of the somewhat bluer UV-continuum slopes β we find using the present ACS+WFC3/IR photometry. Such photometry extends over a wider wavelength baseline than is available using ACS alone (as employed by Bouwens et al. 2009) and therefore allows for much more accurate estimates of β.

The above discussion highlights several clear trends that are present in the properties of star-forming galaxies at z  ∼  4–7 as a function of their UV luminosity and presumably as a function of their stellar mass. The focus of this discussion was the UV-continuum slopes β of galaxies—and by inference—their overall dust extinction.

However, we also noted that the observed β's in the apparent SF sequence showed a modest dependence on the redshift of the sources. At face value, this suggests that the dust extinction in galaxies must increase, from high redshift to low redshift (as per the discussion in Sections 5.2– 5.3). We might expect such an evolution based on the gradual build up of both metals and mass in star-forming galaxies with cosmic time, as seen in the evolution of the mass–metallicity relationship (e.g., Tremonti et al. 2004; Erb et al. 2006a; Maiolino et al. 2008; Mannucci et al. 2009; Laskar et al. 2011), or based on the evolution in the observed correlation between dust extinction and bolometric luminosity from z  ∼  2 to z  ∼  0 (e.g., Reddy et al. 2006b, 2010; Buat et al. 2007). Such a change in the mean dust content of galaxies is also supported by the simulations of Finlator et al. (2011) who predict essentially the same evolution in mean β with redshift (Figure 7), and almost all of the change in β (≳75%) comes from a change in the dust content (K. Finlator 2011, private communication).

Of course, it is always possible that the observed evolution in the UV-continuum slopes β might be due to a change in the overall dust composition or extinction curve with cosmic time, as might occur if the dust composition depended on the age of the stellar population in a galaxy. Indeed, we might expect some change in the dust composition of galaxies as a result of the fact that dust from SNe would be expected to form much earlier in the lifetime of a galaxy than dust from AGB stars (e.g., Maiolino et al. 2004; Maiolino 2006; Reddy et al. 2006b; Gallerani et al. 2010; Finkelstein et al. 2012). While it seems clear that such changes may affect the colors of very young star-forming galaxies, it is not clear how important they are in driving the trends we observe with cosmic time. After all, even without considering such changes in dust composition, the detailed cosmological hydrodynamical simulations of Finlator et al. (2011) are successful in reproducing the approximate evolution in β we observe.

In the previous sections, we saw that the overall shape of the SED for star-forming galaxies—both in the UV-continuum slope and the UV–optical colors (González et al. 2012)—exhibits a very similar dependence on luminosity at all redshifts that we examined z  ∼  4, 5, 6, and 7, suggesting there is a standard star formation sequence for galaxies at high redshift.

In this section, we derive an approximate relationship between the SFR and stellar mass of galaxies that reside on the "star-forming sequence." To do this, we utilize the observed UV luminosities, UV-continuum slopes, and M/L ratios inferred from the observations. We transform the UV luminosities into SFRs by using the canonical Kennicutt (1998) and Madau et al. (1998) UV luminosity-to-SFR conversion factor. A dust correction is made at different UV luminosities using the Meurer et al. (1999) IRX–β relationship and the UV-continuum slope β distribution determined in this paper. As we have seen, these dust corrections are quite significant (∼3×) at higher luminosities around L*_{z = 3}, but are essentially zero at lower luminosities (≲ 0.1L*_{z = 3}: ≳ − 18.5 AB mag).

Stellar masses can be calculated using the luminosity-dependent M/L ratios derived by González et al. (2011). González et al. (2011) derived these M/L ratios utilizing the HST optical + HST near-IR + Spitzer IRAC photometry for a large sample of z  ∼  4 galaxies within the CDF-South ERS field (including those sources that are not individually detected in the IRAC observations). González et al. (2011) find that the M/L ratios scale as (M/L_UV)∝L^0.7_UV. González et al. (2011) find a steeper dependence on luminosity than found by Stark et al. (2009), where the M/L ratio scale as (M/L_UV)∝L^0.175_UV. The z  ∼  4 results should be fairly indicative of the results at higher redshift given the lack of clear evolution in the M/L ratio from z  ∼  7 to z  ∼  4 for galaxies at a fixed UV luminosity (M_UV ∼ −22 to −18 mag; Stark et al. 2009; Labbé et al. 2010a, 2010b; González et al. 2011).

The resulting SFRs and stellar masses for galaxies in various UV luminosity bins are shown in the top panel of Figure 16. Fitting a line to these points in log–log space, we find that the SFR varies as (13⁺⁷_{− 5} M_☉ yr⁻¹)(M_*/10⁹ M_☉)^{0.73 ± 0.32}. Without any dust correction, we find that the SFR varies as (6⁺³_{− 2} M_☉ yr⁻¹)(M_*/10⁹M_☉)^{0.59 ± 0.32}. The equivalent results for the specific star formation rate (SSFR) and dust extinction are shown in the lower two panels of Figure 16 (see also Section 5.6).

Zoom In Zoom Out Reset image size

Interestingly, the SFR versus stellar mass relationship we derive including the effects of dust extinction is much more linear than we would derive without it. Without any dust correction, the results of Stark et al. (2009) and González et al. (2011) imply that SFR∝M^0.85 and SFR∝M^0.59, respectively. However, correcting for dust extinction, these relationships become a much more linear SFR∝M^1.05 and SFR∝M^0.73, respectively. An approximately linear proportionality, i.e., SFR vs. M_*, is exactly what is expected from cosmological hydrodynamical simulations (e.g., Davé et al. 2006; Finlator et al. 2011; Dayal & Ferrara 2012). Indeed, it points to a scenario where galaxies build up exponentially with time along a well-defined star-forming sequence (e.g., Stark et al. 2009; Papovich et al. 2011; Lee et al. 2011).

We emphasize that the SFRs and stellar masses shown in Figure 16 (and also quoted for the "star-forming sequence") are representative values for a luminosity-selected sample. They were derived using the median M/L ratios found in specific bins of UV luminosity (González et al. 2011). In general, one would expect different results for the M/L ratios using mass rather than luminosity-selected samples. Different results are also expected using mean rather than median M/L ratios. In particular, a mass selection would yield higher values for the M/L ratios (though the size of the effect will depend substantially on the scatter in the M/L ratios). Also use of mean rather than median M/L ratios for these calculations would increase the quoted masses. The reason is that medians will not account for the large amounts of mass in the tail of the distribution that extends to high masses. Both effects work in the same sense and would tend to increase the masses of galaxies at each point on the sequence. Overall, by correcting for these effects, we would expect somewhat higher stellar masses for galaxies on the sequence and somewhat lower values of the SSFRs. The precise corrections depend, of course, upon the duty cycle for star formation and scatter in the M/L ratios, but factor of ∼2 corrections would not be surprising. Reddy et al. (2012b) include an extended discussion of the effect of source selection on observed SFR vs. stellar mass relations in their Appendix B and Figure 26.

The results of the previous section allow us to update previous estimates of the SSFR at high redshift to include a correction for the dust extinction. The SSFR, the SFR divided by the stellar mass, has been of considerable interest recently due to the evidence that the SSFR may not evolve very rapidly at high redshift (Stark et al. 2009; González et al. 2010), in significant contrast to that expected from theory (e.g., Bouché et al. 2010; Davé 2010; Dutton et al. 2010; Weinmann et al. 2011; but, see Krumholz & Dekel 2012).

One shortcoming of these early SSFR determinations at high redshift (Stark et al. 2009; González et al. 2010) was that no dust correction was applied (González et al. 2010). Such an approach seemed appropriate at the time, given the very blue UV-continuum slopes β observed for galaxies in the redshift range z ≳ 5 (Bouwens et al. 2009) and large uncertainties on those UV-continuum slopes (e.g., Bouwens et al. 2009). However, now that sufficiently deep, wide-area near-IR data are available we can establish the UV-continuum slopes β and approximate dust corrections more accurately.

Our new estimates of the dust extinction at z  ∼  4–7 allow us to correct previous SSFR estimates at z  ∼  4–7. The results are presented in Figure 17. Typical corrections result in a factor of ∼2–3 increase in the SSFR.

Zoom In Zoom Out Reset image size

We emphasize that the above corrections to the SSFR at z ≳ 4 are very schematic in nature. A proper determination of the SSFR in this regime requires a fairly extensive, self-consistent analysis of the observations. In addition to the issue of dust extinction, other issues that need to be considered are (1) the selection of the samples by mass in contrast to selection by luminosity and (2) emission-line contamination of our broadband flux measurements (e.g., Schaerer & de Barros 2010).

Already there are many theoretical predictions about how the SSFR should evolve with cosmic time (e.g., Bouché et al. 2010; Weinmann et al. 2011; Krumholz & Dekel 2012). How well do these expectations agree with our corrected SSFRs? The evolution of the SSFR in many models follows the specific accretion rate $\dot{M}/M$, which scales as (1 + z)^2.5 (dotted line in Figure 17; Neistein & Dekel 2008; Weinmann et al. 2011). This implies a factor of ∼10 decrease in the SSFR from z  ∼  7 to z  ∼  2. By comparison, our revised estimates of the dust extinction at z  ∼  4–6 imply an SSFR that decreases by a factor of ∼3 to z  ∼  2. While this still does not match the evolution predicted by standard models, the agreement is better. Are there additional ingredients that might lead to further changes? One possibility that has been discussed includes accounting for metallicity dependencies in one's SFR prescription (e.g., Krumholz & Dekel 2012) or changes to the stellar IMF (Davé 2010; Schaye et al. 2010).

The distribution of UV-continuum slopes β shows a remarkably small intrinsic scatter σ_β for a fixed UV luminosity. The typical scatter σ_β observed is just ∼0.34 (Table 4). A similar scatter was found by Labbé et al. (2007) and Bouwens et al. (2009) for galaxy samples at z  ∼  1–2.7 and z  ∼  4, respectively. Such a small scatter allows us to set upper limits on variations in the star formation histories and dust content of galaxies (assuming that changes in one variable do not offset changes in other variables).

Figure 18 illustrates the impact of scatter in dust, age, or metallicity would have on scatter in the β distribution. For the typical galaxy, we assume an E(B − V) of ∼0.15, a metallicity [Z/Z_☉] of −0.7, a Salpeter IMF, and an approximately constant star formation history with some scatter in the SFR. These parameters are fairly representative for what has been found in stellar population models of luminous z  ∼  2–4 galaxies (e.g., Papovich et al. 2001; Shapley et al. 2001, 2005; Erb et al. 2006b; Reddy et al. 2006a).

Zoom In Zoom Out Reset image size

We begin by using the observed scatter in the UV-continuum slope β distribution to set constraints on scatter in the ages (or star formation history) of star-forming galaxies at high redshift. One potentially promising approach is to consider star formation history for galaxies with stochastic variations in the SFR over 50 Myr intervals (a typical timescale over which one might imagine the SFR in a galaxy might be correlated). Assuming similar variations in the SFRs of z > 4 galaxies to that observed at z  ∼  0–1 (where a scatter of ∼0.3 dex in observed the SFR–stellar-mass relationship; Noeske et al. 2007) results in just a 0.14 scatter in β. We computed this scatter by (1) running a simulation with 1000 input galaxies, (2) dividing up the star formation history for each galaxy into ten 50 Myr segments, (3) randomly selecting an SFR for each 50 Myr segment from a log-normal distribution with 0.3 dex scatter, (4) computing the resultant β's for each galaxy based on its star formation history, and (5) calculating the scatter in the derived β distribution for the simulated galaxies. For simplicity, each galaxy in the simulation is taken to have an age of 500 Myr.

The predicted scatter in β, i.e., 0.14, is considerably less than we observe. We can, of course, increase the predicted scatter in β by considering star formation histories with larger variations in the SFRs. For example, a 0.9 dex scatter in the SFRs translates into a 0.34 dex scatter in β,which is a good match to intrinsic scatter in β (σ_β). This is similar to the scatter found by González et al. (2011) in modeling the distribution of M/L ratios for z  ∼  4 star-forming galaxies. In the above modeling, no account is made for changes to the total luminosity of galaxies, as a result of a stochastic star formation history.

Scatter in the dust content can add significantly to the scatter in the UV-continuum slope β, but the magnitude of the scatter will depend directly on how dusty galaxies are in the luminosity range one is considering. If the dust extinction is low, for example, dust has very little impact on the β observed, and therefore small multiplicative changes to the dust extinction factor would have similarly little impact. On the other hand, if the dust extinction is non-negligible, multiplicative changes to the total dust extinction factor would have a big impact on the value of β one observes. For the fiducial luminous galaxy at z  ∼  2.5, with an E(B − V) ∼ 0.15 and a Calzetti et al. (2000) dust law, a ∼0.3 dex scatter in the dust extinction would result in an observed scatter σ_β of ∼0.34 (Figure 18).

Scatter in the metallicity adds very little to scatter in the UV-continuum slopes (Figure 18), and therefore one cannot use the observed scatter in the UV-continuum slope β distribution to set strong limits on scatter in galaxy metallicity.

In summary, the observed scatter in the UV-continuum slopes β distribution allows us to set upper limits on variations in dust content and instantaneous SFR of galaxies (assuming that changes in one variable do not offset changes in other variables). Scatter in the dust extinction of galaxies appear to be ≲0.3 dex and scatter in the instantaneous SFR is ≲0.9 dex.

In this subsection, we determine the SFR density using our current estimates of the dust extinction. As in previous work (Bouwens et al. 2009, 2011b), we base our SFR density determinations on our most recent LF determinations at z  ∼  4–8 (Bouwens et al. 2007, 2011b) and search results at z  ∼  10 (Bouwens et al. 2011a; Oesch et al. 2012a). Luminosity densities are derived by integrating these LFs down to −17.7 AB mag (0.05 L*_{z = 3}) which is the limit to which we probe the LFs at both z  ∼  7 and z  ∼  8. These luminosity densities are then converted into SFR densities using the canonical Madau et al. (1998) and Kennicutt (1998) relation:

$$\begin{equation} L_{\rm UV} = \left(\frac{{\rm SFR}}{M_{\odot } {\rm \,yr}^{-1}} \right) 8.0 \times 10^{27} {\rm erg}\, {\rm s}^{-1}\, {\rm Hz}^{-1}, \end{equation} \tag{ 3 }$$

where a 0.1–125 M_☉ Salpeter IMF and a constant star formation rate of ≳ 100 Myr are assumed. Finally, for our dust extinction estimates, we will use those from Table 6 calculated using the Meurer et al. (1999) IRX–β relationship.

Our latest UV luminosity density and SFR density estimates are summarized in Table 7 and presented in Figure 19. The new estimates are in broad agreement with previous estimates (Bouwens et al. 2009), but we find a lower SFR density at z  ∼  4, as expected, given the lower dust extinction we infer at these redshifts. The change is significant, with the SFR density decreasing by a factor of ∼1.5–2 at this redshift (see also Castellano et al. 2012).

Zoom In Zoom Out Reset image size

Table 7. UV Luminosity Densities and Star Formation Rate Densities to −17.7 AB mag (0.05 L*_{z = 3})^a

Dropout		${\rm log}_{10} \mathcal {L}$	log10 SFR Density	Corrected	Incl. ULIRGb
Sample		(erg s−1	(M☉ Mpc−3 yr−1)
	〈z〉	Hz−1 Mpc−3)	Uncorrected
B	3.8	26.38 ± 0.05	−1.52 ± 0.05	−1.21 ± 0.05	−1.12 ± 0.05
V	5.0	26.08 ± 0.06	−1.82 ± 0.06	−1.54 ± 0.06	−1.51 ± 0.06
i	5.9	26.02 ± 0.08	−1.88 ± 0.08	−1.72 ± 0.08	−1.71 ± 0.08
z	6.8	25.88 ± 0.10	−2.02 ± 0.10	−1.90 ± 0.10	−1.90 ± 0.10
Y	8.0	25.65 ± 0.11	−2.25 ± 0.11	−2.13 ± 0.11	−2.13 ± 0.11
J c	10.3	24.1+0.5− 0.7	−3.8+0.5− 0.7	−3.8+0.5− 0.7	−3.8+0.5− 0.7
J c	10.3	<24.2d	<−3.7d	<−3.7d	<−3.7d

Notes. See Sections 6.1–6.2. ^aIntegrated down to 0.05 L*_{z = 3}. Based upon LF parameters in Table 2 of Bouwens et al. (2011b; see also Bouwens et al. 2007; see Section 6.1). The SFR density estimates assume ≳ 100 Myr constant SFR and a Salpeter IMF (e.g., Madau et al. 1998). Conversion to a Chabrier (2003) IMF would result in a factor of ∼1.8 (0.25 dex) decrease in the SFR density estimates given here. ^bSee Section 6.2. ^cz  ∼  10 determinations and limits are from Oesch et al. (2012a; see also Bouwens et al. 2011a) and assume 0.8 z  ∼  10 candidates in the first case and no z  ∼  10 candidates (i.e., an upper limit) in the second case. ^dUpper limits here are 1σ (68% confidence).

Download table as:  ASCIITypeset image

Dusty, infrared-luminous galaxies also potentially contribute quite meaningfully to the SFR density at high redshift. However, it can be quite challenging to account for this contribution on the basis of optical/near-IR surveys with the Hubble Space Telescope. Not only can it be difficult to identify such sources in these surveys (due to their faintness in the UV or red colors which cause them to be excluded from LBG selections), but it is now well established that the dust extinction for the most dusty, infrared luminous (>10¹² L_☉) galaxies at z  ∼  1–3 cannot be accurately estimated using the observed UV-continuum slopes β and IRX–β relationship (e.g., Reddy et al. 2006b; Elbaz et al. 2007).

To accurately account for the SFR density in this population, a better approach is to include that population explicitly by utilizing a LF in the mid-IR/far-IR and integrating down to 10¹² L_☉ (e.g., Reddy et al. 2008; Reddy & Steidel 2009; Bouwens et al. 2009). We replicate that approach here integrating published mid-IR LFs to 10¹² L_☉, converting the total IR luminosity to SFR using the canonical relation in Kennicutt (1998), and then adding the inferred SFR densities to the dust-corrected UV SFR densities. The IR LFs we utilize are Magnelli et al. (2009, 2011) to z  ∼  2 and Daddi et al. (2009) at z  ∼  4. Magnelli et al. (2011) utilize the full set of deep 24 μm and 70 μm observations over the two GOODS fields from the FIDEL program. Since the 70 μm and 24 μm data that Magnelli et al. (2011) utilize are significantly deeper than that used by Caputi et al. (2007) and allow for a self-consistent correction to the total bolometric luminosity, the Magnelli et al. (2011) SFR density estimate represents a noteworthy improvement on the Caputi et al. (2007) estimates we previously utilized at z  ∼  2 (Bouwens et al. 2009).

Including the contribution from IR bright galaxies, we present our total SFR density estimates in Table 7 and Figure 19. Interestingly, but not surprising, this correction makes very little difference to the SFR density derived at very high redshifts z > 4, but adds modestly to the SFR density at late cosmic times (z ⩽ 3). A small contribution of ULIRGs to the SFR density at high redshifts is expected given their position at the very end of the extended buildup process whereby galaxies gradually acquire higher and higher masses in gas, dust, and stars (e.g., Bouwens et al. 2009).

As in our previous work (Bouwens et al. 2011b), we can compare our SFR density estimates with that implied by recent stellar mass density determinations (e.g., Stark et al. 2009; González et al. 2011; Labbé et al. 2010a, 2010b; González et al. 2011). In doing so, we consider the integrated SFR density and stellar mass density to the same luminosity limits 0.05 L*_{z = 3} for self-consistency.¹⁷

We use the following formula to infer an approximate SFR density (SFRD) at z ≳ 4 from the observed stellar mass density (SMD):

$$\begin{equation} {\rm SFRD}(z_i,z_j) = \frac{{\rm SMD}(z_i) - {\rm SMD}(z_j)}{{\rm time}(z_i)-{\rm time}(z_j)}(1-\epsilon)^{-1} f_{{\rm LE}}, \end{equation} \tag{ 4 }$$

where z_i < z_j are the redshifts of adjacent Lyman-break samples and is the gas recycling factor. The 1 − factor accounts for the recycling of gas mass from high-mass stars back into the interstellar medium through SN explosions. Recycling results in only a fraction of the stars formed being locked up in stellar mass, i.e., dM_*/dt = (1 − )SFR where = 0.3 appropriate for a Salpeter IMF (e.g., Bruzual & Charlot 2003). The f_LE factor in the above equation accounts for the fact that new galaxies enter our magnitude-limited samples at all redshifts simply as a result of galaxy growth. Since these galaxies (and their stellar mass) would not have been included in the magnitude-limited sample just above them in redshift, these sources would cause us to overestimate the SFR density (see Section 7.4 of Bouwens et al. 2011b). Accounting for this latter effect reduces the inferred SFR density by a factor of 1.3; therefore, we take f_LE to be 1/1.3.

The SFR densities implied by several recent stellar mass density determinations (e.g., Stark et al. 2009; González et al. 2010; Labbé et al. 2010a, 2010b; González et al. 2011) are presented in Figure 20, and there is remarkably good agreement over the redshift range z  ∼  4–6. This is a useful consistency check and suggests that the dust extinctions we are inferring at z  ∼  4–7 are reasonable and fit into a consistent picture. The general agreement we observe also points toward no clear evolution in the stellar IMF to high redshift (z ⩾ 4; see also Bouwens et al. 2011b; Papovich et al. 2011) since any changes in the IMF would have an effect on the SFR density or stellar mass densities we compute from the observed light and likely result in a mismatch.

Zoom In Zoom Out Reset image size

The observed agreement also suggests that biases in our stellar mass density estimates at high redshift, e.g., due to the contamination of the rest-frame optical light probed by IRAC with strong emission lines (Schaerer & de Barros 2010), are not huge. Of course, we cannot totally rule out modest levels of contamination by emission lines, particularly to the SEDs of z  ∼  6–7 galaxies, and in fact there may be some evidence in the stacked SEDs of z  ∼  5–7 galaxies that emission lines do have some effect on the IRAC fluxes (González et al. 2012).

We would expect small biases in the mean UV-continuum slope β we measure from our observed samples. These biases arise whenever galaxies with specific properties—intrinsic or observed—are easier to select than others. We discuss two such biases: (1) one based on the fact that galaxies with certain intrinsic colors are easier to select than others (Appendix B.1.1) and (2) one arising from the fact that the information used to measure β is not always independent from that used to select the sources (Appendix B.1.2).

The bias described in Appendix B.1.1 arises because not all values for β are selected with the same efficiency, and the bias described in Appendix B.1.2 arises because the average measured value for β does not always equal the intrinsic value.

The purpose of Appendix B.1 is to provide a description of the relevant biases and their approximate size for the most simplistic cases (i.e., if the β distribution is a delta function). In Appendix B.2, we calculate these biases assuming a more realistic input distribution for β.

B.1.1. "Selection Volume Bias": Biases Related to the Intrinsic Selectability of Sources

In general, we would expect our high-redshift selections to be more effective in identifying galaxies with certain UV-continuum slopes than others, resulting in small biases in the mean UV-continuum slope β derived from our samples. Since bluer galaxies are almost always easier to select than redder galaxies in LBG selections, this bias would be toward bluer colors. This effect is often referred to as "template bias" in the literature.

To quantify the extent to which object selection affects the distribution of UV-continuum slopes recovered from the observations, we run extensive Monte Carlo simulations. In these Monte Carlo simulations, we insert artificial galaxies into the real data with a wide variety of luminosities and colors and then reselect these galaxies using the same method as we used on the real sources. We then use these simulations to determine how the selection efficiency and hence selection volume depends upon the intrinsic UV-continuum slope β. A detailed description of these simulations is given in Appendix C.

Figure 21 provides an illustration of the basic results from these simulations in the magnitude ranges 25.0–27.0 mag and 27.0–28.5 mag. Indeed we find the expected trends in the simulations. The selection efficiency is largest for galaxies with the bluest UV-continuum slopes. As in other Lyman-break selections, the selection volumes show a sizeable decrease toward redder UV-continuum slopes, especially for β ≳ 0.0.

The effect of these biases on the mean UV-continuum slope β derived from the observations, however, is relatively modest for typical samples. For example, assuming the intrinsic distribution of UV-continuum slopes β has a mean β of −1.5 and 1σ scatter of 0.3 (see Table 4), the bias in the mean β is just Δβ ∼ 0.1 for a ∼28 mag galaxy in our z  ∼  4 selection. We estimate the size of this bias by starting with an initial distribution of β's, weighting this distribution by the effective volume at each β (for a given magnitude), and then computing the shift in the biweight mean β.

Figure 22 shows the effect of these selection biases on the mean UV-continuum slope β found in our z  ∼  4–7 selections. The biases shown in Figure 22 assume that the intrinsic β distribution has a mean value of −2.0 and 1σ scatter of 0.3. While we note a bias toward bluer UV-continuum slopes β, these biases are not extraordinarily large. Similarly small selection biases are found in our other data sets, i.e., the two other HUDF09 fields, the ERS field, and the CDF-South CANDELS field.

Zoom In Zoom Out Reset image size

Comparing the number of sources observed with a given β with the computed selection volumes for this β gives us another means of assessing the importance of this bias. As can be seen in Figure 23 for our HUDF selections, we can select sources to much redder β's than generally observed in our selections. As in our direct estimates of the bias (Figure 22), this suggests that the selection effects are only having a modest effect on our derived β's.

Zoom In Zoom Out Reset image size

We emphasize that the biases presented in this subsection and the next are simply intended to be illustrative; more realistic estimates of the bias (and those values of the bias we will use to correct the observations) are presented in Appendix B.2.

B.1.2. Photometric Error Coupling Bias

In Appendix B.1, we estimated the approximate bias in the mean UV-continuum slope β we would determine from the observations. We began by estimating the bias that would result from the fact that galaxies with bluer UV-continuum slopes β are easier to select than galaxies with redder UV-continuum slopes (Appendix B.1.1). We then estimated the bias we would expect in the measured values for β as a result of the slight coupling that occurs between source selection and the measurement process (Appendix B.1.2). Of course, both biases affect the mean UV-continuum slope β, so we need to add the above biases together to obtain the total bias in the UV-continuum slope β. The results are presented in Figure 22 (solid lines).

Remarkably enough, the approximate bias in the mean UV-continuum slope β is almost zero. The bias toward bluer slopes (from the intrinsic selection biases; Appendix B.1.1) largely offsets the bias toward redder slopes (from the slight coupling between source selection and β measurements; Appendix B.1.2), resulting in a very small bias overall. Given the very small size of the biases expected here for our faint samples, it is certainly striking how large these same biases are for the "robust" photometric redshift selection of Dunlop et al. (2012). Figures 10 and 22 provide a rather dramatic illustration of the differences. The contrast in biases is especially noteworthy, especially given the large number of candidates in our faint samples.

While we can calculate the approximate effect of the aforementioned biases on β given some input distribution, for the simulations to be accurate we must use the true underlying β distribution for these inputs. Since we are not able to establish the underlying β distribution without first establishing the biases (so we can apply a correction to the observed distribution of β's), it was necessary to follow an iterative approach. The observed β served as a starting point for the computation of the biases. In each iteration, we used the previous best estimate of the UV-continuum slope β distribution to establish the relevant corrections until we obtained convergence.

Biases in the mean UV-continuum slope were estimated in a very similar way for each of our samples (z  ∼  4, z  ∼  5, z  ∼  6, and z  ∼  7).

Flux measurement errors act to significantly broaden the UV-continuum slope β distribution. The added scatter can be substantial, particularly near the selection limits for our samples. Typical uncertainties on individual UV-continuum slope β measurements range from 0.4 to 1.0, with the largest uncertainties being relevant for z  ∼  7 galaxies where only two passbands (J₁₂₅, H₁₆₀) are available to make the measurement and the wavelength baseline is relatively short. An illustration of the extent to which these errors can increase the apparent scatter in the β distribution can be found in Bouwens et al. (2009; Figure 4).

To quantify the extent to which flux measurement error ("noise") increases the spread in the UV-continuum slope β distribution, we again rely on the simulations described in Appendix B.1 (see also Appendix C). For these simulations, we assume that all galaxies have exactly the same UV-continuum slope β, with no scatter, and then measure the scatter in the recovered distribution. The derived 1σ scatter is shown in Figure 24 for z  ∼  4–7 galaxies identified in the HUDF observations. Scatter is largest for faint sources in each of our z  ∼  4–7 samples and also in our z  ∼  6–7 samples. Scatter in our lower redshift samples is smaller because of the much larger number of passbands and longer wavelength baselines over which we use to measure the UV-continuum slopes β (Table 3 and Section 3.3).

Zoom In Zoom Out Reset image size

Another benefit of having run these simulations is that they allow us to quantify the extent to which we can reduce the overall uncertainties in β using the full flux information in the UV continuum, rather than using just a single color (Section 3.3 and Appendix A). These gains are most dramatic for z  ∼  4 and z  ∼  5 galaxies in our fields where we use the information in four passbands (Table 3) rather than just two (Appendix A). In the z  ∼  4 or z  ∼  5 cases, for example, we find a factor of ∼1.5 reduction in the scatter.

In the previous sections, we outlined the procedure we use to determine the distribution of UV-continuum slopes β for high-redshift galaxy samples, after correcting for various selection and measurement biases.

An important check on these results is to compare the distributions we derive at different depths. How well does the distribution of UV-continuum slopes β obtained from the ultra-deep HUDF09 data set agree with that obtained from the wide-area ERS+CANDELS observations? The results are shown in Figure 25 for all four of our Lyman-break selections. While the agreement is excellent overall, this should perhaps not come as a great surprise given the very small corrections that need to be made to each of our selections (Figure 22).

Zoom In Zoom Out Reset image size