z ∼ 2–9 Galaxies Magnified by the Hubble Frontier Field Clusters. II. Luminosity Functions and Constraints on a Faint-end Turnover

We will base the present deep LF results primarily on the sensitive near-UV, optical, and near-IR observations obtained by the HFF program (Coe et al. 2015; Lotz et al. 2017) and a follow-up GO campaign of the HFF clusters with WFC3/UVIS (Siana 2013, 2015; Alavi et al. 2014, 2016). Over each of the six clusters in the HFF program, at least 16 orbits of WFC3/UVIS time (eight and eight in the UV₂₇₅ and U₃₃₆ bands), 70 orbits of optical ACS time (18, 10, and 42 in the B₄₃₅, V₆₀₆, and I₈₁₄ bands, respectively), and 70 orbits of WFC3/IR time (24, 12, 10, and 26 in the Y₁₀₅, J₁₂₅, JH₁₄₀, and H₁₆₀ bands, respectively) were invested into observations of each cluster, with additional observations coming from the CLASH (Postman et al. 2012) and GLASS (Schmidt et al. 2014) programs. We made use of the v1.0 reductions of these observations made publicly available by the HFF team (Koekemoer et al. 2014). For the WFC3/UVIS observations, we made our own reductions, following similar procedures to what we utilized in reducing the UVIS observations obtained as part of the HDUV program (Oesch et al. 2018b).

In addition, we will make use of the z = 2–9 LF constraints available from the comprehensive set of blank-field HST observations recently utilized by Bouwens et al. (2021a). For that analysis, Bouwens et al. (2021a) not only made use of the extremely sensitive optical, near-IR, and rest-UV observations obtained over the HUDF (Beckwith et al. 2006; Koekemoer et al. 2013; Illingworth et al. 2013; Teplitz et al. 2013), but also made use of the sensitive UV, optical, and near-IR data over the five CANDELS fields and Early Release Science field (Grogin et al. 2011; Koekemoer et al. 2011; Oesch et al. 2018a), Hubble Frontier Field parallels (Coe et al. 2015; Lotz et al. 2017), and pure parallel search fields (Trenti et al. 2011; Yan et al. 2011).

We consider the systematic application of z = 2–9 Lyman-break galaxy selection criteria to the HFF near-UV, optical, and near-IR observations available over the six clusters that make up the HFF program. We describe our application of these criteria to the HFF observations in a companion paper (Bouwens et al. 2022b), and it is very similar to what we previously performed in Bouwens et al. (2015a), Bouwens et al. (2016a), and Bouwens et al. (2017b).

In total, we have identified 765, 1176, 68, 59, 274, 125, 51, and 16 sources as part of our z ∼ 2, z ∼ 3, z ∼ 4, z ∼ 5, z ∼ 6, z ∼ 7, z ∼ 8, and z ∼ 9 selections, respectively. The results are summarized in Table 1. z ∼ 4 galaxies are exclusively selected around the two highest-redshift clusters that make up the HFF program, i.e., MACS0717 and MACS1149. Meanwhile, z ∼ 5 galaxies are selected behind the four lowest-redshift clusters that make up the HFF program, i.e., A2744, MACS0416, AS1063, and A370. z ∼ 4 and z ∼ 5 selections are exclusively made behind these specific clusters to minimize contamination from foreground sources in clusters having spectral breaks at similar wavelengths to the breaks used to select Lyman-break galaxies. See Figures 1–4 from Bouwens et al. (2022b).

An illustration of the distribution of the present lensed sample in redshift and UV magnitude is provided in Figure 1 along with the >24,000 source sample we utilize from the Bouwens et al. (2021a) blank-field selection. As should be apparent, the lensed sample reaches ∼10–40× fainter in UV luminosity than does the blank-field sample, providing us with much greater leverage for probing the faint-end form of the z ≥ 2 LFs.

Zoom In Zoom Out Reset image size

Here we describe our basic procedure for constraining the shape of the UV LF for each of our intermediate- to high-redshift samples at z ∼ 2, z ∼ 3, z ∼ 4, z ∼ 5, z ∼ 6, z ∼ 7, z ∼ 8, and z ∼ 9.

One significant challenge, in the derivation of the UV LF results at high redshift from lensed samples, is the impact of errors in the estimated magnification factors for individual sources. As demonstrated in Bouwens et al. (2017b), errors in the magnification factors can effectively wash out faint-end (>−15 mag) turn-overs in the UV LF, making it difficult to observationally test for the presence of such a turn-over in real-life observations.

We have already demonstrated in Bouwens et al. (2017b) how we can overcome the impact of potential errors in the lensing maps using a forward-modeling procedure (see Atek et al. 2018 for a separate approach using forward modeling). The basic idea is to leverage the availability of the many independent lensing models for each cluster to estimate the uncertainties. Model LFs are evaluated by treating one of the lensing maps as the truth and thus constructing a full catalog of observables with that map and then interpreting the observations using another lensing model. In this way, the expected number of sources per unit luminosity for a model LF could be realistically estimated. Figure 6 of Bouwens et al. (2017b) illustrates the basic procedure.

Here we will follow the same forward-modeling procedure as was introduced in Bouwens et al. (2017b). In evaluating model LFs, we treat one flavor of lensing model as the truth and use it to construct a complete catalog of background sources behind each cluster. Sources are added to a search field in proportion to the product of the model volume density, selection efficiency S(m), and the cosmic volume element—which we take to be the cosmic volume element divided by the magnification factor. The selection efficiency S is, in general, a function of both the apparent magnitude m and the magnification factor μ, but we can largely ignore the impact of magnification in the limit that faint sources have very small sizes. Justification for this is provided in the published results of Bouwens et al. (2017a), Kawamata et al. (2018), Bouwens et al. (2022a), and Bouwens et al. (2022b).

We then use a different lensing model to estimate the magnification factors and UV luminosities for individual sources behind the cluster. To maximize the reliability of our results, we made use of the median magnification from the latest parametric lensing models for our fiducial LF determinations. Not only has use of the median magnification model been shown to provide robust estimates of the magnification to magnifications of ≳50 (Bouwens et al. 2017b; Livermore et al. 2017; Bouwens et al. 2022a), but the parametric models were shown to best reproduce the input magnification models from the HFF comparison project (Meneghetti et al. 2017), with median differences (<0.1 dex differences) to magnification factors of 30. The parametric lensing models available for the HFF clusters and utilizing most of the public multiple image constraints, i.e., v3/v4, include the CATS (Jullo & Kneib 2009; Richard et al. 2014; Jauzac et al. 2015a, 2015b; Lagattuta et al. 2017), Sharon (Johnson et al. 2014), GLAFIC (Oguri 2010; Ishigaki et al. 2015; Kawamata et al. 2016), Zitrin-NFW (Zitrin et al. 2013, 2015), and Keeton (2010) results. When computing the median magnification map used to estimate the magnification factors for our forward-modeling approach, any model we treat as the truth is naturally excluded.

While we base our fiducial LF results on the parametric lensing models available for the HFF clusters, we also derive LF results using the nonparametric lensing models available for the HFF clusters. The nonparametric lensing models have been shown to be a good match (<0.1 dex) to the input models from the HFF comparison project (Meneghetti et al. 2017) to magnification factors of ∼10–20. Given the greater flexibility of the nonparametric models relative to the parametric models, results derived from these models allow us to assess the impact lensing models have on our LF results.

We evaluate the likelihood of an LF model by comparing the observed number of sources in various absolute magnitude bins (0.5 mag width) with the expected number of sources assuming that galaxies are Poisson-distributed:

where

N_obs,i and are the observed and expected number of sources in magnitude interval i. To reduce the impact of sources with complex multicomponent or morphological structure on our analysis—which become common at the bright end of the LF—we only consider luminosity bins fainter than −19 mag. As a result of this choice, this analysis relies entirely on blank-field LF results for constraints brighter than −19 mag. Additionally, our forward-modeling simulations typically include ∼200× as many sources as are present in the actual observations (e.g., see Figure 6 from Bouwens et al. 2017b) to guarantee an accurate calculation of for our likelihood estimates.

To evaluate the likelihood of various LF parameters for our z ∼ 6–9 samples, we must account for our combining five separate selections of z ∼ 6, z ∼ 7, z ∼ 8, and z ∼ 9 galaxies in creating our composite sample of z ∼ 6–9 galaxies. This was done to maximize the utility of our results to constrain the faint end of the z ∼ 6–9 LFs. Selections were constructed by leveraging a separate source catalog, each made with different parameters to handle the background parameters and thresholding (Bouwens et al. 2022b). While this results in a larger number of sources in each bin, many sources occur in our final catalog multiple times and are not independent. We can account for the impact of this by quantifying the fraction of sources f_i,j in each bin i that are counted one, two, three, four, and five times and assuming a similar multiplicity fraction in the modeled statistics. As such, the probability _Pi for measuring a specific number of counts in bin i is as follows:

