Observing AGN feedback with CO intensity mapping

TL;DR

This paper discusses how upcoming CO intensity mapping experiments can provide a less biased way to study quasar feedback by measuring the molecular gas content in a large population of AGN, complementing current limited observations.

Contribution

It proposes using CO intensity mapping combined with spectroscopic AGN surveys to measure average CO luminosity of quasars, enabling large-scale, less biased feedback studies.

Findings

Upcoming experiments like COMAP, CCAT-prime, and CONCERTO can detect the cross-correlation with current data.

These experiments can place upper limits or measure mean CO luminosity of quasars.

Future surveys will significantly improve measurement precision.

Abstract

Current models of galaxy formation require star formation in high-mass galaxies to be limited by poorly understood mechanisms of quasar feedback. Feedback processes can be studied by examining the molecular gas content of AGN hosts through the CO rotational ladder, but the complexity of these observations means that current data are limited to only extremely CO-bright objects. Upcoming CO intensity mapping experiments offer an opportunity for a less biased probe of quasar feedback. By correlating intensity maps with spectroscopic AGN surveys, we can obtain a measurement of the mean CO luminosity of a large population of quasars simultaneously. We show that experiments like COMAP, CCAT-prime, and CONCERTO have enough sensitivity to detect this cross-correlation if existing AGN observations are representative of the whole population, and to place interesting upper limits if they are not. Future surveys will be able to increase the precision of these measurements by orders of magnitude, allowing detailed studies of quasar properties across a wide range of cosmic history.

Tables

Table 1. Table 1: Parameters of COMAP surveys, slightly modified from those quoted in Ihle et al. ( 2019 ) . Observing times and survey areas for Phase II are quoted as the value per field times the number of target fields.

Parameter	COMAP I	COMAP II
System Temperature, $T_{\mathrm{sys}}$ (K)	40	40
Number of Detectors $N_{\det }$	19	95
Beam FWHM, ${\theta }_{\mathrm{FWHM}}$ (arcmin)	4	4
Frequency Range, ${\nu }_{\mathrm{obs}}$ (GHz)	26-34	26-34
Channel Width, $\delta \nu$ (MHz)	2.0	2.0
Observing Time, $t_{\mathrm{obs}}$ (hr)	4500	2250$\times$4
Survey Area, ${\Omega }_{\mathrm{surv}}$ (deg2)	1.7	2.25$\times$4
AGN Source	COSMOS	DESI

Table 2. Table 2: Parameters of C  ii surveys. Stage II parameters are based on Silva et al. ( 2015 ) , the others come from Chung et al. ( 2018 ) .

Parameter	CCATp	CONCERTO	TIME	StageII
${\sigma }_{\mathrm{pix}}$ (MJy/sr s1/2)	0.86	11.0	11.0	0.21
$N_{\det }$	20444The physical CCATp instrument has $N_{\det }=1004$. However, the Fabry-Perot inteferometer used for CCATp only observes a fraction of the total spectral range at any given time. The sensitivity is thus roughly equivalent to that of an effective $N_{\det }\sim 20$ (Chung et al., 2018)	3000	32	16000
${\theta }_{\mathrm{FWHM}}$ (arcsec)	46.0	22.5	22.5	30.1
${\nu }_{\mathrm{obs}}$ (GHz)	200-300	200-300	200-300	200-300
$\delta \nu$ (GHz)	2.5	1.5	1.9	0.4
$t_{\mathrm{obs}}$ (hr)	3400	1200	1000	2000
${\Omega }_{\mathrm{surv}}$ (deg2)	1.7	1.4	0.01	100
AGN Source	COSMOS	COSMOS	COSMOS	DESI

Table 3. Table 3: Target CO transitions, redshift ranges, and available AGN densities.

Transition	$z_{\min }$	$z_{\max }$	AGN/deg2
			XMM-COSMOS	DESI
(1-0)	2.4	3.4	31	27
(3-2)	0.15	0.73	99	4
(4-3)	0.54	1.3	200	44
(5-4)	0.92	1.9	175	78

Table 4. Table 4: 95% Upper Limits on CO luminosity of AGN hosts mapped by different surveys. Values are given as physical luminosities in L ⊙ subscript 𝐿 direct-product L_{\odot} and as effective CO(1-0) luminosities in observer units. Compare to L ¯ ( 1 − 0 ) ′ = 4.1 × 10 10 subscript superscript ¯ 𝐿 ′ 1 0 4.1 superscript 10 10 \overline{L}^{\prime}_{(1-0)}=4.1\times 10^{10} K km/s pc 2 from directly imaged AGN in H18.

IM Survey	AGN Survey	${\overline{L}}_Q$ (107 $L_{\odot }$)				${\overline{L}}_{\left(1-0\right)}^{\prime }$ (1010 K km/s pc2)
		1-0	3-2	4-3	5-4	1-0	3-2	4-3	5-4
COMAP I	zCOSMOS	0.18	–	–	–	3.6	–	–	–
CCATp	zCOSMOS	–	1.1	3.7	9.7	–	0.94	1.4	2.3
CONCERTO	zCOSMOS	–	0.69	2.2	5.7	–	0.60	0.80	1.3
TIME	zCOSMOS	–	10.0	31.4	82.0	–	8.7	11.5	19.4
COMAP II	DESI	0.063	–	–	–	1.2	–	–	–
Stage II	DESI	–	0.003	0.005	0.009	–	0.003	0.002	0.002

Equations

$$P_Q(k)=b_Q^2{P}_m(k)+\frac{1}{{\overline{n}}_Q},$$

$$P_{\mathrm{CO}}(k)={\overline{T}}_{\mathrm{CO}}^2{b}_{\mathrm{CO}}^2{P}_m(k)+P_{\mathrm{CO}}^{\mathrm{shot}},$$

$${\overline{T}}_{\mathrm{CO}}=C_{LT}\int L\frac{dn}{dL}dL,$$

$$C_{LT}=\frac{c^3(1+z{)}^2}{8\pi k_B{\nu }_{\mathrm{CO}}H(z)},$$

$$b_{\mathrm{CO}}=\frac{\int Lb(L)dn/dLdL}{\int Ldn/dLdL},$$

$$P_{\mathrm{CO}}^{\mathrm{shot}}=C_{LT}^2\int L^2\frac{dn}{dL}dL.$$

$$P_{\times }(k)={\overline{T}}_{\mathrm{CO}}{b}_{\mathrm{CO}}{b}_Q{P}_m(k)+C_{LT}{\overline{L}}_Q,$$

$${\overline{L}}_Q=\frac{1}{{\overline{n}}_Q}\int L{\frac{dn}{dL}|}_{Q}dL$$

$$C_{ij}(k)=\left(\begin{array}{cc}{P}_{\mathrm{CO}}(k)+P_N/W^2(k)& P_{\times }(k)\\ P_{\times }(k)& P_Q(k)\end{array}\right),$$

$$W(k)=e^{-k^2{\sigma }_{\perp }^2}{\int }_0^1{e}^{-k^2\left({\sigma }_{\parallel }^2-{\sigma }_{\perp }^2\right){\mu }^2}d\mu ,$$

$${\sigma }_{\perp }=r(z){\sigma }_{\mathrm{beam}}$$

$${\sigma }_{\parallel }=\frac{c}{H(z)}\frac{(1+z)\delta \nu }{{\nu }_{\mathrm{obs}}}$$

