Resolving discs and mergers in z ∼ 2 heavily reddened quasars and their companion galaxies with ALMA Open Access

ABSTRACT

1 INTRODUCTION

In our current paradigm of galaxy formation, accreting supermassive black holes play a fundamental role in shaping the evolution of the most massive galaxies. Gas-rich major mergers have for some time been invoked in theoretical models to fuel both star formation and accretion on to the black hole in massive galaxies (e.g. Di Matteo, Springel & Hernquist 2005; Springel, Di Matteo & Hernquist 2005; Hopkins et al. 2008; Narayanan et al. 2010; Gabor et al. 2016; Steinborn et al. 2018). Idealized as well as cosmological simulations of merging galaxies now have the requisite spatial and temporal resolution to be able to demonstrate that the merger leads to enhanced star formation (e.g. Moreno et al. 2015), triggering of an active galactic nucleus (AGN; e.g. Capelo et al. 2015; Angles-Alcazar et al. 2020), and an alteration in the interstellar medium (ISM) properties of the galaxy (e.g. Moreno et al. 2019). Energy injection from the quasar during a feedback-dominated phase (Fabian 2012) is then thought to be responsible for quenching star formation and curtailing the growth of these massive galaxies. Cosmological simulations routinely implement such feedback in order to explain the observed distributions of massive galaxies (e.g. Sijacki et al. 2015; Beckmann et al. 2017). A major aim of current observational research is to find compelling evidence in favour of this standard paradigm. Recently several aspects of this simple picture have been confronted by new observational data. For example, the existence of massive gas discs in the high-redshift Universe (e.g. Tacconi et al. 2013; Förster Schreiber et al. 2018) provides evidence that processes such as smooth accretion via cold gas flows (Dekel et al. 2009) or minor mergers may play a more prominent role in assembling galaxies even at the high-mass end of the galaxy mass function therefore challenging the major-merger paradigm. The effects of AGN feedback too can be more nuanced than previously thought. For example, theoretical models predict that compression of gas in galactic outflows can indeed trigger rather than inhibit star formation (e.g. Ishibashi & Fabian 2012). There have been some claims of observations of such positive feedback effects (e.g. Cresci et al. 2015; Carniani et al. 2016) with evidence now emerging that stars can in fact be formed in situ within a galactic outflow (Maiolino et al. 2017; Gallagher et al. 2019).

In the context of the merger-driven evolutionary scenario for the formation and co-evolution of massive galaxies and their supermassive black holes, populations of obscured quasars could represent an important phase as a galaxy is transitioning from a dust-obscured starburst to an ultraviolet (UV) luminous quasar (e.g. Hickox & Alexander 2018). The idea has its origins in the seminal work of Sanders et al. (1988). With the advent of wide-area surveys in the X-ray, near- and mid-infrared wavelengths, population studies of the luminous, obscured quasar population are now possible (e.g. Banerji et al. 2012, 2013, 2015; Glikman et al. 2012; Assef et al. 2015; Hamann et al. 2017; LaMassa et al. 2017; Temple et al. 2019). The most luminous obscured quasars detected in these wide-area surveys are often seen to be hosted by merging galaxies (e.g. Urrutia, Lacy & Becker 2008; Glikman et al. 2015; Fan et al. 2016; Bischetti et al. 2018; Díaz-Santos et al. 2018) with significant amounts of ongoing star formation (e.g. Urrutia et al. 2012; Banerji, Fabian & McMahon 2014; Banerji et al. 2017; Wethers et al. 2018). Powerful AGN-driven outflows have also been detected in several of these systems (e.g. Urrutia et al. 2009; Cresci et al. 2015; Perna et al. 2015; Zakamska et al. 2016; Temple et al. 2019; Kakkad et al. 2020), with the observations suggesting that these outflows could have both a positive and negative feedback effect on star formation.

Luminous quasars were among the first high-redshift (z ≳ 2), extragalactic sources to be targeted in molecular gas emission (e.g. Solomon & Vanden Bout 2005, and references therein) and submillimetre and millimetre observations provide us with an opportunity to study host galaxy emission from these quasars over a wavelength range where the emission from the quasar itself is sub-dominant. Molecular gas has been detected in quasars out into the epoch of reionization (e.g. Walter et al. 2003, 2004; Maiolino et al. 2005; Wang et al. 2010, 2011; Venemans et al. 2017) and the new capabilities of the Atacama Large Millimeter Array (ALMA) are now making it possible to spatially resolve both the dust and gas emission in high-redshift quasar hosts (e.g. Salomé et al. 2012; Carilli et al. 2013; Carniani et al. 2013; De Breuck et al. 2014; Wagg et al. 2014; Brusa et al. 2018; Shao et al. 2017; Trakhtenbrot et al. 2017; Venemans et al. 2017; Fogasy et al. 2020). Nevertheless, the number of high-redshift quasars with such resolved measurements of the molecular gas emission are still relatively small. Resolved measurements of the host galaxies in the millimetre wavelengths enable unique constraints on the dynamical masses, gas morphologies and gas fractions in such systems.

Observational evidence is mounting that many of the major physical processes thought to govern massive galaxy formation – e.g. mergers, starbursts, high-accretion black hole growth and outflows – are synchronously occurring during the dust-obscured quasar phase. How are these processes triggered and what effect do they have on the quasar host galaxies? We address these questions here via detailed, spatially resolved observations of the CO(3–2) gas emission in the quasars ULASJ1234+0907 (z ∼ 2.50) and ULASJ2315+0143 (z ∼ 2.56). The two quasars are very similar in terms of the AGN properties. Both quasars are extremely luminous (L_6 µm ∼ 10⁴⁷ erg s⁻¹, L_X ∼ 10⁴⁵ erg s⁻¹; see Lansbury et al. 2020) with massive black holes (M_BH > 10⁹ M_⊙) obscured by significant columns of dust and gas (N_H ∼ 10²² cm⁻², A_V ∼ 5−6 mag; Banerji et al. 2012, 2015, 2017; Wethers et al. 2018; Lansbury et al. 2020). Both quasars are known to reside in highly star-forming host galaxies (Banerji et al. 2014; Wethers et al. 2018, 2020). The first detections of molecular gas via the CO(3–2) emission in the two quasars was presented in Banerji et al. (2017) (B17 hereafter). In Banerji et al. (2018) (Paper I hereafter), we combined the CO(3–2) detections with detections of CO(1–0), CO(7–6), and C i(2–1) to better constrain the ISM properties of the quasar host galaxies. While Paper I focused on the integrated line emission properties of our quasar host galaxies, in this paper we concentrate on studying the resolved emission to gain a better understanding of the physical processes that govern the formation of these galaxies. Throughout the paper we assume a flat, |$\Lambda$| cold dark matter cosmology with H₀ = 70 km s⁻¹ Mpc⁻¹, Ω_M = 0.3, and Ω_Λ = 0.7.

2 ALMA BAND 3 OBSERVATIONS