where the summation runs over all N_obj,i,j where the sum Σ_j=1,5 jN_obj,i,j is equal to N_obj,i . We have verified that Equation (2) reduces to the appropriate Poissonian likelihood distribution in the limit that all sources are present in our catalogs with a fixed multiplicity (e.g., once or five times).

Due to the relatively modest depth of the UV₂₇₅ and U₃₃₆ data available to select our z ∼ 2–3 samples, we only include sources brighter than 28.0 mag and 28.5 mag (V₆₀₆ and I₈₁₄ bands, respectively) to mitigate the impact of the uncertain completeness corrections fainter than these limits. For our z ∼ 4 and z ∼ 5 selections, where contamination from evolved galaxies at the redshift of the cluster is a particular concern, we restrict ourselves to using sources brighter than 27.3 and 27.5 mag, respectively.

We used extensive source recovery simulations to estimate selection efficiencies S(m) for each of our intermediate- to high-redshift samples. For each of the simulations, we first constructed mock catalogs of sources over the general redshift ranges spanned by each of our z ∼ 2–9 samples, i.e., z = 1–4, z = 1–4, z = 2.5–5.0, z = 4–6, z = 5–7, z = 5.5–8, z = 6.5–9.5, and z = 7–10 for our z ∼ 2, 3, 4, 5, 6, 7, 8, and 9 selections, respectively. We then created artificial two-dimensional images for each of the sources in these catalogs in all HST bands used for the selection and detection of the sources and then added these images to the real observations. We then repeated both our catalog creation procedure and source selection procedure in the same way as on the real observations (Bouwens et al. 2017a, 2017b).

Motivated by our earlier findings regarding source sizes for faint z ∼ 2–8 samples from the HFFs (Bouwens et al. 2017a see also Bouwens et al. 2017c, 2022a; Kawamata et al. 2018; Yang et al. 2022), we adopted a point-source size distribution in modeling the completeness of sources over the HFFs for our fiducial LF results. Additionally, we took the UV-continuum slope β distribution to a median value of −2 at z ∼ 2–3, −2.1 at z ∼ 4–5, and −2.3 at z ∼ 6–8, consistent with the Kurczynski et al. (2014) and Bouwens et al. (2014) UV-continuum slope measurements. We have verified that for unlensed sources at the faint end of our HFF selections, we estimate almost identical selection volumes to what we estimate by computing the selection volumes using randomly selected z ∼ 2–4 galaxies as templates in our image simulations. We expect that this is due to the combination of the high surface brightness sensitivity of the HFF observations and the faintest sources in our fields having small sizes.

To quantify the possible systematic uncertainties that could result in our LF determinations from our size modeling, we have also repeated these completeness estimates using the size–luminosity relations of both Shibuya et al. (2015) and Bouwens et al. (2022a) to illustrate how large the systematic uncertainties could be, similar to the exercise we performed in Bouwens et al. (2017b). While we include these estimates as a possible systematic error on the derived faint-end LF results, we emphasize that this is a worst-case scenario, as the results from Bouwens et al. (2017a), Kawamata et al. (2018), and Yang et al. (2022) all point toward substantially smaller source sizes.

As an illustration of the substantial leverage available from the HFF samples to constrain the faint-end slope, we show in Figure 2 the number of sources behind the HFF clusters as a function of UV luminosity and compare that to the predicted numbers for specific values for the faint-end slope of the UV LF at z ∼ 3. From this figure and the substantial numbers of sources fainter than −19 mag, it is clear that the faint-end slope of the z ∼ 3 LF must be fairly close to −1.6, and faint-end slopes of −1.4 and −1.8 can both be excluded at high confidence on the basis of the HFF results.

Zoom In Zoom Out Reset image size

As in Bouwens et al. (2017b), we adopt a standard Schechter functional form for our LF

but modified to allow for curvature in the faint-end slope α at the faint end. We implement this using a new curvature parameter δ and multiply the standard Schechter form with the following expression fainter than −16 mag:

As we demonstrate in Bouwens et al. (2017b), positive values of δ result in a turn-over in the LF at

while negative values for δ result in the LF turning concave upwards.

Given our use of four separate parameters to characterize the overall shape of the UV LF, we utilize a Markov Chain Monte Carlo (MCMC) procedure both to determine the maximum likelihood value and to characterize the observational uncertainties. We begin the MCMC optimization process using the best-fit blank-field LF parameters from Bouwens et al. (2021a) and run ∼10,000 iterations to find the best-fit LF parameters and map out the likelihood space.

3.3.1. Leveraging Only Lensed Sources from the HFF Program

As in Bouwens et al. (2017b), we commence our analysis by exclusively making use of the lensed HFF samples to derive our initial UV LF results. The value in doing this first is that it allows us to test for the presence of any systematic errors in lensing-cluster LF determinations vis-à-vis blank-field determinations. The two approaches have their own strengths and are subject to different sources of systematic error, so it is valuable to first look into this issue before combining the results to arrive at the best constraints on the overall shape of the UV LF.

In Bouwens et al. (2017b), we demonstrated we could obtain accurate constraints on both the faint-end slope α and ϕ* at z ∼ 6 relying only on the lensed z ∼ 6 sources from the HFF lensing clusters. There, a faint-end slope α of −1.92 ± 0.04 and a ϕ* of 0.66 ± 0.06 × 10⁻³ Mpc⁻³ were found, consistent (≲1σ) with the −1.87 ± 0.10 slope and Mpc⁻³ normalization found from blank-field observations (Bouwens et al. 2015a). For these determinations, the characteristic luminosity M* was fixed to −20.94 mag, the value obtained from wide-area blank-field studies, due to their being insufficient volume behind the HFF clusters to achieve strong constraints on this luminosity.

We adopt a similar approach here. We begin by fixing the characteristic luminosities M* for the z ∼ 2–9 LFs to the values obtained by the Bouwens et al. (2021a) blank-field analysis and then use the described MCMC approach to identify the values of ϕ*, α, and δ that maximize the likelihood of recovering our z ∼ 2–9 samples from the full HFF program. Uncertainties on the individual LF parameters can be calculated based their fits to the multidimensional likelihood surface derived from our MCMC simulations.

We present the LF results we derive using our lensed HFF samples and fixed values for the characteristic luminosity in Table 2. In addition, in Figure 3, we present our determinations of the faint-end slope α of the LF versus redshift. A simple linear fit to the faint-end slope α results we derive versus redshift yields the relation

which is presented on Figure 3 as a blue line. For context, Figure 3 also shows a recent determination of the faint-end slope α evolution based on blank-field observations alone (red circles; Bouwens et al. 2021a), with the observed trend with redshift (red line):

Encouragingly enough, the new faint-end slope α results we derive from lensing-cluster observations seem very consistent (within the 1σ errors) both in terms of slope and intercept, with blank-field results over the entire redshift range we are examining z ∼ 2 to z ∼ 9.

Zoom In Zoom Out Reset image size

Table 2. Summary of Our Final Fiducial Constraints on the z ∼ 2–9 UV LFs