$$F_{ij}=\frac{1}{2}\sum _k{N}_m(k)\mathrm{Tr}\left[\frac{\partial \mathrm{C}}{\partial {\mathrm{x}}_{\mathrm{i}}}{\mathrm{C}}^{-1}\frac{\partial \mathrm{C}}{\partial {\mathrm{x}}_{\mathrm{j}}}{\mathrm{C}}^{-1}\right],$$

$$N_m=\frac{k^2\Delta k{V}_{\mathrm{surv}}}{4{\pi }^2}$$

$$P_N={\sigma }_{\mathrm{vox}}^2{V}_{\mathrm{vox}},$$

$${\sigma }_{\mathrm{vox}}=\frac{T_{\mathrm{sys}}}{\sqrt{N_{\det }\delta \nu t_{\mathrm{pix}}}},$$

$$P_N=\frac{{\sigma }_{\mathrm{pix}}^2}{t_{\mathrm{pix}}^{\mathrm{CII}}}{V}_{\mathrm{vox}}^{\mathrm{CII}}.$$

$$t_{\mathrm{pix}}^{\mathrm{CII}}=\frac{N_{\det }{t}_{\mathrm{obs}}}{{\Omega }_{\mathrm{surv}}/{\Omega }_{\mathrm{pix}}^{\mathrm{CII}}}$$

$${\sigma }_{\times }(k)=\frac{1}{\sqrt{N_m(k)}}\sqrt{P_{\times }^2(k)+\left(P_{\mathrm{CO}}(k)+\frac{P_N}{W^2(k)}\right)P_Q(k)},$$

$${\sigma }_{\times }(k)\approx \sqrt{\frac{P_{\mathrm{CO}}(k)W^2(k)+P_N}{{\overline{n}}_Q{W}^2(k)N_m(k)}}$$

$$F_{{\overline{L}}_Q{\overline{L}}_Q}=\sum _k\frac{1}{{\sigma }_x^2(k)}{\left(\frac{\partial P_{\times }}{\partial {\overline{L}}_Q}\right)}^2=\frac{1}{{\sigma }_{{\overline{L}}_Q}^2}.$$

$${\sigma }_{{\overline{L}}_Q}=\frac{1}{\sqrt{f_L{\overline{n}}_Q}},$$

$$f_L\equiv C_{LT}^2\sum _k\frac{W^2(k)N_m(k)}{P_{\mathrm{CO}}(k)W^2(k)+P_N}$$

$$\frac{L_{(1-0)}^{\prime }}{\mathrm{K}km/s\mathrm{p}{\mathrm{c}}^2}\equiv 2.04\times 10^5{R}_J{\left(\frac{{\nu }_J}{{\nu }_{(1-0)}}\right)}^{-3}\frac{{\overline{L}}_Q}{L_{\odot }},$$

$$f_{L^{\prime }}\equiv {\left(4.9\times 10^{-5}\frac{L_{\odot }}{\mathrm{K}\mathrm{km}/\mathrm{s}{\mathrm{pc}}^2}\right)}^2{R}_J^2{\left(\frac{{\nu }_J}{{\nu }_{(1-0)}}\right)}^6{f}_L$$

$$P_Q(k)=V⟨{\left|{\tilde{\delta }}_Q(\mathbf{k})\right|}^2⟩,$$

$${\tilde{\delta }}_Q(\mathbf{k})=\frac{1}{V}\int {\delta }_Q(\mathbf{x})e^{\mathrm{i}\mathbf{k}\cdot \mathbf{x}}{d}^3\mathbf{x}=\frac{1}{N_Q}\sum _{i=1}^{N_c}{N}_i{e}^{\mathrm{i}\mathbf{k}\cdot {\mathbf{x}}_i},$$

$$⟨{\left|{\tilde{\delta }}_Q(\mathbf{k})\right|}^2⟩=\frac{1}{N_Q^2}\sum _{i,j=1}^{N_c}⟨N_i{N}_j⟩e^{\mathrm{i}\mathbf{k}\cdot ({\mathbf{x}}_i-{\mathbf{x}}_j)}$$

$$\begin{gathered} {⟨{\left|{\tilde{\delta }}_Q(\mathbf{k})\right|}^2⟩}_{i\ne j}=\frac{1}{N_Q^2}\sum _{i\ne j}{\overline{n}}_Q^2\delta V^2\left[1+{\xi }_Q({\mathbf{x}}_j-{\mathbf{x}}_i)\right] \\ \times e^{\mathrm{i}\mathbf{k}\cdot ({\mathbf{x}}_i-{\mathbf{x}}_j)}, \end{gathered}$$

$$\begin{gathered} {⟨{\left|{\tilde{\delta }}_Q(\mathbf{k})\right|}^2⟩}_{i\ne j}=\frac{1}{V^2}\int e^{\mathrm{i}\mathbf{k}\cdot (\mathbf{x}-{\mathbf{x}}^{\prime })}{d}^3\mathbf{x}{d}^3{\mathbf{x}}^{\prime } \\ +\frac{1}{V^2}{b}_Q^2\int {\xi }_m(\mathbf{x}-{\mathbf{x}}^{\prime })e^{\mathrm{i}\mathbf{k}\cdot (\mathbf{x}-{\mathbf{x}}^{\prime })}{d}^3\mathbf{x}{d}^3{\mathbf{x}}^{\prime }. \end{gathered}$$

Full text

1 Introduction

In order for analytic simulations of galaxy formation to match the observed galactic number density over cosmic time, it is necessary to invoke a process to limit the maximal stellar mass of galaxies in the universe (Springel et al., 2005; Bower et al., 2006; Croton et al., 2006; Silk & Mamon, 2012). At the most massive end of the galaxy mass function, this process is believed to be feedback from the active galactic nuclei (AGN) (Hopkins et al., 2006; Fabian, 2012) which should shut off and prevent further star formation in the host galaxy (Sanders et al., 1988; Hopkins et al., 2006). One tracer of the star formation in galaxies is their molecular gas content. Molecular gas provides the material for star formation and is directly correlated with a galaxy’s star formation rate and its ability to turn cold gas in to stars. Thus, it is crucial to understand how AGN-driven outflows are affecting the molecular gas in their host galaxies.

Tracing molecular gas is difficult, however, because the bulk of the molecular gas in the interstellar medium (ISM) is in the form of molecular hydrogen (H2), which has no emission or absorption features at the relatively cool temperatures expected in the ISM (though quasar-driven winds might enhance the emission of warm H2, Hill & Zakamska 2014). Instead, carbon monoxide (CO) is a common and useful tracer of the molecular gas content, with several rotational transitions at the correct temperature.

Detections of individual AGN in CO remain relatively scarce and an unbiased study of the molecular gas content of AGN hosts is even more elusive. Often individual objects are targeted because of previous submillimeter detections, their intrinsic bolometric luminosity (e.g. Hill et al., 2018) or because of the presence of molecular outflows in other lines (e.g. Cicone et al., 2014). Even with the most sensitive instruments available today, sampling a large area remains difficult. In addition, observations of the molecular gas content of powerful quasars using sub-millimeter continuum and molecular line emission provide some conflicting results. Some objects appear gas-rich and strongly star-forming (Beelen et al., 2004; Coppin et al., 2008; Husemann et al., 2017; Banerji et al., 2017) and others depleted of gas (Cicone et al., 2014; Kakkad et al., 2016; Fiore et al., 2017; Carniani et al., 2017).

