Deep Realistic Extragalactic Model (DREaM) Galaxy Catalogs: Predictions for a Roman Ultra-deep Field

We use a Planck 2018 cosmology (Planck Collaboration et al. 2020; Ω_b = 0.04893, Ω₀ = 0.3111, Ω_Λ = 0.6889, H₀ =67.66 km/s/Mpc, σ₈ = 0.8102, and n_s = 0.9665), with a box size of 115 h⁻¹Mpc and N = 2048³ particles. This corresponds to a particle mass of 1.5 × 10⁷ M_⊙ h⁻¹. The initial conditions are created using Music (Hahn & Abel 2011), and the simulations are run in Gadget-2 (Springel 2005), with a softening length of 1.13 co-moving h⁻¹kpc. The simulation outputs are shown in Figure 3. We generate halo catalogs and merger trees using Rockstar (Behroozi et al. 2013a) and Consistent Trees (Behroozi et al. 2013b). We define halo masses using the Bryan & Norman (1998) virial definition. Halos are required to have at least 20 particles, which corresponds to a minimum mass of 8.50 M_⊙ h⁻¹.

Zoom In Zoom Out Reset image size

Given the halo lightcone catalog constructed in Section 3, we assign a galaxy stellar mass, M_gal, to each halo. The SHAM procedure we use can be expressed mathematically as:

$$\begin{eqnarray}&&{\int }_{x}^{\infty }n(x^{\prime} ,z){\rm{d}}x^{\prime} ={\int }_{{M}_{\mathrm{gal}}}^{\infty }\phi ({M}_{\mathrm{gal}}^{{\prime} },z){\rm{d}}{M}_{\mathrm{gal}}^{{\prime} }\end{eqnarray} \tag{ 1 }$$

where ϕ(M_gal, z) is the SMF (i.e., the co-moving number density of galaxies per unit galaxy mass, per unit redshift), and n(x) is the number density of halos as a function of some mass proxy, x.

In the simplest formulation of SHAM, galaxy mass is directly mapped to halo mass. However, subhalo properties at the time of accretion are a better indicator of galaxy mass, since as subhalos undergo tidal stripping, they lose a significant amount of dark matter before the galaxy is disrupted (e.g., Conroy et al. 2006; Vale & Ostriker 2006). There are many different choices of mass proxies, and each of these choices produces slightly different spatial distributions of galaxies. We adopt the peak of the maximum circular velocity over the entire merger history, V_peak, as the halo mass proxy, because SHAM with V_peak as the mass proxy is known to reproduce the small-scale clustering of galaxies (Hearin et al. 2013; Reddick et al. 2013; Lehmann et al. 2017; Campbell et al. 2018).

We use the SMFs from Williams et al. (2018), which are continuously evolving double-Schechter functions for both star-forming and QGs out to redshifts z > 10. By construction, these SMFs match observational constraints at low redshifts (z ≤ 4), and reproduce known UVLFs when convolved with the distribution of rest-frame UV magnitudes, M_UV, i.e.,

$$\begin{eqnarray}&&\begin{array}{l}{{\rm{\Phi }}}_{\mathrm{UV}}({M}_{\mathrm{UV}},z)\\ \quad ={\displaystyle \int }_{0}^{\infty }\phi ({M}_{\mathrm{gal}}^{{\prime} }){ \mathcal N }[{M}_{\mathrm{UV}},{\bar{M}}_{\mathrm{UV}}({M}_{\mathrm{gal}}^{{\prime} },z),{\sigma }_{\mathrm{UV}}]{\rm{d}}{M}_{\mathrm{gal}}^{{\prime} },\end{array}\end{eqnarray} \tag{ 2 }$$

where ${ \mathcal N }[{M}_{\mathrm{UV}},{\bar{M}}_{\mathrm{UV}}({M}_{\mathrm{gal}},z),{\sigma }_{\mathrm{UV}}]$ is a normal distribution centered on the average UV magnitude ${\bar{M}}_{\mathrm{UV}}$, and σ_UV is the scatter in the M_gal–${\bar{M}}_{\mathrm{UV}}$ relation. Williams et al. (2018) provided a relation to describe the M_gal–${\bar{M}}_{\mathrm{UV}}$ relation and scatter, which we use to model our galaxies (see Section 6.4).

For redshifts z ≤ 4, Williams et al. (2018) fit the double-Schechter function to the data from Tomczak et al. (2014). For QGs, above z ∼ 4, the Schechter parameters are extrapolated to higher redshifts. This extrapolation is in agreement with the few constraints for QGs at z > 3.5. For star-forming galaxies (SFGs), the high-redshift SMFs are inferred from UVLFs. Specifically, Williams et al. (2018) used data from Bouwens et al. (2015) for 4 ≤ z ≤ 8 and from Oesch et al. (2018) at z = 10 to fit the SMF in Equation (2). Beyond z ∼ 10, the SMFs for the SFGs are extrapolated.

We perform SHAM on the halo catalogs using the SMFs described above. Specifically, given N halos, we first find the minimum galaxy mass, ${M}_{\min }$, such that ${N}_{\mathrm{gal}}(\gt {M}_{\min })=N$, where

$$\begin{eqnarray}&&{N}_{\mathrm{gal}}(\gt {M}_{\min })=V{\int }_{{M}_{\min }}^{\infty }\phi ({M}_{\mathrm{gal}}^{{\prime} }){\rm{d}}{M}_{\mathrm{gal}}^{{\prime} },\end{eqnarray} \tag{ 3 }$$

and V is the co-moving volume of the simulation. We then sample N galaxy masses from the SMF above ${M}_{\min }$. Finally, we rank-order the galaxy masses and assign galaxy masses to halos (such that the halo with the largest V_peak value gets the largest galaxy mass).

Scatter is commonly introduced in the relation between galaxy and halo properties by either deconvolving the SMF (e.g., Behroozi et al. 2010) or by directly adding scatter to the stellar masses, re-ranking, and iteratively solving for the galaxy mass (Hearin et al. 2013). Since our main goals are not to reproduce the observed scatter in the stellar-mass–halo-mass relation (SHMR), we did not include scatter in the SHAM procedure. As illustrated in Section 4.3, we are able to produce realistic statistical galaxy properties, such as galaxy clustering with our approach.

We classify each galaxy as either an SFG or a QG by randomly generating a number and comparing it to the probability a galaxy of that mass is star-forming, as calculated from the SMFs. This method does ignore a possible correlation between star formation rate (SFR) and mass accretion histories (e.g., Behroozi et al. 2019), which could potentially impact galaxy clustering with SFR. In particular, QGs are more likely to be found in denser environments (e.g., Kauffmann et al. 2004; Kimm et al. 2009). Finally we remove all galaxies from the catalog with M_gal < 10⁵ M_⊙. As shown in Section 8, the fraction of detectable galaxies below this mass is zero for redshifts z > 1. ¹²

For this work, we consider a survey with a 1 deg² FOV. To reach redshifts z ≳ 10, we tile 60 boxes in the line-of-sight direction. This corresponds to a maximum redshift of 13 and a co-moving distance of 6900 Mpc. As in, e.g., Bernyk et al. (2016), to avoid the replication of structures viewed by the observer, we randomly translate, reflect, and permute the axes of every tiled simulation box. We construct the halo lightcone by finding where each halo first crossed the observer's past lightcone (i.e., the location at which light has had just enough time to reach the observer). Technical aspects associated with this procedure are outlined extensively in the literature (e.g., Evrard et al. 2002; Blaizot et al. 2005; Kitzbichler & White 2007; Merson et al. 2013; Bernyk et al. 2016; Smith et al. 2017; Korytov et al. 2019). In this section, we closely follow the procedure used to generate the CosmoDC2 sky catalog for the Vera C. Rubin Observatory (Hollowed 2019; Korytov et al. 2019).

We assume that the observer is at the origin of the coordinate system (x, y, z) = (0, 0, 0). For each tiled box, we begin by placing all of the host halos on the lightcone. For every host halo in snapshot j, we determine its position in the subsequent snapshot, r_j+1. We extrapolate the halo position from j by assuming constant velocity to find its extrapolated position in snapshot j + 1, r_j+1,extrap. If the halo did not have a descendant located in snapshot j + 1, we set r_j+1 = r_j+1,extrap. Otherwise, we use:

$$\begin{eqnarray}&&{r}_{j+1}={r}_{j+1,\mathrm{extrap}}+\mathrm{int}\left(\displaystyle \frac{{r}_{j+1,\mathrm{extrap}}-{r}_{j+1,\mathrm{desc}}}{L}\right)\times L,\end{eqnarray} \tag{ 4 }$$

where L is the length of the simulation box, "int" represents integer, and r_j+1,desc is the position of the descendant in snapshot j + 1. This approach allows the halo to cross the edge of the simulation box between snapshots j and j + 1, and thus r_j+1 might be outside the domain of the tiled box. In these cases, we also consider the scenario where the halo crosses the lightcone on the other side of the box; i.e., we apply periodic boundary conditions, such that

$$\begin{eqnarray}\begin{array}{rcl}{r}_{j+1}^{p} & = & {r}_{j+1}\,\mathrm{mod}\,L\\ {r}_{j}^{p} & = & {r}_{j}+({r}_{j+1}^{p}-{r}_{j+1}),\end{array}\end{eqnarray} \tag{ 5 }$$

to allow the halo to cross the observer lightcone between positions r_j ^p and ${r}_{j+1}^{p}$.