Data Sets + Method		ϕ* (10−3 Mpc−3)	α	δa	MT d
z ∼ 2
HFF + Blank-field (Parametric, Fiducial)	−20.30 ± 0.08	3.4 ± 0.4	−1.53 ± 0.03	0.09 ± 0.11	>−15.2
HFF + Blank-field (Nonparametric) b	−20.31 ± 0.09	4.0 ± 0.7	−1.53 ± 0.03	0.20 ± 0.13	>−15.3
HFF (Parametric) + M* Fixed	−20.28 (fixed)	3.0 ± 0.6	−1.56 ± 0.06	0.15 ± 0.13
HFF (Parametric) + M* + ϕ* Fixed	−20.28 (fixed)	4.0 (fixed)	−1.49 ± 0.03	0.05 ± 0.11
Σ(μ) Fit c			≤ − 1.61 ± 0.06
Blank-field (B21a)	−20.28 ± 0.09		−1.52 ± 0.03	⋯
z ∼ 3
HFF + Blank-field (Parametric, Fiducial)	−20.84 ± 0.07	2.3 ± 0.4	−1.60 ± 0.03	−0.06 ± 0.05	>−13.1
HFF + Blank-field (Nonparametric) b	−20.96 ± 0.09	1.9 ± 0.8	−1.65 ± 0.04	0.05 ± 0.07	>−12.3
HFF (Parametric) + M* Fixed	−20.87 (fixed)	2.0 ± 1.5	−1.62 ± 0.09	−0.03 ± 0.06
HFF (Parametric) + M* + ϕ* Fixed	−20.87 (fixed)	2.1 (fixed)	−1.60 ± 0.03	−0.06 ± 0.04
Σ(μ) Fit c			≤ − 1.71 ± 0.04
Blank-field (B21a)	−20.87 ± 0.09		−1.61 ± 0.03	⋯
z ∼ 4
HFF + Blank-field (Parametric, Fiducial)	−20.93 ± 0.07	1.6 ± 0.3	−1.69 ± 0.03	−0.18 ± 0.18	>−15.3
HFF + Blank-field (Nonparametric) b	−20.95 ± 0.07	1.8 ± 0.3	−1.71 ± 0.03	−0.06 ± 0.12	>−15.6
HFF (Parametric) + M* Fixed	−20.93 (fixed)	1.5 ± 0.6	−1.71 ± 0.09	−0.12 ± 0.22
HFF (Parametric) + M* + ϕ* Fixed	−20.93 (fixed)	1.69 (fixed)	−1.68 ± 0.05	−0.19 ± 0.17
Σ(μ) Fit c			≤ − 1.87 ± 0.15
Blank-field (B21a)	−20.93 ± 0.08		−1.69 ± 0.03	⋯
Bouwens et al. (2007)	−20.98 ± 0.10	1.3 ± 0.2	−1.73 ± 0.05	⋯
Bouwens et al. (2015a)	−20.88 ± 0.08		−1.64 ± 0.04
z ∼ 5
HFF + Blank-field (Parametric, Fiducial)	−21.13 ± 0.09	0.73 ± 0.13	−1.78 ± 0.04	−0.07 ± 0.20	>−14.7
HFF + Blank-field (Nonparametric) b	−21.15 ± 0.10	0.74 ± 0.14	−1.80 ± 0.05	−0.02 ± 0.24	>−14.8
HFF (Parametric) + M* Fixed	−21.10 (fixed)	0.66 ± 0.15	−1.82 ± 0.08	−0.03 ± 0.21
HFF (Parametric) + M* + ϕ* Fixed	−21.10 (fixed)	0.79 (fixed)	−1.77 ± 0.06	−0.07 ± 0.20
Σ(μ) Fit c			≤ − 1.94 ± 0.16
Blank-field (B21a)	−21.10 ± 0.11		−1.74 ± 0.06	⋯
Bouwens et al. (2015a)	−21.17 ± 0.12		−1.76 ± 0.05
z ∼ 6
HFF + Blank-field (Parametric, Fiducial)	−20.87 ± 0.07	0.57 ± 0.11	−1.87 ± 0.04	0.05 ± 0.10	>−14.3
HFF + Blank-field (Nonparametric) b	−20.98 ± 0.08	0.45 ± 0.10	−1.98 ± 0.06	0.36 ± 0.18	>−15.2
HFF (Parametric) + M* Fixed	−20.93 (fixed)	0.61 ± 0.29	−1.85 ± 0.05	0.00 ± 0.14
HFF (Parametric) + M* + ϕ* Fixed	−20.93 (fixed)	0.51 (fixed)	−1.87 ± 0.02	0.09 ± 0.09
Σ(μ) Fit c			≤ − 1.90 ± 0.09
Blank-field (B21a)	−20.93 ± 0.09		−1.93 ± 0.08	⋯
Bouwens et al. (2015a)	−20.94 ± 0.20		−1.87 ± 0.10
z ∼ 7
HFF + Blank-field (Parametric, Fiducial)	−21.13 ± 0.08	0.20 ± 0.05	−2.05 ± 0.06	0.24 ± 0.20	>−15.2
HFF + Blank-field (Nonparametric) b	−21.21 ± 0.09	0.17 ± 0.05	−2.13 ± 0.07	0.56 ± 0.26	>−15.4
HFF (Parametric) + M* Fixed	−21.15 (fixed)	0.17 ± 0.12	−2.09 ± 0.13	0.36 ± 0.28
HFF (Parametric) + M* + ϕ* Fixed	−21.15 (fixed)	0.19 (fixed)	−2.06 ± 0.04	0.25 ± 0.16
Σ(μ) Fit c			≤ − 1.84 ± 0.12
Blank-field (B21a)	−21.15 ± 0.13		−2.06 ± 0.11	⋯
Bouwens et al. (2015a)	−20.87 ± 0.26		−2.06 ± 0.13
z ∼ 8
HFF + Blank-field (Parametric, Fiducial)	−20.90 ± 0.19	0.096 ± 0.065	−2.20 ± 0.09	0.23 ± 0.28	>−15.2
HFF + Blank-field (Nonparametric) b	−20.96 ± 0.19	0.094 ± 0.051	−2.27 ± 0.09	0.50 ± 0.30	>−15.3
HFF (Parametric) + M* Fixed	−20.93 (fixed)	0.11 ± 0.09	−2.17 ± 0.15	0.15 ± 0.31
HFF (Parametric) + M* + ϕ* Fixed	−20.93 (fixed)	0.09 (fixed)	−2.21 ± 0.05	0.32 ± 0.17
Blank-field (B21a)	−20.93 ± 0.28		−2.23 ± 0.20	⋯
Bouwens et al. (2015a)	−20.63 ± 0.36		−2.02 ± 0.23
z ∼ 9
HFF + Blank-field (Parametric, Fiducial)	−21.15 ± 0.12	0.018 ± 0.009	−2.28 ± 0.10	0.53 ± 0.41	>−15.6
HFF + Blank-field (Nonparametric) b	−21.15 ± 0.11	0.021 ± 0.011	−2.34 ± 0.11	0.52 ± 0.41	>−15.6
HFF (Parametric) + M* Fixed	−21.15 (fixed)	0.019 ± 0.012	−2.27 ± 0.12	0.53 ± 0.29
HFF (Parametric) + M* + ϕ* Fixed	−21.15 (fixed)	0.021 (fixed)	−2.21 ± 0.05	0.33 ± 0.19
Blank-field (B21a)	−21.15 (fixed)		−2.33 ± 0.19	⋯

Notes.

^a Best-fit curvature in the shape of the UV LF fainter than −16 mag (Section 3.2). ^b The nonparametric LF results presented in this table represent a mean of our derived LFs treating the nonparametric lensing models, i.e., grale and Diego, as the truth. ^c Constraints based on the observed trend in source surface density versus model magnification factor μ (Bouwens et al. 2022b) ^d Brightest luminosity at which the current constraints from the HFF permit a turn-over in the UV LF (95% confidence).

Download table as:  ASCIITypeset images: 1 2

The observed agreement between the two results is remarkable given the enormous differences between the two approaches and their different challenges. For example, while blank-field probes provide us with less leverage in luminosity to constrain the faint-end slope α of the UV LF, they are much less sensitive to assumptions about source size near the detection limits (since almost all sources at these limits are unresolved) and sample sufficient volume to obtain much better constraints on the bright end of the LF. Meanwhile, for lensing-cluster probes, despite the sensitivity of this approach to an accurate modeling of the sizes (e.g., Bouwens et al. 2017a), the additional leverage provided by lensing allows us to probe substantially fainter in luminosity, providing us with much greater leverage to probe the faint-end slope. Figure 2 illustrates this leverage quite well.

We emphasize that our use of very small sizes is critical for achieving consistent faint-end slope α results to blank-field studies. As demonstrated in both Bouwens et al. (2017a) and Bouwens et al. (2022a), if we had instead assumed that galaxies followed a more standard size–luminosity relation (e.g., Huang et al. 2013; Shibuya et al. 2015), completeness would be a decreasing function of the magnification of the source. This would result in much steeper faint-end slopes α and significantly higher volume densities of sources >−15 mag (Bouwens et al. 2017a, 2022a).

This is the first time such consistent faint-end slope results have been found between blank-field and lensing-cluster analyses over such an extended range in redshift. The independence of the two approaches and their consistency strongly suggests that we can combine the two approaches to arrive at even more robust determinations of the overall shape and evolution of the UV LF.

3.3.2. Leveraging Both Blank-field and HFF Results

Having demonstrated we can use lensed samples of z = 2–9 galaxies to obtain constraints on the faint-end slope α consistent with blank-field studies, we now leverage both data sets to arrive at our best overall estimate on the shape of the UV LFs at z = 2–9.

For our blank-field constraints on the z ∼ 2–9 LF results, we rely on the likelihood constraints Bouwens et al. (2021a) derived from the comprehensive set of HST fields. Following the treatment we provided in Bouwens et al. (2017b), we again allow for a 20% uncertainty in the relative normalization of the value of blank-field and cluster LF results at z ∼ 2, 3, 6, 7, 8, and 9. This 20% uncertainty includes a ∼15% contribution from cosmic variance (Robertson et al. 2014) and a ∼10% uncertainty in the selection volume calculation for both the bright and faint ends of the LFs. For our z ∼ 4 and z ∼ 5 results, we allow for a 26% and 18% cosmic variance uncertainty (and 28% and 22% relative uncertainty in total) owing to our consideration of lensed sources behind only two and four clusters, respectively, due to the challenge in identifying z ∼ 4–5 galaxies behind those clusters without substantial contamination.

Additionally and motivated by the modest uncertainties in the size distribution of the faintest sources, we allow for modest uncertainties in the selection efficiency of galaxies in the faintest magnitude intervals. In particular, we consider both the impact of increasing the selection efficiencies in the 27.5–28.0 mag, 28.0–28.5 mag, 28.5–29.0 mag, and >29.0 mag intervals by 2%, 5%, 20%, and 30% and decreasing the selection efficiencies by the same amount and then marginalize the results over both the regularly calculated efficiencies and the two others. These percentage changes indicate the approximate impact of ∼10% differences between the model and derived S/N for sources of a given magnitude. For our derived LF parameters at z ∼ 2 and 3, similar adjustments were made to selection volumes to determine the impact of possible systematics in the estimated volumes, but 1.2 and 0.7 mag brighter (corresponding to the shallower UV₂₇₅ and U₃₃₆ coverage over the HFF clusters and also shallower z ∼ 2 and z ∼ 3 selections).

