Distinguishing Mergers and Disks in High Redshift Observations of Galaxy Kinematics

TL;DR

This study uses simulations to assess how often galaxy mergers at high redshift are misclassified as disks in observations, finding that the misclassification rate is small and does not significantly bias disk fraction estimates.

Contribution

It provides a systematic analysis of the likelihood that merging galaxies are mistaken for disks in high-redshift observations, quantifying the impact on measured disk fractions.

Findings

Merger galaxies can pass as disks in observations with probabilities from 0 to 100%.

The misclassification rate slightly overestimates true disk fractions by up to 5-15%.

The effect of mergers on disk identification is weak and does not significantly bias results.

Abstract

The majority of massive star-forming galaxies at $z\sim2$ have velocity gradients suggestive of rotation, in addition to large amounts of disordered motions. In this paper, we demonstrate that it is challenging to distinguish the regular rotation of a disk galaxy from the orbital motions of merging galaxies with seeing-limited data. However, the merger fractions at $z\sim2$ are likely too low for this to have a large effect on measurements of disk fractions. To determine how often mergers pass for disks, we look to galaxy formation simulations. We analyze $\sim$24000 synthetic images and kinematic maps of 31 high-resolution simulations of isolated galaxies and mergers at $z\sim2$. We determine if the synthetic observations pass criteria commonly used to identify disk galaxies, and whether the results are consistent with their intrinsic dynamical states. Galaxies that are intrinsically mergers pass the disk criteria for anywhere from 0 to 100$\%$ of sightlines. The exact percentage depends strongly on the specific disk criteria adopted, and weakly on the separation of the merging galaxies. Therefore, one cannot tell with certainty whether observations of an individual galaxy indicate a merger or a disk. To estimate the fraction of mergers passing as disks in current kinematics samples, we combine the probability that a merger will pass as a disk with theoretical merger fractions from a cosmological simulation. Taking the latter at face-value, the observed disk fractions are overestimated by small amounts: at most by $5\%$ at high stellar mass ($10^{10-11}$ M$_{\odot}$) and $15\%$ at low stellar mass ($10^{9-10}$ M$_{\odot}$).

Equations

$$V_{rot}=\frac{V_{max}-V_{min}}{2\times \sin i}$$

$$f_{obs}=(1-f_{int})\times f_{pair}\times X_{contamination}+f_{int}$$

$$F_{late}=F_{pair}\times \left(\frac{{\tau }_{late}}{{\tau }_{late}+{\tau }_{pair}}\right)$$

$$\vec{r}$$

$$=\rho \sin \varphi \cos \theta \hat{x}+\rho \sin \varphi \sin \theta \hat{y}+\rho \cos \varphi \hat{z}$$

$$R$$

$$=\sqrt{(\rho \sin \varphi \cos \theta {)}^2+(\rho \sin \varphi \sin \theta {)}^2}$$

$$=\rho \sin \varphi$$

$$<R>$$

$$=\frac{{\int }_0^{2\pi }{\int }_0^{\pi }(\rho \sin \varphi )(\sin \varphi d\varphi d\theta )}{{\int }_0^{2\pi }{\int }_0^{\pi }\sin \varphi d\varphi d\theta }$$

$$=\rho \frac{{\int }_0^{\pi }{\sin }^2\varphi d\varphi }{{\int }_0^{\pi }\sin \varphi d\varphi }$$

$$=\rho \left(\frac{\pi }{4}\right)$$

Full text

1. Introduction

The cosmic star-formation rate density peaks between $1\,<\,z\,<\,3$ and marks a critical period for galaxy assembly (Madau & Dickinson, 2014). The processes thought to be most relevant for regulating the mass growth and structural transformation of galaxies at this time, e.g., stellar and AGN feedback, violent disk instabilities, cold mode accretion, minor/major mergers (see e.g., review by Somerville & Davé 2015, and references therein), may also be dynamically disruptive, destroying ordered disk rotation on timescales comparable to the dynamical time of the galaxy.

The observed fraction of high redshift galaxies that are disk-like thus provides important insight into the frequency and relevance of these processes. If the fraction is high, it implies that early galaxy assembly is governed by dynamically calm processes, ones that are hospitable for disk formation and survival. If the fraction is low, it implies that early galaxy assembly is dominated by dynamically disruptive processes, those that tend to stall the formation and settling of well-ordered disks.

The internal kinematics of galaxies offer the most relevant observable for determining if a galaxy is a disk. A disk galaxy will at least have more ordered motions than random disordered motions, as quantified by the rotation velocity $V_{rot}$ and the gas velocity dispersion $\sigma_{g}$, respectively.

Emission-line kinematics of the ionized gas in galaxies have now been measured for several hundred star-forming galaxies at $1\lesssim z\lesssim 3$ (e.g., Erb et al. 2003; Genzel et al. 2006; Förster Schreiber et al. 2006; Erb et al. 2006; Law et al. 2007; Kassin et al. 2007; Wright et al. 2007; Genzel et al. 2008; Yang et al. 2008; Puech et al. 2008; Shapiro et al. 2008; Law et al. 2009; Förster Schreiber et al. 2009; Lemoine-Busserolle & Lamareille 2010; Jones et al. 2010; Gnerucci et al. 2011; Swinbank et al. 2012; Law et al. 2012; Kassin et al. 2012; Contini et al. 2012; Newman et al. 2013; Glazebrook 2013; Wisnioski et al. 2015; Mieda et al. 2016; Price et al. 2016; Olivares et al. 2016; Simons et al. 2016; Contini et al. 2016; Stott et al. 2016; Rodrigues et al. 2017; Mason et al. 2017; Straatman et al. 2017; Turner et al. 2017; Übler et al. 2017; Straatman et al. 2018) — with samples that are large and representative enough to derive meaningful population statistics.

At $z\,=\,2$, approximately 50% of low mass galaxies (M${}_{*}=10^{9-10}$ M*⊙) and 70% of high mass galaxies (M${}_{*}=10^{10-11}$ M⊙) have $V_{rot}/\sigma_{g}\,>\,1$, i.e., more rotational support than dispersion support (Simons et al. 2017). Integral field spectroscopy (IFS) allows one to check for additional disk signatures in two-dimensional velocity and velocity dispersion maps. For the KMOS3D* IFS survey, Wisnioski et al. (2015) define a series of increasingly strict criteria that are used to identify disk-like systems. These disk criteria include: a continuous single gradient in the velocity map, a rotation velocity that exceeds the velocity dispersion, a maximum in the velocity dispersion map that is coincident with the dynamical center, an alignment between the photometric and kinematic major axes, and a spatial coincidence of the continuum and dynamical centers. Using these, they conclude that $\sim$50 - 70% of massive galaxies at $z\sim 2$ are disk-like, with the exact fraction depending on the number of criteria used (Wisnioski et al., 2015).