Our previous work (B17) presented detections of the CO(3–2) line in a sample of four heavily reddened quasars (HRQs) at z ∼ 2.5 that were initially discovered in Banerji et al. (2012, 2015). These four quasars are among the reddest quasars found in our near-infrared search for luminous, obscured quasars that are not present in wide-field optically selected spectroscopic quasar samples such as those from the Sloan Digital Sky Survey (SDSS) and BOSS. Briefly, the quasars have extremely red near-infrared colours of (H − K)_Vega > 2, corresponding to dust extinctions towards the quasar continuum of E(B − V) ≳ 1 (A_V ≳ 3 mag for R_V = 3.1 mag). The bolometric luminosities and black hole masses are comparable to the most luminous quasars selected in optical surveys such as the SDSS. As a detection experiment, the Cycle 3 ALMA observations in B17 made use of the most compact available configuration of the telescope leading to an angular resolution of ∼2.5–3 arcsec. Two of the quasars from the sample – ULASJ1234+0907 (z ∼ 2.50; J1234 hereafter) and ULASJ2315+0143 (z ∼ 2.56; J2315 hereafter) – were marginally resolved in the Cycle 3 data and were therefore re-observed during ALMA Cycle 4 at ∼3× better angular resolution with the aim of spatially resolving the gas emission in both sources and constraining the gas dynamics. J1234 has two nearby gas-rich galaxies that were first discovered via their CO(3–2) detections in B17. Both companion galaxies (named G1234N and G1234S) were also imaged at higher resolution as part of the observations presented in this work. Further details of the new ALMA Cycle 4 observations can be found in Table 1.

The spectroscopic observations of both quasar systems made use of a correlator configuration with four dual polarization bands, each of bandwidth 1.875 GHz. In the case of J1234, observations made use of the time division mode (TDM) of the correlator, providing a channel width of 15.625 MHz. In the case of J2315, the frequency division mode (FDM) of the correlator was instead used, providing a much finer channel width of ∼488 kHz across the entire bandwidth, which we then binned up to 15.625 MHz resolution.

The raw data were processed using casa (v4.7.0), by executing the appropriate calibration scripts for each set of observations as provided by the ALMA Observatory. Line cubes were produced after u–v plane continuum subtraction employing the line-free channels. We then used the clean algorithm with a natural weighting of the visibilities to transform to the image plane. Continuum images were also produced by running clean in mfs mode on the line-free spectral windows, and once again using a natural weighting. All images have had a primary beam correction applied to them using the casa task pbcor. Images have been generated with a pixel-scale of 0.1 arcsec pixel⁻¹. The beam sizes as well as RMS noise levels reached by the data for each of the sources is summarized in Table 1. The spectra have all had Hanning smoothing applied to them as part of the pipeline processing.

3 OBSERVED MOLECULAR GAS AND DUST CONTINUUM PROPERTIES

In Figs 1 and 2, we present the CO(3–2) images as well as integrated spectra for all galaxies in the J1234 and J2315 systems, respectively, obtained via our ALMA Cycle 4 program. We infer initial source sizes by fitting a two-dimensional Gaussian to the zeroth moment image using casa. The spectra have been extracted using apertures of variable sizes (diameter 1–2 arcsec) centred on the positions of each galaxy stated in Table 2. The aperture sizes chosen depend on the size of the CO emitting region for each galaxy, which is also quoted in Table 2 and discussed further in Sections 5 and 6. These new observations provide higher signal-to-noise ratio (S/N) spectra for all sources compared to the Cycle 3 data presented in B17, as well as higher S/N measurements of the continuum flux densities at observed frequencies of ∼86–103 GHz. These measurements are summarized in Table 2.

Figure 1.

Maps of the CO(3–2) emission in the three galaxies making up the J1234 system (left). North is up and East is to the right. The quasar position is marked with the star in the top-left image and the synthesized beam size and position angle are indicated by the blue ellipse at the bottom right of each image. The right-hand panels show the integrated CO(3–2) spectra for the three galaxies together with the best-fitting double Gaussian fit to these line profiles. Double-peaked profiles are evident for all three galaxies. Zero velocity, indicated by the dotted vertical line corresponds to the systemic redshifts of these galaxies obtained from our ALMA Cycle 3 observations in B17.

Open in new tabDownload slide

Figure 2.

Top: Map of the CO(3–2) emission in the J2315 system. North is up and East is to the right. The source is clearly resolved into two spatial components and the quasar position is indicated by the star. The synthesized beam size is indicated by the blue ellipse at the bottom-right of the image. Bottom: CO(3–2) spectra extracted over the entire system (left), at the quasar position (middle), and at the position of the northern companion galaxy (right) together with the best-fitting single Gaussian used to fit these line profiles. Zero velocity, indicated by the dotted vertical line corresponds to the systemic redshift obtained from our ALMA Cycle 3 observations in B17.

Open in new tabDownload slide

All three galaxies in the J1234 system are spatially resolved in our data based on the fact that casa is able to deconvolve the intrinsic emission from the convolved image. We have also checked that the radial profiles of all three galaxies are more extended than would be expected for a point source. The spectra exhibit clear double-peaked CO(3–2) emission line profiles. Such clear double-peaked profiles were not evident in the original Cycle 3 spectra in B17, and indicate that we have now resolved multiple kinematic components in the gas emission. The J2315 system, which was already spatially resolved in CO emission in B17, now breaks up into two distinct galaxies in a close-separation major merger in our new data. We are therefore able to extract CO(3–2) spectra for both the quasar host and the companion galaxy lying to the North of the quasar. Table 2 presents the line and continuum flux densities for both our quasars and their companion galaxies. These measurements have been analysed in more detail in Paper I, where we combine them with observations of the CO(1–0), CO(7–6), and C i(2–1) lines as well as the far-infrared to radio continuum in order to infer the properties of the ISM in these galaxies. In Paper I, we find that the quasar J1234 as well as its companion galaxy G1234S are comprised of highly excited, dense gas consistent with what is seen in luminous AGNs. J2315 on the other hand has excitation properties that are more similar to that observed in star-forming galaxies. Both the quasars J1234 and J2315 as well as their companion galaxies were inferred to be highly star forming based on their dust continuum detections in Paper I.

The focus of the current paper is on a detailed kinematic analysis of the CO(3–2) emission and inferred constraints on the gas dynamics and gas morphologies in our two quasar systems. Before turning to such an analysis we first present a description of the dynamical models we will employ.

4 DYNAMICAL MODELLING

4.1 ^3DBarolo

We model the CO(3–2) traced gas dynamics in our quasars and their companion galaxies using the tilted ring model fitting code ^3DBarolo – 3D-Based Analysis of Rotating Objects from Line Observations (Di Teodoro & Fraternali 2015). The approach is similar to our previous work (B17) employing the GIPSY rotcur task, but with the key improvement that ^3DBarolo fits models to the entire data cube rather than just the velocity field. By convolving the model with the restoring beam and comparing data and model on a channel-by-channel basis, the effects of beam smearing are reduced, and additional dynamical aspects may be examined. Such an approach has been used previously to study the rotation curves of nearby galaxies (e.g. de Blok et al. 2008; Iorio et al. 2017; Shelest & Lelli 2020) and is now being routinely applied at high redshifts (e.g. Jones et al. 2017; Xue et al. 2018; Tadaki et al. 2020).