In Figure 4, we show the 68% and 95% confidence intervals on our LF results marginalizing over the results for the four families of parametric models, i.e., GLAFIC, CATS, Sharon/Johnson, and Keeton, while comparing against independent constraints on the bright end of the UV LFs from various blank-field probes (Bouwens et al. 2007, 2015a; 2019, 2020, 2021a; Bowler et al. 2015, 2017, Harikane et al. 2022). Along with the formal confidence intervals, we also show the allowed LF results (light shaded region) if the mean sizes of lower-luminosity galaxies are more consistent with the Shibuya et al. (2015) scalings. We stress that the observational results of Bouwens et al. (2017a), Kawamata et al. (2018), Bouwens et al. (2022a), and Yang et al. (2022) strongly suggest the true size distribution is significantly smaller than this, but we show these allowed regions to present the possible systematic error.

Zoom In Zoom Out Reset image size

For convenience, in Table 3, we include a compilation of our derived 68% confidence intervals for our z = 2–9 LF results from z ∼ 9 to z ∼ 2. The formal uncertainties provided in Table 3 include not only the formal 68% confidence results from our MCMC fitting results using the parametric models, but include systematic uncertainties on the volume densities if the true sizes of lower-luminosity galaxies are much larger than implied by the Bouwens et al. (2017a, 2017c), Kawamata et al. (2018), and Bouwens et al. (2022a) results.

Table 3. Sixty-eight Percent Likelihood Intervals Derived for the z ∼ 2–9 UV LFs Using Forward Modeling (Section 3.3.2) ^a

M1600,AB b	ϕk	M1600,AB b	ϕk	M1600,AB b	ϕk	M1600,AB b	ϕk
	(Mpc−3 mag−1)		(Mpc−3 mag−1)		(Mpc−3 mag−1)		(Mpc−3 mag−1)
z ∼ 2 galaxies		z ∼ 4 galaxies		z ∼ 6 galaxies		z ∼ 8 galaxies
−18.75		−18.75		−18.75		−18.75
−18.25		−18.25		−18.25		−18.25
−17.75		−17.75		−17.75		−17.75
−17.25		−17.25		−17.25		−17.25
−16.75		−16.75		−16.75		−16.75
−16.25		−16.25		−16.25		−16.25
−15.75		−15.75		−15.75		−15.75
−15.25		−15.25		−15.25		−15.25
−14.75		−14.75		−14.75		−14.75
−14.25		−14.25		−14.25		−14.25
−13.75		−13.75		−13.75		−13.75
−13.25		−13.25		−13.25		−13.25
−12.75		−12.8		−12.75		−12.75
−12.25		−12.2		−12.25		−12.25
z ∼ 3 galaxies		z ∼ 5 galaxies		z ∼ 7 galaxies		z ∼ 9 galaxies
−18.75		−18.75		−18.75		−18.75
−18.25		−18.25		−18.25		−18.25
−17.75		−17.75		−17.75		−17.75
−17.25		−17.25		−17.25		−17.25
−16.75		−16.75		−16.75		−16.75
−16.25		−16.25		−16.25		−16.25
−15.75		−15.75		−15.75		−15.75
−15.25		−15.25		−15.25		−15.25
−14.75		−14.75		−14.75		−14.75
−14.25		−14.25		−14.25		−14.25
−13.75		−13.75		−13.75		−13.75
−13.25		−13.25		−13.25		−13.25
−12.75		−12.75		−12.75		−12.75
−12.25		−12.2		−12.25		−12.25

Notes.

^a Sixty-eight percent confidence intervals on the z ∼ 6 UV LF we achieve (Section 3.3.2) by applying our forward-modeling formalism to observations of all six HFF clusters in Section 3.3.2. The quoted constraints give the geometric mean of the forward-modeling results using the GLAFIC, CATS, Sharon/Johnson, and Keeton parametric models as inputs. Examples of these constraints for individual lensing models are provided in Figures 14 and 15. The 1σ upper limits indicate the upper limits if one adopts the Bouwens et al. (2022a) size–luminosity scalings (implemented in a similar manner to Section 5.4 of Bouwens et al. 2017b). As the intervals provided here are based on our parameterized modeling, the error bars on individual bins are not independent. ^b Derived at a rest-frame wavelength of 1600 Å.

Download table as:  ASCIITypeset image

We also include a simple binned version of our LF results in Figure 4 using the relation

where ϕ_m is the derived volume density of sources in magnitude bin m, N_m is the number of sources in magnitude bin m, and V_m is the estimated selection volume in magnitude bin m. We calculate V_m as follows:

where dA is a differential area in the image plane, C(z, m, μ) is the selection completeness as a function of redshift z, apparent magnitude m, an the magnification factor μ, and is a differential volume element. The plotted uncertainties in Figure 4 include not only the formal Poissonian uncertainties but also the systematic uncertainties on the volume densities if the true sizes of lower-luminosity galaxies are as given by the Bouwens et al. (2022a) size–luminosity relation instead of adopting point-source sizes. No account is made for the impact of uncertainties in magnification factors on the binned constraints. For convenience, we list the binned constraints shown in Figure 4 in Table 4.

Table 4. Binned Determinations of the Rest-frame UV LF at z ∼ 2–9 (Section 3.3.2) ^a

M1600,AB b	ϕk c	M1600,AB b	ϕk c	M1600,AB b	ϕk c	M1600,AB b	ϕk c
	(Mpc−3 mag−1)		(Mpc−3 mag−1)		(Mpc−3 mag−1)		(Mpc−3 mag−1)
z ∼ 2 galaxies		z ∼ 4 galaxies		z ∼ 6 galaxies		z ∼ 8 galaxies
−18.75		−18.75		−18.75		−18.75
−18.25		−18.25		−18.25		−18.25
−17.75		−17.75		−17.75		−17.75
−17.25		−17.25		−17.25		−17.25
−16.75		−16.75		−16.75		−16.75
−16.25		−16.25		−16.25		−16.25
−15.75		−15.75		−15.75		−15.75
−15.25		−15.25		−15.25		−15.25
−14.75		−14.75	<0.088	−14.75		−14.75	<0.068
−14.25		−14.25	<0.35	−14.25		−14.25
−13.75		−13.75		−13.75		−13.75
−13.25				−13.25		−13.25
−12.75	<0.11			−12.75	<0.97
−12.25
z ∼ 3 galaxies		z ∼ 5 galaxies		z ∼ 7 galaxies		z ∼ 9 galaxies
−18.75		−18.75		−18.75		−18.75
−18.25		−18.25		−18.25		−18.25
−17.75		−17.75		−17.75		−17.75
−17.25		−17.25		−17.25		−17.25
−16.75		−16.75		−16.75		−16.75
−16.25		−16.25		−16.25		−16.25
−15.75		−15.75		−15.75		−15.75
−15.25		−15.25	<0.051	−15.25		−15.25	<0.037
−14.75		−14.75		−14.75		−14.75	<0.097
−14.25		−14.25		−14.25		−14.25	<0.27
−13.75		−13.75	<3.0	−13.75		−13.75	<0.98
−13.25				−13.25	<0.56	−13.25	<5.8
−12.75				−12.75	<3.7
−12.25

Notes.

^a The LF results presented in this table are derived using Equation (6) from Section 3.3.2. One advantage of the results presented in this table is that each binned value is independent of the other bins. On the other hand, no account is made for uncertainties in the magnification map in deriving constraints on the overall shape of the UV LF, so this is one disadvantage. ^b Derived at a rest-frame wavelength of 1600 Å. ^c The 1σ upper limits indicate the upper limits if one adopts the Bouwens et al. (2022a) size–luminosity scalings (implemented in a similar manner to Bouwens et al. 2017b).

Download table as:  ASCIITypeset image

Results derived using the individual parametric and nonparametric magnification models are presented in Table 7 of Appendix A. Figures 14 and 15 from Appendix A show the 68% and 95% confidence intervals on the z ∼ 2 and z ∼ 6 LF results, respectively, using the same set of parametric and nonparametric models. In general, very similar LF results are obtained utilizing different lensing models, with derived parameters typically differing by much less than the formal uncertainties derived for a single model. In general, LFs derived using the nonparametric lensing models showed higher values for the curvature parameter δ. This appears to be a consequence of the larger differences seen between the magnification factors from these models and those in the parametric lensing models and the impact this has in washing out potential turn-overs at the faint end of the UV LFs.