Given the positions r_j and r_j+1 for each host halo, we calculate the time the halo would cross the past lightcone, t_e (Equations (27)–(29) in Korytov et al. 2019). If the halo crosses between the snapshot times t_j and t_j+1, we calculate the position on the lightcone, r_e from:

$$\begin{eqnarray}\begin{array}{rcl}{v}_{\mathrm{lin}} & = & ({r}_{j+1}-{r}_{j})/({t}_{j+1}-{t}_{j})\\ {r}_{e} & = & {r}_{j}+{v}_{\mathrm{lin}}({t}_{e}-{t}_{j})\,.\end{array}\end{eqnarray} \tag{ 6 }$$

We do not include any halos where r_e is beyond the domain of the tiled box.

The final step to generating the lightcone is to assign halo properties (e.g., mass and substructure) to each object. One approach to assigning halo properties is to allow the merger of halos to happen at a time randomly between the two snapshot (e.g., Smith et al. 2017). However, this approach will result in halos being double counted, and the lightcone needs to be carefully pruned. Since our simulations have very fine time resolution (500 snapshots between redshifts z = 0 and z = 20), we can avoid the technical difficulties and assumptions introduced by needing to prune the catalog, and follow the same approach from Korytov et al. (2019). Specifically, every halo that crosses the lightcone between snapshots j and j + 1 is assigned properties from snapshot j. When including the substructure, we ensure subhalos have the same position and velocity offset from the host as in snapshot j.

Given the lightcone catalog consisting of halo properties, galaxy masses, and positions, we cut out a wedge in the survey volume corresponding to 1 deg². Following Bernyk et al. (2016), we convert the co-moving (x, y, z) positions of each galaxy to an angular position as follows:

$$\begin{eqnarray}\begin{array}{rcl}d & = & \sqrt{{x}^{2}+{y}^{2}+{z}^{2}}\\ {\rm{R}}.{\rm{A}}. & = & \displaystyle \frac{360}{2\pi }\arctan (y/x)\\ \mathrm{decl}. & = & \displaystyle \frac{360}{2\pi }\arcsin (z/d).\end{array}\end{eqnarray} \tag{ 7 }$$

We only consider galaxies with R.A.< 1 deg and decl. < 1 deg, and then center the galaxies on (R.A., decl. ) = (0, 0). ¹³

Figure 4 shows the resulting galaxy mass density as a function of position. The lightcone procedure results in a realistic distribution of galaxies that traces the underlying cosmic web. Since the survey will reach very deep distances, the survey is tiled in six rows, with each row showing 10 tiled boxes. The survey wedge is complete to a maximum distance of 6588 Mpc/h (which corresponds to a cosmological redshift of ∼10.5); at higher redshifts, the angle of the survey is larger than the width of the simulation box. There are discontinuities in the galaxy density where the simulation boxes are tiled, which is a common feature for lightcones. As addressed in Bernyk et al. (2016), statistical properties of the galaxy catalog will be accurate on scales smaller than the box size.

Zoom In Zoom Out Reset image size

To verify the lightcone pipeline results in a realistic galaxy population, we calculate the halo mass functions (HMFs), SMFs, SHMRs, and galaxy clustering of the synthetic galaxy catalog.

4.3.1. Halo Mass Functions

To demonstrate that the HMF of the halo lightcone catalog is correct, we plot the recovered HMF of host halos in Figure 5. The theoretical curve is from the parameterization from Despali et al. (2016), and calculated using the HMF routine in Colossus (Diemer 2018). The theoretical curve and simulation results agree very well. The catalog was constructed to be complete above galaxy masses M_gal = 10⁵ M_⊙, which corresponds to a halo mass limit of approximately M_halo ∼ 10^9.5 h⁻¹ M_⊙

Zoom In Zoom Out Reset image size

4.3.2. Stellar Mass Functions and the Stellar-to-halo Mass Relations

Additionally, we show the SMF of the lightcone catalog in Figure 6 for both SFGs and QGs. For comparison, we show the volume-weighted average SMFs from Williams et al. (2018; dotted lines), averaged over the redshift bin. The catalog is complete at all redshifts for galaxy masses greater than 10⁵ M_⊙. The galaxies in the lightcone agree with the SMFs from Williams et al. (2018), which verifies that the lightcone pipeline and abundance matching procedure reproduce realistic galaxy counts. As discussed in Section 3.2, producing the desired SMF is important for reproducing observed luminosity functions.

Zoom In Zoom Out Reset image size

4.3.3. Galaxy Clustering

Galaxy clustering is commonly described using the two-point correlation function (2PCF), ξ(r), defined as

$$\begin{eqnarray}&&{\rm{d}}P={n}_{g}(r)[1+\xi (r)]{\rm{d}}V,\end{eqnarray} \tag{ 8 }$$

where dP is the excess probability above the Poisson noise, n_g (r) is the mean density of galaxies at separation r, and dV is the differential volume. In practice, the projected 2PCF, w_p , is often used. The 2PCF depends on the projected separation, r_p , and the line-of-sight separation, π:

$$\begin{eqnarray}&&{w}_{p}={\int }_{-\infty }^{\infty }\xi ({r}_{p},\pi ){\rm{d}}\pi \approx 2\times {\int }_{0}^{{\pi }_{\max }}\xi ({r}_{p},\pi ){\rm{d}}\pi ,\end{eqnarray} \tag{ 9 }$$

where the second equality assumes isotropy. The upper limit, ${\pi }_{\max }$, needs to be chosen large enough to nullify the effect of peculiar velocities, but not so large that it creates artificial edge effects from the survey boundary. The optimal value for ${\pi }_{\max }$ depends on the underlying survey volume, but is typically in the range 40–80 Mpc (Zehavi et al. 2005; van den Bosch et al.2013).

We calculate ξ based on the commonly used Landy-Szalay estimator (Landy & Szalay 1993),

$$\begin{eqnarray}&&\xi \approx \displaystyle \frac{{DD}-2{DR}+{RR}}{{RR}},\end{eqnarray} \tag{ 10 }$$

where DD, DR, and RR are the normalized number of counts from data–data, data–random, and random–random pairs. The random catalog contains 10³ points per ${\mathrm{arcmin}}^{2}$, with distances assigned such that the random catalog has the same redshift distribution as the synthetic catalog, i.e.,

$$\begin{eqnarray}&&\begin{array}{l}N({z}_{\min }\lt z\lt {z}_{\max })\\ \quad =\theta \sin (\theta ){\displaystyle \int }_{\chi ({z}_{\min })}^{\chi ({z}_{\max })}{\chi }^{2}{\displaystyle \int }_{{M}_{1}}^{{M}_{2}}\phi (M,z(\chi )){\rm{d}}M{\rm{d}}\chi \end{array}\end{eqnarray} \tag{ 11 }$$

where θ is the survey angle, and ϕ is the Williams et al. (2018) SMF. Values for R.A. and decl. coordinates for each random galaxy were selected so that galaxies are distributed isotropically:

$$\begin{eqnarray}\begin{array}{rcl}{\rm{R}}.{\rm{A}}. & = & {r}_{1}\theta \\ \mathrm{decl}. & = & \arcsin \left[2{r}_{2}\sin \left(\displaystyle \frac{\theta }{2}\right)\right]\,,\end{array}\end{eqnarray} \tag{ 12 }$$

where r₁, r₂ are random numbers selected uniformly between 0 and 1. The random catalog is then centered on (R.A., decl.) = (0,0), in the same way as the synthetic catalog.

To measure the projected 2PCF, we use the package Corrfunc (Sinha & Garrison 2020), with ${\pi }_{\max }\,=60\,{h}^{-1}\,\mathrm{Mpc}$. We estimate errors in the measured 2PCF by bootstrapping the galaxy catalog data with 200 subsamples. We present this measured 2PCF in Figure 7, along with Sloan Digital Sky Survey (SDSS) data for redshift z ≈ 0.1 galaxies from Yang et al. (2012). To compare our synthetic catalog 2PCF to the data, we use galaxies in the redshift range 0.1 < z < 0.2 within a 10 deg × 10 deg survey. The low-redshift clustering of the synthetic galaxies agrees with observations within the 1σ error bars, indicating that the lightcone halo clustering is realistic.

Zoom In Zoom Out Reset image size

Galaxies sizes are known to decrease with increasing redshift at a fixed stellar mass (e.g., van der Wel et al. 2014; Shibuya et al. 2015; Curtis-Lake et al. 2016), typically evolving as

$$\begin{eqnarray}&&{R}_{e}\propto {\left(1+z\right)}^{-\alpha },\end{eqnarray} \tag{ 13 }$$

where R_e is the half-light radius. Further, SFGs and QGs evolve differently, with SFGs having lower values of α (e.g., van der Wel et al. 2014; Ma et al. 2018). ¹⁴ In addition to these observed trends, theoretical models commonly show a correlation between the size of a galaxy and that of its host dark matter halo (e.g., Mo et al. 1998; Kravtsov 2013; Jiang et al. 2019). Given that our synthetic catalog begins with the underlying dark matter structure, we assign galaxy sizes based on the dark matter halo sizes, and then test that the size–redshift–mass relation is consistent with current data.

Specifically, we characterize the size of each galaxy as the half-light radius in the semimajor axis R_eff. ¹⁵ To assign R_eff to each galaxy, we use the relation between halo and galaxy size (e.g., Kravtsov 2013; Zanisi et al. 2020):