The three-dimensional model fitting considers some of the same variables that are included in two-dimensional fits: inclination, i, position angle, PA, rotational velocity, v_rot, systemic redshift, v_sys, and the spatial centroid [x0, y0]. By expanding to the spectral dimension, we can also consider velocity dispersion σ_v, radial brightness profile I(r), and scale height Z_o as free parameters. While ^3DBarolo has been well tested on low-S/N and low-resolution data sets (e.g. Di Teodoro & Fraternali 2015), the best-fitting results are strongly dependent on both the model geometry (i.e. number of rings, ring width) and initial estimates for each parameter. To reduce any implicit bias, we automate the fitting process by using signal isolation and moment maps.

We first utilize the ^3DBarolo SEARCH utility (based on DUCHAMP; Whiting 2012) to create a signal mask. This utility automatically finds the RMS noise level of an input data cube and identifies signal peaks above a user-provided S/N threshold (S/N_upper = 3.5). A search is then conducted around these peaks for emission above a second user-provided S/N_lower = 3.0, creating a three-dimensional mask of significant emission. This approach includes fewer non-signal pixels than a 2D spaxel mask or a mask based on a single S/N threshold.

Next, we use the casa task immoments to create moment 0 (integrated intensity) maps using all channels containing significant line emission. By fitting this map with a two-dimensional Gaussian (casa task imfit), we find a best-fitting beam-deconvolved full width at half-maximum (FWHM) of the major and minor axes, position angle, and morphological central position. We also apply immoments to all pixels identified as signal to create moment 1 (line-of-sight velocity) and moment 2 (velocity dispersion) maps. These maps can be seen in the left-most panels of Figs 3, 4, 5, and 8.

Figure 3.

Line intensity (top; M0), velocity (middle; M1), and velocity dispersion (bottom; M2) maps for ULASJ1234. The observed data are shown in the left-most panels and the best-fitting models from ^3DBarolo together with the fit residuals are shown in the middle and right-hand panels, respectively.

Open in new tabDownload slide

Figure 4.

Line intensity (top; M0), velocity (middle; M1), and velocity dispersion (bottom; M2) maps for G1234N. The observed data are shown in the left-most panels and the best-fitting models from ^3DBarolo together with the fit residuals are shown in the middle and right-hand panels, respectively.

Open in new tabDownload slide

Figure 5.

Line intensity (top; M0), velocity (middle; M1), and velocity dispersion (bottom; M2) maps for G1234S. The observed data are shown in the left-most panels and the best-fitting models from ^3DBarolo together with the fit residuals are shown in the middle and right-hand panels, respectively.

Open in new tabDownload slide

With this mask and these maps in hand, we are able to provide physically motivated estimates for most model parameters. We base the maximum model radius on the 2D Gaussian fit to the moment zero map: 1.0× (the deconvolved FWHM of the major axis). In order to not over-resolve the data, we adopt a minimum ring width of (FWHM of the minor axis of the restoring beam)/2.5 – see also Talia et al. (2018) and Fan et al. (2019). The number of model rings is then found by dividing the maximum model radius by the minimum ring width and rounding down to the nearest integer, while the ring width is taken to be the maximum model radius divided by this number. For all galaxies modelled in this work, this process returns models with three rings. The central position of the model is fixed to be the centroid of the 2D Gaussian fit to the moment zero map. Since each emission line is well detected, the systemic velocity is expected to be 0 km s⁻¹, with respect to the redshifts listed in Table 3. ^3DBarolo provides an estimate of the rotational velocity based on the moment 1 map, while we provide the value of the moment 2 map at the morphological centre as an estimate of the velocity dispersion. The inclination estimate is set to 45°, and allowed to vary between 10° and 80°, while the position angle is provided by the user. We assume a thin disc with a scale height of ∼100 pc.

The uncertainties on the dynamical mass only include the output rotational model uncertainties from ^3DBarolo (δv_rot,max) and the width of the model ring (δr_max) and the significant uncertainty on the inclination is not included. Since the fitting process is three-dimensional and multistaged, the inclusion of the inclination uncertainty in the error budget is non-trivial and we simply note that the formal uncertainties on the dynamical masses from ^3DBarolo are likely to be underestimated

4.2 Exponential disc modelling

^3DBarolo makes no assumption regarding the underlying density distribution and mass profile of the galaxy and simply fits the CO emitting material across the extent over which it is seen in our ALMA data. We therefore use an independent method to check the mass estimates by approximating the kinematics of the galaxies as a thin rotating disc (Carniani et al. 2013; Pensabene et al. 2020). Such a model assumes that the surface brightness of the disc has an exponential profile |$\Sigma \propto \mathrm{ e}^{-r/R_d}$|⁠, where R_d is the scale radius that corresponds to a half-light radius of R_1/2 = 1.67 × R_d. By assuming that the mass surface density has the same profile as the gas emission, the circular velocity associated with the exponential profile can be calculated following equation 2.157 of Binney & Tremaine (2008). Both intensity and velocity model maps are then matched to the spatial resolution of the observations by taking into account the observed synthesized beam. The free parameters of this dynamical model are the spatial coordinates of the disc centre, the half-light radius (or scale radius), inclination, position angle, dynamical mass, and systemic velocity. The dynamical mass is the disc mass enclosed within a radius of 10 kpc. Because this radius is much larger than the typical disc scale radius at z ∼ 2 (Cresci et al. 2009), this mass corresponds to the total baryonic mass in the disc. Note also that this radius is much larger than the radius over which the CO emission is visible, and therefore the radius over which ^3DBarolo fits are conducted. The full parameter space is explored by using the Markov chain Monte Carlo algorithm emcee (Foreman-Mackey et al. 2013)), which allows one to obtain a probability distribution of all free-parameters along with marginalized uncertainties. The disc centre coordinates and half-light radius are derived from modelling the intensity maps, while the probability distributions of the other parameters are driven from the fitting of the velocity maps.

5 THE J1234 SYSTEM: A QUASAR WITH TWO MASSIVE COMPANION GALAXIES

5.1 Dynamical masses, sizes, and gas fractions

The quasar host galaxy J1234 (z = 2.501) and its two CO-luminous companions – G1234N (z = 2.514; ∼170 kpc projected distance from the quasar) and G1234S (z = 2.497; ∼90 kpc projected distance from the quasar) – are from hereon referred to as the J1234 system. The CO images shown in Fig. 1 have been constructed by collapsing all the spectral channels containing the line in each source. Initial source sizes for the gas emission were derived by fitting a single two-dimensional Gaussian to the CO images using the casa task imfit. Using these constraints as an initial guess, the fitting was then re-done in the u–v plane using uvmodelfit. The results from the two fitting approaches are consistent for all the galaxies and the FWHM of the major axis is 0.48 ± 0.14, 0.57 ± 0.18, and 0.53 ± 0.24 arcsec for ULASJ1234, G1234N, and G1234S, respectively, after deconvolving from the beam. The FWHM values can be converted to a half-light radius using the equation, R_1/2 = 0.5 × FWHM for a Gaussian profile, and correspond to physical scales of ∼2 kpc at the redshifts of these galaxies. We summarize the half-light radius constraints for all three galaxies in Table 3. A second estimate of the half-light radius can be estimated from fitting the exponential disc models (Section 4.2) to the intensity maps. The half-light radii derived from the exponential disc models are also specified in Table 3 and are generally larger than those derived from the Gaussian fits – ∼3–4 kpc at the redshifts of these galaxies. This is to be expected given the 2D Gaussian profiles are not able to reproduce the full spatial extent of the CO emission seen from these galaxies.