Upcoming line intensity mapping observations offer a possible solution to this lack of statistics. Rather than aiming to image individual known galaxies, intensity maps measure the total emission in a particular spectral line with relatively coarse angular resolution (see Kovetz et al., 2017, for a review). By observing at many closely-spaced frequency bands, intensity mapping surveys can probe structure in three dimensions with sub-percent redshift accuracy. An intensity map thus probes the combined emission of all sources within a given volume, including faint objects which would be difficult or impossible to target individually. While the technique cannot tell us the line intensity of any individual object, it is sensitive to the statistical properties of all emitting galaxies, and can be used to make many useful inferences about the nature of high-redshift galaxies (see, e.g., Breysse et al., 2016; Comaschi & Ferrara, 2016; Serra et al., 2016; Lagache, 2018).

Though intensity mapping is a relatively young field, it has attracted significant experimental investment. Intensity mapping of the 21 cm spin-flip transition in neutral hydrogen has long been seen as a useful technique for studying cosmology (Pritchard & Loeb, 2012, and references therein), with several surveys completed or in progress across a wide range of redshifts (Tingay et al., 2013; van Haarlem et al., 2013; Bandura et al., 2014; Ali et al., 2015; Xu et al., 2015; Newburgh et al., 2016; DeBoer et al., 2017), and a pair of detections in cross-correlation (Masui et al., 2013; Anderson et al., 2018). In addition to this, there has been a surge of recent interest in intensity mapping of other lines, including fine-structure lines of C ii (Gong et al., 2012; Yue et al., 2015) and other hydrogen transitions such as H$\alpha$ and Ly$\alpha$ (Doré et al., 2014; Pullen et al., 2014; Silva et al., 2018; Cooray et al., 2016) among others. Most interesting for our purposes, however are experiments targeting rotational CO transitions.

The CO intensity mapping signal was first modeled in Righi et al. (2008), who were primarily interested in it as a possible foreground for cosmic microwave background observations. It was quickly realized that CO fluctuation maps were interesting in their own right as tracers of molecular gas and star formation in high-redshift galaxies (Lidz et al., 2011; Pullen et al., 2013; Breysse et al., 2014). A tentative detection of small-scale CO(1-0) intensity fluctuations at $z\sim 3$ was reported by the CO Power Spectrum Survey (COPSS) using inteferometric observations (Keating et al., 2015; Keating et al., 2016). Currently, the CO Mapping Array Pathfinder (COMAP) experiment is underway which will have sufficient sensitivity for a high-significance detection over a similar redshift range and a much larger volume. “Phase II” of COMAP will target the same redshift at higher sensitivity (Li et al., 2016) and a planned “Phase III” will add a second frequency channel to target z=6-8 at the end of the Epoch of Reionization (EoR).

Surveys which do not have CO as their primary target can also be useful for studying molecular gas in AGN. Specifically, several higher-$J$ CO transitions between $z\sim 0-2$ fall into the frequency ranges of experiments targeting the C ii line at the EoR. These lines are typically treated as foreground contaminants for reionization science, and substantial effort has gone into studying how they can be subtracted from C ii maps (Breysse et al., 2015; Cheng et al., 2016; Lidz & Taylor, 2016; Sun et al., 2018). However, these so-called “foregrounds” carry important astrophysical information in their own right. We will consider here three C ii experiments which are planned to take data over the next few years: the Tomographic Ionized-Carbon Intensity Mapping Experiment (TIME, Crites et al. 2014), the Cerro Chajnantor Atacama Telescope-prime (CCATp, Stacey et al. 2018), and the Carbon C ii line in post-reionization epoch (CONCERTO, Lagache 2018). We will also consider a hypothetical futuristic C ii experiment which has been termed C ii Stage II (abbreviated StageII hereafter) defined in Silva et al. (2015). An additional C ii experiment, the Experiment for Cryogenic Large-aperture Intensity Mapping (EXCLAIM) has recently been announced (Padmanabhan, 2018), but they target much lower redshift than the other experiments we consider here ($z\sim 3$ for C ii rather than $z\sim 6-8$), and as such the CO lines they are sensitive to will be either higher-$J$ or lower redshift. We therefore leave detailed forecats of EXCLAIM sensitivities to future work.

As stated previously, intensity maps are sensitive to the total emission from all emitters in a surveyed volume. In order to specifically study AGN, we need a way to isolate the contribution from AGN to the total emission. This can be done using cross-correlation. Wolz et al. (2017) (W17 hereafter) showed that, on small scales, the cross-power spectrum111This same statement can be expressed mathematically equivalently in terms of a stacking analysis, or in terms of the real-space two point function. However, as the intensity mapping literature predominantly prefers to work with power spectra, we choose to follow W17 and express our forecasts in Fourier space. between a line intensity map and a galaxy survey depends on the mean line luminosity of only those sources contained in the galaxy catalog. This means that by cross-correlating a CO intensity map and a map of quasar positions we can directly target the molecular gas contained within the mapped AGN. In doing so, we will be able to simultaneously observe dozens or hundreds of AGNs without many of the selection biases that plague targeted molecular gas observations. It is this powerful technique that will allow for a systematic study of the molecular gas content of AGN hosts over large portions of cosmic time.

This paper is organized as follows: Section 2 describes our cross-correlation technique and how it can be used to forecast constraints on molecular gas in AGN hosts, then proceeds to give the parameters of the fiducial surveys we consider. Section 3 gives the results of our forecasts, and we describe the science implications of these forecasts in Section 4. Further discussion is found in Section 5, and we conclude in Section 6.

2 Method

Cosmological density fields, whether they be maps of galaxy density or line intensity, are often described in terms of their power spectra. For a map of the number density $n_{Q}(\mathbf{x})$ of AGN, we typically consider the power spectrum of the dimensionless density contrast $\delta_{Q}\equiv(n_{Q}(\mathbf{x})-\overline{n}_{Q})/\overline{n}_{Q}$. As quasars are a biased, discrete tracer of the underlying dark matter density field, we can write this power spectrum as

[TABLE]

where $b_{Q}$ is the quasar bias and $1/\overline{n}_{Q}$ is the shot noise caused by Poisson randomness in the galaxy positions. For the convenience of the reader, detailed derivations of this and the other power spectra discussed in this section are provided in Appendix A.

Qualitatively, the first term in Eq. (1) comes from the correlation between pairs of galaxies, while the shot-noise term comes from the correlation of each galaxy with itself. In the intensity mapping regime, we weight each galaxy by its intensity in a chosen line, and take the power spectrum of the resulting intensity field, which we define here in brightness temperature units $T(\mathbf{x})$. The auto-spectrum of a map of CO is given by

[TABLE]

(see, for example Lidz et al. (2011); Breysse et al. (2014)). The CO emitting galaxies have their own bias $b_{\rm{CO}}$, and the large-scale clustering term is weighted by the sky-averaged mean intensity of the line $\overline{T}_{\rm{CO}}$. There is again a scale-independent shot noise component, as the CO emission is still coming from a population of discrete emitters.

If we assume a luminosity function $dn/dL$ for the CO emitters, we can write the mean intensity as

[TABLE]

For compactness, we have defined the conversion factor between CO luminosity density and observed brightness temperature at redshift $z$ as

[TABLE]

where $c$ is the speed of light, $k_{B}$ is Boltzmann’s constant, $H(z)$ is the Hubble parameter, and $\nu_{\rm{CO}}$ is the rest frequency of the target CO transition. The CO bias is conventionally defined as

[TABLE]

where $b(L)$ is the bias of a galaxy with CO luminosity $L$. The shot noise again comes from the correlation of a given CO emitter with itself, and so is proportional to the mean-square CO luminosity. In full,

[TABLE]

Equations (3), (5), and (6) can be modified to be written in terms of the halo mass function $dn/dM$ by assuming a mass-luminosity relation L(M).