However, these observed kinematic signatures of disks may not be unique to disks. The main challenge comes from mergers. The orbital motions of merging galaxies, once convolved with the typical observational seeing-limit at $z=2$ ($\sim 0.\arcsec 6$ or 5 kpc), can mimic the regular ordered rotation of a disk. While techniques have been developed for distinguishing mergers and disks in kinematics observations (e.g., Kinemetry, Krajnović et al. 2006), they are typically restricted to a small subset of observations at $z\sim 2$ with sufficient resolution and signal-to-noise (e.g., Shapiro et al. 2008). Using a set of artificially-redshifted local mergers and synthetic observations of idealized binary merger simulations, Hung et al. (2015, 2016) demonstrate that the ability to distinguish between disks and mergers in kinematic data depends strongly on the interaction phase, becoming progressively more difficult during later stages.

Merger fractions are found to increase from $z=0$ to $z=2$ in both observations and simulations (Lotz et al., 2011; Rodriguez-Gomez et al., 2015; Mantha et al., 2018). This means that the confusion of mergers as disks is likely more important to the observed disk fractions at higher redshifts. Rodrigues et al. (2017) re-analyzed the Wisnioski et al. (2015) sample at $z\sim 1$ and conclude that 58$\%$ are involved in a merger (from first approach to post-collision) and that only one-third are isolated and virialized disks. It is clear that mergers present a potentially important problem for determining disk fractions in seeing-limited kinematics observations at high redshift, but how significant this problem is remains an open question.

In this paper, we quantify the significance of merger contamination in current high redshift seeing-limited kinematic surveys using state of the art galaxy formation simulations. First, we create and analyze synthetic observations of a suite of high-resolution zoom-in hydrodynamic simulations (Ceverino et al., 2014) to determine the probability that a merger in a given stage will pass as a disk in the observations. These simulations are run in a cosmological context and include realistic cosmic accretion and galaxy-galaxy mergers. We then use theoretically-derived merger fractions from the large Illustris cosmological simulation (Genel et al., 2014; Vogelsberger et al., 2014) to determine if there are enough galaxies in mergers (at all stages) to significantly affect global disk fractions at $z\sim 2$.

The paper is structured as follows. In §2 we describe the simulation suite. In §3 we characterize the intrinsic dynamical properties of the simulated galaxies and identify mergers and close pairs. In §4 we use the SUNRISE dust-radiative transfer program to generate synthetic observations of the simulations. In §5, we measure photometric properties from the synthetic imaging and kinematic properties from the synthetic spectra. In §6 we compare the interpretation from the synthetic data with the intrinsic dynamical state of each galaxy and determine the probability that a merger will pass as a disk as a function of projected separation. In §7 we use merger fractions from the Illustris simulation to determine if merger contamination plays a significant role in current global disk fraction statistics at $z\,=\,2$. In §8 we summarize our conclusions. We adopt a $\Lambda$CDM cosmology defined with ($h$, $\Omega_{m}$, $\Omega_{\Lambda}$) = (0.7, 0.27, 0.73).

2. VELA Simulations

We analyze a set of 31 cosmological zoom-in galaxy formation simulations from the VELA simulation suite (Ceverino et al., 2014; Zolotov et al., 2015). The simulations are briefly reviewed here and we refer to Ceverino et al. (2014) for further details.

The simulations were run with the N-body and Eulerian gas dynamics Adaptive Refinement Tree code (ART, Kravtsov et al. 1997), including subgrid recipes for: gas and metal cooling, UV-background photoionization, stochastic star formation, gas recycling and metal enrichment from stellar winds, thermal feedback from supernovae and feedback from young stars through radiation pressure and radiative heating (Ceverino et al., 2010, 2012, 2014). The adaptive mesh grid has a maximum resolution between 17 and 35 physical parsecs and this is typically reached for gas densities between $\sim 10^{-2}$ – $100$ cm*-3*. The mass of the star particles ranges from $10^{3}$ to $10^{6}$ M*⊙, depending on the gas mass of the parent cell, with a typical value of $10^{4}$ M⊙*.

The simulated galaxies have stellar masses spanning $9.3<\log$ M*∗*/M${}_{\odot}<10.7$ at $z=2$ — comparable to the typical stellar mass range of seeing-limited kinematics samples in the literature at this redshift (e.g., Wisnioski et al. 2015; Simons et al. 2016). The VELA simulations have been used in the past to study the assembly and transformation of high redshift galaxies, in particular the origin and evolution of stellar elongation (Ceverino et al., 2015b; Tomassetti et al., 2016), clumps (Moody et al., 2014; Inoue et al., 2016; Mandelker et al., 2017), morphologies (Snyder et al., 2015), inflows and outflows (Ceverino et al., 2016a, b), galaxy-halo spin correlations (Jiang et al., 2018), and compact galaxy structure (Zolotov et al., 2015; Tacchella et al., 2016a, b; Huertas-Company et al., 2018).

The subset of independent simulations analyzed in this paper span $1\,<\,z\,<\,3$ and are separated in equal intervals of cosmological scale factor $\Delta a$ = 0.01. This corresponds to time intervals of 110 Myr at $z=3$, 125 Myr at $z=2.0$ and 160 Myr at $z=1.0$ and totals over 1000 simulation snapshots.

3. Intrinsic States of the Simulated Galaxies

The objective of this paper is to compare the intrinsic properties of the simulations with the interpretation one would reach from their synthetic observations. The intrinsic properties of the galaxies are measured directly from the simulation data. For each simulation snapshot (i.e, time step), we use the three-dimensional positions and velocities of the star particles to (i.) determine if the most massive galaxy in the simulation box, hereafter referred to as the central galaxy, has a rotationally-supported star-forming disk and (ii.) identify all galaxies that are nearby the central galaxy, hereafter referred to as companion galaxies. As follows, the central galaxy in each simulation snapshot is categorized into one of the following: isolated disk galaxy, disk galaxy with a companion, isolated non-disk galaxy, or a non-disk galaxy with a companion.

3.1. Using Intrinsic Properties of the Simulated Galaxies to Identify Disks

In this section, we determine whether the central galaxy of each simulation snapshot has an intrinsic disk. Galaxy kinematics at high redshift are almost exclusively measured from emission lines, tracing the ionized gas in the star-forming regions of galaxies. For an appropriate comparison to the observations, we use the young star particles ($<20$ Myr) to search for disks in the simulated galaxies.