We have already noted the double-peaked line profiles for all three sources seen in Fig. 1, which provide evidence for kinematically distinct components of the gas emission. Two Gaussian components are required to adequately model the one-dimensional spectral line profiles and these best spectral-fit models are overplotted on the observed spectra in Fig. 1. The velocity separation between the two Gaussian components is 510 ± 160, 340 ± 50, and 350 ± 40 km s⁻¹ for the quasar, G1234N, and G1234S, respectively. Such double-peaked profiles could arise from the approaching and receding sides of a rotating gas disc, a coalescing merger of two distinct galaxies, or large-scale inflows and outflows affecting the molecular gas.

As can be seen in the moment 1 maps in Figs 3, 4, and 5, a velocity gradient is evident across the full spatial extent of all three galaxies that comprise the J1234 system. The best-fitting ^3DBarolo models together with the fit residuals can also be seen in the same figures. In general, the fit residuals are small, at least in the central portions of these galaxies, suggesting that a smooth velocity gradient provides a reasonable approximation to these data. The resulting dynamical masses are summarized in Table 3. With dynamical masses of ∼10¹¹ M_⊙, all three galaxies in the J1234 system can be considered among the most massive systems seen at these epochs.

We also derive dynamical masses for all three galaxies using the exponential disc modelling. The resulting dynamical masses summarized in Table 3 are larger than those obtained from ^3DBarolo by 0.9, 0.1, and 0.3 dex for the three galaxies. As stated earlier, the exponential disc models are sensitive to mass over a much larger physical radius than the ^3DBarolo models. Furthermore, the inclination of these galaxies is poorly constrained and adds significant additional uncertainty to the mass estimates. Given these factors, the discrepancies in the dynamical mass estimates from the two methods is not surprising.

In Fig. 6, we show the dynamical mass versus the half-light radius of the CO emission from the exponential disc models for all three galaxies in the J1234 system and compare to both field and cluster galaxy samples from the literature at z ∼ 1−2.5. The field galaxies from the PHIBSS sample correspond to CO(3–2) measurements (Tacconi et al. 2013), whereas the z ∼ 1.6 cluster galaxy measurements are taken from the CO(2–1) observations presented in Noble et al. (2019). We take the galaxy mass as the sum of the stellar mass and molecular gas mass in the case of the literature samples and as the dynamical mass in the case of our quasar host galaxy and its companions. In general, all three galaxies in the J1234 system are at the more massive and more compact end of the distribution of z ∼ 2 galaxies from the literature. We also show the disc half-light radii for z ∼ 4−6 quasars from Pensabene et al. (2020), which uses the same disc modelling as we have employed in this work. The higher redshift quasars are generally found to be more compact in size but this likely reflects the differing spatial distribution of the [C ii] gas in the case of the z ∼ 4−6 QSOs and the relatively low-J CO gas traced by the z ∼ 2 samples.

Figure 6.

Host galaxy mass versus disc half-light radius for the three galaxies in the J1234 system and compared to various samples from the literature. The PHIBSS field galaxy sample is taken from Tacconi et al. (2013), the cluster sample measurements are taken from Noble et al. (2019), and the high-redshift quasar sample is from Pensabene et al. (2020).

Open in new tabDownload slide

In Paper I, we derived gas masses for all galaxies in the J1234 system using several different estimators, thus mitigating against the uncertainties implicit in using a single molecular tracer of the total gas mass. The resulting gas masses were found to be ∼3 × 10¹⁰, ∼2 × 10¹⁰, and ∼1−7 × 10¹⁰ M_⊙ for ULASJ1234, G1234N, and G1234S, respectively. In the case of G1234S, the large range results from a much larger C i-derived gas mass compared to the CO-derived gas mass (Paper I). Combining these gas mass estimates with the dynamical masses and assuming that the inner parts of the discs traced by our CO observations are indeed baryon-dominated, we can infer the molecular gas fractions. The gas fraction is traditionally defined as the ratio of the gas mass to the total stellar mass. In the absence of stellar mass estimates for our high-luminosity quasars, we infer the stellar masses from the dynamical masses via M_* = M_dyn − M_gas −M_BH, as has been done in previous studies of high-luminosity quasars (e.g. Bischetti et al. 2021). In the case of the quasar, very different gas fractions are implied depending on whether we consider the dynamical mass from ^3DBarolo, which represents the mass in the innermost regions of the galaxy, which are seen in CO emission, or whether we consider the exponential disc model mass, which extends out to 10 kpc. In the former case, we obtain gas fractions of ∼85 per cent whereas the gas fraction relative to the extended disc model mass is more like ∼6 per cent. It is reasonable to assume that the galaxy is baryon dominated even out to 10 kpc, in which case there must be a substantial stellar component already in place with M_* ∼ 10¹¹ M_⊙ in this quasar. G1234N has a gas fraction of ∼14–18 per cent whereas in G1234S the gas fraction is more uncertain owing to the discrepant gas masses from the CO and C i observations, and could be in the range 7–90 per cent.

The star formation rates for the two companion galaxies are poorly constrained due to a lack of photometric data points sampling the peak of the galaxy dust SED. However, for the quasar host galaxy we find an SFR of >1000 M_⊙ yr⁻¹ based on Herschel and ALMA observations (Paper I; Banerji et al. 2014). The well-established star-forming ‘main sequence’ of galaxies at high redshifts (Daddi et al. 2007) can be used to predict how much star formation a typical star-forming galaxy of a certain stellar mass is expected to have at these redshifts. We find that an M_* ∼ 10¹¹ M_⊙ galaxy at z ∼ 2 would have a star formation rate of ∼200 M_⊙ yr⁻¹ were it to lie on the main sequence and the quasar host galaxy is therefore clearly experiencing a prolific starburst. Both starburst galaxies as well as the most massive main-sequence galaxies in the high-redshift Universe have molecular gas fractions of ∼30–50 per cent Tacconi et al. (2013), Casey, Narayanan & Cooray (2014) with some suggestions that the gas fractions might be even higher in richer protocluster and cluster environments (Noble et al. 2017; Hayashi et al. 2018).