$$\begin{eqnarray}&&{R}_{\mathrm{eff}}=A\,{R}_{\mathrm{vir}},\end{eqnarray} \tag{ 14 }$$

where R_vir is the radius of the halo in physical units. We assign R_eff values from this mean relation, with log-normal scatter, σ_R . We use the coefficients A and scatter σ_R from Zanisi et al. (2020), as summarized in Table 1. Zanisi et al. (2020) used size distributions of central galaxies from the SDSS DR7 sample (Abazajian et al. 2009; Meert et al. 2016). While these coefficients capture the dependence of R_eff on mass and on whether a galaxy is star-forming or quiescent, they are calibrated for z = 0 galaxies.

Table 1. Coefficients, A, for the R_eff–R_vir Relation (Equation (14)) and the Scatter σ_R ; Values Are from Zanisi et al. (2020)

${\mathrm{log}}_{10}({M}_{\mathrm{gal}}/{M}_{\odot })$	A (SFGs)	A (QGs)	σR
≤9.5	0.018	0.006	0.2
(9.5,10]	0.019	0.007	0.2
(10,10.5]	0.019	0.010	0.15
(10.5,11]	0.019	0.011	0.15
(11,11.5]	0.019	0.015	0.15
>11.5	0.024	0.016	0.1

Download table as:  ASCIITypeset image

To test whether this relation extends reasonably to higher redshifts, we compare the median R_eff–z relation from our catalog to data from van der Wel et al. (2014) and Shibuya et al. (2015), as shown in Figure 8. The data sets shown in Figure 8 consist of galaxies from the 3D-HST+CANDELS catalog. Redshifts and masses are from Skelton et al. (2014), and galaxies are classified as either star-forming or quiescent using UVJ diagram criteria. The van der Wel et al. (2014) R_eff measurements are slightly higher than the Shibuya et al. (2015) measurements (see Figure 3 in Shibuya et al. 2015), due to the differences in the method used to fit R_eff. Our results are more consistent with Shibuya et al. (2015).

Zoom In Zoom Out Reset image size

Figure 8 demonstrates that our simple assumption that Equation (14) holds for all redshifts agrees remarkably well with the data. We capture the decrease in R_eff with both decreasing mass and increasing redshift. We also capture the effect that QGs are smaller than SFGs on average. Our ability to capture redshift trends, despite only using known relations between R_eff and R_vir at z = 0, likely reflects the fact that virial radii evolve as R_vir ∝1/(1 + z), and therefore higher-redshift galaxies have smaller virial radii. Our finding that the redshift evolution in the mass–size relation can be captured through the redshift evolution R_vir has also been found by Mowla et al. (2019), who show that the ratio r₈₀/R_vir is constant for redshifts z < 3.

The distribution of galaxy axis ratios and Sérsic indices should differ between SFGs and QGs (e.g., Franx et al. 2008; Bell et al. 2012; Mortlock et al. 2013) and with redshift (e.g., van der Wel et al. 2011; Guo et al. 2015). To assign shape (defined as the projected axis ratio q = b/a, where a is the semimajor half-light size and b is the semiminor half-light size) and Sérsic indices, n_s , we use the method directly from Williams et al. (2018). Specifically, using data compiled from van der Wel et al. (2012) and Skelton et al. (2014), Williams et al. (2018) found the distribution of q and n_s as a function of redshift. Due to limited data, morphologies for QGs with z > 4 are drawn from the 3 ≤ z ≤ 4 distribution, and morphologies for SFGs with z > 6 are drawn from the 5 ≤ z ≤ 6 distributions. We use these resulting distributions to draw q and n_s values to assign to each galaxy.

We assign PAs uniformly between 0 and 2π. Synthetic galaxy catalogs often assume that galaxies are oriented isotropically (e.g., Williams et al. 2018), and this assumption will not affect the galaxy clustering, number counts, or cosmic SFH of our synthetic catalog.

One of the main observables included in the synthetic catalog is the photometry in the Roman filters, which requires accurate modeling of the galaxy SEDs. We model the SEDs for each galaxy using the software FSPS (Conroy et al. 2009; Conroy & Gunn 2010). This section serves as an overview of the model used to generate the galaxy SEDs.

We use a Chabrier (2003) initial mass function (IMF), and include the IGM absorption model from Madau (1995). The SFH model is described in Section 6.2.1. Additionally, for SFGs, we include the FSPS nebular emission model (Byler et al. 2017), which is controlled by the gas ionization parameter and the gas metallicity. As in Williams et al. (2018), we approximate the stellar and ISM metallicities as a single metallicity value, Z_met.

The dust modeling includes dust emission and dust absorption. Dust emission follows Draine & Li (2007), which is a silicate-graphite-PAH grain model. For dust absorption, we use the Calzetti et al. (2000) attenuation curve, where dust attenuation is applied to all starlight equally, and therefore depends on one parameter, ${\bar{\tau }}_{v}$. This parameter gives the opacity at 5500 Å, and will set the normalization of the Calzetti dust attenuation curve. We assume dust attenuation is the same for emission line and continuum. We also include the asymptotic giant branch circumstellar dust models from Villaume et al. (2015). Including circumstellar dust can make a significant contribution to the IR emission from galaxies with little diffuse gas, and will need to be included in stellar modeling to accurately interpret data from upcoming IR facilities, including JWST.

Overall, this leaves seven free parameters to describe the galaxy SEDs:

• 1.  
  M_gal: the stellar mass of the galaxy.
• 2.  
  z: the redshift of the galaxy (or alternatively, t_age, the age of the universe at that redshift).
• 3.  
  Z_met: the galaxy metallicity. The gas-phase and stellar metallicities are assumed to be equal.
• 4.  
  t_start: the age of the universe at the start of star formation.
• 5.  
  τ: the e-folding time for star formation.
• 6.  
  ${\bar{\tau }}_{v}$: the dust attenuation parameter, defined as the opacity at 5500 Å.
• 7.  
  U_S : the gas ionization parameter.

In addition to these seven parameters, other relevant quantities included in the synthetic catalog include the SFR, the UV magnitude, M_UV, the slope of the UV continuum, β, and the rest-frame colors of the galaxies. Section 6.2 describes the calculation of each of these quantities.

Each galaxy begins with a redshift z (from the lightcone pipeline), a mass, M_gal (from abundance matching), and is labeled as either an SFG or a QG, as sampled from the SMF (see Section 3.2). Sections 6.4 and 6.5 below describe the SED pipeline for assigning the remaining FSPS parameters (Z_met, t_start, τ, ${\bar{\tau }}_{v}$, and U_S ) for SFGs and QGs, respectively. The pipeline is summarized in Figure 9.

Zoom In Zoom Out Reset image size

In brief, we first generate "parent catalogs" for both SFGs and QGs, spanning a realistic range of parameters. These parent catalogs serve as a lookup table for realistic galaxy SEDs. Details on how the parent catalogs are calculated are given in Section 6.3. We propose parameters for each galaxy, and find the 10 nearest neighbors in the parent catalog for each set of galaxy parameters using k-d trees. ¹⁶ For SFGs, the proposed parameters are M_gal, z, M_UV, and β while for QGs, the proposed parameters are M_gal and z. The distance metrics used to find the nearest neighbors are detailed in Sections 6.4 and 6.5.

We assign the five free FSPS parameters (Z_met, t_start, τ, ${\bar{\tau }}_{v}$, and U_S ) to each galaxy by taking the weighted average of the nearest neighbors. Given these five parameters, along with the mass and redshift of the galaxy, we use FSPS to generate an SED. We then calculate the UV properties, SFRs and rest-frame colors directly from the spectrum. These calculated properties will not be identical to the proposed parameters; however, we still retain the imposed scaling relations in the parent catalogs (see Appendix B).

For the SED pipeline outlined above, we need to accurately measure the SFR, rest-frame colors, and UV properties of the galaxy SEDs. This section details the methods we use to measure these quantities. Since FSPS normalizes the stellar modeling, such that one stellar mass is created over the formation history, we scale all of the spectra. Specifically, we scale all fluxes and the SFR by M_gal/f, where f is the surviving mass fraction in the stellar population, not including stellar remnants.

6.2.1. Star Formation History

Galaxy SFRs dictate the CSFRD of the universe, which in turn constrains the amount of ionizing photons produced by galaxies in the EoR. Given the importance of SFRs to galaxy evolution and the sources of ionizing photons during the EoR, we need to ensure our synthetic galaxies have SFRs consistent with observations.

The SFR depends on the age, e-folding time, and mass of the galaxy. We model the SFH, ψ(t), using a "delayed-tau model," in which:

$$\begin{eqnarray}&&\psi (t)\propto t\exp -t/\tau ,\end{eqnarray} \tag{ 15 }$$

where t is the time since the start of star formation, t_start. As discussed in Williams et al. (2018), this parameterization achieves the expectation from simulations that high-redshift galaxies have rising SFHs (Finlator et al. 2011), and accurately reproduces the colors and mass-to-light ratios of galaxies in smoothed particle hydrodynamics simulations (e.g., Simha et al. 2014).

To determine the SFR for each galaxy, we use the FSPS calculated SFR, and scale it by M_gal/f, as described above.

6.2.2. Rest-frame Colors

6.2.3. UV Properties

Roman will constrain the number of ionizing photons that are produced by galaxies in the EOR by precisely measuring the UVLF at high redshift. The UVLF depends on the SMF and the UV magnitude, M_UV (see Equation (2)).