In a disk galaxy, stars will tend to have orbital motions aligned with the net galaxy rotation, while in a dispersion-supported galaxy they will have no preferred orbital direction.

To assess the alignment of each star particle’s orbit with the net galaxy rotation, we adopt the widely-used circularity parameter $j_{z}/j_{circ}$ (e.g., Abadi et al. 2003; Okamoto et al. 2005; Governato et al. 2007; Scannapieco et al. 2008, 2009; Sales et al. 2012; Ceverino et al. 2015a). The quantity $j_{z}$ is the component of a particle’s specific angular momentum aligned with the net angular momentum of the galaxy. The quantity $j_{circ}$ is the specific angular momentum of a theoretical circular orbit at the particle’s position, irrespective of the orbital direction, and is calculated from the total mass inside of a shell at that position, $j_{circ}(r)\,=\,\,r\,v_{circ}(r)\,=\,\sqrt{r\,M(<r)G}$, where spherical symmetry is assumed for simplicity.

The direction of the net angular momentum of the central galaxy ($\hat{z}=\vec{J}_{gal}/J_{gal}$) is defined using an iterative fit of a cylindrical disk to both the cold gas ($T<1.5\times 10^{4}$ K) and the stars with ages less than 100 Myr (see §3.2 of Mandelker et al. 2017 for details).

The circularity parameter is calculated for every star particle in the simulation box. Those particles associated with a disk lie near $j_{z}/j_{circ}\,=\,1$, i.e., on circular orbits aligned with the net rotation of the galaxy. In a dispersion-supported system, star particles are distributed symmetrically around $j_{z}/j_{circ}\,=\,0$. Particles with $j_{z}/j_{circ}\,=\,-1$ are traveling retrograde to the net galaxy rotation.

The joint distribution of stars in $j_{z}/j_{circ}$ versus distance from the central galaxy is shown for three example snapshots in Figure 1. External galaxies appear as distinct clusters of points in this plane, i.e., their stars are close to one another in space and are traveling on similar trajectories.

For each simulation snapshot, we create a segmentation map of the joint distribution of $j_{z}/j_{circ}$ versus distance. This is used to detect significant clusters of star particles. The segmentation procedure iterates until the number of detected objects does not change with an increasing clipping level. To prevent spurious detections, a threshold is set such that a companion galaxy must have a minimum total stellar mass of 106 M*⊙*. Once this procedure is complete, each star particle in the simulation is associated with a unique galaxy. We record the mass and mass-weighted position of all galaxies detected within 100 kpc of the central galaxy.

For each timestep of each simulation, we calculate the mean circularity of the young star particles that are associated with the central galaxy.

Each central galaxy of each simulation snapshot is determined to host a disk if the average circularity of its young stars $\left<j_{z}/j_{circ}\right>_{young}$ is greater than $1/\sqrt{3}$, i.e., the young stars are strongly rotationally-supported, and to not host a disk if $\left<j_{z}/j_{circ}\right>_{young}$ is less than 1/3, i.e., the young stars are strongly dispersion-supported (Figure 2). This classification applies only to their young stellar component — the older stars are typically more dispersion-supported (Figure 1).

These thresholds in $\left<j_{z}/j_{circ}\right>_{young}$ are motivated by the following argument. The Jeans equation relates the circular velocity, the rotation velocity and the 1D velocity dispersion as $V_{circ}^{2}=V_{rot}^{2}+\alpha\,\sigma_{r}^{2}$, where $\alpha=2-3$ depending on the potential. For an isothermal sphere, $\alpha=2$ and $j_{z}/j_{circ}>1/\sqrt{3}$ and $<1/3$ are equivalent to $V_{rot}/\sigma_{r}\,>1$ and $<$ 0.5, respectively. The former is the standard observational threshold for classifying a galaxy as a disk — and is one of the observational thresholds used later in this paper — and the latter is sufficiently low for the galaxy to be considered dispersion-dominated. For simplicity, we ignore galaxies with intermediate $\left<j_{z}/j_{circ}\right>_{young}$ to these two cases.

Over $1\,<\,z\,<\,3$, the majority of the central galaxies in VELA have $\left<j_{z}/j_{circ}\right>_{young}\,>\,1/\sqrt{3}$ — i.e., they have a young stellar disk according to our classification (Figure 2). We stress that the goals of this paper do not rely on the true disk fraction of the VELA simulations. Our focus is to test how well one can recover the intrinsic properties of a galaxy from its synthetic observations — an exercise that is impartial to the rarity of those intrinsic properties in the simulation set.

3.2. Using Intrinsic Properties of the Simulated Galaxies to Identify Companion Galaxies

As described above, we use the joint distribution of stars in radial position and $j_{z}/j_{circ}$ to identify sufficiently massive ($>10^{6}$ M*⊙*) external galaxies in the simulation box.

Companion galaxies are defined to be those nearby the central galaxy (3D separation $\,<35$ kpc) with sufficiently high relative stellar mass (stellar mass ratio with the central less than 30). As such, we include major (1:1 - 1:4 mass ratios), minor (1:4 - 1:10), and very minor (1:10-1:30) companions. The latter two are referred to as simply ‘minor’ for the rest of this paper. If more than one companion galaxy is detected in a snapshot, we only consider the most massive one for this paper.

A central galaxy is considered isolated if it meets the following two criteria: there are no galaxies identified within 75 kpc in the current snapshot and it did not have a companion galaxy in the previous snapshot. The second criteria ensures that we exclude galaxies that experienced a merger within the past 150 Myrs.

4. Synthetic Observations of the Simulations: images and spectra

To generate synthetic observations, each simulation snapshot is post-processed with the dust radiative transfer software SUNRISE (Jonsson, 2006; Jonsson et al., 2010; Jonsson & Primack, 2010). Synthetic images and spectra are created from the SUNRISE output.

The synthetic HST images used in this paper are available as a high-level science product from MAST111https://archive.stsci.edu/prepds/vela: https://doi.org/10.17909/t9-ge0b-jm58 (catalog https://doi.org/10.17909/t9-ge0b-jm58).

4.1. Dust Radiative Transfer with SUNRISE

SUNRISE contains two primary steps. It first generates a spectral energy distribution (SED) for each source of radiation in the simulation. It then propagates the associated polychromatic rays through the 3D gas grid of the simulation, taking into account dust absorption and scattering, towards a user-defined camera(s).