If we consider the inner region of the J1234 quasar, the gas fraction is very high and comparable to what is seen in the most gas-rich starburst galaxies. However, over more extended scales the galaxy would appear to be substantially depleted of gas. In the case of the quasar host galaxy, it is plausible that feedback from the AGN has begun to deplete the molecular gas reservoir, e.g. as suggested by Brusa et al. (2018). However, the two companion galaxies also exhibit somewhat depleted gas reservoirs with currently no evidence for an actively accreting AGN within them. Tadaki et al. (2019) suggest that the most massive galaxies with masses of >10¹¹ M_⊙ in the highest density environments at these redshifts could have gas accretion suppressed due to inefficient cooling of massive haloes, which could explain the fact that both G1234N and G1234S have gas fractions that are at the lower end of what is seen in typical star-forming galaxies at these epochs. We caution however that there are considerable uncertainties in both the dynamical and gas masses as illustrated by the discrepancies in these estimates using different dynamical modelling methods and different gas mass tracers.

5.2 CO line profiles

In Fig. 7, we compare the line profiles for the various different CO transitions presented here and in Paper I for all three galaxies in the J1234 system. Progressively higher J transitions of CO trace warmer, denser gas and a comparison of the line profiles therefore helps illustrate any potential differences in the kinematics of these different gas components. For the quasar J1234 and G1234N all the CO lines appear to be well matched in terms of their velocity and line shape suggesting that the cooler lower density, and warmer, denser gas are experiencing similar large-scale motions in these systems. A double-peaked line profile is evident in the CO(7–6) emission from G1234S, which matches what is seen in the CO(3–2) emission. However, for this particular galaxy we also see a very broad wing in the CO(7–6) emission extending out to velocities of >1000 km s⁻¹. Note, the blue side of the broad-wing is not traced by our observations as the C i(2–1) emission line affects this portion of the spectrum. Curiously, no such broad wings are seen in the CO(3–2) line profile of this galaxy. If the broad wing is ascribed to a molecular outflow, it appears to only be affecting the dense gas in the galaxy. Molecular outflows with similar observational properties to what is seen for G1234S, have previously been observed in luminous, high-redshift quasars (Maiolino et al. 2012; Cicone et al. 2015), where the AGN is responsible for driving the gas to velocities of ∼1000 km s⁻¹. Currently we have no evidence for an AGN being present in G1234S. However, no optical/infrared counterpart to the galaxy exists in our imaging data of this field down to i-band magnitudes of <25.1 (Banerji et al. 2014), which suggests that the galaxy is highly obscured. We therefore speculate that G1234S with its broad wings seen in CO(7–6) emission and high excitation CO SLED (Paper I) may well host an obscured AGN although high-resolution, deep X-ray or mid-infrared data would be needed to confirm or refute this claim. It is also curious that only the dense gas is being affected by any putative outflow in this galaxy. Such an observation would be consistent with recent claims that star formation can happen in situ within galactic outflows due to the dense outflowing gas cooling and fragmenting. Indeed recent observations of the prototypical local obscured quasar Mrk231 have shown that the fast outflowing gas in this galaxy shows an unusually high fraction of dense gas, which is conducive to forming stars (Cicone et al. 2020)

Figure 7.

Comparison of the CO line profiles in the three galaxies in the J1234 system. The CO(1–0) line profiles have been scaled up by a factor of 10 for ease of visualization.

Open in new tabDownload slide

6 THE J2315 SYSTEM: A MAJOR MERGER

6.1 Dynamical mass and gas fraction

Our compact configuration ALMA observations for J2315 in B17 had already resolved the CO(3–2) emission in this quasar over physical scales of >20 kpc. In those observations we found that the velocity field could be satisfactorily modelled as a single (albeit very large) rotating gas disc. As can be seen in Fig. 2, the source breaks up into two distinct components in the higher spatial resolution observations presented in this paper, which clearly highlights that many high-redshift, massive galaxies that only have single component gas emission seen in low angular resolution observations, could well be merging systems. The two components are spatially separated by 15.0 ± 0.5 kpc and very close together in velocity, with the northern companion galaxy blueshifted relative to the quasar host by only 170 ± 20 km s⁻¹. We fit two 2D Gaussians to the zeroth moment CO image shown in Fig. 2, finding major axis FWHM of 1.0 ± 0.2 and 1.4 ± 0.2 arcsec for the quasar host galaxy and the northern companion, respectively, after deconvolving from the beam. This corresponds to half-light radii of ∼4–6 kpc at the redshifts of these galaxies, which is more than 2× the source sizes inferred for all three galaxies in the J1234 system based on the 2D Gaussian fitting (Table 3).

As was done for the galaxies in the J1234 system, we construct intensity, velocity, and dispersion maps for both galaxies using the casa task immoments. These maps can be seen in Fig. 8. The two galaxies of the merger are clearly resolved in the maps with comparable CO line intensities in each merging component. No velocity structure is evident in the individual galaxies at the current resolution, although we clearly see that the northern companion galaxy is blueshifted relative to the quasar host.

Figure 8.

Line intensity (top), velocity (middle), and line dispersion (bottom) maps for ULASJ2315. The source is clearly resolved into a major merger of two gas-rich galaxies with comparable line intensities. The northern companion galaxy is blueshifted relative to the quasar host. The ^3DBarolo best-fitting dynamical models and residuals shown for completeness in the middle and right-most panels do not provide a good description of these data as expected for a merging system.

Open in new tabDownload slide

For completeness, we also show the best-fitting ^3DBarolo models fit to the J2315 system as well as the fit residuals in the middle and right-most panels of Fig. 8. The fits are poor as expected as the system cannot be well modelled as a rotating disc and we do not consider the ^3DBarolo fitting results in detail in the paper as we do not expect them to be robust. A simple estimate of the total mass of the system can instead be obtained using the virial theorem assuming a projected distance between the galaxies of 15 kpc and a velocity separation of 170 km s⁻¹. This leads to M_tot = 1.0 × 10¹¹ M_⊙ for both galaxies assuming that the galaxies are moving along the line of sight.

Alternatively, we can assume that both the galaxies in the J2315 system are individual rotating discs and use the FWHM of the CO line profiles shown in Fig. 2 together with the rotating disc estimator (Neri et al. 2003) to calculate dynamical masses for each galaxy. Inclination angles of 52° and 61° are estimated for the quasar host galaxy and the northern companion, respectively, by considering the inverse cosine of the minor-to-major axial ratio derived from the CO image. Correcting for these inclination values we find that the quasar host galaxy has M_dyn = 2.3 × 10¹⁰ M_⊙ and the northern companion has M_dyn = 1.0 × 10¹⁰ M_⊙, corresponding to a total mass of 3.3 × 10¹⁰ M_⊙. Thus, the two galaxies appear to be undergoing a major merger with a ∼1:2 mass ratio.

In Paper I, we used various different estimators to derive gas masses for this system. The total gas mass is ∼3–4 × 10¹⁰ M_⊙. Combining this with the dynamical mass estimates above and estimating the stellar mass as M_* = M_dyn − M_gas −M_BH as was done in the case of J1234, we find very high gas fractions of >45 per cent. The lower limit assumes that the merging galaxies are moving along the line of sight and the gas fraction could be substantially higher if the inclination angle is sufficiently different from face-on.

6.2 CO line profiles