Following Robertson et al. (2013), we define the UV magnitude, M_UV, as the average magnitude at rest-frame wavelength in a flat filter, in the range 1450–1550 Å. We measure M_UV directly from the rest-frame galaxy SEDs by calculating the average flux density:

$$\begin{eqnarray}&&\langle {f}_{\lambda }\rangle =\displaystyle \frac{{\int }_{{\lambda }_{1}}^{{\lambda }_{2}}{f}_{\lambda }d\lambda }{{\lambda }_{2}-{\lambda }_{1}},\end{eqnarray} \tag{ 16 }$$

with λ₁ = 1450Å and λ₂ = 1550Å.

In addition to the UV magnitude, the slope of the UV continuum, β, helps determine the role of galaxies in ionizing the universe. The slope β is defined as f_λ ∝ λ^β (e.g., Meurer et al. 1999). The presence of very blue UV continuum slopes in Roman imaging may provide a signpost of very high escape fractions in early galaxies. We use the method from Dunlop et al. (2012) to determine β. Specifically, we take the SEDs shifted to redshift z = 7, and calculate

$$\begin{eqnarray}&&\beta =4.43(J125-H160)-2\end{eqnarray} \tag{ 17 }$$

using the FSPS included filters F125W and F160W for WFC3.

We create SFG and QG parent catalogs that serve as a lookup table between the assigned FSPS parameters (M_gal, z, Z_met, t_start, τ, ${\bar{\tau }}_{v}$, and U_S ) and the derived spectra quantities (SFR, M_UV, β, and rest-frame colors). We create the parent catalogs by sampling from known scaling relations. By only populating regions of parameter space corresponding to the desired galaxy population, we can fully sample the relevant space, while also limiting the size of the parent catalog. The distributions we use to construct the parent catalogs are summarized in Tables 2 and 3 for SFGs and QGs, respectively, and outlined in detail in this section.

Table 2. Parameters for SEDs of SFGs

Parameter	Proposed Parameter	Parent Catalog	Distance Metric	Assigned Parameter
Mgal	fixed	uniform in ${\mathrm{log}}_{10}({M}_{\mathrm{gal}}/{M}_{\odot })\in [5,12]$	${\mathrm{log}}_{10}({M}_{\mathrm{gal}}/10{M}_{\odot })$	fixed
z	fixed	uniform in a = 1/(1 + z) ∈ [0.07, 1]	a = 10/(1 + z)	fixed
Zmet	...	FMR	...	nearest neighbors
US	...	US –Zmet relation	...	nearest neighbors
${\bar{\tau }}_{V}$	...	${ \mathcal N }(\mu =0,\sigma =0.5)\in [0,4]$	...	nearest neighbors
tstart	...	uniform in [1 Myr, tage]	...	nearest neighbors
τ	...	uniform in [0.1, 100] Gyr	...	nearest neighbors
ψ	...	calculated; cut ${\mathrm{log}}_{10}(\psi /{M}_{\odot }{\mathrm{yr}}^{-1})\in [-5,4]$	...	calculated
MUV	MUV–Mgal	calculated; cut in [−25, −5]	MUV + 20	calculated
β	β–MUV	calculated; cut in [−3, −1]	β	calculated
U-V	...	calculated	...	calculated
V-J	...	calculated	...	calculated

Note. The synthetic galaxies have fixed parameters M_gal and z and "proposed" parameters M_UV and β. We identify the 10 closest neighbors in the parent catalog (defined by the distance metric), and assign the final catalog properties Z_met, U_S , ${\bar{\tau }}_{V}$, t_start, and τ from the nearest neighbors.

Download table as:  ASCIITypeset image

6.3.1. Star-forming Galaxy Parent Catalog

For the SFG parent catalog, we sample the seven FSPS parameters for 10⁶ galaxies. The mass, redshift, and SFH parameters t_start and τ are sampled from uniform distributions, as summarized in Table 2. We sample the dust parameter, ${\bar{\tau }}_{V}$, uniformly from a normal distribution centered on zero with a standard deviation of 0.5, truncated between 0 and 4.

Since the metallicity of galaxies depends on the metal production from stars, SFR is related to metallicity. More massive galaxies have higher metallicities on average (the mass–metallicity relation; e.g Lequeux et al. 1979; Tremonti et al. 2004), and the scatter in this relation depends strongly on SFR. The relation between M_gal–Z_met–ψ is known as the fundamental metallicity relation (FMR; e.g., Mannucci et al. 2010; Hunt et al. 2016), and does not display any significant redshift evolution up to at least redshift z > 3 (e.g., Cresci et al. 2019, and references therein). Given these observed trends, we first assign realistic SFRs for galaxies, before assigning metallicities from the FMR.

To assign realistic SFRs, we use the parameterization from Schreiber et al. (2017), constructed to match observations in the redshift range 0 < z < 7. For SFGs,

$$\begin{eqnarray}&&\begin{array}{l}{\mathrm{log}}_{10}[\psi /({M}_{\odot }{\mathrm{yr}}^{-1})]\\ \quad ={\mathrm{log}}_{10}({M}_{\mathrm{gal}}/{M}_{\odot })-9.5+1.5{\mathrm{log}}_{10}(1+z)\\ \qquad 0.3{\left[\max (0,{\mathrm{log}}_{10}({M}_{\mathrm{gal}}/{M}_{\odot })-9.36-2.5{\mathrm{log}}_{10}(1+z))\right]}^{2}\end{array}\end{eqnarray} \tag{ 18 }$$

with a log-normal scatter of 0.3 dex.

We then propose Z_met by drawing from the FMR presented in Williams et al. (2018; based off of data from Hunt et al. 2016)

$$\begin{eqnarray}&&\begin{array}{l}{\mathrm{log}}_{10}({Z}_{\mathrm{met}}/{Z}_{\odot })+8.7\approx 12+{\mathrm{log}}_{10}({\rm{O}}/{\rm{H}})\\ \quad =\,-0.14\mathrm{log}+10(\psi /{M}_{\odot }{\mathrm{yr}}^{-1})\\ \qquad +\,0.37{\mathrm{log}}_{10}({M}_{\mathrm{gal}}/{M}_{\odot })+4.82\end{array}\end{eqnarray} \tag{ 19 }$$

with scatter from a Student's t-distribution with 3° of freedom, a standard deviation of 0.3 in ${\mathrm{log}}_{10}({Z}_{\mathrm{met}}/{Z}_{\odot })$, and truncated between − 2.2 < Z_met < 0.24. A Student's t-distribution is more heavily tailed than a Gaussian, and agrees better with observational data. As discussed in Williams et al. (2018), imposing a maximum metallicity, ${\mathrm{log}}_{10}({Z}_{\mathrm{met}}/{Z}_{\odot })=0.24$, in the synthetic galaxies reproduces the observed turnover in the mass–metallicity relation, despite the fact that Equation (19) is linear in mass. We use the same FMR to assign metallicities to both SFGs and QGs, as we find that existing measurements of QGs are consistent with the FMR used to create this synthetic catalog (e.g., Peng et al. 2015; Leethochawalit et al. 2018).

Gas-phase metallicity is also an indirect tracer of gas flows in galaxies. Metallicity and stellar mass both correlate with SFH, but gas flows will introduce scatter in the metallicity and stellar mass relation. Thus, the gas ionization parameter, U_S (which represents the ratio of the number density of ionizing photons to the number density of hydrogen atoms) correlates with metallicity. To assign the gas ionization parameter, we use the low-redshift relation from Carton et al. (2017):

$$\begin{eqnarray}&&{\mathrm{log}}_{10}{U}_{S}=-0.8{\mathrm{log}}_{10}({Z}_{\mathrm{met}}/{Z}_{\odot })-3.58\end{eqnarray} \tag{ 20 }$$

with a scatter of 0.3 dex sampled from a Student's t-distribution with three degrees of freedom, truncated between $-4\lt {\mathrm{log}}_{10}{U}_{S}\,\lt -1$. As in Williams et al. (2018), we assume the relation holds at all redshifts.

Given the seven FSPS parameters, we then calculate M_UV, β, and ψ as outlined in Section 6.2. We discard galaxies that have unrealistic ψ, M_UV, or β values, as summarized in Table 2. The resulting SFG parent catalog has ∼ 9 × 10⁵ galaxies.

6.3.2. Quiescent Galaxy Parent Catalog

We generate the parent catalog for the QGs in a similar way to the parent catalog for the SFGs. For 10⁷ galaxies, we assign galaxy mass, redshifts, and SFH parameters t_start and τ from uniform distributions, as summarized in Table 3. We sample the dust parameter, ${\bar{\tau }}_{V}$ uniformly from a normal distribution centered on zero with a standard deviation of 0.5, truncated between 0 and 4.

Table 3. Parameters for SEDs of QGs