For each snapshot, the ART output is first loaded into the YT analyses software (Turk et al., 2011) and converted into the appropriate SUNRISE format. This is carried out using an updated version of the pipeline used in Moody et al. (2014).

An SED is assigned to each star particle based on its age, mass, and metallicity using STARBURST99 (Leitherer et al., 1999) assuming a Kroupa (2001) initial mass function. Emission lines are generated from young star particles, less than 10 Myrs old, using a starburst model from the photoionization and shock modeling code MAPPINGS III (Dopita et al., 2005; Groves et al., 2008). This model includes a central massive young stellar cluster, whose SED is set with STARBURST99, a surrounding HII region, and a photodissociation region with a covering fraction of 0.2 (Groves et al., 2008). To account for $10^{4}$ K thermal broadening, we convolve the original MAPPINGS III H$\alpha$ emission line with a 10 km s*-1* 1D Gaussian kernel.

The dust density is assumed proportional to the gas metal mass using a dust-to-metals mass ratio of 0.4 (Dwek 98). Consequently, the dust geometry is resolved to the typical gas cell size of 17-35 pc. We adopt the dust grain size distribution from Weingartner & Draine (2001) and Draine & Li (2007) and set the slope of the dust extinction law, $R_{v}$, to the Milky Way value of $3.1$ (e.g., Gordon et al. 2003).

SUNRISE cameras are placed at 19 positions around the galaxy. The camera parameters are user-selected and not default to SUNRISE. Of the cameras, 6 are random but remain fixed in comoving space over time, 8 are random but change from timestep to timestep, and the remaining 5 are fixed with respect to the net angular momentum of the galaxy and, as such, change from timestep to timestep. Those fixed to the galaxy include one which is face-on, i.e., along the gas angular momentum axis, one which is reverse face-on, two which are edge-on, and one which is placed at a 45*∘* angle from face-on. Multiple sightlines allow us to overcome the intrinsic randomness associated with a single viewing angle. To avoid highly-uncertain inclination-corrections, we do not use the 2 face-on cameras or any of the random cameras that are oriented within $20^{\circ}$ of the face-on direction for a given snapshot.

We run the radiative transfer step of SUNRISE twice, once to create broadband images and once to create high-resolution spectra. A sufficient number of polychromatic rays are used for both runs, 107, to overcome Monte Carlo noise.

In the first mode, SUNRISE is run at low spectral resolution but with coverage over the NUV–NIR (0.02 – 5 $\micron$). Each camera contains 800 x 800 pixels and covers a physical FOV of either 50 or 100 kpc. The larger FOV is used for 6 of the largest galaxies in the simulation suite to ensure sufficient coverage in their outskirts at low redshift. The output spectral cube is integrated over the spectral response function of the far-UV through IR filters available with HST/ACS–WFC3, taking into account the redshift of the simulation and the corresponding cosmological Doppler shift and surface brightness dimming of the SED.

In the second mode, SUNRISE is run at high spectral resolution around the H$\alpha$ emission line (0.65 – 0.66 $\micron$). The camera contains 400 x 400 pixels and covers the same physical FOV (50 or 100 kpc) as in the first mode. The pixel-by-pixel surface brightness in the resulting SUNRISE cube is scaled with redshift to account for cosmological surface brightness dimming.

This process is repeated for each snapshot of each galaxy simulation. The final suite contains thousands of idealized UV–IR band images and idealized spectral cubes around the H$\alpha$ emission line.

4.2. Synthetic HST/WFC3 imaging

Following Snyder et al. (2015), the idealized images are degraded in noise and spatial resolution to simulate the typical imaging quality of modern deep-field galaxy surveys with Hubble/WFC3 (Grogin et al., 2011; Koekemoer et al., 2011). The angular sizes and fluxes of the pixels are scaled with the angular size distance and luminosity distance of the simulation redshift. The images are convolved with the typical point-spread function of the instrument/filter (Krist et al., 2011) and the pixel size is rebinned to the mosaic pixel scale of the CANDELS survey (0.″06 pixel*-1* for HST/WFC3-IR; Koekemoer et al. 2011). Finally, shot noise is added to reach a specific surface brightness limit. We generate images for two surface brightness limits — 25 mag arcsec*-2* and 27 mag arcsec*-2*.

In this paper, we use the synthetic HST/WFC3 F160W image (hereafter referred to as H-band) at the deeper surface brightness limit of 27 mag arcsec*-2*. The observed H-band traces the rest $\sim$V-band at z = 1 and the $\sim$I-band at z = 3.

4.3. Synthetic VLT/KMOS Integral Field Spectroscopy

The idealized high spectral resolution cube from SUNRISE is degraded to the surface brightness depth, spectral and spatial resolution of typical ground-based seeing-limited data. The majority of the large IFS kinematics surveys at $1\,\lesssim\,z\,\lesssim\,3$ are performed with the KMOS near-infrared spectrograph on the VLT, e.g., KMOS-3D (Wisnioski et al., 2015), KROSS (Stott et al., 2016), KDS (Turner et al., 2017). For an apt comparison with the literature, we use the specifications of this instrument to generate the synthetic spectra. We emphasize, however, that our conclusions are general to all seeing-limited instruments/surveys.

The physical pixel size of the SUNRISE cube is 0.25 kpc px*-1* and the corresponding angular scale is determined from the simulation redshift. Each spectral slice of the cube is blurred with a 2D Gaussian of 0.″6 FWHM to simulate the typical atmospheric seeing in the NIR. The pixels are rebinned to 0.″2, the angular pixel scale of KMOS, and the spectral dimension is convolved with a 1D Gaussian of $R=\Delta\lambda/\lambda=3800$ to simulate the KMOS spectral resolution. Finally, shot noise is added to each pixel to reach a typical 5$\sigma$ point source depth of 24 AB magnitudes, which is equivalent to a typical 8 hour exposure with KMOS.

Although we aim to simulate realistic ground-based observations, we do not include sky lines and ignore that the sky transmission is significantly lower over the observed wavelengths of H$\alpha$ at $1.8\,<z\,<2.2$. These assumptions do not affect the conclusions of this paper.

5. Measurements from the synthetic observations

5.1. Inclination, Photometric Position Angles, and Continuum Centroid

Photometric properties — inclinations and position angles — are measured from the synthetic HST H-band images. Unique sources in each image are identified and characterized using the Photutils package of Astropy (Astropy Collaboration et al., 2013), with the requirement that at least 5 connected pixels have a flux that exceeds the background rms.