In Fig. 9, we compare the CO line profiles for the various different CO rotational transitions detected in the J2315 galaxies. The CO(1–0) and CO(7–6) profiles are integrated over the entire merging system as the individual components of the merger are not resolved in these data. We do however show separate CO(3–2) profiles extracted for the quasar host galaxy and the northern companion galaxy. These profiles are considerably less well matched compared to the J1234 quasar host galaxy (Fig. 7) and the kinematics of the gas in this major-merger system is clearly messier as already indicated by the poor ^3DBarolo fits. The CO(7–6) emission is significantly broader than the CO(3–2) emission with an FWHM of twice that of the CO(3–2) line (see Table 3 in Paper I). However, as can be seen in Fig. 9, the broadening almost certainly arises because there isn’t sufficient angular resolution in the ALMA Band 6 data to disentangle the companion galaxy from the quasar host. The separate CO(3–2) line profiles extracted for the quasar and the companion galaxy, straddle the broad CO(7–6) emission seen in the total system. As already noted in Paper I, the CO(1–0) emission is significantly redshifted by ∼500 km s⁻¹ relative to the other CO lines. Thus, we have clear evidence for different gas dynamics affecting the cold/diffuse and warm/dense gas components in this major merger. In Paper I, we also saw tentative evidence for the CO(1 0) emission originating from closer to the northern companion galaxy rather than the quasar, which could explain the velocity differences observed. However, it is interesting to note that the CO(3–2) emission in the companion galaxy is actually blueshifted rather than redshifted relative to the CO(3 –2) in the quasar host. Higher angular resolution observations of all the CO lines are necessary to be able to locate co-spatial components of the gas emission and therefore disentangle the origin of the complicated velocity shifts observed in Fig. 9.

Figure 9.

Comparison of the CO line profiles in J2315. The CO(1–0) and CO(7–6) line profiles are derived over the entire merging system as the individual components of the merger are not resolved in these data. Separate CO(3–2) line profiles for the QSO host and the northern companion galaxy are however shown. The CO(1–0) line profile has been scaled up by a factor of 5 for ease of visualization.

Open in new tabDownload slide

7 MULTIWAVELENGTH PROPERTIES OF THE HRQS

In this section, we briefly review existing observational constraints on the multiwavelength properties of the two HRQs and combine the information provided by these multiwavelength data with the new constraints on the properties of the quasars from this work, in order to place the quasars and their host galaxies in context within an overall picture of galaxy black hole co-evolution. Both sources were selected as red, infrared luminous quasar candidates with broad emission lines present in the near-infrared spectra (Banerji et al. 2012, 2015). X-ray observations (Lansbury et al. 2020) place constraints on the accretion properties of the black hole and the line-of-sight gas column densities. Both quasars have similarly high X-ray luminosities – log₁₀|$(L_{2{\!-\!}10\, \rm {keV}})\gt 45$| – and Eddington ratios of >0.1, indicating that we are seeing the black holes during a phase of rapid accretion. The gas column densities are also significant – ∼10²² cm⁻², consistent with the quasars being in an obscured phase. In terms of the properties of the central engine, both HRQs are therefore very similar.

Deep optical photometry from current wide-field imaging surveys allows us to trace the rest-frame UV emission from the HRQs. As the quasar is highly obscured, in the majority of HRQs the rest-frame UV emission likely arises from unobscured star-forming regions in the quasar host galaxy (Wethers et al. 2018). J1234 and its companion galaxies are undetected in the i band down to a 5σ magnitude limit of ∼25.1 (Banerji et al. 2014). Far-infrared to millimetre observations (Banerji et al. 2014, 2018) however suggest significant levels of obscured star formation in the host galaxy. J2315 is detected in data from the Dark Energy Survey with a measured i-band magnitude of 22.5 (Wethers et al. 2018). The source is morphologically extended in the Dark Energy Survey imaging confirming that the UV emission arises from unobscured regions within the quasar host galaxy rather than the quasar itself. The UV-derived star formation rate for the galaxy is estimated to be ∼100–200 M_⊙ yr⁻¹ (Wethers et al. 2018). Recent, deeper and higher spatial resolution observations from the HyperSuprimeCam Survey Data Release 2 (Aihara et al. 2019) clearly shows resolved emission spatially co-incident with the location of the quasar (Fig. 10; left-hand panel) as well as more extended emission in the direction of the northern companion galaxy. The rest-frame UV observations of the quasar host galaxies therefore indicate very different levels of obscuration towards the star-forming regions despite similar levels of obscuration towards the black holes, with the host galaxy of J1234 being completely obscured in the UV whereas J2315 shows evidence for unobscured star-forming regions, albeit the main bulk of star formation activity in this galaxy is also likely obscured (Banerji et al. 2018; Wethers et al. 2020).

Figure 10.

Left: A 10 × 10 arcsec multicolour grz band composite image of the J2315 system from HSC Data Release 2. North is up and East is to the right. The white contours denote the CO(3–2) emission from this paper and the yellow dot–dashed contours denote the JVLA CO(1–0) emission from Paper I. The red star denotes the position of the quasar from near-infrared imaging. Morphologically ‘fuzzy’ blue emission is visible both at the location of the quasar and extending to the North co-incident with the location of one of the CO(1–0) peaks. Right: A portion of the XShooter H-band spectrum (paper in preparation) centred on the position of the quasar marked with the red star on the left. The red vertical lines mark the expected positions of H β and the [O iii] doublet based on the ALMA redshift from this paper. The [O iii] 5008.239 emission is clearly blueshifted relative to systemic signifying the presence of AGN-driven outflows.

Open in new tabDownload slide

In Paper I, we presented constraints on the ISM properties of the two HRQs using detections of multiple molecular emission lines of CO and C i. The excitation conditions in J1234 are consistent with those seen in other luminous AGNs with the CO spectral line energy distribution (SLED) peaking at high J values. Photodissociation region modelling suggests high radiation field strengths in the ISM of the galaxy. J2315 on the other hand appears to be dominated by lower excitation gas. The CO SLED peaks around J_upper = 3 and the radiation field strength is ∼1.5–2 dex lower than that inferred for J1234. The ISM conditions in J2315 are therefore more reminiscent of star-forming galaxies rather than AGNs.

Temple et al. (2019) study narrow-line region outflows in HRQs and find that [O iii] is undetected in J1234, which, even after correcting for the significant levels of obscuration in the galaxy, suggests that the narrow-line region has either been overionized or physically removed from the galaxy. This is also the case for ∼10 per cent of very high luminosity optically selected unobscured quasars (Coatman et al. 2019) and consistent with AGN feedback having had a significant effect on the narrow line region gas. In the case of J2315, broad blue wings are seen in the [O iii] emission line profile extending to velocities of ∼2000 km s⁻¹, providing direct evidence for high-velocity outflowing gas being powered by the luminous AGN. The [O iii] and H β region from a recent XShooter spectrum of J2315 are shown in the right-hand panel of Fig. 10 (paper in preparation).