Finally, we can write down the cross-spectrum we would obtain by correlating these two maps:

[TABLE]

where

[TABLE]

is the mean CO luminosity of the quasars, which are assumed to have a CO luminosity function $\left.dn/dL\right|_{Q}$. The origin of the first term is clear, each pair of galaxies contributes one factor of $\overline{T}_{\rm{CO}}b_{\rm{CO}}$ and one factor of $b_{Q}$. The cross-shot power comes from galaxies which simultaneously emit CO and host AGN, and is thus proportional to the mean CO luminosity $\overline{L}_{Q}$ of these and only these galaxies. In other words, a measurement of the shot power in a CO-AGN cross-spectrum provides a measure of the mean CO luminosity of the targeted quasars. As with the prior auto-spectra, a derivation of Eq. (7) can be found in Appendix A (also see Eqs. (13) and (14) of W17).

In all of our power spectra thus far we have assumed a linear, scale-independent bias for all tracers. This approximation is valid on large scales, but deviations from this simple biasing scheme will likely be important on the scales where the power spectrum transitions from clustering to shot-noise dominated (see discussions in W17, Chung et al. 2019). Wolz et al. (2018) demonstrate that halo models can be developed which can fully account for the nonlinear power in simulated intensity maps. For simplicity here, however, we will neglect these nonlinear contributions, as errors due to this assumption are likely small compared to the order-of-magnitude uncertainties on the modeling of the CO signal. To be conservative, we will limit our forecasts here to upper limits on $\overline{L}_{Q}$, with the understanding that a full measurement would require subtracting a nonlinear component from the cross-spectrum.

We can use the Fisher matrix formalism to estimate how well a given instrument can measure $\overline{L}_{Q}$. Any given survey may not have sufficient resolution to access fully shot-noise dominated scales, so we consider possible degeneracies between $\overline{L}_{Q}$ and the amplitude of the large-scale clustering term. We can write down the covariance matrix of our two density fields as

[TABLE]

where $P_{N}$ is the noise power spectrum of the CO survey assuming a noise $\sigma_{N}$ in each voxel. The only noise in the AGN power spectrum is the shot noise, which we have included in our definition of $P_{Q}$.

The quantity $W(k)$ defines the sensitivity cutoff at the resolution limit of the telescope, given by

[TABLE]

where

[TABLE]

is the cutoff in the plane of the sky for a beam with width $\sigma_{\rm{beam}}$222$\sigma_{\rm{beam}}=\theta_{\rm{FWHM}}/\sqrt{8\ln 2}$ for beam full width at half maximum $\theta_{\rm{FWHM}}$. and a survey at comoving distance $r(z)$. The cutoff in the line-of-sight is set by the frequency resolution $\delta\nu$, and is given by

[TABLE]

for observed frequency $\nu_{\rm{obs}}$. Eqs. (10-12) are derived in Appendix C of Li et al. (2016).

For free parameters $x_{i}=\left(\overline{T}_{\rm{CO}}b_{\rm{CO}},b_{Q},\overline{L}_{Q}\right)$, the Fisher matrix is then

[TABLE]

(Fonseca et al., 2018; Tegmark et al., 1997), where the sum is over bins in $k$ space with width $\Delta k$, each of which contains

[TABLE]

modes for a survey covering volume $V_{\rm{surv}}$. Including both the auto- and cross-spectra in the Fisher matrix allows us to use information in the auto-spectra of the two tracers to break degeneracies in the cross-spectra.

2.1 Fiducial Surveys

Here we describe the available AGN maps and the parameters of our planned intensity mapping surveys.

2.1.1 Quasar Reference Samples

In order to harness the predictive power of a CO line intensity mapping survey, we require a large sample of quasars for cross-correlation. While large photometric surveys provide at least an order of magnitude more objects than spectroscopic surveys, the increased redshift uncertainty makes them unsuited for our task (see Chung et al., 2019). Modeling our result with photometric errors is beyond the scope of this paper and will be treated in future work. Instead we focus on current and upcoming deep spectroscopic surveys with sufficient numbers of quasars.

In order to obtain a sufficient number density of current spectroscopically-confirmed AGN, we rely on data from deep field spectroscopic surveys. In particular, the Cosmic Evolution Survey (COSMOS; Scoville et al., 2007) field has a wealth of spectroscopic follow-up of X-ray sources, the majority of which are AGN. The entirety of the COSMOS field has been observed by XMM-Newton for a total of 1.55 Ms (Hasinger et al., 2007; Cappelluti et al., 2007, 2009) with the final catalog containing 1887 (point-like) sources detected in either the soft (0.5-2.0 keV), hard (2-10 keV), or ultrahard (5-10 keV) bands. Various surveys (Trump et al., 2009; Lilly et al., 2009) have conducted spectroscopic follow-up of AGN candidates in this field using IMACS on Magellan (Trump et al., 2009) or VIMOS on the VLT (zCOSMOS bright and faint; Lilly et al., 2009) which together have covered about half of all AGN-candidates in COSMOS. A catalog of the spectroscopic counterparts to XMM-COSMOS sources is presented in Brusa et al. (2010) which includes a total of 890 unique spectroscopic objects. To construct our sample we include only those objects from Brusa et al. (2010) which are part of the XMM-COSMOS flux limited sample, have only one reliable counterpart with a secure redshift, and are classified by their optical spectroscopy and x-ray detections as being an AGN, not a star-forming galaxy leaving us with a total sample of 730 AGN from XMM-COSMOS.

Looking towards the future, large, deep, all-sky spectroscopic surveys will allow for a cross-correlation with much larger samples of AGN. For example, the Dark Energy Spectroscopic Instrument (DESI) survey (DESI Collaboration, 2016) will operate from the Mayall 4-meter telescope on Kitt Peak over the course of 5 years and observe 14,000 deg2 of the sky. As part of its mandate to study baryon acoustic oscillations (BAO) and redshift-space distortions, DESI will target bright emission line galaxies out to z=1.7, and quasars at higher redshifts ($2.1<z<3.5$) with a final sample including approximately 17.1 and 1.7 million objects in each sample respectively. For our purposes, the expected number density of spectroscopic quasars from DESI has been taken from the DESI Science Motivations and Requirements white paper333https://www.desi.lbl.gov/final-design-report/ Tables 2.3 & 2.7, depending on the redshift range considered. To derive the total number of quasars expected for a specific CO transition, we interpolate between the numbers provided in the DESI Science Motivations and Requirements document. While the number density of quasars is obviously not going to be as dense as a deep field, DESI will cover a large portion of the sky and therefore large-area C ii experiments will be able to conduct cross-correlations with a very large number of quasars, which will improve their sensitivity.

2.1.2 CO Surveys

The surveys we consider here can be broken into two categories. First we have Phase I and II of the COMAP survey, which target the CO(1-0) transition at $z\sim 3$. Both phases will have some sensitivity to CO(2-1) at $z\sim 7$, but as there are few spectroscopic galaxies at such high redshifts we limit ourselves to the $z\sim 3$ observations here. Second, we have the CCATp, CONCERTO, TIME, and C ii-Stage II surveys, which primarily target the C ii transition during reionization. For these experiments, which we will group together under the general name “C ii surveys”, we are interested here in higher-order CO rotation transitions which are foregrounds to the target C ii line. As the CO(1-0) and C ii experiments observe in very different frequency ranges, there are different conventions in the literature for computing $P_{N}$.

For COMAP-type surveys, we typically write

[TABLE]

where

[TABLE]