Each image is centered on the central galaxy and so we select the central source as the main target. The ratio of the isophotal semi-major and semi-minor axis, $a$ and $b$, respectively, is used to determine the inclination, as $\cos^{2}i=((b/a)^{2}-q_{0}^{2})/(1-q_{0}^{2})$. We adopt a value of $q_{0}=0.25$, the standard literature assumption for a thick disk. The photometric position angle is defined as the direction of the semi-major axis. We determine the continuum centroid by calculating the light-weighted center of the central source in the synthetic H-band image.

5.2. Kinematics

The line-of-sight velocity and velocity dispersion are measured for each spaxel of the synthetic KMOS cube. To do this, we fit the H$\alpha$ emission line in the spaxel with a 1D Gaussian profile. The width, amplitude and center are left as free parameters in the fit. The center of the best-fit profile is taken as the mean line-of-sight velocity of that spaxel. The velocity dispersion of the spaxel is taken from the RMS width of the best-fit after subtracting the spectral resolution of the instrument in quadrature ($\sigma(R)=\sqrt{\sigma^{2}_{measured}-\sigma^{2}_{instrument}}$). This is repeated for all spaxels. The result of this procedure is a 2D map of the velocity and velocity dispersion. This is repeated for all of the synthetic KMOS cubes. Example kinematic maps are shown in Figure 3.

The kinematic major axis is defined as the line intersecting the the centroid of the pixels with the maximal and minimal 10% velocities, $V_{max}$ and $V_{min}$, respectively. A 3-spaxel wide slit is placed along the kinematic major axis. The rotation curve and velocity dispersion profile is measured along this slit by averaging the velocity and velocity dispersion of the 3 spaxels perpendicular to the slit axis.

The rotation velocity uncorrected for inclination, $V_{rot}\times\sin i$, is taken as the difference between the maximum and minimum velocities of the rotation curve. This quantity is then inclination-corrected using the inclination derived from the synthetic HST/WFC3 H-band image:

[TABLE]

The seeing will tend to smear unresolved velocity gradients and artificially elevate the velocity dispersion in the centers of galaxies, where the rotation curve is steepest (see Fig. 9 of Weiner et al. 2006 and the Appendix of Simons et al. 2015). The severity of beam smearing depends on several factors, including the intrinsic shape of the rotation curve and the ratio of the seeing to the galaxy size (Burkert et al., 2016). To minimize the effects of beam smearing, we adopt the technique used in Wisnioski et al. (2015). We assume that the velocity dispersion is constant across the face of the galaxy and measure its intrinsic value using the average uncertainty-weighted mean of the velocity dispersion profile on the galaxy outskirts, i.e., those at and beyond the maximal and minimal 10% velocities (typically near $\sim$1.5 $\times$ effective radius).

The dynamical center is defined along the kinematic major axis as the location where $v_{sys}=0$, i.e., the midpoint of the maximal and minimal velocities. In practice, the position is measured by interpolating between the locations of the 4 pixels that are contiguous to the pixel that is nearest to $v_{sys}=0$ (2 on either side). The location of the peak velocity dispersion is measured from the 2D velocity dispersion map. The uncertainties on both of these quantities are calculated using a Monte Carlo technique: the values of the velocity and velocity dispersion maps are resampled by their uncertainties 1000 times and, for each iteration, the dynamical center and the location of the velocity dispersion peak are re-measured. The width of the resulting distributions of measurements is taken as the uncertainty.

5.3. Criteria for classifying a galaxy as a disk

To determine whether the synthetic observations are consistent with a disk-like system, we adopt the disk criteria outlined in Wisnioski et al. (2015). The first three criteria are:

The velocity map exhibits a continuous velocity gradient along a single axis. 2. 2.

The inclination-corrected rotation velocity exceeds the ionized gas velocity dispersion, i.e., $V_{rot}$/$\sigma_{{\mathrm{H}\alpha}}\,>\,1$. 3. 3.

The position of the steepest velocity gradient is coincident with the peak in the 2D velocity dispersion map, within the uncertainties.

Aside from the inclination-correction for $V_{rot}$, these first three criteria rely only on information from the spatially-resolved kinematic maps. Wisnioski et al. (2015) report that 68% of star-forming galaxies at $z\sim 2$ satisfy these criteria. The last two criteria include:

The photometric and kinematic axes are aligned within 30 degrees. 2. 5.

The centroid of the continuum center is coincident, within the uncertainties, with the position of the steepest velocity gradient.

We will also explore one additional criterion:

The HST H-band image contains only a single nucleus.

6. Confusing Mergers for Disks in Resolved Kinematic Maps

6.1. Synthetic Observations of a Single Simulated Merger

We first examine the synthetic observations of a single simulation, namely VELA 20. The central galaxy in the simulation box merges with a slightly less massive companion galaxy at z = 1.6. Prior to their encounter, the stellar mass of the central and companion are 4.1 $\times 10^{10}$ M*⊙* and 1.2 $\times 10^{10}$ M*⊙*, respectively, qualifying this as a major merger. We follow the merger over 4 timesteps, each separated by $150$ Myr. The angular momentum profile, synthetic HST/WFC3 H-band image, synthetic VLT/KMOS H$\alpha$ flux map, and synthetic KMOS H$\alpha$ velocity and velocity dispersion maps are shown for each timestep in Figure 3.

Over the first two snapshots, $1.6<z<1.7$, the central galaxy has a rotationally-supported star-forming disk (Figure 3; top left) — the young star particles are distributed near $j_{z}/j_{circ}=1$ and extend out to $15$ kpc from the galaxy center. As expected, the H$\alpha$ velocity and velocity dispersion maps reflect this ordered rotation (Figure 3; bottom left) — the primary velocity gradient is smooth and continuous, the kinematic and photometric axes are aligned, the velocity dispersion map peaks near the center of the steepest velocity gradient, and the inclination-corrected rotation velocity (345 $\pm$ 11 km s*-1*) exceeds the typical velocity dispersion (61 $\pm$ 7 km s*-1*). This galaxy satisfies the observational disk criteria during these time steps, illustrating their success when applied to disk galaxies with no nearby companions.

The observations are significantly more challenging to interpret during the first interaction at $z=1.56$. Before and during the merger, the central galaxy retains its rotationally-supported disk and the orbital direction of the merger is aligned with the rotational direction of the disk (top panel). At the time of collision, the maximum velocity of the velocity map, i.e., the most redshifted component, traces the rotational signature of one side of the disk, while the minimum velocity of the velocity map, i.e., the most blueshifted component, traces the orbital motion of the companion galaxy. When taken together, these two define the kinematic major axis, i.e., the direction between the maximum and minimum velocities. Instead of tracing the rotational motion of the disk — which is still present in the map — the kinematic major axis follows the velocity gradient defined by one edge of the disk and the companion galaxy. This measurement has significantly less physical meaning than either the rotation velocity of the disk or the relative velocity between the two galaxies.