Parameter	Proposed Parameter	Parent Catalog	Distance Metric	Assigned Parameter
Mgal	fixed	uniform in ${\mathrm{log}}_{10}({M}_{\mathrm{gal}}/{M}_{\odot })\in [5,12]$	${\mathrm{log}}_{10}({M}_{\mathrm{gal}}/10{M}_{\odot })$	fixed
z	fixed	uniform in [0, 13]	a = 10/(1 + z)	fixed
Zmet	...	FMR	...	nearest neighbors
${\bar{\tau }}_{V}$	...	${ \mathcal N }(\mu =0,\sigma =0.5)\in [0,4]$	...	nearest neighbors
tstart	...	uniform in [1 Myr, tage]	...	nearest neighbors
τ	...	uniform in [0.01, 10] Gyr	...	nearest neighbors
ψ	...	calculated; cut ${\mathrm{log}}_{10}(\psi /{M}_{\odot }\,{\mathrm{yr}}^{-1})\in [-4,1]$	...	calculated
U-V	...	calculated; cut in UVJ diagram	...	calculated
V-J	...	calculated; cut in UVJ diagram	...	calculated

Note. The synthetic galaxies have fixed parameters M_gal and z. We identify the 10 closest neighbors in the parent catalog (defined by the distance metric), and assign the final catalog properties Z_met, ${\bar{\tau }}_{V}$, t_start, and τ from the nearest neighbors.

Download table as:  ASCIITypeset image

We assume QG metallicities, Z_met, follow the same FMR as SFGs, as given in Equation (19), which is consistent with existing measurements of QGs (e.g., Peng et al. 2015; Leethochawalit et al. 2018). We propose SFRs from Schreiber et al. (2017) as before. For QGs,

$$\begin{eqnarray}&&\,{\mathrm{log}}_{10}[\psi /({M}_{\odot }{\mathrm{yr}}^{-1})]=0.5{\mathrm{log}}_{10}(M/{M}_{\odot })+{\mathrm{log}}_{10}(1+z)-6.1,\end{eqnarray} \tag{ 21 }$$

with a log-normal scatter of 0.45 dex.

As before, given the seven FSPS parameters, we then calculate ψ and the UVJ colors, as outlined in Section 6.2. We remove any galaxies in the QG parent catalog that do not fall in the proper space in the UVJ diagram, to avoid including any galaxies with unrealistic colors. Specifically, we use the criteria that a galaxy can be classified as SF if it satisfies any of the following conditions:

$$\begin{eqnarray}&&\begin{array}{l}U-V\lt 1.3\\ V-J\gt 1.6\\ U-V\lt 0.88(V-J)+0.49\end{array}\end{eqnarray} \tag{ 22 }$$

(Williams et al. 2009). After cuts in SFR and UVJ, the resulting parent catalog has 5 × 10⁴ galaxies. Though SFR is not set explicitly for QGs, these cuts ensure that QGs have a low SFR, differentiating them from the SFG population (see Appendix B).

This section outlines our procedure for assigning realistic SEDs for SFGs. To accurately capture the observed CSFRD and UVLFs, galaxies must have realistic SFRs and UV magnitudes. We also aim to produce realistic UV continuum slopes, β, the FMR, and the relationship between metallicity and gas ionization. Overall, we propose five parameters to define the SFG SEDs, ψ, M_UV, β, Z_met, and U_S . Details of how these proposed parameters are generated, and how the FSPS parameters are determined for the SFGs are given in this section, and summarized in Table 2.

To reproduce observed UVLFs, the synthetic galaxies need to follow the same ${\bar{M}}_{\mathrm{UV}}$–M_gal relation that Williams et al. (2018) used to generate the high-redshift SMFs used in this work:

$$\begin{eqnarray}\begin{array}{rcl}{\bar{M}}_{\mathrm{UV}} & = & -1.66({M}_{\mathrm{gal}}-{b}_{{M}_{\mathrm{UV}}})-20\\ {b}_{{M}_{\mathrm{UV}}} & = & \left\{\begin{array}{ll}0.12{z}^{2}-0.98z+11.4175 & \mathrm{for}\ z\leqslant 3.75\\ -0.12z+9.88 & \mathrm{for}\ 3.75\lt z\lt 8\\ 8.92 & \mathrm{for}\ \,z\geqslant 8\end{array}\right..\end{array}\end{eqnarray} \tag{ 23 }$$

Additionally, Williams et al. (2018) found that scatter observed in 3D-HST data (Skelton et al. 2014) is constant in both stellar mass and redshift, with an average value of σ_UV = 0.7. Therefore, we propose M_UV values from Equation (23) with a Gaussian scatter of σ_UV = 0.7.

We also use the relations from Williams et al. (2018) to propose UV continuum slopes, β,

$$\begin{eqnarray}\begin{array}{rcl}\beta & = & {m}_{\beta }({M}_{\mathrm{UV}}+19.5)+{b}_{\beta }\\ {b}_{\beta } & = & \left\{\begin{array}{ll}-0.09z-1.49 & \mathrm{for}\ \ z\leqslant 8\\ -0.77 & \mathrm{for}\ \ z\gt 8\end{array}\right.\\ {m}_{\beta } & = & \left\{\begin{array}{ll}-0.007z-0.09 & \mathrm{for}\ \ z\leqslant 8\\ -0.146 & \mathrm{for}\ \ z\gt 8\,.\end{array}\right.\end{array}\end{eqnarray} \tag{ 24 }$$

Given the proposed galaxy properties, we next determine which FSPS parameters to assign. To generate SEDs that match closely to the proposed parameters, we find the 10 nearest neighbors to the SFG parent catalog. Table 2 outlines the distance metric we use to find the nearest neighbors. We average the five free FSPS parameters (Z_met, t_start, τ, ${\bar{\tau }}_{v}$, and U_S ) of these neighbors, with a weight based on their distance to the proposed parameters. In cases where the assigned t_start value is greater than t_age, we set t_start = t_age. This method allows us to "fit" the target (proposed) parameters, even with the very large number of galaxies, something that would not be computationally realistic to achieve iteratively.

Given the large degeneracies between different galaxy properties (e.g., age, metallicity, dust, SFH), we find that our solution depends sensitively on how realistically we populate the parent catalogs, and also how we define the distance metric. In particular, we find that to ensure the final parameters were close to the proposed parameters, it is especially important to find neighbors that are close in redshift.

In addition to SFGs, we include a QG population. We turn off the nebular emission model for QGs, as low-redshift galaxies are known to have very low gas content. Though galaxy modeling often assumes QGs do not contain any gas or dust (e.g., Williams et al. 2018), recent results from Gobat et al. (2018) suggest that higher-redshift galaxies should have significant dust content. We therefore include dust for the QGs. We assign QG properties in a similar manner as to how we assign SFG properties. Given the QG parent catalog described in Section 6.3, we determine FSPS parameters using the 10 nearest neighbors in mass and redshift. We summarize the details of this procedure in Table 3.

We present the resulting SEDs for two example galaxies at redshift z ∼ 8 in Figure 10. The Roman WFI filters include seven photometric filters: R062, Z087, Y106, J129, H158, F184, and F213 ¹⁷ , a wide-field filter W146, a grism, and a prism. We plot the flux calculated in the seven photometric filters R062, Z087, Y106, J129, H158, F184, and F213. Both galaxies show a clear Lyman break in Y106.

Zoom In Zoom Out Reset image size

The SFG (top panel) should be selected as a Lyman-break galaxy (LBG) in Roman using the common Lyman-break method. This example also displays strong Lyα emission. Lyα emitters (LAEs) are a common way to study the EoR, and Roman will detect many such galaxies. Quantifying the science returns of LAEs in a Roman UDF will be the subject of future work.

Due to the SMFs we used to construct our catalog, we model a small population of high-redshift QGs, such as the one shown in the bottom panel of Figure 10. It is unclear whether redshift z = 8 galaxies exist, or what the mechanism of their formation would be. Due to the large area a deep survey with Roman will cover, a Roman UDF could help identify rare high-redshift QGs if they exist. This will be discussed further in Section 9.2.

Placing strong constraints on the faint end of the UVLF during the EoR is a main goal of a Roman UDF. To reproduce observed UVLFs, the synthetic galaxies must follow the correct distribution of M_UV, as discussed in Section 6.2.3. On average, the UV magnitude, M_UV, decreases (becomes brighter) with galaxy mass, M_gal (e.g., Stark et al. 2009). For galaxies with ${\mathrm{log}}_{10}({M}_{\mathrm{gal}}/{M}_{\odot })\gt 10$, the relationship flattens (e.g., Spitler et al. 2014). Apart from this flattening at high masses, the observed M_UV–M_gal relation has a constant slope, with the normalization evolving with redshift (e.g., Duncan et al. 2014; Stefanon et al. 2017). Similarly, the UV continuum slope, β, correlates with M_UV. Though the exact relation between β and M_UV is still being studied, β appears bluer (decreases) with increasing (fainter) M_UV and increasing redshift (e.g., Bouwens et al. 2014).

We calculate the UV properties, M_UV and β, of the synthetic SFGs as described in Section 6.2.3, and presented in Figure 11. We reproduce trends for both M_UV and β imposed from Williams et al. (2018). The M_UV–mass relation enables us to reproduce observed UVLFs, as discussed below. We also closely reproduce the β–M_UV relation for β > − 2.5. Below β < − 2.5, our galaxy catalog agrees with studies that indicate that the β–M_UV relationship flattens for M_UV ≳ − 19 galaxies (e.g., Bouwens et al. 2014).

Zoom In Zoom Out Reset image size

We show the UVLFs of our synthetic SFGs in Figure 12, compared to the UVLFs that are expected from convolving the SMF with a normal distribution, centered on M_UV (Equation (2)). The galaxy catalog agrees very well with the imposed relation, and is consistent with current observations, as shown in Section 8.3. Our catalog has slightly fewer galaxies at faint magnitudes (M_UV > − 17) than expected due to the fact that we do not perfectly recreate the M_UV–M_gal relationship. These magnitudes are fainter than the detection limits of a 30 m_AB survey, and our synthetic catalog UVLF is consistent with current observations, as discussed in Section 8.3.

Zoom In Zoom Out Reset image size

As discussed in Section 6.2.1, the CSFRD describes the ionizing photon contribution of galaxies as a function of cosmic time (for a review, see, e.g., Madau & Dickinson 2014). We display the evolution of the CSFRD in the top panel of Figure 13. We calculate the CSFRD by taking the sum of galaxy SFRs, ψ (averaged over the last 100 Myr), in each bin divided by the co-moving volume of the bin. As in Williams et al. (2018), to compare to Madau & Dickinson (2014), we imposed a luminosity limit of 0.03 L_*. ¹⁸

Zoom In Zoom Out Reset image size

Similar to Williams et al. (2018), we find that the CSFRD agrees roughly with the model presented in Madau & Dickinson (2014), but more closely to data (e.g., Cucciati et al. 2012; Oesch et al. 2018). Most notably, our synthetic galaxies agree with high-redshift (z > 8) measurements of the CSFRD. As in Williams et al. (2018), our synthetic galaxies have a slight deficit at redshifts 1 < z < 3, which reflect the underlying SMFs. As discussed in Williams et al. (2018), the SMFs we used do not include the dusty SFGs currently missed in the UV-selected samples.

We show the evolution of the specific star formation rate (sSFR) in the bottom panel of Figure 13, for galaxies with $8.8\lt {\mathrm{log}}_{10}({M}_{\mathrm{gal}}/{M}_{\odot })\lt 10$. We calculate the median sSFR by dividing the SFR averaged over the last 100 Myr by the galaxy mass. For comparison, we plot data from Tasca et al. (2015). The synthetic galaxies and observations agree very well. The sSFR of galaxies at redshifts z > 4 is a matter of active research, but we find an increase in sSFR, consistent with current studies (Stark 2016). Given that our galaxy catalog is consistent with observations for both galaxy clustering and cosmic SFHs, we now turn our attention to making preliminary predictions for the science returns of a Roman UDF.

Understanding of galaxy formation increased immensely when techniques to select high-redshift (z > 5) galaxies were developed (for a review, see, e.g., Dunlop 2013). The most commonly used method is the Lyman-break technique, which selects LBGs based on a step in the blue UV continuum emission at λ_rest = 1216 Å. To determine whether a UDF with Roman would select each galaxy in the synthetic catalog, we use a magnitude limit based on the Lyman break, in the same manner as Williams et al. (2018). Specifically, for each filter, we assign a nonoverlapping wavelength range (the dividing point was taken to be halfway between the centers of two adjacent bands). Each galaxy is assigned a "dropout" filter, according to the wavelength that corresponds to its Lyman break, λ = 1216(1 + z)Å. If a galaxy is brighter than the magnitude limit in the band immediately redder to its dropout filter, we consider it to be detectable. We use the magnitude limit m_AB = 30 for a 5σ detection.

We show the fraction of detectable galaxies in Figure 14. The mass at which ∼50% of galaxies ranges from ∼10⁶ at z ≈ 1–10⁷ at z = 10. We find that nearly 100% of galaxies are detectable above masses of M_gal = 10⁸ M_⊙, which indicates that there should be very tight constraints on galaxy properties (such as galaxy clustering and UVLFs) above this mass. Given this, we can measure the clustering of faint galaxies (M_UV ≲− 18) at the EoR. With completeness corrections, we can measure quantities such as the UVLFs for galaxies brighter than M_UV ≈ − 17 out to redshifts z ≈ 10. Section 8.3 will explore these potential constraints in more detail.

Zoom In Zoom Out Reset image size

In the following sections, we assume that any detectable galaxy is selected. This approach assumes that there is no scatter in the photometric redshifts, and no galaxies are obscured by forefront galaxies. Future work will examine Roman UDF constraints in more detail, taking into account these effects. The catalogs and images presented in this work will be instrumental to explore these sources of systematic errors.

Here we provide number counts of the detectable galaxies (5σ), to understand the science returns of a Roman UDF and measure the completeness of the survey. Figure 15 shows both the cumulative number counts (top), and the differential number counts (bottom). We find that a 1 deg² Roman UDF could detect more than 10⁴ galaxies at redshifts z > 7. For comparison, HST has detected on the order of 10³ galaxies at redshifts z ∼ 4–10 (Koekemoer et al. 2013; Lotz et al. 2017), and Cycle 1 JWST programs are expected to detect ∼5000 galaxies at redshifts z > 6 (Williams et al. 2018).

Zoom In Zoom Out Reset image size

Though a UDF catalog will consist largely of SFGs, there will also be a number of QGs, as demonstrated in Figure 16. These predictions indicate that a 1 deg² will contain a handful of detectable QGs above redshift z = 7 if they exist. In combination with, e.g., JWST spectroscopy, this could identify z > 7 QGs. Given that we currently have only detected QGs out to redshifts of z ∼ 4–5 (Merlin et al. 2019; Valentino et al. 2020), this will greatly enhance our understanding of the QG population. Determining the number and mass distribution of high-redshift QGs is fundamental to understanding galaxy evolution, and for testing different models for the emergence of QGs (see, e.g., Cortese et al. 2021, for a recent review on quenching).

Zoom In Zoom Out Reset image size

Due to the nature of these kinds of predictions, the underlying galaxy properties in the catalog are extrapolated beyond our current knowledge (in particular, the number of low-mass objects is dependent on extrapolations of the SMFs). Therefore, the number counts presented in this section represent a reasonable estimate of what a UDF might see, but there is much to be learned in the future. Regardless, the large number of galaxies that Roman will detect will provide the best census to date of faint high-redshift galaxies, including those thought responsible for reionizing the universe.

Some of the quantities that a Roman UDF will aim to constrain include the SMF, the UV luminosity function, galaxy clustering and the galaxy–halo connection at high redshifts. In this section we provide preliminary projected constraints for a 1 deg² Roman UDF, based on the detectable galaxies identified in Section 8.1. We have not included uncertainties from cosmic variance or from errors on photometric redshifts in our predictions, as that is beyond the scope this paper. However, we note that our synthetic catalog provides the basis for rigorous studies of systematic errors in a Roman UDF, and will be the focus of future work.

8.3.1. Stellar Mass Functions

The evolution of the SMF of galaxies over cosmic time provides information on how galaxy populations have evolved, the SFH of the universe, and the galaxy–halo connection. Figure 17 shows the SMF for the detectable galaxies in the synthetic catalog, with the underlying SMFs shown as dotted lines for both SFGs (top) and QGs (bottom). We perform a completeness correction by dividing the number of detectable galaxies by the fraction of detectable galaxies at the corresponding mass and redshift. As shown in Figure 6, our synthetic catalog matches the SMF nearly exactly by construction.

Zoom In Zoom Out Reset image size

For our catalog, the SMF for the detectable SFGs is measurable to masses of M_gal > 10^6.5 M_⊙ for all redshifts. For the detectable QGs, the SMF is measurable to redshifts z ∼ 7 (masses M_gal > 10^9.5 M_⊙). In comparison, observations have only recently placed any constraints on the SFG SMF at redshift z = 8–10. Further, these constraints only exist for galaxies with M_gal < 10⁸ M_⊙, and have very large uncertainties (about 100%; Stefanon et al. 2021). Measurements of QG SMFs only extend to redshifts z ≈ 4 (e.g., Girelli et al. 2019).

The exact constraints on an SMF from Roman will depend on stellar mass measurements and the accuracy of completeness corrections. This work provides a first calculation of the completeness of a Roman UDF (Section 8.1), and highly realistic simulated data to enable future rigorous calculations. In practice, stellar mass estimates of high-redshift galaxies will require photometry redshifts beyond the reddest Roman filters. However, JWST will provide tight constraints on the M_UV–M_gal relation, which can be used to estimate stellar masses of galaxies in a Roman UDF.

8.3.2. UV Luminosity Functions

To clearly answer whether there were enough galaxies during the EoR to reionize the universe, we need strong constraints on the faint end of the UVLF during the EoR. Bouwens et al. (2021) provide the most recent, comprehensive constraints on the UVLF, and are able to constrain the faint end (M_UV < − 17) out to redshift z ∼ 10. At redshift z ∼ 10, only eight sources are used to constrain the UVLF, so the uncertainties are very large (∼50%–100%). In comparison, a Roman UDF could detect ∼10³ sources at z ∼ 10, and thus greatly reduce the uncertainties in the UVLF.

We show the UVLF of the detectable synthetic galaxies in Figure 18, compared to the data from Bouwens et al. (2021). As in the previous section, we perform a completeness correction by dividing the number of selected galaxies by the fraction of selected galaxies at the corresponding mass and redshift. We predict a Roman UDF will be able to constrain the UVLF beyond redshift z ≈ 10, with galaxies M_UV < − 17. Due to the vast number of galaxies in a Roman UDF, it could constrain the faint-end UVLF to within ∼1 percent to redshifts z ≳ 10, far tighter than existing limits.

Zoom In Zoom Out Reset image size

We note that the bright end of the synthetic galaxy catalog UVLF is higher than the Bouwens et al. (2021) data at low redshifts. As shown in Figure 12, our galaxy catalog agrees with the underlying UVLF model we used at the bright end, and therefore this difference can be attributed to the SMFs used in this work. This excess of low-redshift (z < 4), bright galaxies is consistent with the possibility of inefficient mass quenching, low dust obscuration, or hidden AGN activity, as suggested by Harikane et al. (2021).

8.3.3. Galaxy Clustering

Roman's wide, contiguous FOV will allow for measurements of the 2PCF of faint galaxies at the time of reionization and thus their underlying dark matter halo masses. To date, halo masses have only been strongly constrained to redshifts z ∼ 6, and are limited to bright (M_UV ∼ − 20) galaxies at redshifts z = 4–6 (e.g., Harikane et al. 2018). JWST will potentially provide the first measurement of the clustering of high-redshift galaxies (Endsley et al. 2020), but due to its smaller FOV, it will likely suffer from cosmic variance, and will not have the same capability as Roman to study clustering in different environments.

We measure the 2PCF of the detectable galaxies in the 1 deg² survey, as outlined in Section 4.3. We included all galaxies brighter than M_UV < − 18. Figure 19 shows the measured 2PCF. The 2PCF for these faint galaxies is measured to very high accuracy (within 1%) out to redshifts z ≈ 7, with constraints out to z ≈ 10 (to within 10%).

Zoom In Zoom Out Reset image size

8.3.4. Galaxy–Halo Connection

In addition to constraining the number and spatial distribution of faint, high-redshift galaxies, Roman will also provide constraints on the galaxy–halo connection. Galaxy clustering can be used to infer dark matter halo mass, as the underlying dark matter will dictate the gravitational field (Mo & White 1996). Together with galaxy mass measurements, galaxy clustering gives a direct measure of the SHMR (which summarizes the connection between galaxy masses and their host dark matter halos).

As discussed above, halo masses have only been measured directly for redshifts z ≲ 6 (e.g., Harikane et al. 2018), while at high redshifts, SHMR constraints have been achieved with abundance matching between observed stellar masses and dark matter only simulations (e.g., Stefanon et al. 2021). There is a potential disagreement between these two techniques—the SHMR measured from abundance matching is typically three to four times higher than that measured from clustering (Stefanon et al. 2021). The difference between SHMRs derived from clustering and AM indicates that direct halo mass measurements are needed to understand the galaxy–halo connection at high redshifts. A Roman UDF will provide extraordinary data, which may elucidate the origin of this discrepancy.

We show the SHMR for our synthetic galaxy catalog in Figure 20, and recent constraints from Harikane et al. (2018) and Stefanon et al. (2021). The Stefanon et al. (2021) data were scaled down by a factor of 1.7 to convert from a Salpeter (1955) IMF to a Chabrier (2003) IMF. We use all of the 5σ Lyman-break selected galaxies in the catalog. Our synthetic catalog clustering measurements agree well with the Stefanon et al. (2021) data (within ∼1σ), and are slightly higher than the Harikane et al. (2018) data.

Zoom In Zoom Out Reset image size

In this work we do not attempt to predict how well a Roman UDF will be able to constrain the SHMR. As discussed in Section 8.3.1, stellar mass measurements for high-redshift galaxies likely require photometry at wavelengths longer than what Roman will provide, but they may be possible to estimate with scaling relations. Since our synthetic galaxy catalogs contain information regarding host dark matter halo properties, this work can form the basis of future analyses of the galaxy–halo connection in UDFs.

In addition to the galaxy catalog, we present synthetic images of a Roman UDF. These images are intended for developing analysis tools, and studying systematics (see Section 9.5). We include FITS files for each Roman filter, along with full-resolution versions of the RGB images presented in this section in our data release (see Appendix A for details).

We create the synthetic images using GalSim (Rowe et al. 2015). GalSim contains a module specifically for Roman observations (Troxel et al. 2021), which includes five of the Roman filters: Z087, Y106, J129, H158, and F184. We model each galaxy as a Sérsic profile, with an index n_s , axis ratio q, and PA, as described in Section 5. We truncate the distribution of Sérsic indices between 0.3 and 6.2, to avoid numerical inaccuracies at more extreme values. The scale radius of the Sérsic profile, r₀, is directly related to the assigned half-light radius, R_eff:

$$\begin{eqnarray}&&{r}_{0}=\displaystyle \frac{{R}_{\mathrm{eff}}}{{b}^{{n}_{s}}\sqrt{q}},\end{eqnarray} \tag{ 25 }$$

where b ≈ 2n_s − 1/3 (Moriondo et al. 1998). Figure 21 shows a composite of the full synthetic galaxy catalog using Z087 (blue), Y106 (green), and H158 (red) filters. We note that we do not include simulated stars in any of the released images.

Zoom In Zoom Out Reset image size

To simulate a Roman 30 m_AB survey, we convolve the image with the Roman point-spread function (PSF), and add noise to the image. We calculate the image noise assuming proposed exposure times from Koekemoer et al. (2019), to reach a 5σ limit of m_AB ≈ 30 (see Table 4). We add noise to the image in the same manner as Troxel et al. (2021). First, we generate the sky background accounting for the stray light and thermal emission from the telescope, which is added to the image. Given this, we add errors associated with Poisson noise, reciprocity failure, dark current, the calibration of the Roman detectors, inter-pixel capacitance, and instrument read noise.

We present an RGB visualization of the full Roman catalog, with the PSF and noise included in Figure 22, with one Roman footprint overlaid on top. Each Roman pointing will have 18 detectors. This visualization shows the rich amount of structure that will be contained in a Roman UDF. As described in Koekemoer et al. (2019), a 1 deg² UDF would consist of three Roman pointings. Though a UDF based on three tiled WFI pointings would not be perfectly square, and the exposure time would not be perfectly uniform, these synthetic images demonstrate the richness of a 1 deg² UDF and will be incredibly useful in designing and preparing for wide, deep surveys.

Zoom In Zoom Out Reset image size

We include FITS maps of the galaxy catalog for all five filters currently included in GalSim in our data release. The top row of Figure 23 shows the flux in one WFI detector (each detector has an area of ∼$7\buildrel{\,\prime}\over{.} 3\,\times \,7\buildrel{\,\prime}\over{.} 3$). The bottom row shows an example z = 9.4 galaxy (with the background subtracted in each panel). This galaxy would be selected as a Y106 dropout. The dropout galaxy has a size of R_eff = 0.38 kpc, and therefore is not resolved with the WFI (which has an angular resolution of $0\buildrel{\,\prime}\over{.} 11$ pixel⁻¹).

Zoom In Zoom Out Reset image size

The synthetic images we release are at the resolution of the WFI detector, $0\buildrel{\,\prime}\over{.} 11$ pixel⁻¹. An actual Roman UDF would likely be dithered, which would sample the sky on a subsampled pixel scale, which would improve the PSF sampling and angular resolution. Since this process would create correlations between adjacent pixels, we did not include this in our data release. Future work will build upon DREaM to examine dithering strategies of a Roman UDF.

The EoR is the final frontier for galaxy surveys. Given the difficulty in measuring galaxies at high redshifts, this period in the universe's history is remarkably unconstrained. High-redshift low-mass galaxies were likely the major source of the ionizing photons in the EoR (e.g., Finkelstein et al. 2019), and observations indicate that reionization was a "patchy" process (e.g., Furlanetto & Oh 2005; Villasenor et al. 2021). To fully understand the EoR, we need a complete census of galaxies and their ionizing photon contribution.

To determine whether there are enough high-redshift, faint galaxies to cause reionization, we must accurately model the CSFRD, which depends on the UVLF. Though there are some constraints on the UVLF to z ∼ 10 (e.g., Bouwens et al. 2021), there is still much improvement to be made on constraining both the bright steep end of the UVLF and the faint end (e.g., Bowler et al. 2014, 2015) at high redshifts. JWST will make improvements on this front. For example, Kauffmann et al. (2020) predicted that the $100\,{\mathrm{arcmin}}^{2}$ Cosmic Evolution Early Release Science Survey will constrain the faint end of the UVLF to a precision of 0.25 at z ≥ 8, but that a survey would need to be at least $300\,{\mathrm{arcmin}}^{2}$ to constrain the bright end up to z = 8. We predict a Roman UDF will capture enough high-redshift galaxies to constrain the UVLF to within 1% at redshifts beyond z ∼ 10. This will either confirm that there are enough faint galaxies to account for the ionizing photons needed to cause reionization, or indicate that another source, such as AGNs (e.g., Madau & Haardt 2015), must contribute.

In addition to the abundance of galaxies, the total ionizing photon budget depends on the Lyman-continuum (LyC) production efficiency, ξ_ion, and the LyC escape fraction, f_esc. To account for reionization, galaxies need higher escape fractions at high redshifts than have been observed at low redshifts (Davies et al. 2021). Since the neutral IGM absorbs LyC photons, f_esc is very difficult to constrain at the EoR. However, f_esc could possibly be measured during the EoR using indirect methods (e.g., Leethochawalit et al. 2016; Zackrisson et al. 2017; Chisholm et al. 2018, 2020). A Roman UDF will detect tens of thousands of galaxies during the epoch of reionization, including rare, bright sources that are ideal for spectroscopic follow-up to measure f_esc. Further, since Roman can map out the density around each galaxy, it will allow for the measurement of the escape fraction in different environments.

There exists a clear bimodality in SFRs, indicating two distinct populations (star-forming and quiescent). The decline in the CSFRD since redshift z ≈ 2 was likely caused by the quenching of galaxies (e.g., Renzini 2016). However, the mechanisms that transform galaxies from star-forming to quiescent are not fully understood. Possible mechanisms that have been proposed include feedback from stars and AGNs, or the removal of gas through tidal or ram pressure stripping. Quenching is likely caused by a combination of these mechanisms, and the dominant processes may be redshift dependent (Kalita et al. 2021).

Advances in the understanding of the origin of QGs is being greatly improved by surveys such as the spectroscopic survey Gemini Observations of Galaxies in Rich Early ENvironments (Balogh et al. 2021), which targets QGs in clusters around redshift z ∼ 1. Roman has the potential to identify QG populations out to much higher redshifts. In particular, the synthetic catalog presented in this work predicts that a UDF with Roman will contain ∼10⁵ detectable QGs, including QGs beyond redshift z = 7. Though the exact number of QGs relies heavily on the underlying assumed SMF, Roman will likely detect the highest redshift QG to date if they exist, allowing for the study of quenching mechanisms at high redshifts.

A Roman UDF will likely provide the first strong clustering measurements of the faint, high-redshift galaxies responsible for reionization, which will enable a measurement of the underlying dark matter mass. In our preliminary predictions, we have shown that the coverage from a UDF will constrain the 2PCF for faint galaxies out to redshifts z = 10 (at the 10% level). These measurements will likely provide the first direct measurement of halo masses at z = 10, placing constraints the galaxy–halo connection at high redshift. Since the galaxy–halo connection depends on both galaxy formation physics and the underlying cosmological model, this will provide important tests for current models of galaxy formation.

Roman's extensive FOV will allow for the contiguous, deep imaging of galaxies, reducing cosmic variance, and probe the environment around individual galaxies. The large number of galaxies detected by a Roman UDF will reduce Poisson noise, and increase the number of rare objects that will be detected. Since a single Roman pointing is wide enough to capture several reionization bubbles at the height of reionization, a Roman UDF will allow for the study of differences between galaxies in ionized and neutral regions. For instance, there may be possible variations in the faint-end slope of the UVLF with environment, which can be measured with a deep, wide, galaxy survey.

Before the launch of Roman, JWST will begin to address questions regarding the EoR, galaxy formation, and the galaxy–halo connection. In particular, the JWST Advanced Deep Extragalactic Survey (JADES) is a $236{\mathrm{arcmin}}^{2}$ planned imaging and spectroscopy survey. ¹⁹ Though JADES will provide revolutionary data, it will detect less than ∼10 × the number of objects at all redshifts compared to a Roman UDF. In addition, cosmic variance may dominate over Poisson noise for future high-z surveys (Trapp & Furlanetto 2020), and JADES will cover a much smaller area.

Though Roman will cover a huge area on the sky, it will be difficult to obtain photometry at wavelengths greater than 2 μm. At high redshifts, accurate measurements of stellar masses and QG identification require supplementary imaging from another observing facility. JWST will probe very far into the infrared, which will potentially allow for Balmer-break selection of galaxies, improving redshift measurements and stellar mass estimates (though JWST is not sensitive to light blueward of 1 μm). Therefore, JWST may provide valuable spectroscopic follow-up to rare detections from a Roman UDF. For example, Roman will detect many galaxies on the bright end of the UVLF. Bright galaxies (e.g., starburst galaxies, with M_UV < − 22) are ideal for spectroscopic follow-up with JWST, and can be used to observe nebular lines and estimate f_esc.

In addition to the simulated galaxy catalog, we provide synthetic images of a Roman UDF. These images can be used to determine the impact of source blending, line confusion, and potential problems with SED fitting (e.g., Borlaff et al. 2019; Kauffmann et al. 2020; Massara et al. 2021). Having very deep synthetic images will be useful for studying WFI systematics, processing issues (e.g., low surface brightness issues), and secondary analysis (e.g., photo z studies). Quantifying these issues will also be beneficial for many other Roman surveys, such as the High Latitude Survey. For instance, in Section 8.1, we assumed that all galaxies that are detectable will be selected. However, since a UDF will be so richly populated with structure, a fraction of the high-redshift galaxies will be obscured by forefront galaxies. We intend to use the synthetic images and catalogs in this work to quantify this effect.

In this paper we have only included photometric predictions, but the realistically modeled galaxy SEDs can also be used to generate grism predictions. In particular, the MOSFIRE Deep Evolution Field (MOSDEF) survey (Kriek et al. 2015) measured the detailed rest-frame optical emission-line SEDs in galaxies z = 2–3. More recent measurements with Keck/MOSFIRE have presented spectroscopic measurements out to z ∼ 8 (Topping et al. 2021). We can combine these observations with our synthetic catalog, to guide the Roman grism data reduction pipeline at high redshifts.

The synthetic catalog can also be used to study LAEs. LAEs produce a large amount of ionizing photons, and both UVLFs and clustering of LAEs are important to characterize EoR. The Lyα UVLF decreases toward the early stage of EoR, since neutral hydrogen absorbs Lyα photons (e.g., Ouchi et al. 2018). The clustering of LAEs emitters constrains the ionized fraction and topology of reionization.

As discussed in Section 6.2.2, UVJ diagrams differentiate between star-forming and QG populations. We display the U−V and V−J colors for our two galaxy populations in Figure 24, along with the UVJ selection box from Williams et al. (2009). The two populations occupy two distinct regions in this parameter space, demonstrating that the synthetic galaxy catalog does capture the bimodal population. Observational constraints for QGs do not currently exist for galaxies past redshifts z ≈ 4; the UVJ distribution of the synthetic galaxies at low redshifts does agree with observations (e.g., Schreiber et al.2015).

Zoom In Zoom Out Reset image size

We do not explicitly assign SFRs to the galaxies in the synthetic catalog, but we reproduce realistic trends in SFRs. By assigning realistic M_UV values to the SFGs, and constraining QGs to the appropriate place in the UVJ color diagrams, we accurately model a bimodal population, with SFGs having higher SFRs than QGs. Additionally we reproduce the observed cosmic star formation rate density (CSFRD) as demonstrated in Section 7.2.

Figure 25 shows the SFR–mass relation of the synthetic catalog compared to the relations from Schreiber et al. (2017). We closely match the Schreiber et al. (2017) SFG SFR–mass relations. The SFRs for QGs are lower than the Schreiber et al. (2017) relation. However, the SFRs–mass relations from Schreiber et al. (2015) assumed that all IR emission from QGs originated from residual star formation. Alternative explanations are that the IR emission originates from AGN torus emission, dust heating, or incorrect classification of SFGs. The QG SFRs still agree with what current observations can predict.

Zoom In Zoom Out Reset image size

Galaxy ages are a measure of SFHs. SFHs of individual galaxies give a direct measurement of the evolution of different galaxy populations. Trends in galaxy ages have been well established; for instance, QGs are older than SFGs, and low-redshift galaxies are older than high-redshift galaxies (e.g., Webb et al. 2020). However, galaxy ages for individual galaxies are difficult to measure due to degeneracies between other parameters (metallicity, dust), sensitivity to priors in fitting (e.g., Leja et al. 2019), and the similarity between SEDs in galaxies older than ∼5 Gyr (e.g., Gallazzi et al. 2005). Our synthetic galaxy catalogs provide age estimates, and realistic SEDs to examine this further.

We calculate the mass-weighted age of the synthetic galaxies as

$$\begin{eqnarray}&&{a}_{\mathrm{gal}\ }=\displaystyle \frac{{\int }_{0}^{{t}_{\mathrm{age}}}({t}_{\mathrm{age}}-t)\psi {dt}}{{\int }_{0}^{{t}_{\mathrm{age}}}\psi {dt}}\,.\end{eqnarray} \tag{ B1 }$$

For the delayed-tau model,

$$\begin{eqnarray}&&{a}_{\mathrm{gal}\ }={t}_{\mathrm{tot}}-2\tau +\displaystyle \frac{{t}_{\mathrm{tot}}^{2}}{\tau {e}^{{t}_{\mathrm{tot}}/\tau }-(\tau +{t}_{\mathrm{tot}})},\end{eqnarray} \tag{ B2 }$$

where t_tot = t_age − t_start is the total time of star formation.

Figure 26 shows the average galaxy ages as a function of redshift for both SFGs (blue) and QGs (red). QGs are older than SFGs at all redshifts, and age decreases with increasing redshift, as expected. At low redshifts, QGs are ≈ 3 Gyr older than SFGs.

Zoom In Zoom Out Reset image size

Galaxy metallicity can greatly affect quantities such as galaxy color, and therefore must be accurately modeled. We assigned metallicities from the FMR (Equation (19)), to the parent catalogs used in the SED pipeline. Figure 27 shows the FMR for the synthetic galaxies. We show for comparison the FMR (dotted lines), where we have used the assumption $12+{\mathrm{log}}_{10}({\rm{O}}/{\rm{H}})\approx {\mathrm{log}}_{10}({Z}_{\mathrm{met}}/{Z}_{\mathrm{sol}})$.

Zoom In Zoom Out Reset image size

Our synthetic galaxies follow the expected FMR, where metallicity increases with mass, with a turnover at high masses, as seen in observations (see discussion in Section 6.3). The low-mass galaxies have lower metallicities than the high-mass galaxies as expected, but do not vary much with SFR. However, at these low masses, very few galaxies will have SFRs ψ > −1.