is the noise in a single survey voxel, $T_{\rm{sys}}$ is the system temperature, $N_{\rm{det}}$ is the number of detector feeds, and $t_{\rm{pix}}$ is the time spent observing each pixel (Li et al., 2016). A “pixel” here is defined to cover a solid angle $\Omega_{\rm{pix}}=\sigma^{2}_{\rm{beam}}$, so for a survey area $\Omega_{\rm{surv}}$, $t_{\rm{pix}}=t_{\rm{obs}}(\Omega_{\rm{pix}}/\Omega_{\rm{surv}})$. The volume $V_{\rm{vox}}$ of a single voxel can then be set by separating each 2D pixel into voxels with line-of-sight width set by $\delta\nu$.

For COMAP Phase I, we consider correlations with the 1.7 deg2 XMM-COSMOS AGN survey, while we consider correlating Phase II with the much larger DESI survey. Table 1 gives the parameters we assume here for COMAP Phases I and II. These are slightly modified from those given in Ihle et al. (2019). For Phase I, we reduce the observing time by $1.7/2.25$, the ratio between the full COMAP survey area and that of the COSMOS field. For Phase II, the full survey plans to target four separate fields, and since we are correlating with DESI we assume that there will be quasar data available in each of them. We also use a better frequency resolution for both phases than that quoted in Ihle et al. (2019) (K. Cleary, private communication).

For the noise power spectrum for C ii surveys, we write following Chung et al. (2018):

[TABLE]

The quantities $t_{\rm{pix}}^{\rm{CII}}$ and $V_{\rm{vox}}^{\rm{CII}}$ are defined by convention slightly differently than their equivalents in the COMAP band, as the solid angle spanned by a single pixel is set to be $\Omega_{\rm{pix}}^{\rm{CII}}=2\pi\sigma_{\rm{beam}}^{2}$. For CCATp, CONCERTO, and TIME we use the $\sigma_{\rm{pix}}$ values directly from Table 1 of Chung et al. (2018). For Stage II we take the Noise Equivalent Flux Density from Table 8 of Silva et al. (2015) and convert it to $\sigma_{\rm{pix}}$ by dividing by $\Omega_{\rm{beam}}$. We can combine these with $V_{\rm{vox}}^{\rm{CII}}$ computed using the modified $\Omega_{\rm{pix}}$ and the frequency bandwidth, and

[TABLE]

to get the noise power.

Table 2 shows the parameters we assume for the C ii-type surveys. Though each has a somewhat different true frequency coverage, we consider a 200-300 GHz frequency band with each instrument for ease of comparison between them. Stage II parameters are based on Silva et al. (2015), the others come from Chung et al. (2018). As with COMAP , we assume that CCATp can only observe AGN in the 1.7 deg2 COSMOS field, so we reduce its survey area and observing time accordingly. TIME and CONCERTO plan to observe smaller fields and thus require no modification. Stage II is a futuristic experiment, so we assume it can cross-correlate with DESI, and thus will see AGN throughout its 100 deg2 footprint.

2.2 CO Emission Models

The suite of experiments we consider here are sensitive to several different CO rotational transitions at a variety of redshifts. Here we summarize the literature models we use to forecast the non-AGN CO emission. The COMAP surveys primarily target the CO(1-0) transition, while the C ii experiments have the CO(3-2) and higher transitions as “foregrounds” to their reionization-era C ii signal. The COMAP surveys are also in principle sensitive to CO(2-1) from the EoR, but we will not consider that line here as large enough spectroscopic quasar samples do not exist at those redshifts. Table 3 shows the transitions we consider and the redshift ranges over which we can access them, along with the number density of quasars accessible in the XMM-COSMOS and DESI surveys. Note that though DESI has a lower density of quasars, it spans a much larger area of the sky, so the total number of quasars available for correlation is much larger than in COSMOS.

In order to maintain consistency with past predictions in the literature, we model the non-AGN luminosities of the CO(1-0) line and the higher frequency lines somewhat differently. For CO(1-0), we use the model from Li et al. (2016), which assumes a power law relation between the luminosity of a halo and its star formation rate (SFR), and uses the SFR($M,z$) predicted in Behroozi et al. (2013) to get a mass-luminosity relation. This model also assumes some scatter between the mean $L_{\rm{CO}}(M)$ function and the luminosity of any given halo. For CO(3-2) and above, we use the model from Silva et al. (2015), which provides in their Table 4 a direct parameterization of $L_{\rm{CO}}(M)$ for each line. We neglect the evolution in CO luminosity with redshift, and simply compute power spectra at the center of each frequency band. For all of these models, we use the Tinker mass function and bias (Tinker et al., 2008; Tinker et al., 2010), and we compute matter power spectra using CAMB (Lewis et al., 2000).

Figure 1 shows cross-spectra between CO and AGN for different values of $\overline{L}_{Q}$. The left panel shows spectra for CO(1-0) observed with COMAP, the right CO(3-2) observations with CCATp. When forecasting the detectability of $\overline{L}_{Q}$ with these experiments, we will take $\overline{L}_{Q}=0$. Forecasted errors on $\overline{L}_{Q}$ will in this case give the upper limits a given experiment could place on $\overline{L}_{Q}$ without a detection. This is the most conservative assumption we can make given our lack of knowledge of the true $\overline{L}_{Q}$ and the fact that we are not modeling the nonlinear clustering power.

2.3 Scaling with Survey Parameters

To see how sensitivity to $\overline{L}_{Q}$ scales with instrument parameters, we now make some approximations to examine the output of Eq. (13) more closely. Consider a case where we use the cross-spectrum alone and neglect uncertainties on $\overline{T}_{\rm{CO}}b_{\rm{CO}}$ and $b_{Q}$. This is equivalent to assuming that the auto-spectra provide sufficient constraints on these quantities enough to avoid cross-spectrum degeneracies. This assumption turns out to be valid at the  few percent level for the surveys we consider here, though care must be taken when applying it in general.

Under this assumption, the error on the cross-power can be written as

[TABLE]

(Villaescusa-Navarro et al., 2015, W17). Eq. (19) contains contributions from the sample variance of the cross-spectrum and both auto-spectra, along with the shot noise in the CO and AGN maps (implicit in our definitions of $P_{\rm{CO}}$ and $P_{Q}$), and the instrumental noise in the intensity map. The black lines in Figure 1 show the different contributions to this error, with the full $\sigma_{\times}$ shown as a solid curve, and errors neglecting instrumental noise and/or quasar shot noise shown as broken curves.

It can clearly be seen from Figure 1 that the dominant source of error in both of these measurements is the shot noise from the limited number of AGN. We can simplify Eq. (19) by taking only terms proportional to $1/\overline{n}_{Q}$, yielding

[TABLE]

As we are neglecting errors on the clustering terms, we only need the Fisher matrix element for $\overline{L}_{Q}$ to compute our upper limit:

[TABLE]

Combining Eqs. (20) and (21) with our definitions of the power spectra yields

[TABLE]

where we have defined

[TABLE]

This $f_{L}$ quantity we have defined encapsulates all of the relevant sensitivity information in the CO survey into a single quantity, including sample variance, shot noise, instrument noise, and resolution effects. The information in the AGN map is in turn contained entirely in the value of $\overline{n}_{Q}$. Eq. (22) now provides a simple way to compare the power of different surveys. Despite all of our simplifying assumptions, it remains accurate at the $\sim$few percent level in comparison to forecasts using the full Fisher matrix from Eq. (13).

3 Results