A summary of the mean LF parameters M*, ϕ*, α, and δ we derive on the basis of the four parametric magnification models (GLAFIC, CATS, Sharon/Johnson, Keeton) and the two nonparametric models (grale and Diego) is provided in Table 2. They are indicated by the descriptors "Parametric" and "Nonparametric," respectively. Finally, Figure 5 shows our best-fit z ∼ 2, z ∼ 3, z ∼ 4, z ∼ 5, z ∼ 6, z ∼ 7, z ∼ 8, and z ∼ 9 LFs, along with the z ∼ 10 results from Oesch et al. (2018a). The z ∼ 10 LF results from Oesch et al. (2018a) rely on both blank-field and lensing-field results. From the plotted constraints, it is clear that the UV-bright galaxies undergo a much more rapid evolution in volume density than galaxies at the faint end of the LF. In the next section, we parameterize the evolution of the UV LF in terms of convenient fitting formulas.

Zoom In Zoom Out Reset image size

The availability of deep lensed samples of z = 2–9 galaxies behind the HFF clusters has made it possible to significantly improve our constraints on the faint-end slope of the UV LF and therefore the evolution of the UV LF from z ∼ 10 to z ∼ 2.

Our estimates of the Schechter parameters in Table 2 provide a good illustration of this. Comparing the uncertainties on the faint-end slope α from our most recent blank-field determinations, i.e., Bouwens et al. (2021a), and those combining blank-field constraints with the lensing-cluster constraints, we have been able to reduce the uncertainties on the faint-end slope α by factors of ∼1.5, 2, 2, 2, and 2 at z ∼ 5, 6, 7, 8, and 9, respectively.

Taking advantage of these new constraints, we are in a position to further refine our characterization of the evolution in each Schechter parameter. Following Bouwens et al. (2021a), we fit the evolution in α as a linear function of redshift, log₁₀ ϕ* as a quadratic function of redshift, and M* as a linear function of redshift, with a break at z ∼ 2.5. Simultaneously fitting the present LF constraints over the redshift range z ∼ 2–9 as well as the Oesch et al. (2018a) results at z ∼ 10, we find the following best-fit relation:

where z_t = 2.42 ± 0.09. Figure 6 compares the observed evolution of α, M*, and ϕ* with the above relation.

Zoom In Zoom Out Reset image size

As we previously noted (Bouwens et al. 2015a, 2021a), the evolution in the faint-end slope α agrees remarkably well with the change in slope expected based on the evolution of the halo mass function over the same range in redshift, i.e., d α/dz ∼ − 0.12 (Bouwens et al. 2015a). This is fairly similar to the faint-end slope evolution we recovered in Bouwens et al. (2015a), i.e., d α/dz ∼ − 0.10 ± 0.03, the d α/dz ∼ − 0.11 ± 0.01 trends derived by Parsa et al. (2016) and Finkelstein (2016) fitting to the available z ∼ 0–10 and z ∼ 4–10 LF results, respectively, and the d α/dz ∼ − 0.09 ± 0.05 trend recently derived by Bowler et al. (2020; however, they express the LF evolution they derive using a double power-law form).

The slow evolution in the characteristic luminosity M* at z > z_t seems very likely to be a consequence of the fact that UV luminosity reaches a maximum value of ∼−21 to −23 mag due to the increasing importance of dust extinction in the highest stellar mass and SFR sources (Bouwens et al. 2009a; Reddy et al. 2010). Finally, as Bouwens et al. (2015a, 2021a) have argued, the evolution in ϕ* can readily be explained by evolution in the halo mass function and no significant evolution in the SFE. As in Bouwens et al. (2021a), we find that log₁₀ ϕ* depends on redshift with a clear quadratic dependence (but here significant at 6σ instead of 4σ).

Before discussing in more detail the implications of the present LF determinations, it makes sense to compare our new results with the many previous determinations of these LFs that exist in the literature. Comparing with previous work is very valuable for improving LF determinations in general, as it allows us to identify differences in the results and ascertain the best path to improve LF results in the future.

To this end, we present our new z ∼ 2, 3, 4, 5, 6, 7, 8, and 9 LF results in Figures 7 and 8 using the 68% and 95% confidence intervals we have derived (light-gray and dark-gray shaded regions, respectively). For the comparisons we provide, we will focus on previous LF results that provided results on the faint end of the z = 2–9 LFs. In the preparatory work to this using blank-field data (Bouwens et al. 2021a), we already provided a significant discussion of previous LF results, which are more relevant for the bright end of the UV LFs.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

To help with the discussion of various faint-end LF results, we also include in Figures 7 and 8 the approximate regions in parameter space allowed if we take galaxy sizes to follow the Bouwens et al. (2022a) size–luminosity relation (light shaded regions above the nominal 68% and 95% contours).

We will structure the comparisons we make to previous work as an increasing function of redshift:

z ∼ 2–5: At present, there have been only two studies that have reported LF results on the extreme faint end of the UV LF at z ∼ 2–5, i.e., at >−16 mag, for Lyman-break galaxies: one by Alavi et al. (2016) based on lensed samples of z ∼ 1–3 galaxies identified behind two HFF clusters A2744 and MACS0717 and A1689, and a second by Parsa et al. (2016) using a deep (∼30.5 mag) photometric redshift selection over the HUDF. The faint constraints reported by Alavi et al. (2016) and Parsa et al. (2016) lie above and below our own constraints and reach very different conclusions. Alavi et al. (2016) report a faint-end slope from their LF results of −1.72 ± 0.04 and −1.94 ± 0.06 at z ∼ 2 and z ∼ 3, respectively, while Parsa found −1.31 ± 0.04 and −1.31 ± 0.04 at approximately the same redshifts.

As illustrated in Figure 2, the leverage available from our lensed HFF samples should allow us to quantify the faint-end slope for the UV LF z ∼ 2 and z ∼ 3 with great precision, if an accurate account can be made for various sources of systematic errors. Given the consistency of our own blank-field and lensed measurements of the faint end of the z ∼ 2 and z ∼ 3 LF, what might drive the lower and higher volume density results obtained by Parsa et al. (2016) and Alavi et al. (2016)? For the faint end of the Parsa et al. (2016) probe, the answer is not entirely clear, as their LF determinations are in good agreement with both the Bouwens et al. (2021a) blank-field LF results and our lensed results brighter than −16.5 mag and differences in our faint-end slope inferences may follow from our different M* determinations. ⁸ Nevertheless, fainter than −16.5 mag, differences between the Parsa et al. (2016) LF measurements and our own are more significant. One potential concern for the LF determinations of Parsa et al. (2016) is the inclusion of sources to ∼30.5–31.0 mag, i.e., at essentially the detection limit of the HUDF. At such faint magnitudes, segregating sources into different redshift bins is more challenging, especially given that the UV observations Parsa et al. (2016) report utilizing (essential for probing the position of the Lyman break) have a 5σ depth of 28 mag, i.e., ∼10× brighter than the z ∼ 1–4 sources being included in their LFs.

Alavi et al.'s (2016) z ∼ 2–3 LF results are higher than our own for the faintest sources (>−15 mag). This almost certainly results from the completeness corrections Alavi et al. (2016) calculated extrapolating the Shibuya et al. (2015) size–luminosity relation to lower luminosities. In fact, making use of a similar size–luminosity relation, we find a three times higher volume density for the faintest sources, largely reproducing Alavi et al.'s (2016) results. It is worth emphasizing that use of the size–luminosity relation from Shibuya et al. (2015) for completeness calculations was very reasonable, given the lack of information that existed 5 yr ago regarding the sizes of the faintest sources at z ∼ 2–3.

Finally, we include some constraints on the z ∼ 2–5 LFs of galaxies by Weisz et al. (2014) that rely on resolved stellar population analyses of nearby dwarf galaxies, abundance matching, and evolving the luminosities of these dwarf galaxies backwards in time (see also Boylan-Kolchin et al. 2015). While such analyses are obviously very interesting to pursue and should be useful in providing indicative constraints on the volume density of faint sources during the first few billion years of the universe, they also necessarily involve a number of significant assumptions regarding the precise history of the stars that make up these analyses (e.g., mergers and galaxy disruption). Given the uncertainties, it is encouraging how consistent the Weisz et al. (2014) LF inferences and our own results are.

z ∼ 6–9: At lower luminosities, our new z ∼ 6–7 LF results from the HFFs are in excellent agreement with the z ∼ 6–7 results from Atek et al. (2015b) using the first three HFF clusters, the z ∼ 5–10 results from Castellano et al. (2016) based on the first two HFF clusters, our previous results obtained from the first four HFF clusters (Bouwens et al. 2017b), the z ∼ 6–7 results from Atek et al. (2018) based on the full HFF program, the Yue et al. (2018) results based on the first four HFF clusters, and the LF results obtained by both McLure et al. (2013) and Oesch et al. (2013) at z ∼ 9. The z = 6–9 LF results of Bhatawdekar et al. (2019) also appear to be in reasonable agreement with our own LF results, particularly if the comparison is made against their results modeling the sizes of faint galaxies as "disk galaxies" (with a mean size of 015). Given the very small sizes adopted in our analysis, we would have expected the "point-source" results from Bhatawdekar et al. (2019) to be most consistent with our own, but the point-source results from Bhatawdekar et al. (2019) are ∼1.5× lower; it is not clear why this would be the case.