While we should be cautious about drawing conclusions from a sample of only two quasars studied here in detail, the comprehensive multiwavelength observations that have been assembled for the two HRQs paint a consistent picture of them being seen in different phases during their evolution. J2315 is conceivably an earlier phase, where two galaxies are viewed in the process of undergoing a major merger. While the quasar itself and a significant amount of the ongoing star formation is obscured by dust, regions of unobscured star formation remain visible in the rest-frame UV with their location broadly consistent with the extended reservoir of molecular gas. Evidence for quasar-driven outflows is provided by the broad blue wings seen in the [O iii] emission line profiles. Nevertheless, the molecular gas fraction is high and the excitation conditions of the ISM are more consistent with those of star-forming galaxies rather than luminous quasars, which might suggest that the effects of AGN feedback are yet to have a dramatic impact on the properties of the host galaxy.

J1234 can putatively be assumed to be seen at a more evolved stage. The quasar together with both companion galaxies are completely obscured at rest-frame UV and optical wavelengths. A significant starburst appears to be occurring in the quasar host galaxy based on the ALMA dust continuum detections (Paper I), which may or may not have been triggered by a merger. While the velocity fields of the CO(3–2) molecular gas observed in this work are consistent with rotational dynamics, the relatively modest resolution of the observations do allow the possibility of a merging system seen close to or soon after final coalescence when the gas is perhaps settling into a rotation pattern. Unlike in the case of J2315, the star formation is entirely obscured in the rest-frame UV. No [O iii] emission is detected in the galaxy, suggesting that powerful AGN-driven winds have already overionized or cleared the galaxy of its [O iii] gas. Further evidence for efficient feedback is provided by the low molecular gas fraction observed in this system. The ISM properties are very similar to those observed in other high-luminosity quasars with highly excited warm and dense phases of the molecular gas being dominant.

8 THE M_BH–M_GAL RELATION

Our spatially resolved observations of the HRQ host galaxies have allowed us to estimate their dynamical masses. Such dynamical mass estimates allow luminous quasars to be placed on the M_BH–M_gal relation (Willott, Bergeron & Omont 2015, 2017; Shao et al. 2017). The majority of optically selected quasars with M_BH ≳ 10⁹ M_⊙ have black holes that are significantly overmassive relative to their host galaxies, and compared to the local relation (Kormendy & Ho 2013). While this is largely a selection effect of optical surveys targeting the most massive black holes at the rapidly declining end of the mass function (e.g. Lauer et al. 2007), it is also symptomatic of a much larger scatter in the M_BH–M_gal relation at high redshifts relative to the local Universe. On the other hand, the population of SCUBA selected submillimetre galaxies at z ∼ 2 were inferred by Alexander et al. (2008) to have black holes that are undermassive relative to the local relation. Dust-obscured quasars are a plausible transition population between the starburst and optical quasar phases. The emergence of wide-area surveys at mid- and far-infrared wavelengths such as the WISE All-Sky Survey and the Herschel Extragalactic Surveys in the last decade, has allowed a wider population of obscured galaxies and AGNs to be assembled, with higher dust temperatures than the initial SCUBA detected galaxies. At least some of these also have black holes that are overmassive relative to their hosts (e.g. Assef et al. 2015; Matsuoka et al. 2018) although the number of obscured galaxies with both host and black hole mass measurements remains limited. The location of our HRQs on the M_BH–M_gal plane could therefore place interesting constraints on the relative time-scales for black hole and galaxy growth as these massive galaxies are being assembled.

Figure 11.

Black hole mass versus galaxy mass. In the case of observational data, the galaxy masses all refer to dynamical masses inferred using resolved observations of molecular emission lines arising from the ISM in the quasar host galaxy. In the case of our HRQs, the uncertainties on the black hole masses are shown to be 0.5 dex, which encodes the level of systematic uncertainty expected in these masses due to various factors. The resolved dynamical mass estimates for the z ∼ 4–6 quasars are taken from Pensabene et al. (2020) and the black hole masses come from De Rosa et al. (2011), De Rosa et al. (2014), Trakhtenbrot et al. (2011), and Venemans et al. (2015). The grey points show the total baryonic mass contained within 2× the stellar half-mass radius for Illustris galaxies (Sijacki et al. 2015) at z ∼ 2. The dashed line denotes the Kormendy & Ho (2013) relation.

Open in new tabDownload slide

In Fig. 11, we also show the total baryonic mass contained within 2× the stellar half-mass radius, versus black hole mass for z ∼ 2 galaxies in the Illustris simulations (Sijacki et al. 2015), which is described in more detail in the next section. It is interesting to note that J1234, which is the more ‘evolved’ of our two systems appears to be broadly consistent with an extrapolation of the Illustris galaxies to the highest masses, whereas J2315, the major merger, lies further away from the relation suggesting some further evolution of the total galaxy mass is necessary in this system in order to bring it on to the M_BH–M_gal relation predicted at these redshifts.

8.1 Comparison to Illustris

It is interesting to use the outputs of the Illustris simulations (Vogelsberger et al. 2014; Nelson et al. 2015; Sijacki et al. 2015) to consider both the progenitors of our z ∼ 2.5 quasars as well as the kinds of systems they would evolve into at the present day. The Illustris-1 simulations used here simulate a total comoving volume of (106.5 cMpc)³ and self-consistently evolve five different types of resolution elements from z = 127 to the present day. The public data release paper (Nelson et al. 2015) contains more details regarding the simulation outputs that we make use of. At a z ∼ 2.5, corresponding to the redshift of our quasars, the simulated volume contains 22 079 black holes, 345 of which have M_BH > 10⁸ M_⊙ and 29 of which have M_BH > 10⁹ M_⊙. While the simulation volume is not large enough to contain >10⁸ M_⊙ black holes at the highest redshifts of z ∼ 6, the simulations can nevertheless be used to track back the z ∼ 2.5 population to understand the progenitors of these massive black holes at higher redshifts.

Figure 12.

The red points and lines track the evolution of galaxy and black hole mass in the Illustris simulation for the most massive black holes at z = 2.5 tracked forward to z = 0.5 (left) and for the most massive black holes at z = 6 tracked forward to z = 2.5 (right). The background set of grey points trace the z = 0.5 (left) and z = 2.5 (right) galaxy populations in Illustris.

Open in new tabDownload slide

Given the location of the z ∼ 4−6 quasars in Fig. 11 relative to our z ∼ 2.5 quasars, it is also interesting to ask whether the evolution from z ∼ 6 to z ∼ 2.5 also occurs in a predominantly horizontal direction in the M_BH–M_gal plane. As can be seen in the right-hand panel of Fig. 12, this is not the case. The Illustris volume is too small to contain a typical high-luminosity z ∼ 6 quasar, but the most massive black holes in Illustris at z ∼ 6 with M_BH ∼ 10⁷ M_⊙ could putatively evolve to M_BH ∼ 10⁹ M_⊙ quasars at z ∼ 2.5. The evolution at these epochs happens in a diagonal direction, i.e. between z ∼ 6 and z ∼ 2.5 black holes tend to grow more rapidly than their host galaxies. Selecting the 10 most massive black holes at z1 = 6 and tracking them forward to z2 = 2.5 results in a mean slope in the M_BH–M_gal plane of 1.88 ± 0.4, considerably steeper than the slope of the evolution from z = 2.5 to z = 0.5. The highest redshift quasars shown in Fig. 11 are therefore unlikely to be the progenitors of our z ∼ 2.5 quasar population. Indeed, in Illustris the z ∼ 6 progenitors of the most massive z ∼ 2.5 black holes lie ∼0.13 ± 0.18 dex below the best-fitting scaling relation, suggesting that although the z ∼ 6 quasars will continue to grow along with their host galaxies, the z ∼ 2.5 quasars likely originated from z ∼ 6 black holes which were undermassive relative to their host galaxies.