Table 4 shows our forecasted 95% upper limits on $\overline{L}_{Q}$ for the surveys we describe above. In addition to physical $L_{\odot}$ units, we also quote our forecasts in observer-based (K km/s pc2) units. To make the observer-unit values more readily comparable between different lines, we follow e.g. Hill et al. (2018) (hereafter referred to as H18) and use the CO line ratios measured for QSO’s in Carilli & Walter (2013) to convert our upper limits to an effective CO(1-0) luminosity $L^{\prime}_{1-0}$. This effective luminosity is given by

[TABLE]

(Li et al., 2016), where $\nu_{J}$ is the frequency of the target transition and $R_{J}$ is the $L`_{(1-0)}/L`_{J}$ ratio. We can modify our $f_{L}$ parameter from Eq. (23) to work in observer units by writing

[TABLE]

To judge the usefulness of our forecasted measurements, we can compare the upper limits in Table 4 to values obtained through direct observations of CO in quasars. For example, we can compare our results to the 13 QSOs observed in CO(3-2) by H18. As currently observations of individual quasars in molecular gas can be very expensive, many objects chosen for observation have been previously detected in the sub-mm or are in some way extreme. This holds for the objects observed in H18 which are considered to be rare, hyper-luminous and UV-bright QSOs at $2.5<z<2.9$ which, additionally, were all previously detected in the sub-mm with SCUBA-2. Thus, we might expect the actual value of $\overline{L}_{Q}$ to be substantially lower than this value.

H18 find an average CO(3-2) luminosity of $\overline{L}^{\prime}_{(1-0)}=4.1\times 10^{10}$ K km/s pc2. As shown in H18 figure 5, using the average value of $\overline{L}^{\prime}_{(1-0)}$ from H18 is a good indicator of the typical value currently observed in high-z ($\gtrsim 1.4$) quasar samples. The observed values for low-z AGN ($\lesssim 0.4$) could be an order of magnitude lower (Scoville et al., 2003). Nevertheless, we see little redshift evolution in our cross-correlation sensitivity and most of our results probe samples at $z>0.4$ so we find H18 to be an appropriate comparison. Thus we can use this value of $\overline{L}^{\prime}_{(1-0)}$ as a benchmark for comparing the quality of the surveys we consider here. In other words, if all quasars had CO luminosities similar to the H18 sample, we would expect to detect cross-shot power at reasonable significance in CCATp, CONCERTO, and COMAP II, along with an extremely precise measurement in C ii StageII. The COMAP Phase I survey in our forecast has enough sensitivity to marginally detect the H18 quasars, and the TIME survey appears to be not well-optimized for this type of observation.

To illustrate more clearly the dependence of our forecast on survey design, Figure 2 shows the upper limit a survey could obtain on $\overline{L}^{\prime}_{(1-0)}$ as a function of $\overline{n}_{Q}$ and the $f_{L}$ parameter defined in Eq. (23). Though this only approximates the values used in Table 4, as we stated previously the approximation should be valid at the $\sim$few percent level. Plotting the forecasts in this way allows us to completely separate out the dependence on the AGN survey and the CO survey.

Of the current-generation surveys, we find that CCATp and CONCERTO can make the most sensitive measurements of $\overline{L}_{Q}$. CONCERTO does slightly better than CCATp, largely because it has higher frequency resolution, making it more sensitive to the small spatial scales we care about here. The TIME survey does not do as well as the other C ii experiments as it covers a much smaller volume, limiting it to a small number of Fourier modes and relatively few quasars. COMAP Phase I, because it observes at higher redshift than the C ii experiments, also has fewer AGN to work with and less sensitivity to small scales. The latter is compounded by the fact that COMAP targets the lowest-frequency CO line. However, the added sensitivity of COMAP Phase II/DESI is enough to significantly improve the Phase I constraints. The C ii Stage II experiment is by far the most ambitious experiment we consider here. Because its high sensitivity allows it to map an order of magnitude more sky than the other experiments, it has a much lower $f_{L}$ and can measure $\overline{L}_{Q}$ at extremely high precision.

The forecasts for the three current-generation C ii experiments evolve similarly with the targeted CO transition. In the COSMOS field, the total number of AGN per sky area is highest in the CO(4-3) band, while the total comoving volume of the surveys decreases monotonically to higher $J$’s. This causes the comoving quasar density to decrease for higher-order lines, and worsens the forecasted limits. The $f_{L}$ parameter for these surveys also gets slightly worse with increasing $J$, though at slightly different rates for the three surveys as they have very different parameters. For the Stage II survey, the constraints improve significantly between CO(4-3) and CO(5-4). This is because, unlike in COSMOS, the number of AGN per square degree actually increases from CO(4-3) to CO(5-4). The comoving number density remains roughly constant between these two lines, leading to the observed behavior in Figure 2.

4 Science Applications

The CO cross-power spectrum has strong implications for the study of molecular gas in quasars and we envision several possible science applications. Of course, our ability to extract scientifically-interesting information from this sample is influenced by our ability to convert from $\overline{L}_{\rm{Q}}$ to a quantity of scientific interest, specifically the gas mass, or $M_{\rm{gas}}$. There is disagreement over estimates for this conversion factor, $\alpha_{\rm{CO}}$ for AGN hosts (Bolatto et al., 2013) with some papers adopting the Milky Way value, $\alpha_{CO(1-0)}\approx 4$M*⨀(K km s-1* pc2)-1 (Feruglio et al., 2010; Herrera-Camus et al., 2019) while others use the value typical of ULIRGs $\alpha_{CO(1-0)}\approx 0.8-1.0$M*⨀(K km s-1* pc2)-1 (Cicone et al., 2014; Fiore et al., 2017; Veilleux et al., 2017; Hill et al., 2018; Fluetsch et al., 2019). Measurements to confirm this value for a relevant sample of AGN, will be necessary for proper interpretation of our results.

If the mean CO luminosity of DESI quasars is around the H18 level, then under our assumptions a survey like C ii StageII will provide sub-percent level measurements of the mean CO luminosity in AGN hosts. At that point, our understanding of the details of AGN feedback would likely be limited by the precision of our models, as detailed computation of the physics which determines CO emission is a highly complicated task (see, e.g. Popping et al. 2019). However, this level of precision (and the excellent redshift resolution offered by intensity mapping) means that one can break the full frequency bandwidth of the map into several shells. This will sacrifice sensitivity in each band, but could pay substantial dividends by allowing one to map the detailed redshift evolution of AGN feedback processes and search for changes in the average gas mass of AGN hosts as a function of redshift.

Several of the proposed C ii experiments will, in fact, give access to several different CO transitions, and thus allow us to explore the CO spectral line energy distribution (SLED) for quasars (assuming a proper accounting for the redshift evolution of the average CO luminosity of quasars as mentioned above). There is emerging evidence to suggest that quasars may have more highly excited molecular gas than submillimter galaxies (SMGs) or main sequence galaxies. Riechers et al. (2011) find evidence for a single, highly excited gas component in a sample of five redshift $z\sim 2.5$ lensed quasars. Evidence suggests that all of the molecular gas is associated with star formation. In addition, it is possible that either wind- or jet-mediated quasar feedback would introduce a second more excited component of molecular gas. For example, Dasyra et al. (2016) found that the ratio of CO(4-3)/CO(2-1) in IC 5063 was $\approx 1$ for the disk but closer to 5 in regions near the AGN-driven jet (necessitating that the gas be optically thin). Our technique would allow us to probe the average quasar SLED which will facilitate the comparison of future observations made at different CO transitions and reveal details of the excitation state of the molecular gas in AGN hosts.