Our LF results at z = 6–8 also show good overall agreement with the results from Livermore et al. (2017) brighter than −17 mag. At lower luminosities, the differences are larger. We refer interested readers to Bouwens et al. (2017a) and Bouwens et al. (2017b) for discussions of the differences. Our results also show a broad similarity to the z ∼ 6–9 LF results of Ishigaki et al. (2018), but we note a slight excess in the volume density of sources they found in their z ∼ 6–7 LF results and z ∼ 9 results (upper- and lower-right panels of Figure 8) at ∼−19 mag and a slight deficit in their z ∼ 8 LF results at ∼−18 to −17 mag. The slightly higher volume densities that Ishigaki et al. (2018) found in their LF results at z ∼ 9 may be derived from their probing a slightly lower redshift range with their z ∼ 9 selection than we do in our own determinations. Use of the photometric redshift estimates from Ishigaki et al. (2018; their Table 8) substantiates this assertion, as Ishigaki et al. (2018) probed a mean redshift of z ∼ 8.5 with their z ∼ 9 selection, and we probe a mean redshift of z ∼ 8.9. The slight deficit Ishigaki et al. (2018) report at −18 and −17 mag in their z ∼ 8 LF results appears to be derived from the limited number of sources Ishigaki et al. (2018) have in the faintest bins, i.e., 1, 1, and 3, respectively.

Given the impressively deep UV LF results we have over the redshift range z = 2–9, it is clearly interesting to use these measurements to map out the evolution of the UV luminosity density of galaxies to very faint luminosities from z ∼ 9 to z ∼ 2.

To maximize the utility of this exercise, we have elected to compute the luminosity density to the faint-end limits −17, −15, −13, and −10 mag. We adopted those faint-end limits due to their frequent use in both blank-field LF studies and reionization calculations (e.g., Robertson et al. 2013; Bouwens et al. 2015b; Robertson et al. 2015; Ishigaki et al. 2018).

To compute the luminosity densities implied by our new LF results, we simply compute the luminosity density implied by various parameterizations and marginalize across the likelihood distribution. We compute both 68% and 95% confidence intervals on the luminosity densities and have tabulated our results in Table 5. The results are also presented in Figure 9. From the results presented in both the table and figures, it is clear that probes to −13 mag contribute meaningfully to the UV luminosity density vis-à-vis probes to −17 mag, increasing the total UV luminosity by 0.1 dex (st z ∼ 2) and by 0.4 dex (at z ∼ 8).

Zoom In Zoom Out Reset image size

Table 5. Sixty-eight percent and 95% Confidence Intervals on the UV Luminosity Density at z ∼ 2–9 to Various Limiting Luminosities ^{a , b}

	(UV Luminosity Density)
	(ergs s−1 Hz−1 Mpc−3)
	Lower Bound		Upper Bound
Faint-end Limit	95%	68%	68%	95%
	z ∼ 2
MUV < − 17	26.38	26.41	26.46	26.49
MUV < − 15	26.46	26.49	26.54	26.57
MUV < − 13	26.48	26.51	26.57	26.60
MUV < − 10	26.48	26.51	26.59	26.68
	z ∼ 3
MUV < − 17	26.45	26.51	26.64	26.70
MUV < − 15	26.52	26.58	26.72	26.77
MUV < − 13	26.56	26.62	26.76	26.81
MUV < − 10	26.58	26.66	26.88	26.94
	z ∼ 4
MUV < − 17	26.42	26.47	26.56	26.60
MUV < − 15	26.52	26.56	26.66	26.70
MUV < − 13	26.56	26.62	26.82	27.01
MUV < − 10	26.58	26.71	29.18	31.17
	z ∼ 5
MUV < − 17	26.24	26.27	26.35	26.38
MUV < − 15	26.35	26.38	26.46	26.50
MUV < − 13	26.39	26.44	26.62	26.87
MUV < − 10	26.39	26.45	28.52	30.84
	z ∼ 6
MUV < − 17	26.08	26.12	26.20	26.24
MUV < − 15	26.23	26.26	26.34	26.38
MUV < − 13	26.28	26.32	26.43	26.50
MUV < − 10	26.29	26.34	26.58	27.04
	z ∼ 7
MUV < − 17	25.92	25.95	26.01	26.05
MUV < − 15	26.11	26.14	26.23	26.27
MUV < − 13	26.15	26.20	26.33	26.41
MUV < − 10	26.15	26.20	26.41	26.68
	z ∼ 8
MUV < − 17	25.55	25.61	25.72	25.78
MUV < − 15	25.85	25.89	26.01	26.07
MUV < − 13	25.90	25.96	26.20	26.41
MUV < − 10	25.90	25.96	26.43	28.01
	z ∼ 9
MUV < − 17	24.99	25.07	25.26	25.34
MUV < − 15	25.33	25.39	25.54	25.60
MUV < − 13	25.36	25.43	25.65	25.75
MUV < − 10	25.36	25.43	25.67	25.82

Notes.

^a The derived luminosity densities represent a geometric mean of the LF results derived treating the parametric models as the truth and recovering the LF results using a median of the other parametric models. ^b See also Figure 9 to see these luminosity densities (integrated to −13 mag) presented as the equivalent SFR densities as a function of redshift.

Download table as:  ASCIITypeset image

Remarkably, the present observational results allow us to constrain the luminosity densities to an uncertainty of <0.05 dex at z ∼ 2–3 and z ∼ 6–7 and <0.12 dex at z ∼ 4–5 and z ∼ 8–9, equivalent to <13% and <30% relative uncertainties, respectively. The UV luminosity densities we derive brighter than −10 mag are much less well constrained. While the 95% confidence intervals we derive from our z ∼ 2–3, z ∼ 6–7, and z ∼ 9 results span a ≲0.5 dex range, these same intervals span a ≳2 dex range at z ∼ 4–5 and z ∼ 8.

It is interesting to convert our new determinations of the UV luminosity density into equivalent SFR densities using the conversion factors in Madau & Dickinson (2014). Assuming a Chabrier (2003) initial mass function, a constant SFR, and metallicity Z = 0.002 Z_⊙, the conversion factor _FUV is 0.7 × 10⁻²⁹ M_⊙ yr⁻¹ erg⁻¹s Hz. The specific value adopted for the metallicity does not have a huge impact on this conversion factor (<1.5%; Madau & Dickinson 2014). The equivalent SFR densities to the integrated UV luminosity densities to −13 mag are also shown on the left vertical axis of Figure 9.

We have compared the implied unobscured SFR density from our UV LF results in Figure 9 to the obscured SFR density inferred from the Bouwens et al. (2020) ASPECS HUDF study, the obscured SFR density results at z = 0–2 from Magnelli et al. (2013), and the unobscured SFR density results from z = 0 to z = 1.5 from Wyder et al. (2005) and Moutard et al. (2020). We also present the unobscured SFR density inferred at z ≥ 2 integrated down to a brighter limit, −17 mag, similar to what Madau & Dickinson (2014) utilized in their SFR density figures.

For convenience, Table 6 presents our new SFR density estimates along with estimates of the obscured SFR densities Bouwens et al. (2020) derived combining their own inferences of the obscured star formation from their ASPECS HUDF results with the ULIRG results of Wang et al. (2019), Franco et al. (2020a), and Dudzevičiūtė et al. (2020). We adopt a fiducial 0.1 dex uncertainty in the obscured SFR densities at z ≥ 2 given the considerable challenges in deriving SFRs from far-IR spectral energy distributions (SEDs) resulting from the uncertain SED shapes, an uncertain contribution from active galactic nuclei (AGNs), and uncertainties in the selection volume. Figure 10 compares the unobscured SFR densities to the total SFR densities.

Zoom In Zoom Out Reset image size

Table 6.  UV Luminosity Densities and Star Formation Rate Densities to −13.0 AB mag (0.0006 : Section 4.2)

Lyman			SFR Density
Break		(ergs s−1	(M⊙ Mpc−3 yr−1)
Sample	< z >	Hz−1 Mpc−3) a	Unobscured	Obscured b , c	Total
U275	2.1	26.54 ± 0.03	−1.61 ± 0.03	−1.07 ± 0.10	−0.96 ± 0.08
U336	3.0	26.69 ± 0.07	−1.46 ± 0.07	−1.36 ± 0.10	−1.11 ± 0.07
B	3.8	26.74 ± 0.12	−1.41 ± 0.12	−1.47 ± 0.10	−1.14 ± 0.08
V	4.9	26.54 ± 0.09	−1.62 ± 0.09	−1.81 ± 0.10	−1.40 ± 0.07
i	5.9	26.38 ± 0.05	−1.78 ± 0.05	−2.28 ± 0.10	−1.66 ± 0.05
z	6.8	26.27 ± 0.06	−1.89 ± 0.06	−2.93 ± 0.10	−1.85 ± 0.06
Y	7.9	26.08 ± 0.12	−2.07 ± 0.12	−3.31 ± 0.10	−2.05 ± 0.11
J	8.9	25.54 ± 0.11	−2.61 ± 0.11	⋯	−2.61 ± 0.11