9 CONCLUSIONS

We have presented spatially resolved observations of the CO(3–2) emission from ALMA in two high-luminosity, HRQ systems seen at the peak epoch of galaxy formation at z ∼ 2.5. Both quasars have very similar luminosities, black hole masses, and dust obscurations and the new ALMA observations clearly demonstrate that there is a range in gas morphologies, gas fractions, and dynamical masses seen in these systems.

• The quasar J1234 has two massive, CO-luminous companion galaxies. The CO(3–2) emission shows a clear velocity gradient in the quasar host as well as the two companion galaxies located ∼ 170 and ∼ 90 kpc (projected distance) from the quasar. By modelling these velocity gradients using the tilted-ring model ^3DBarolo, we infer dynamical masses of ∼10¹¹ M_⊙ in all three galaxies. Dynamical masses estimated using an exponential disc model are larger in all cases, presumably because they trace mass over a much larger physical extent of ∼10 kpc.

• The gas fractions in all three galaxies in the J1234 system are somewhat lower than expected from typical star-forming galaxies at these redshifts, particularly when deriving these gas fractions using the full dynamical mass encompassed within ∼10 kpc. This is consistent with the process of AGN feedback and inefficient cooling in the most massive haloes suppressing the accretion of gas into the system.

• All three galaxies in the J1234 system are at the more massive and more compact end of the distribution of both field and cluster galaxies at z ∼ 2 on a mass–size relation.

• We compare the velocity profiles of the CO(3–2) emission to the line profiles of other CO lines detected in the J1234 system galaxies. For the galaxy G1234S, we uncover very broad wings in the CO(7–6) gas emission that are not seen in the CO(3–2) line profiles. This could indicate the presence of a molecular outflow with a much higher fraction of warm, dense gas in the outflow. The implied outflow velocities are high and lead us to speculate about the presence of a ‘hidden’, highly obscured AGN within G1234S, although this needs to be confirmed with deep, high-resolution X-ray and/or mid-infrared data.

• We resolve the quasar J2315 into a close-separation major merger with a ∼1:2 mass ratio. The two galaxies in the merger have a physical separation of ∼15 kpc and a velocity separation of ∼170 km s⁻¹ providing evidence that major mergers are associated with luminous accretion and prolific star formation in high-redshift, obscured quasar systems. We also present new rest-frame UV and optical images of the system from HyperSuprimeCam where we see unobscured UV emitting regions spatially co-incident with the extended molecular gas reservoir.

• The gas fraction in J2315 is ≳45 per cent, which is consistent with or higher than other star-forming galaxies and much higher than in the J1234 system. This could suggest that AGN feedback is yet to fully kick in and that this system is viewed at an earlier evolutionary stage compared to J1234.

• The various CO line profiles in the J2315 merging system show significant velocity shifts relative to one another suggesting that the warm/dense and cool/diffuse gas phases are undergoing different motions and are likely not co-spatial. In particular, the CO(7–6) emission is much broader than the CO(3–2) emission most likely as the two merging galaxies are not resolved in these data, and the CO(1–0) emission is redshifted by ∼500 km s⁻¹ relative to the other CO lines.

• We use the dynamical mass estimates from our resolved observations to place our HRQs on the M_BH–M_gal relation. J1234 is consistent with the predictions from the Illustris simulations for z ∼ 2 galaxies. J2315 on the other hand, has a black hole mass that is larger than would be implied by the host galaxy mass. This is consistent with J1234 representing a later phase in the evolution of these massive galaxies compared to J2315.

• We use the outputs from the Illustris simulations to understand both the progenitors and likely descendants of our z ∼ 2.5 quasars. We find that most of the mass in the black hole for such massive quasars is expected to be assembled by z ∼ 2.5 consistent with very massive black holes already being in place at these redshifts. From this epoch the galaxies mainly grow in stellar mass until the present day. We also demonstrate that the high-luminosity z ∼ 6 quasar population is unlikely to evolve into our z ∼ 2.5 quasars. Instead our massive z ∼ 2.5 quasars likely evolved from black holes that were undermassive relative to their host galaxies at z ∼ 6.

Overall these observations highlight the diversity in gas morphologies, environments and gas fractions even within a uniformly selected quasar sample with very similar black hole properties and AGN obscuration. The differences could be consistent with the two quasars being seen in different evolutionary stages along a major-merger sequence. During final approach, when the galaxies are separated by a few tens of kpc (as is the case for the J2315 major-merger system), the gas fractions are still elevated as the effects of AGN feedback are yet to fully kick in. These feedback processes likely rapidly become more efficient following final coalescence, when the AGN begins to deplete the gas supply leading to the low gas fractions we observe in J1234. It is also important to point out however that there are numerous examples of unobscured, UV-luminous quasars with similar host galaxy properties to both J1234 (e.g. Coppin et al. 2008) and J2315 (e.g. Wagg et al. 2014) suggesting that AGN continuum obscuration is not necessarily a pre-requisite for gas-rich, starburst/merger host galaxies. With a sample of only two obscured quasars analysed in this work, it is difficult to draw any firm conclusions regarding how they fit into a standard, merger-driven paradigm of massive galaxy formation. Crucially, matched control samples of unobscured, UV-luminous, high-redshift quasars with similar resolved observations of the gas emission are still missing. Assembling these large samples is costly, but nevertheless essential if we are to put the obscured quasar population in context within galaxy formation scenarios.

ACKNOWLEDGEMENTS

We thank the MNRAS referee for their comments, which have helped improve the paper. MB acknowledges support from The Royal Society via a University Research Fellowship. GCJ acknowledges ERC Advanced Grant 695671 ‘QUENCH’ and support by the Science and Technology Facilities Council (STFC). CD acknowledges support from ERC starting grant 638707 and support from the STFC. MB would like to thank Matt Auger, Chris Carilli, Paul Hewett, Debora Sijacki, and Matthew Temple for many useful discussions. This paper makes use of the following ALMA data: ADS/JAO.ALMA 2016.1.01308.S. ALMA is a partnership of ESO (representing its member states), NSF (USA), and NINS (Japan), together with NRC (Canada), NSC and ASIAA (Taiwan), and KASI (Republic of Korea), in cooperation with the Republic of Chile. The Joint ALMA Observatory is operated by ESO, AUI/NRAO, and NAOJ.

DATA AVAILABILITY

The ALMA data presented in this paper can be downloaded from the ALMA archive using the Program ID 2016.1.01308.S (PI: Manda Banerji).

REFERENCES