While this paper explores the possible result of cross correlating CO intensity mapping experiments with AGN samples, this technique is applicable for any number of galaxy samples. Leaving aside for now the science questions to be answered by different galaxy samples, it will also be possible to compare the results of the CO cross-power spectrum for quasar and other galaxy samples. This is where we see the most power to shed insight on the question of quasar feedback. Comparing the CO cross-power spectrum of AGN hosts to a mass-matched sample of star-forming and quenched galaxies will allow us to search for a difference in the mean CO luminosity. As the samples are independently-selected on small scales where the shot noise term dominates the power spectrum, these results will be dominated by the measurement error on the smaller quantity, expected to be the quasar sample due to their relative lack of abundance. Thus, the measurement errors on our ability to calculate the difference in $\overline{L}_{\rm{CO}}$ between a sample of quasars and a complimentary sample of galaxies will be similar to those expressed in Table 4.

If AGN are in fact responsible for clearing out the molecular gas from their hosts, the depletion timescale is expected to be on the order of a few hundred Myr (Cicone et al., 2014; Fiore et al., 2017). Thus, we would expect the average gas mass in the AGN hosts to be comparable to that of a mass-matched sample of passively evolving galaxies and significantly gas-deficient compared to a mass-matched sample of star-forming galaxies (alternatively, one could compare to a population of galaxies selected to have the same large-scale linear bias as the AGN sample). Such results have not been confirmed with current AGN samples (see, for example Shangguan et al., 2018) though there does appear to be a dearth of molecular gas in quasars as compared to star-forming galaxies as a function of stellar mass if accounting for the redshift evolution of such trends for star-forming galaxies (Genzel et al., 2015; Fiore et al., 2017). Alternatively, some observations suggest that swept-up gas from AGN outflows may produce additional molecular material (Richings & Faucher-Giguère, 2018) in which case quasar samples may have higher gas masses than any comparison sample. Our calculations show that upcoming intensity mapping surveys will allow us to distinguish between these two scenarios.

Finally, the search for feedback from AGN also rests on the presence of molecular outflows as measured by an excess of blueshifted emission in the CO transitions of AGN hosts (Cicone et al., 2014; Veilleux et al., 2017). In an instrument with sufficient frequency resolution, outflows would create an extra source of small-scale anisotropy in our power spectra. This signature comes from the fact that the CO emission in quasars would be extended in frequency space while remaining compact in the plane of the sky. Searches for quasar outflows tend to look for gas with outflow speeds $\gtrsim 500$ km/s (Cicone et al., 2014). The frequency resolution of COMAP or a Stage II C ii experiment would be sufficient to resolve gas at such speeds and thus we would expect to be able to detect such a signature in CO-AGN cross correlations. A determination and detectability forecast of the average quasar outflow speed is beyond the scope of this paper but will be discussed in future work. As with the other probes we can compare cross-correlations of CO with AGN with correlations between CO and star forming galaxies to search for specifically quasar-driven outflows.

5 Discussion

We described here the potential of correlating CO and AGN surveys to study AGN feedback. This is far from the only useful science which can be gleaned from cross-correlations using the surveys we discuss here. We discussed correlating star forming and quiescent galaxies in COSMOS and DESI with our CO intensity maps as a means of comparing molecular gas properties between galaxies which do and do not host AGN. However, these additional correlations are interesting in their own right, and could be used to study how molecular gas mass scales with a variety of different galaxy properties in a similar manner to the work done in W17. One could also cross-correlate with other surveys. For example, the correlation between COMAP and the Hobby-Eberly Telescope Dark Energy Experiment (HETDEX, Hill et al. 2008) discussed in Chung et al. (2019) can be used to study the CO emission from Lyman $\alpha$ emitters. Beyond that, the cross-spectra we consider here only use a small fraction of the vast DESI survey area. All-sky intensity maps like those produced by the upcoming Spectro-Photometer for the History of the Universe, Epoch of Reionization, and Ices Explorer satellite (SPHEREx, Doré et al. 2014) will be able to correlate with DESI to produce results integrating over orders of magnitude more galaxies than we consider here, opening an unprecedented window into the high-redshift universe.

One way this cross-correlation method can potentially be improved would be if one could make effective use of photometric galaxy surveys. Typical number densities of AGN with photo-$z$ measurements are around an order of magnitude higher than spectroscopic quasars. The exception here may be COSMOS due to the extensive spectroscopic follow-up of AGN candidates. Only about 50% of photometrically-classified AGN have not been confirmed with follow-up spectroscopy in this field. On the other hand, the upcoming Large Synoptic Survey Telescope (LSST; Ivezić et al., 2008) will identify approximately $10^{7}$ AGN (LSST Science Collaboration et al., 2009), with a significant fraction at $z>2.4$ but will have only photometric redshifts for the vast majority of these objects. As quasar shot noise is the dominant source of error in our forecasts, using photo-$z$ AGN would naïvely seem to dramatically improve the power of our forecasts, as can be seen in Figure 2. However, Chung et al. (2019) demonstrated that decreasing the redshift accuracy of a galaxy survey to photometric levels significantly reduces the strength of the correlation, implying that some of this gain in sensitivity would be lost. Beyond raw statistical significance, this also would add a potentially difficult systematic to the measurement, as one would need to know the true redshift distribution of the photo-$z$ AGNs to reconstruct what the true cross-shot noise would be absent this signal loss. We leave for future work the question of whether the increased quasar number counts are worth dealing with these challenges.

We consider COSMOS here as a cross-correlation field for COMAP Phase I, but Chung et al. (2019) found that cross-correlating with HETDEX would be a better option for overall cross-correlation signal to noise. Specifically, it was found that the best signal-to-noise would come from targeting the Spitzer/HETDEX Exploratory Large-Area field (SHELA), which covers $\approx$24 degrees of the Sloan Digital Survey “Stripe 82” (Annis et al., 2014). HETDEX is designed to pick out Lyman $\alpha$ emitters between $1.9<z<3.5$ using wide-field integral field unit (IFU) spectrographs and so will undoubtedly contain a large number of AGN with spectroscopic redshifts for which a similar cross-correlation can be conducted. Estimates for the exact number density currently vary though assuming a similar detection rate from their test field (Hill et al., 2008) of one AGN per 40 arcmin with a filling factor of 1/3 implies a number density of $\sim$five AGN/deg2 in the SHELA field where the filling factor will be of order unity (Papovich et al., 2016).

We have neglected in this paper any discussion of continuum foregrounds in our CO maps. This is perhaps the most significant challenge facing 21 cm intensity mapping surveys, as foreground contamination at low frequencies is several orders of magnitude brighter than the cosmological signal (see for example Liu & Tegmark 2011). Continuum foregrounds are less of a challenge in maps of other lines, but they can be dealt with, as they are much less dominant over the signal and largely remain contained to the lowest line-of-sight Fourier modes (Keating et al., 2015). Even if residuals remain after foreground cleaning, cross-spectra like the ones we propose here are very robust against foreground contamination, as the cosmological signal will be spatially correlated with the AGN locations and the foregrounds generally will not be, as demonstrated in Switzer et al. (2019) and Cunnington et al. (2019).

6 Conclusion