Notes.

^a From Table 5. ^b From Table 8 from ASPECS HUDF analysis of the infrared excess (Bouwens et al. 2020). The obscured SFR density from Bouwens et al. (2020) explicitly includes the ULIRG results from Wang et al. (2019), Franco et al. (2020a, 2020b), and Dudzevičiūtė et al. (2020). ^c In light of the considerable challenges in deriving SFRs from far-IR SEDs resulting from the uncertain SED shapes, an uncertain contribution from AGN, and uncertainties in the selection volume, we assume a fiducial 0.1 dex uncertainty in the obscured SFR densities at z ≥ 2.

Download table as:  ASCIITypeset image

It is clear from these results that unobscured star formation dominates the SFR density in the high-redshift universe, and obscured star formation dominates the SFR density at intermediate and low redshifts. As has been found previously (e.g., Bouwens et al. 2009a, 2016b; Dunlop et al. 2017; Zavala et al. 2021), we find that the cross-over point between these two regimes is at z ∼ 4. Additionally, it is interesting to note the impact that the faint-end limit can have on the transition redshift between the two regimes. In the Bouwens et al. (2020) ASPECS study, the transition redshift was z ∼ 5 using a brighter faint-end limit of −17 mag, but using a faint-end limit of −13 mag, enabled by our new HFF results, we find a transition redshift z ∼ 4.

Thanks to the faintness of lensed HFF samples, our new constraints put us in position to set constraints where the UV LF might turn over at the faint end. This question is relevant both because it provides insight into the efficiency of star formation in lower-mass galaxies (Behroozi et al. 2013) and because it allows for a more accurate accounting for the contribution of ionizing photons from especially faint star-forming galaxies.

As in Bouwens et al. (2017b), we can constrain the brightest magnitude where a turn-over in the LF can possibly occur using Equation (3) and the overall likelihood distribution we have derived on the three parameters M*, α, and δ. By marginalizing over the results we have obtained on the UV LF at z = 2–9 based on the four different families of parameterized lensing models, we can derive 68% and 95% confidence intervals on the UV luminosity of the turn-over.

We present in the upper panel of Figure 11 the constraints we obtain on the position of the turn-over in the UV LF for star-forming galaxies from z ∼ 9 to z ∼ 2. We find that our results rule out the presence of a turn-over brighter than −15.5 mag (95% confidence) for all z = 2–9 samples we consider. We obtain our tightest constraints on the luminosity of a possible turn-over at z ∼ 3, where our results rule out the presence of a turn-over brighter than −13.1 mag (95% confidence). We note that Alavi et al. (2014) previously presented evidence for the UV LF at z ∼ 2 extending so faint. At z ∼ 6, our results rule out the presence of a turn-over brighter than −14.3 mag.

Zoom In Zoom Out Reset image size

Interestingly enough, our z ∼ 9 LF results seem consistent with a turn-over at −15 mag. While indeed this would be interesting if this were the case, it is challenging to establish the robust presence of a turn-over in the LF at z ∼ 9 due to the small number of sources expected fainter than −16 mag, and therefore our results cannot even establish the presence of a turn-over at 2σ significance.

To make matters even more challenging, another complicated factor is the impact that incompleteness in faint samples could have on the results. If faint sources have larger sizes than adopted here based on several recent observational probes (Bouwens et al. 2017a, 2022a; Kawamata et al. 2018; Yang et al. 2022), this would result in faint sources being much less complete in our selections than in our simulations and cause us to systematically underestimate the volume density of faint galaxies. This would cause the actual luminosity of a turn-over in the UV LF to be substantially fainter than what we infer.

In general, the present constraints on the luminosity of a turn-over in the LF are consistent with most previous results in the literature. Atek et al. (2015b), Atek et al. (2018), Castellano et al. (2016), Livermore et al. (2017), Bouwens et al. (2017b), and Yue et al. (2018) all agree that the HFF results provide strong evidence against the UV LF showing a turn-over brighter than −15 mag. At z ≥ 6 and fainter than −15 mag, there has been a wide variety of differing conclusions drawn about whether firm constraints can be set regarding the existence of a turn-over and how faint those constraints extend (e.g., Castellano et al. 2016; Bouwens et al. 2017b; Livermore et al. 2017; Yue et al. 2018; Atek et al. 2018).

Finally, it is useful for us to compare the current constraints on the evolution of the UV LF with those available from a number of recent theoretical models and cosmological hydrodynamical simulations. While we had previously looked at this in Bouwens et al. (2017b), here we have the advantage that we can compare the simulation/theory results with our new UV LF results over a much more extended baseline in both redshift and cosmic time, reaching from z ∼ 9 to z ∼ 4. While no comparisons are made at z ∼ 2 and z ∼ 3 due to the lack of published model LF results in these two redshift intervals, it would be presumably possible to derive such LF results in the near future based on a number of ongoing simulation efforts, e.g., the IllustrisTNG (Naiman et al. 2018; Nelson et al. 2018; Pillepich et al. 2018; Springel et al. 2018) and NewHorizons (Dubois et al. 2021) simulations.

We consider the following theoretical models:

DRAGONS [Liu et al. 2016]: The LF results from Liu et al. (2016) rely on the Dark-ages Reionization And Galaxy-formation Observables from Numerical Simulations (DRAGONS) ⁹ project, which builds semi-numerical models of galaxy formation on top of halo trees derived from N-body simulations done over different box sizes to probe a large dynamical range. The semi-numerical models include gas cooling physics, star formation prescriptions, and feedback and merging prescriptions, among other components of the model. The turn-over in the LF results of Liu et al. (2016) at ∼−11.5 mag correspond to the approximate halo masses ∼10⁸ M_⊙ where the gas temperature is 10⁴K. Above this temperature, atomic cooling processes become efficient.

CROC [Gnedin 2014, Gnedin 2016]: The model LF results for the Cosmic Reionization On Computers (CROC) were computed using gravity + hydrodynamical simulations executed with the Adaptive Refinement Treement (ART) code (Kravtsov 1999; Kravtsov et al. 2002; Rudd et al. 2008). A wide variety of physical processes, including gas cooling and heating, molecular hydrogen chemistry, star formation, stellar feedback, and radiative transfer of ionizing and UV light from stars is included in these simulations and done 20 h⁻¹ Mpc boxes at a variety of resolutions. There is a flattening in the effective slope of CROC LFs to fainter magnitudes, with a peak at ∼−12 mag. However, the peak at ∼−11.5 mag is reported to depend on the minimum particle size in the simulations and thus not to be a robust result of the simulation.

Finlator et al. (2015), Finlator et al. (2016), Finlator et al. (2017) [F17]: The Finlator et al. (2015, 2016, 2017) LF results are derived from a cosmological simulation of galaxy formation in a Mpc³ volume of the universe including both gravity and hydrodynamics. It is implemented in the GADGET-3 code (Springel 2005). Gas cooling has been added to this code through collisional excitation of hydrogen and helium (Katz et al. 1996). Metal line cooling is implemented using the collisional ionization equilibrium tables from Sutherland & Dopita (1993). Star formation is included using the Kennicutt–Schmidt law, with supernovae (SNe) feedback implemented following the "ezw" prescription from Davé et al. (2013) and metal enrichment from SNe as in Oppenheimer & Davé (2008). Less efficient gas cooling at lower halo masses results in a flattening of the Finlator et al. (2015, 2016, 2017) LF results at the faint end, with the turn-over in the LF occurring at ∼−11.5 mag.

Park et al. (2019) [P19]: The Park et al. (2019) LF results are based on a flexible, physically motivated modeling of star formation in galaxy halos. In their model, Park et al. (2019) took the SFE, the SFE scaling with halo mass, and a turn-over mass to be free parameters, which they then fit to the LF constraints from Bouwens et al. (2015a, 2017b) and Oesch et al. (2018a). In their fits, Park et al. (2019) allowed for the turn-over mass to be between 10⁸ M_⊙ and 10¹⁰ M_⊙. Given that the tuning of the Park et al. (2019) LF model to match the observations of Bouwens et al. (2015a) and Bouwens et al. (2017b), it is not especially surprising that their model fits our new observational constraints quite well. The approximate turn-over luminosity in the Park et al. (2019) results occurs at ∼−11.3 mag.

Yung et al. (2019) [Y19]: The Yung et al. (2019) LF results are based on a recent version of the Santa Cruz semi-analytic model (Somerville et al. 2015), which includes not only merger trees constructed by a standard Press-Schechter formalism (Lacey & Cole 1993), but also gas cooling, star formation, chemical evolution, and SNe-driven winds, photoionization feedback, and a critical molecular hydrogen surface density necessary for star formation. Yung et al. (2019) produced their results to provide semi-analytical model forecasts for JWST and rely on halos with circular velocities V_vir ≈ 20–500 km s⁻¹. Yung et al. (2019) report that SNe feedback plays the dominant role in flattening the LF at the faint end. The turn-over at M_UV,AB ∼ − 9 mag is imposed as a result of the atomic cooling limit in halos with V_vir ≈ 20 km s⁻¹ and is thus not a resolution effect.