At face-value, the velocity map is consistent with rotational motions — the velocity gradient between the maximum and minimum velocity is smooth and continuous and velocity dispersion peaks along the steepest observed gradient.

The line-of-sight velocity difference between the edge of the disk and the merger is 546 km s*-1*, while the line-of-sight velocity difference of the disk is 358 km s*-1*. As such, a measurement of $V_{rot}\times\sin i$ along the kinematic major axis overestimates the true rotation velocity of the disk by a factor of 1.5 or 0.18 dex.

Figure 4 shows the distribution of the observed inclination-corrected rotation velocities, as measured from the kinematic major axis, for 16 sightlines around VELA 20 at $z=1.56$. For the majority of sightlines, 13/16, the derived rotation velocity along the kinematic major axis is higher than both the intrinsic rotation velocity of the central disk and the maximum circular velocity — the latter being unphysical.

In this example, the galaxy is sufficiently large, $20$ kpc or $2.^{\prime\prime}5$, compared to the size of the ground-based seeing, $0.^{\prime\prime}6$, that the complicated kinematic signatures of the merger are apparent in the velocity maps, specifically that there are two distinct velocity gradients in the map. As follows, we examine cases where the galaxy is less resolved in the seeing-limit.

6.2. Mergers Passing as Disks

We now examine the full suite of simulations and consider only two cases: those snapshots in which the central galaxy has no disk and a nearby companion and those snapshots in which the central galaxy has a disk and is isolated. For the purposes of analyzing clean and distinct samples, other scenarios (i.e., those snapshots with a disk and a companion or with a non-disk that is isolated) are ignored.

We show synthetic data for example isolated disk galaxies from the VELA simulations in the left panel of Figure 5. The central galaxies in all of these examples exhibit regular morphological and kinematic properties. More specifically, they satisfy the observational disk criteria outlined in §5.3 and in Wisnioski et al. (2015).

We consider all the isolated galaxy snapshots with a rotationally-supported young stellar disk in the VELA simulations, and find that the observational disk criteria work with high confidence. In Figure 6, we show the distribution of the fraction of sightlines in which an isolated intrinsic disk is classified as a disk. This distribution is highly skewed towards a value of 1, indicating that, in a majority of sightlines, isolated disks are correctly classified as disks in the synthetic observations.

There is, however, a shallow tail in the distribution, indicating that a non-negligible fraction of isolated intrinsic disks are misclassified as non-disks. This is most often due to the misalignment of photometric and kinematic axes and/or the offset between the continuum center and the dynamical center (i.e., failure of criteria 4 and 5). In these cases, the young star-forming regions are structurally dissimilar to the bulk of the stars — the former determining the kinematic morphology and the latter determining the continuum morphology.

We present synthetic data for examples in which galaxies are dispersion-supported and have a nearby companion in the right panel of Figure 5. Although these galaxies do not intrinsically have disks, the synthetic data exhibit many of the same characteristics that we expect for disks. In all of the examples presented, two of the first three criteria outlined in §5.3 are satisfied — there is a smooth and continuous velocity gradient and the observed inclination-corrected “rotation-velocity” (in actuality, the orbital motion of the merger) exceeds the velocity dispersion (as measured at the galaxy outskirts). The third criterion — that the dynamical center is coincident with a peak in the velocity dispersion — is satisfied in all but the second row. The fourth criterion — that the photometric major axis is aligned with the kinematic major axis — is satisfied in all rows. The fifth criterion — that the center of the continuum is coincident within the uncertainties of the dynamical center — is satisfied in all but the fifth row.

However, it is also apparent from the high-resolution HST H-band image that these galaxies are either highly disturbed or, as we discuss below, clearly have two nuclei.

These examples illustrate the challenge in distinguishing disks and mergers in seeing-limited data. Isolated disk galaxies are typically identified as disks using the observational disk criteria, but so too are mergers and sufficiently close pairs.

6.3. Contamination as a Function of 3D Pair Separation

In Figure 7, we show the fraction of sightlines meeting different subsets of the disk criteria for each snapshot of each galaxy. We consider only those galaxies that are dispersion-supported ($\left<j_{z}/j_{circ}\right>_{young}<1/3$) and with a nearby companion. The fraction is shown as a function of the 3D separation between the two galaxies.

In the top panel, we show the fraction of galaxies meeting the first 3 disk criteria outlined in §5.3. Between $10-100\%$ of the sightlines in these galaxies — with no intrinsic disk — are classified as disks from these criteria. We note that there is only a weak apparent dependance on the projected separation of the central and companion galaxies.

These first 3 criteria only use information from the kinematic maps. When a merger is misclassified as a disk, it is primarily due to the orbital motion of the merging galaxies. This orbital motion will often determine the direction of the kinematic major axis, while still satisfying the criteria $V_{rot}/\sigma_{g}>1$ — where $V_{rot}$ is, in reality, measuring the velocity difference of the two galaxies. The central peak in the velocity dispersion results from the steep velocity step-function midway between the two galaxies, i.e., beam smearing of the orbital motions.

The situation improves when we consider the photometric criteria. In the middle panel, we include the 4th and 5th disk criteria — the photometric and kinematic axes are aligned within 30 degrees and the steepest velocity gradient is coincident within the errors of the centroid of the continuum center. Adding these two criteria reduces the fractions to between 0 and 50%.

Finally, in the bottom panel, we explore an additional 6th criterion §5.3 — that the HST image contains only a single nucleus. We note that when high-resolution ancillary imaging and spectra are available, several studies adopt a similar criterion to select against merging galaxies. This is typically carried out through visual inspection or through quantitative and objective measures such as projected distances, redshift separations, and mass ratios (e.g., Übler et al. 2017). This final criterion dramatically improves upon the previous 5. The majority of non-disk galaxies are correctly identified as non-disks, except in two cases where the galaxies are $<5$ kpc average projected separation.

7. Quantifying the Effect of Mergers on Measurements of the Global Disk Fraction at High Redshift

In this next section, we use theoretically-derived pair fractions to estimate the contamination of mergers to the observed fraction of disks at high redshift ($1\lesssim z\lesssim 3$).