In this paper, we have demonstrated that a cross-correlation between a CO intensity map and a spectroscopic quasar survey can provide valuable new insights into the process of AGN feedback. Current or near-future intensity mapping experiments such as COMAPII, CCATp, and CONCERTO will have sufficient sensitivity to study the molecular gas content of AGN hosts with far less bias than current directed observations. This will provide crucial insight in to the effect of “quasar mode” feedback on the molecular reservoirs of a uniformly-selected sample of AGN hosts, where the selection process can be done independently of knowledge of the molecular content or star formation rate of the host but solely as a function of AGN luminosity. Such a measurement will allow us to determine, with robust statistics, the effect of this feedback on the gas mass of the host galaxy and shed insight on to contradictory evidence for negative and positive feedback in AGN hosts. Comparing quasar samples with similar samples of star forming and quiescent galaxies will allow us to determine whether the molecular reservoirs of AGN hosts more closely resemble those of actively star-forming galaxies or quiescent galaxies with exhausted molecular reservoirs.

Future experiments like C ii StageII and DESI will dramatically improve on these near-term measurements, offering deep maps over large enough volumes to study the molecular gas content of AGN hosts in exquisite detail. In particular, larger samples of quasars paired with large area CO intensity mapping surveys will enable us to chart the evolution in the gas mass of quasars as a function of redshift or other relevant parameters, and to search for evidence of molecular outflows with speeds $>500$ km/s. Such measurements will allow us to directly view the impact of the AGN population on their host galaxies and search for evidence of the depletion of molecular material in their hosts due to molecular outflows. Molecular outflows specifically can be a “smoking gun” signal of AGN feedback impacting the host galaxy, making future intensity mapping experiments with higher frequency resolution highly desirable for cross-correlation.

We have highlighted here one of many exciting applications for cross-correlations between galaxy surveys and intensity maps. As increasing amounts of data become available, this and other similar analyses will open fascinating new windows into the physics of distant galaxies and expand the research possibilities available because of intensity mapping surveys.

Acknowledgements

We thank Alberto Bolatto, Kieran Cleary, Dongwoo Chung, Neal Dalal, Hamsa Padmanabhan, Tina Peters, Anthony Pullen, and Eric Switzer for useful discussions. R.M.A. would like to acknowledge funding from the National Science and Engineering Research Council of Canada.

Appendix A Power Spectra Derivations

Here we will derive the relevant auto- and cross-spectra between our CO and AGN maps, loosely following the derivation from W17. As these derivations are rarely found in the literature, we present them here in full detail for the convenience of the reader.

A.1 AGN Auto Power

First consider a map containing the locations of $N_{Q}$ quasars in volume $V$ (this derivation is true for any discrete galaxy survey, but we focus on quasars here as that is what we use above). To derive the power spectrum of such a map, we first split the volume into $N_{c}$ infinitesmial cells with small volume $\delta V$ such that cell $i$ contains $N_{i}=0$ or $1$ quasars. We can then take advantage of the useful property that $N_{i}^{2}=N_{i}$. The power spectrum $P_{Q}(k)$ of the quasar density contrast $\delta_{Q}$ is defined by the estimator

[TABLE]

where $k$ is the magnitude of Fourier wavenumber $\mathbf{k}$. The tilde indicates a Fourier transform

[TABLE]

where in the second equality we have replaced the integral with a sum over our infinitesimal cells (note that the second equality is only valid for $\mathbf{k}\neq 0$).

We can then write

[TABLE]

This sum can be divided into two parts, one where $i\neq j$, i.e. the contribution from pairs of sources, and one where $i=j$, the contribution from each source correlated with itself. In the first case, we have

[TABLE]

where $\xi_{Q}$ is the two-point correlation function of the quasar population. If we assume a simple linear bias $b_{Q}$, this is related to the correlation function of the underlying dark matter field $\xi_{m}$ by $\xi_{Q}=b_{Q}^{2}\xi_{m}$. We can rewrite Eq. (29) as

[TABLE]

The first integral vanishes due to symmetry, and we can simplify the second by recognizing that statistical isotropy and homogeneity requires that $\xi_{m}$ depend only on $\mathbf{r}\equiv\mathbf{x}-\mathbf{x}^{\prime}$. This yields

[TABLE]

which is simply the definition of the matter power spectrum $P_{m}(k)$. Thus we have

[TABLE]

For the $i=j$ component, we can instead rewrite Eq. (28) as

[TABLE]

Since $\left<N_{i}^{2}\right>=\left<N_{i}\right>=\overline{n}_{Q}\delta V$, we have

[TABLE]

The sum of the volumes of all the cells by definition must add up to $V$, so we have

[TABLE]

Combining Eq’s (32) and (35) with Eq. (26) yields the full form of the quasar power spectrum:

[TABLE]

with clustering and shot-noise components identical to those given in Eq. (1).

A.2 CO Auto Power

Next we consider the power spectrum of our CO intensity map. We again divide our map into infinitesmial cells, but we now will start by considering the intensity contributed only by halos with CO luminosities between $L$ and $L+dL$. The intensity in cell $i$ from these halos is

[TABLE]

where $N_{i}(L)=0$ or 1 is the number of emitters with the given luminosity and $C_{LT}$ is the conversion factor defined in Eq. (4). The power spectrum of an intensity map is defined similarly to that of a galaxy map, through

[TABLE]

Now we need to break our estimator again into contributions from pairs of sources and from single sources, and also down into sources with different luminosities. Starting with the two-source part, we have

[TABLE]

Similar to before, the expectation value $\left<N_{i}(L)\right>=dn/dL\delta_{V}$, where $dn/dL$ is the luminosity function of the CO line. We can also define a bias $b(L)$ between the dark matter density and the CO emitter density. We then can write

[TABLE]

This can be simplified following the same logic used to derive Eq. (32), yielding

[TABLE]

Finally, integrating over all source luminosities gives

[TABLE]

where we have assumed that $L_{1}$ and $L_{2}$ are uncorrelated, an assumption which should hold on large scales.

For the single-source term, we have

[TABLE]

Taking advantage again of the properties of $N_{i}$, we can write that

[TABLE]

where in this case the Dirac delta function enforces that $L_{1}=L_{2}$ since both describe the same source. We thus have

[TABLE]

which simplifies to

[TABLE]

Integrating over luminosities gives

[TABLE]

Combining Eq.’s (42) and (47) with (38) gives the full form of the intensity mapping power spectrum from Eq. (2):

[TABLE]

where $\overline{T}_{\rm{CO}}$, $b_{\rm{CO}}$, and $P^{\rm{shot}}_{\rm{CO}}$ have their definitions from Section 2.

A.3 CO/AGN Cross Power

Finally, we want to compute the cross-power spectrum $P_{\times}(k)$ between our two observables, defined via

[TABLE]

For our purposes, we will assume that the biases for both tracers have no scale dependence. We follow the same procedure as before. Consider first the two-source component, for CO emitters with luminosities between $L$ and $L+dL$,

[TABLE]

The volume has still been split into infinitesimal cells, $N_{i}(L)$ in this case refers to the number of CO emitters with luminosity between $L$ and $L+dL$ in cell $i$, $N_{j}$ refers to the number of quasars in cell $j$. Both $N_{i}(L)$ and $N_{j}$ can still only be 0 or 1. From here, we proceed similarly to the previous two cases:

[TABLE]

[TABLE]

with the final average over CO luminosities yielding

[TABLE]

The calculation of the contribution from single sources also begins similarly to before:

[TABLE]

Terms in the sum are now only nonzero if $N_{i}(L)$ and $N_{i}$ are both 1, i.e. if the source in a cell is both a quasar and a CO emitter. Thus we have

[TABLE]

where $dn/dL|_{Q}$ is the CO luminosity function of only the directly-observed quasars. Then,

[TABLE]

The full cross-spectrum is then, after integrating over luminosity,

[TABLE]

which is the form given in Eq. (7).