CoDa2 [Ocvirk et al. 2016, 2020–O20]: The Cosmic Dawn (CoDa) simulations use the RAMSES-CUDATON code (Ocvirk et al. 2016) to execute a full modeling of both gravity + hydrodynamics + radiative transfer for a large ∼(100 Mpc)³ volume of the universe. The simulations include standard prescriptions for star formation and SNe explosions following standard recipes (Ocvirk et al. 2008; Governato et al. 2009, 2010). One key feature of the CoDa simulations is the inclusion of radiative transfer into the simulations through the ATUN code (Aubert & Teyssier 2008), in the sense that hydrodynamics and radiative transfer are now fully coupled. As a result, the effects of photoionization heating on low-mass galaxies are fully included in the CoDa simulations. Ocvirk et al. (2016, 2020) report that radiative feedback plays a big role in suppressing star formation in low-mass galaxies and modulating the very faint end (M_AB > −11) of the LF, resulting in a faint-end turn-over to the UV LF at ≈−11 mag.

FirstLight [Ceverino et al. 2017]: The model LF results from the FirstLight project are based on zoomed-in simulations of galaxies with circular velocity between 50 km s⁻¹ and 250 km s⁻¹. The galaxy simulation results are executed using the ART gravity + hydrodynamics code (Kravtsov et al. 1997; Kravtsov 2003). This code also includes gas cooling (atomic hydrogen, helium, metal, and molecular hydrogen), photoionization heating, star formation, radiative feedback, and SNe feedback. Ceverino et al. (2017) report that stellar feedback drives a flattening of their LF results at the faint end, i.e., M_UV,AB > − 14 mag, with an approximate turn-over luminosity ≈−11.5 mag.

Renaissance [O'Shea et al. 2015]: The "Renaissance" simulations (O'Shea et al. 2015) are zoomed-in simulations of a (28.4Mpc/h)³ volume of the universe, powered by the Enzo code (Bryan et al. 2014). This code self-consistently follows the evolution of gas and dark matter, includes H₂ formation and destruction from photodissociation, and includes star formation and SNe physics. Ionizing and UV radiation are produced as given by Starburst99 (Leitherer et al. 1999). Individual dark matter particles in the simulations have masses of 2.9 × 10⁴ M_⊙, and thus the smallest resolved halos in the simulation have masses of 2 × 10⁶ M_odot (∼70 particles/halo). A detailed description of the implementation of the physics and subgrid recipes is provided in Chen et al. (2014) and Xu et al. (2013, 2014). In the "Renaissance" simulations, flattening in the UV LF directly results from the decreasing fraction of baryons converted to stars in the lowest-mass halos, as a result of radiative feedback and less efficient gas cooling, with an approximate turn-over luminosity at ≈−8.5 mag. The presented results from O'Shea et al. (2015) are at z ∼ 12 where results are available and compared with our z ∼ 9 LF results.

Delphi [Dayal et al. (2014, 2022)]: Delphi uses a binary merger tree approach (Parkinson et al. 2008) to jointly track the buildup of dark matter halos and their baryonic components (gas, stellar, metal, and dust mass). This model follows the assembly histories of z ∼ 4.5 galaxies with halo masses 8–14 up to z ∼ 40. The SFE in any halo is the minimum between that required to eject the rest of the gas and an upper maximum threshold value. The flattening of the UV LF with decreasing redshift is driven by a flattening of the halo mass function coupled with a decrease in the gas mass as a result of the SN feedback experienced by a galaxy over its entire assembly history. The approximate turn-over luminosity ranges from −12 mag to −10 mag.

ASTRAEUS [Hutter et al. 2021] Jeans Mass + Early Heating: The Astraeus framework couples an N-body simulation (160 comoving Mpc; dark matter mass resolution of 6.2 × 10⁶ h⁻¹ M_⊙) with a modified version of the Delphi model for galaxy formation and a semi-numerical scheme for reionization. The authors introduce a filtering mass below which baryonic fluctuations can be suppressed due to reionization heating and explore six models for such reionization feedback. Here, we explore two models: (i) the early heating model where reionization feedback is time-delayed and has a weak to intermediate impact; and (ii) the Jeans mass model, which results in an instantaneous and maximum radiative feedback. The flattening at the faint end of the UV LF is a result of the impact of feedback (both SN and radiative) and the simulation resolution, with a turn-over luminosity occurring between ∼−12 and −14 mag for the early heating simulation and between −8 and −12 mag for the Jeans mass simulation.

THESAN PROJECT [Smith et al. 2022; Garaldi et al. 2022; Kannan et al. 2022]: The galaxy UV LF results from the THESAN project are based on a radiation-magneto-hydrodynamics simulation of a large volume of the universe (95.5 cMpc on a side) that models both the large-scale statistical properties of the intergalactic medium and the galaxies that lie within the volume. The flagship simulation resolves baryonic and dark matter masses down to 5.8 × 10⁵ M_⊙ and 3.1 × 10⁶ M_⊙, respectively. The simulations are executed using efficient code AREPO-RT (Kannan et al. 2019), a radiation hydrodynamics extension to the hydrodynamics code AREPO (Springel 2010). Star formation, black accretion, SNe winds, and other subgrid physics are implemented using the recipes developed and tested as part of the IllustrisTNG simulations (Vogelsberger et al. 2014). The UV LF predicted by the highest-mass component of the THESAN project shows a faint-end turn-over at ∼−12 mag, which is largely set by the resolution limit of the lowest-mass galaxies in the simulation.

Yue et al. (2016) [Y16]: Yue et al. (2016) made use of a semi-analytic formalism to predict the evolution of the UV LF. Yue et al. (2016) started with the halo mass function, broke up the star formation history of each halo into segments according to which the halo grows in mass by a factor of two, and then assumed that the SFR must be such to maintain a constant stellar mass–halo mass relation, which they calibrated to the z ∼ 5 LF of Bouwens et al. (2015a) and Bouwens (2016c). This approach is very similar to what Mason et al. (2015) employed in predicting the evolution of the UV LF (see also Trenti et al. 2010 and Tacchella et al. 2013). Yue et al. (2016) then looked into the impact that radiative feedback would have during the reionization era. Yue et al. (2016) then considered the star formation to be quenched in galaxies below some fixed circular velocity. Here we show their results where the quenching occurs below a circular velocity of 30 km s⁻¹, where the turn-over in the UV LF occurs at ∼−11.5 mag.

Of course, the aforementioned simulation efforts of the galaxy UV LF at z ≥ 4 are not exhaustive and are restricted to those studies that probe the faint-end form of the UV LF, i.e., ≪−16 mag, and thus do not include other very sophisticated modeling efforts such as conducted by ASTRID (Bird et al. 2022), Bluetides (Wilkins et al. 2017), and FLARES (Vijayan et al. 2021).

Finally, and as in Bouwens et al. (2017b), we again compare with the empirical constraints on the form of the UV LFs at z ∼ 7 from Boylan-Kolchin et al. (2015).

Boylan-Kolchin et al. (2015): Boylan-Kolchin et al. (2015) constrained the faint end of the z ∼ 7 LF by leveraging sensitive probes of the color–magnitude relationship of nearby dwarf galaxies to estimate their luminosity at z ∼ 7. By comparing the inferred luminosity distribution for these dwarfs with their implied numbers extrapolating z ∼ 7 LFs to −10 mag, Boylan-Kolchin et al. (2015) inferred a break in the LF at z ∼ −13 mag and transitioned from a faint-end slope of ∼−2 to ∼−1.2.

We compare the predicted UV LFs from these theoretical and empirical models with our new observational constraints on the faint-end form of the UV LFs in Figures 12 and 13. We also present the approximate turn-over luminosities and faint-end slopes derived from the theory LFs in Figure 11 as a function of redshift. The dark- and light-gray shaded regions show the 68% and 95% confidence constraints, respectively, by marginalizing over our LF fit results.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In general, we find good overall agreement between our observational findings and the general form of the model and empirical UV LFs. The faint-end slopes for all models show roughly the same trend with redshift as we see in the observations, although several models lie above and below the trends we derive here.

Additionally, the approximate UV luminosity of the turn-over in the theory LFs, which lie in the range ∼−13.5 to ∼−9 mag, is consistent with our observational constraints (Figure 11), which generally constrain the turn-over to be fainter than −15.5 or −14.3 mag depending on the redshift. At z ∼ 3, we obtain the tightest constraints on the UV luminosity of the turn-over. We find the turn-over to be fainter than −13.1 mag (95% confidence). Unfortunately, none of the theoretical models provide predicted LFs at z ∼ 3, but those that do so at z ∼ 4, i.e., Finlator et al. (2017) and Yung et al. (2019), are consistent with showing a turn-over at ∼−13.1 or fainter.