We have demonstrated that, due to seeing and small galaxy sizes, galaxies at $z\sim 2$ with a nearby companion and no disk can display the same kinematic signatures as those expected from disks (i.e., those outlined in §5.3). The fraction of sightlines that an observer will incorrectly identify a merger as a disk depends on both the observables available and the physical separation of the merging galaxies — ranging between 10 – 100% if only kinematic information is available and 0 – 50% if both kinematic and photometric information is available (Figure 7).

To estimate the degree of this effect on real observations, one must know the fraction of galaxies with companions or late-stage mergers — i.e., those at a projected separation $<\,$25 kpc — as a function of mass and redshift.

Through real observations of high redshift galaxy pairs, Mantha et al. (2018) find that the fraction of massive galaxies (M∗/M${}_{\odot}\,>2\times 10^{10}$) with major nearby companions ($4:1$ mass ratio; $5-50$ kpc projected separation) increases with redshift, to \sim$$15\% at $z\sim 1$, before beginning to decline. This is consistent with the close pair fractions computed by Rodrigues et al. (2017) ($\sim 15-17\%\,\pm\,6\%$, depending on the stringency of selection choice) at $z\sim 1$, once the close-pair selection cuts employed by both studies are made consistent.

As an important aside — given that the kinematic maps are luminosity-weighted, the appropriate statistic is the H$\alpha$-flux pair fraction— i.e., companion galaxies that have sufficiently high H$\alpha$ flux to register velocity values in otherwise unregistered spaxels. It is not straightforward to convert between a mass ratio-selected pair fraction and an H$\alpha$ ratio-selected pair fraction. For instance, Mantha et al. (2018) find that \sim$$20\% of galaxies qualify as a major pair at $z\sim 2.5$ when selecting by H-band flux ratio (4:1 flux ratio), while only $5\%$ qualify when selecting by stellar mass ratio (4:1 mass ratio).

7.1. The Fraction of Central Galaxies with a Companion in the Illustris Simulation

In lieu of using observed pair fractions, which are subject to incompleteness at low stellar masses and/or high mass ratios, we use results from the (106.5 Mpc)3 Illustris cosmological galaxy formation simulation box (Genel et al., 2014; Vogelsberger et al., 2014). Illustris mergers have been studied in detail by Rodriguez-Gomez et al. (2015). Pair fractions are calculated by Snyder et al. (2017) using lightcones, a technique for generating realistic synthetic surveys from cosmological simulations. Specifically, they followed Kitzbichler & White (2007) to generate three 140 square arcmin lightcone catalogs from which to select pairs using common observational criteria. Using simulation results here has its own limitations owing to their uncertain galaxy formation physics, and the challenge of properly assigning mass at small galaxy-galaxy separation (see extensive discussion in Rodriguez-Gomez et al. 2015), but the simulated pair statistics results agree reasonably well with data in the regions of parameter space where they overlap (Snyder et al., 2017).

For consistency with the rest of the paper, we first select pairs using the 3D separation of galaxies in the synthetic Illustris lightcone catalog. Then, using the results of the Appendix, we convert the 3D pair fraction statistics to the more observationally-accessible 2D statistics — pair fraction as a function of average projected 2D separation.

Figure 8 shows cumulative pair fractions as a function of projected separation on the sky at $1.5\,<\,z\,<\,3$ from the synthetic Illustris lightcone catalog. These are broken into two bins of stellar mass — $10^{9}\,<\,$M*∗/M${}_{\odot}<10^{10}$ (left panel) and $10^{10}<\,$M∗*/M${}_{\odot}<10^{11}$ (right panel). We further split these bins into two bins in stellar mass ratio (i.e., stellar mass of the companion: stellar mass of the central) — 1:1 – 1:4 (i.e., major pairs), 1:4 – 1:30 (i.e., minor and very minor pairs).

7.2. Estimating the Contamination to Observed Disk Fractions from Mergers in the Pair Stage

The observed fraction of galaxies with disks ($f_{obs}$) is related to the true fraction of galaxies with disks ($f_{int}$), the fraction of central galaxies with a sufficiently massive companion ($f_{pair}$), and the probability that, given a random sightline, a merger will pass as a disk ($X_{contamination}$):

[TABLE]

We determined that $X_{contamination}$ spans the full range between 0.0 and 1.0 for projected separations less than 25 kpc, depending strongly on the specific criteria used and weakly on the separation of the merging galaxies (Figure 7). For simplicity, we ignore that an isolated disk can be observed as a non-disk in a small fraction of sightlines (as found in 6). We use the Illustris lightcone pair fractions to determine $f_{pair}$ (Figure 8).

In Figure 9, we use the analytic model of Eq. 2 to translate from observed disk fractions to true disk fraction at $z\sim 2$. We adopt three values for $X_{contamination}$ that span the scatter in Figure 7: 10% (dashed-dotted line), 40% (solid line), and 70% (dashed line). The correction is smaller for high mass galaxies (M∗ = $10^{10-11}$ M*⊙) than it is for low mass galaxies (M∗ = $10^{9-10}$ M⊙*).

To estimate how large the correction is for real observations, we consider observed disk fractions from both slit-based (Simons et al. 2017) and IFS-based surveys (Wisnioski et al., 2015). We note that the slit-based measurement identify disks using only criteria 2 and 3 (and need to assume 1 and 4), but are in general agreement with the IFS-derived fractions (Simons et al., 2017). The observed value of $f_{obs}$ at $z\sim 2$ is $\,\simeq\,70\%$ (Wisnioski et al., 2015; Simons et al., 2017) for high mass galaxies and the difference between $f_{obs}$ and $f_{int}$ is $<\,5\%$. For low mass galaxies, the observed disk fraction at $z\sim 2$ is $f_{obs}\,\simeq\,50\%$ (Simons et al., 2017), and the difference is $\,<15\%$.

In Figure 10, we estimate the true disk fraction given observed disk fractions from Simons et al. (2017) with the most extreme assumption one can make about merger contamination — that 100% of pairs inside some projected separation ($D_{contamination}$) masquerade as disks. We explicitly show the observed value of the disk fraction, with their uncertainties (dashed line and shaded region). The true disk fraction derived from Eq. 2, i.e., $f_{int}$, is shown as a solid line.

Even with this extreme assumption, the low mass galaxies only need a significant correction factor ($25-50\%$) to $f_{obs}$ for $D_{contamination}$ = 15 - 25 kpc. However, these cases almost always show two nuclei in the HST imaging and should identified in samples if imaging is available (Figure 7). At high mass, the correction factor is still minimal ($<5\%$), even for $D_{contamination}=25$ kpc.

The discrepancy between the true disk fraction and the observed disk fraction due to merger contamination is minimal for both low-mass galaxies ($10^{9}\,<\,\,$M*∗/M${}_{\odot}\,<\,10^{10}$) and high-mass ($10^{10}\,<\,\,$M∗*/M${}_{\odot}\,<\,10^{11}$) galaxies. Taking the Illustris pair fractions at face-value, the difference between the observed disk fraction and the true disk fraction is estimated to be $\lesssim 5\%$ and $\lesssim 15\%$ at high and low mass, respectively.

7.3. Estimating the Contamination to Observed Disk Fractions from Mergers in the Late-Stage

In the two subsections above, we focus on the pair and close pair stages — where the two galaxies have a projected separation on the sky of at least 5 kpc. This owes to the uncertainty of the Illustris-derived pair fractions when the simulated galaxies are near-coalescence and overlapping (Rodriguez-Gomez et al., 2015).

We now estimate the fraction of central galaxies in the late-stage/pre-coalescence stage, i.e., the mergers that are too close to be included in the pair fractions but are not yet relaxed. The fraction of galaxies in the pair stage $F_{pair}$ ($5<d<50$ kpc) and the fraction of galaxies in the late stage $F_{late}$ ($d<5$ kpc) are related by their relative timescales:

[TABLE]

where $\tau_{pair}$ and $\tau_{late}$ are the time a single merger will spend in the pair stage and the late stage, respectively.

The time for a merger to pass from the pair stage to the late-stage/coalescence stage, i.e., $\tau_{pair}$, was measured in Illustris as $\sim 700$ Myr at $z\sim 2$ (Snyder et al., 2017).

The time for a merger to dynamically relax after entering the late-stage is an unknown, but we use the following simple argument to estimate its value. The kinematics observations trace the ionized gas in galaxies. As gas is collisional, it should dynamically relax within a $\sim$ few crossing times ($\sim$100 Myr) after the mergers enter the late-stage.

Let us consider $\tau_{late}=350$ Myr and $\tau_{pair}=700$ Myr. To associate the timescale $\tau_{pair}$ with the correct pair fraction in Eq. 2, we must use the pair fraction as measured at the average separation of all pairs in the pair stage. The Illustris-derived pair fractions at $z\sim 2$ ($F_{pair}$) — the fraction of galaxies with a pair within 25 kpc projected separation — are 20% at high mass and 35% at low mass (Fig. 8). Given these values, $F_{late}$ as calculated from Eq. 3 is only 7% and 12%, respectively. We then use Eq. 2 to estimate the correction to the observed disk fraction. Even if we assume that 100% of these late-stage mergers masquerade as disks (i.e., the most extreme assumption), the correction to the observed disk fractions is only $\sim 1.8\%$ and 5.5$\%$ for high and low mass galaxies, respectively (i.e., negligible).

8. Conclusions

A majority ($\sim 70\%$) of massive star-forming galaxies at $z\sim 2$ are observed to have the following disk-like qualities: continous velocity gradients, rotation velocities that exceed their average local ionized gas velocity dispersions, central peaks in their velocity dispersion maps, and aligned kinematic and photometric major axes (Wisnioski et al., 2015). It is generally inferred that these galaxies are disk-like, albeit with large amounts of disordered motions. However, with typical ground-based seeing, the orbital motions of merging galaxies can appear regular and disk-like in resolved kinematic maps (Hung et al., 2015, 2016). Given the increased merger activity at $z\sim 2$ compared to the local Universe, it is important to quantify the degree to which mergers are misclassified as disks in current high redshift samples.

This paper has two main conclusions. First, from seeing-limited kinematic data and Hubble images, one cannot tell with certainty whether any individual galaxy at $z\sim 2$ is a merger or a disk. Second, the estimated merger rates at $z\sim 2$ are not high enough for mergers to have a significant effect on measured disk fractions.

To come to these conclusions, we created synthetic Hubble and KMOS IFU observations of cosmological galaxy formation simulations at $z\sim 2$, and measured the probability that a merger will pass as a disk as a function of merger stage. We found that merging galaxies pass the disk criteria in a non-negligible fraction of sightlines. This fraction ranges from 0 to 100$\%$ and depends strongly on the specific disk criteria adopted, and weakly on the separation of the merging galaxies (Fig. 7).

To estimate how often mergers pass as disks in current observational samples, we combined these fractions with estimates for the number of galaxies in a merger at $z=2$ from the Illustris simulation. For galaxies with stellar masses above $10^{9}$ M*⊙, the merger rates are not high enough for this effect to have a significant impact on measurements of disk fractions (Fig. 9). Taking the simulations at face-value, the observed disk fractions are overestimated by at most $15\%$ and $5\%$ at low ($10^{9-10}$ M⊙) and high ($10^{10-11}$ M⊙*) stellar mass, respectively.

Acknowledgements

RCS would like to acknowledge support for this project through a STScI/DDRF grant to SAK. GFS appreciates support from a Giacconi Fellowship at the Space Telescope Science Institute, which is operated by the Association of Universities for Research in Astronomy, Inc., under NASA contract NAS 5-26555. KBM acknowledges support from the HST archival research grant HST-AR-15040. This work was partly supported by the US-Israel BSF grant 2014-273 (AD, JRP) and by the NSF AST-1405962 (AD). The Flatiron Institute is supported by the Simons Foundation. DC has been partly funded by the ERC Advanced Grant, STARLIGHT: Formation of the First Stars (project number 339177). NM acknowledges support from the Klauss Tschira Foundation through the HITS-Yale Program in Astrophysics (HYPA). We thank the referee for a useful report. The simulations were performed at NASA advanced Supercomputing (NAS) at NASA Ames Research Center and at the National Energy Research Scientific Computing Center (NERSC), Lawrence Berkeley National Laboratory. This work has made us of the publicly available software YT (Turk et al., 2011). This research made use of Astropy, a community-developed core Python package for Astronomy (Astropy Collaboration et al., 2013).

Appendix A From 3D separation to average 2D projected separation

Throughout this paper, we substitute a physical quantity measured from the simulations — the 3D separation between a central galaxy and a companion galaxy — with a quantity more relevant for comparing with observations — their average projected separation on the sky. As shown below, these quantities are related by a factor $\pi/4$. Let us define an arbitrary 3D point ${\vec{r}}=(x,y,z)$ with a 3D distance from the origin $\rho$:

[TABLE]

The projected separation $R$ of this point and the origin for an observer normal to $\hat{z}$ is:

[TABLE]

The average of $R$ over a sphere of fixed radius $\rho$ is:

[TABLE]

That is, the average 2D projected separation of two points is $\pi/4$ times their 3D separation.