Abstract

Using the most recent releases of WISE and Planck data, we perform updated measurements of the bias and typical dark matter halo mass of infrared (IR)-selected obscured and unobscured quasars, using the angular autocorrelation function and cosmic microwave background lensing cross-correlations. Since our recent work of this kind, the WISE ALLWISE catalogue was released with improved photometry, and the Planck mission was completed and released improved products. These new data provide a more reliable measurement of the quasar bias and provide an opportunity to explore the role of changing survey pipelines in results downstream. We present a comparison of IR colour-selected quasars, split into obscured and unobscured populations based on optical–IR colours, selected from two versions of the WISE data. Which combination of data is used impacts the final results, particularly for obscured quasars, both because of mitigation of some systematics and because the newer catalogue provides a slightly different sample. We show that ALLWISE data is superior in several ways, though there may be some systematic trends with Moon contamination that were not present in the previous catalogue. We opt currently for the most conservative sample that meet our selection criteria in both the previous and new WISE catalogues. We measure a higher bias and halo mass for obscured quasars (bobsc ∼ 2.1, bunob ∼ 1.8) – at odds with simple orientation models – but at a reduced significance (∼1.5σ) as compared to our work with previous survey data.

1 INTRODUCTION

Large astronomical surveys over a wide range of wavelengths have led to a dramatic increase in public data that is mined by the community for studies of all kinds, especially systematic studies of large samples. Such surveys are generally multiyear efforts, involving multiple data releases as more observations are carried out, sky coverage expanded, and reduction pipelines improved. These changes can propagate through to the scientific results, and impact findings that need to be revisited. Further, the way that survey data are handled – samples selected, regions discarded/weighted, etc. – can shift results in important ways (e.g. Scranton et al. 2002; Huterer et al. 2006; Myers et al. 2006; Ross et al. 2011; Ho et al. 2012, 2015; Agarwal, Ho & Shandera 2014; Leistedt & Peiris 2014)

Quasars, the luminous accreting supermassive black holes in the nuclei of massive galaxies, are relatively rare objects seen primarily in the early Universe. Because of this rarity, their study has seen rapid improvements with the dramatic increase in sample sizes from astronomical surveys such as (to name a few) the Large Bright Quasar Survey (Hewett, Foltz & Chaffee 1995), Faint Images of the Radio Sky at Twenty Centimeters (Becker, White & Helfand 1995; Helfand, White & Becker 2015), the Sloan Digital Sky Survey (SDSS; York et al. 2000), and the 2dF QSO Redshift Survey (Croom et al. 2004). Quasars are a key component to studying the growth of black holes over cosmic time (e.g. Elvis et al. 1994; Richards et al. 2006; Alexander & Hickox 2012), as well as probing the potential links between those periods of growth, quasar host galaxies, and parent dark matter (DM) haloes (e.g. Hopkins et al. 2008; Booth & Schaye 2010; Volonteri, Natarajan & Gültekin 2011; Sabra et al. 2015).

Large samples of both spectroscopic and photometric quasars have permitted analysis of their large-scale distribution and clustering, which reflects their distribution in the underlying DM density field (the quasar ‘bias’, bq) and provides insight into the masses of their DM haloes. Studies of optically bright (unobscured, or type 1)1 quasars have revealed that they tend to reside in haloes of similar mass (∼3 × 1012h−1 M) across all redshifts (Porciani, Magliocchetti & Norberg 2004; Croom et al. 2005; Coil et al. 2007; Myers et al. 2007; da Ângela et al. 2008; Padmanabhan et al. 2009; Ross et al. 2009; Krumpe, Miyaji & Coil 2010; White et al. 2012; Shen et al. 2013; Eftekharzadeh et al. 2015). This suggests a link between black hole fuelling and the growth of large-scale structure.

The DM haloes of quasars not only impact their spacial distribution, but deflect the photons of the cosmic microwave background (CMB, which backlights the whole sky) travelling past them via gravitational lensing. Full-sky maps of the CMB have been steadily improving in depth and resolution, with the current state-of-the-art data being provided by the Planck satellite (Planck Collaboration I 2011). The lensing signature of large-scale structure has now been detected in numerous studies (Das et al. 2011; van Engelen et al. 2012; Planck Collaboration XVII 2014; Planck Collaboration 2015b). Combined with estimates of the intrinsic CMB power spectrum (Seljak & Zaldarriaga 1999; Hu 2001), lensing maps of the CMB can trace the projected mass along a given line of sight back to the surface of last scattering at z ∼ 1100.

CMB lensing measurements are a particularly powerful tool for studying the haloes of quasars, which peak in number density at z ∼ 2 (Croom et al. 2004; Richards et al. 2005; Fan et al. 2006), coinciding with the peak of the CMB lensing kernel. Additionally, CMB lensing measurements are subject to different systematics than clustering, providing independent follow-up to such studies. The first significant detection of a cross-correlation of unobscured quasars and the CMB lensing convergence found a typical halo mass in agreement with clustering results (Sherwin et al. 2012).

However, a subset of the quasar population has remained largely hidden from study because their optical (and in more extreme cases even X-ray) light is dramatically diminished by intervening gas and dust. The existence of these obscured (type 2) quasars has been known for some time (Setti & Woltjer 1989; Comastri et al. 1995), but only recently have large infrared (IR) data sets from Spitzer (Werner et al. 2004) and the Wide-Field Infrared Survey Explorer (WISE; Wright et al. 2010) allowed detailed study of their demographics (Lacy et al. 2004, 2013, 2015; Stern et al. 2005, 2012; Hickox et al. 2007; Assef et al. 2013, 2015; Mateos et al. 2013). However, the nature of the obscuration in these sources is still unclear (e.g. Netzer 2015), with two prominent models being geometric obscuration by a dusty torus (‘unification by orientation’, well supported at low-L and low-z; e.g. Antonucci 1993) or larger, galaxy-scale obscuration (e.g. Goulding et al. 2012) that may be a product of an evolutionary sequence (Sanders et al. 1988; Hopkins et al. 2008; Croton 2009; Booth & Schaye 2010).

A simple test of orientation models for quasars is to compare their DM haloes. If obscured quasars are simply unobscured quasars seen from a dustier line of sight, such as through a torus, then they should have the same halo mass, on average. Some evolutionary scenarios, however, predict a difference in halo mass between the subclasses as the halo and black hole grow as a product of major galaxy mergers. Hickox et al. (2011), Donoso et al. (2014), and DiPompeo et al. (2014, hereafter D14) measured the halo masses of IR-selected quasar samples split into obscured and unobscured populations via their optical–IR colours (Hickox et al. 2007). All found that obscured quasars seem to cluster more strongly, and thus reside in higher mass haloes, though the levels of significance varied considerably. DiPompeo et al. (2015a, hereafter D15) followed up on D14 by cross-correlating the WISE-selected quasar samples with a Planck CMB lensing map, and found excellent agreement with the clustering results.

However, there have been other recent studies that suggest no difference in the bias and halo masses of obscured and unobscured quasars. Geach et al. (2013) cross-correlated WISE-selected quasars with a CMB map from the South Pole Telescope (and Planck as well), and found that the bias of obscured and unobscured quasars was roughly consistent. Mendez et al. (2015) used a spectroscopic sample over a reduced area (∼10°) that benefits from individual source redshifts (instead of relying on an estimate of the ensemble average; see Section 2.2) compiled from several fields and find no significant difference between obscured and unobscured halo masses. While the additional redshift information assures that the evolution of the bias with z does not skew the results, D15 illustrated that their mean redshift estimates would need to be offset by an unreasonably large amount to fully account for the difference they measured.

Clearly the relative halo masses of obscured and unobscured quasars is not yet conclusively determined, and so we follow up here on the work of D14 and D15 using updated data from both WISE and Planck. We find that the measurements using the new and old data in various combinations can produce some significant variation in the results. Our goal here is to provide both a quantitative analysis of the difference between the samples selected from the original and new catalogues, as well as to identify the most reliable set of measurements. We then provide an updated measurement of the obscured and unobscured quasar bias based on the best possible sample.

All models and measured properties use a cosmology of H0 = 70.2 km s−1 Mpc−1, ΩM = ΩCDM + Ωb = 0.229 + 0.046 = 0.275, ΩΛ = 0.725, and σ8 = 0.82 (Komatsu et al. 2011).

2 DATA

2.1  WISE

2.1.1 ALLSKY and ALLWISE catalogues

The WISE mission mapped the entire sky at 3.4, 4.6, 12, and 22 μm (W1, W2, W3, and W4) with angular resolutions of 6.1, 6.4, 6.5, and 12 arcsec, respectively (Wright et al. 2010). The survey reached at least 0.08, 0.11, 1, and 6 mJy 5σ point-source sensitivities in each band (in unconfused regions), with this depth increasing towards the ecliptic poles due to the observing strategy. The WISE full cryogenic mission phase in 2010 resulted in the ALLSKY (AS) data release. The AS source catalogue includes objects with signal-to-noise ratio (SNR) >5 in any band, at least five good measurements, and not flagged as spurious in at least one band.2

After both cryogen tanks were exhausted, The NEOWISE Post-Cryogenic Mission surveyed the entire sky again using the two shorter-wavelength bands, W1 and W2 (Mainzer et al. 2011). Combining data from the original WISE survey with the NEOWISE data, along with improved reduction and calibration pipelines, led to the updated ALLWISE (AW) catalogue release, with improved photometric sensitivity and accuracy, better astrometric precision, and new information on source motions and variability.3

In this work, our goal is to update the analyses of D14 and D15, which used the WISE AS catalogue to select quasars and measure their bias and host halo masses, by using the improved AW catalogue (Sections 4.1 and 4.2). We will also provide a detailed comparison of the objects selected as quasars from the two catalogues (Section 4.4).

2.1.2 Quasar selection

WISE is ideal for identifying quasars via their characteristic hot dust emission, which causes a rising power-law spectrum in the mid-IR while stellar populations tend to peak around 1.5 μm (Lacy et al. 2004; Stern et al. 2005, 2012; Donley et al. 2007; Assef et al. 2013; Mateos et al. 2013). Critically, mid-IR selection is efficient at identifying both optically luminous unobscured and optically faint obscured quasars, the latter of which may make up around half of the full quasar population (Hickox et al. 2007; Assef et al. 2015).

Various methods and criteria using mid-IR data have been developed to select quasars (and lower luminosity active galactic nuclei, e.g. Lacy et al. 2004; Stern et al. 2005; Donley et al. 2012; Mateos et al. 2012; DiPompeo et al. 2015b; Myers et al. 2015), and used in an array of studies of the IR-selected quasar population (e.g. Hickox et al. 2011; Goulding et al. 2014; Satyapal et al. 2014; Smith, Koss & Mushotzky 2014; Ellison, Patton & Hickox 2015; Mendez et al. 2015). However, even a simple colour-cut of W1 − W2 > 0.8 (along with W2 < 15.05, the 10σ flux limit in this band, which helps reduce contamination from high-redshift star-forming galaxies as well as faint stars) identifies quasars at 80 per cent completeness and a contamination rate of only 5 per cent (Stern et al. 2005, 2012). These cuts have also been used to successfully study unobscured and obscured quasars (D14; D15; Donoso et al. 2014), and we adopt them here. WISE photometry is not corrected for Galactic extinction, as the rapidly dropping near-IR extinction curves of Fitzpatrick & Massa (2009) show that this will affect WISE minimally, and we avoid the Galactic plane where extinction is more prevalent (see below).

To provide a direct comparison between quasar selection with AW and AS and the resulting effect on bias measurements, we restrict our sample to the same region as Donoso et al. (2014), D14, and D15: 135° < RA < 226° and 1° < Dec < 54°. This area is sufficiently far from the Galactic plane to limit high stellar densities and the majority of Galactic reddening, but is also not affected by depth changes and source confusion in WISE. An extension beyond this footprint to utilize the full potential of millions of quasars in WISE (e.g. Secrest et al. 2015) while properly handling the selection function across the whole sky is reserved for a future paper.

In the AW catalogue, 225 303 objects satisfy the selection criteria within this footprint. This is a reduction of about 10 per cent from the AS catalogue, which had 250 163 objects satisfying these cuts.

2.1.3 Cleaning the data

An accurate data mask is necessary to properly handle the normalization of the angular autocorrelation function by comparison with a randomly distributed sample, as well as to remove pixels from the Planck maps where quasar data are discarded. D14 highlighted the importance of proper masking of the WISE data, and by more conservatively removing regions around flagged WISE data found a notable drop in the IR-selected obscured quasar bias compared to Donoso et al. (2014). We develop a similar mask for the AW data here,4 using the spherical cap utility mangle (Hamilton & Tegmark 2004; Swanson et al. 2008). Full details of these components can be found in D14, unless a change is mentioned here:

  • Regions with high Galactic extinction in the g band (Ag > 0.18).

  • WISE Atlas tiles (the main imaging product of WISE) with contamination from the Moon. We mask tiles with moon_lev >1 in W4. Using W4 makes the mask more conservative, as the longer wavelength bands can be affected by scattered light as far away from the Moon as ∼30°, while the shorter bands used for selection can be affected to ∼10°. The new values of moon_lev in the updated Atlas tiles are used (and are a major source of difference between the two masks; see Fig. 1).

  • Regions around highly grouped sources in the AW flagged data (cc_flags ≠ 0 in W1 or W2, ext_flag ≠ 0, or n_b >2; see fig. 1 of D14).

  • Regions with poor photometric quality based on the WISEph_qual flags, which is not included in the previous AS mask. We pixelize the sky using healpix5 (Górski et al. 2005) with nside = 64 (pixel areas of ∼0.8 deg2), and mask any pixel that has more than one object with photometric quality not set to ‘A’ (SNR ≥10) in W1 and W2. This procedure is tuned based on the WISE data over the full sky, primarily to remove prominent strips of low-quality WISE data, but does remove some regions in the footprint considered here. We do not discard other objects with lower quality, but these only make up a small fraction of our sample (<1 per cent), due to the W2 < 15.05 cut.

  • The SDSS bright star mask, which masks circular regions around bright stars from the Yale and Tycho-2 bright star catalogues (Warren & Hoffleit 1987; Høg et al. 2000). The vast majority of bright stars are already masked well by the WISE flagged data mask in (iii).

Top: the distribution on the sky of ALLWISE-selected quasars using the ALLWISE mask. Center: ALLWISE-selected quasars that fall within the mask generated for the ALLSKY data. Bottom: the distribution on the sky of ALLWISE-selected quasars after applying both masks.
Figure 1.

Top: the distribution on the sky of ALLWISE-selected quasars using the ALLWISE mask. Center: ALLWISE-selected quasars that fall within the mask generated for the ALLSKY data. Bottom: the distribution on the sky of ALLWISE-selected quasars after applying both masks.

The final region, after all masking is completed, has an area of 3422 deg2, and contains 175 911 quasars from the AW catalogue (Table 1). The sample distribution on the sky using this mask is shown in the top panel of Fig. 1. The usable area and sample size is similar to what was found based on the AS sample (3338 deg2, 180 606 quasars), but is distributed somewhat differently on the sky – the middle panel of Fig. 1 shows the objects that are not masked by the AW mask, but fall within the AS mask of D14. The majority of the difference is due to the Moon level between the two catalogues, as it is clear that the strips contaminated by the Moon in AS are wider than in AW.

Table 1.

Summary of samples with various masks applied. Summary of the areas (in deg2) and number N of each subsample using the two WISE catalogues and various masks. Note that when samples are split to obscured/unobscured, the additional SDSS bad fields mask is applied, reducing the sample size and area further. An asterisk indicates a sample with the masks from both ALLSKY and ALLWISE applied, and ‘both’ indicates only objects that meet the Stern et al. (2012) quasar selection criteria in both ALLSKY and ALLWISE. The second half gives the same information but after also applying Planck mask and discarding any partially masked healpixnside = 2048 pixels, to limit errors in estimating the used acmbrea of partial pixels in the quasar density calculation. Obscured/unobscured percentages are given in the relevant rows – the fractions are very consistent across the samples.

ALLSKYALLWISEALLSKY*ALLWISE*Both
Area (WISE)33733422325032503250
N all IR180 860175 911174 597167 364154 862
Area (WISE+SDSS)33393387321632163216
N obscured74 972 (42.1 per cent)72 587 (41.9 per cent)72 717 (42.1 per cent)69 031 (41.7 per cent)62 715 (41.4 per cent)
N unobscured102 858 (57.9 per cent)101 247 (58.1 per cent)99 894 (57.9 per cent)96 334 (58.3 per cent)88 834 (58.6 per cent)
+CMB lensing mask
Area (WISE)33003189303230323032
N all IR176 243165 044163 356157 036145 238
Area (WISE+SDSS)30783150299429942994
N obscured70 093 (42.1 per cent)67 935 (41.7 per cent)67 839 (42.1 per cent)64 637 (41.7 per cent)58 688 (41.4 per cent)
N unobscured96 401(57.9 per cent)94 885 (58.3 per cent)93 375 (57.9 per cent)90 256 (58.3 per cent)83 187 (58.6 per cent)
ALLSKYALLWISEALLSKY*ALLWISE*Both
Area (WISE)33733422325032503250
N all IR180 860175 911174 597167 364154 862
Area (WISE+SDSS)33393387321632163216
N obscured74 972 (42.1 per cent)72 587 (41.9 per cent)72 717 (42.1 per cent)69 031 (41.7 per cent)62 715 (41.4 per cent)
N unobscured102 858 (57.9 per cent)101 247 (58.1 per cent)99 894 (57.9 per cent)96 334 (58.3 per cent)88 834 (58.6 per cent)
+CMB lensing mask
Area (WISE)33003189303230323032
N all IR176 243165 044163 356157 036145 238
Area (WISE+SDSS)30783150299429942994
N obscured70 093 (42.1 per cent)67 935 (41.7 per cent)67 839 (42.1 per cent)64 637 (41.7 per cent)58 688 (41.4 per cent)
N unobscured96 401(57.9 per cent)94 885 (58.3 per cent)93 375 (57.9 per cent)90 256 (58.3 per cent)83 187 (58.6 per cent)
Table 1.

Summary of samples with various masks applied. Summary of the areas (in deg2) and number N of each subsample using the two WISE catalogues and various masks. Note that when samples are split to obscured/unobscured, the additional SDSS bad fields mask is applied, reducing the sample size and area further. An asterisk indicates a sample with the masks from both ALLSKY and ALLWISE applied, and ‘both’ indicates only objects that meet the Stern et al. (2012) quasar selection criteria in both ALLSKY and ALLWISE. The second half gives the same information but after also applying Planck mask and discarding any partially masked healpixnside = 2048 pixels, to limit errors in estimating the used acmbrea of partial pixels in the quasar density calculation. Obscured/unobscured percentages are given in the relevant rows – the fractions are very consistent across the samples.

ALLSKYALLWISEALLSKY*ALLWISE*Both
Area (WISE)33733422325032503250
N all IR180 860175 911174 597167 364154 862
Area (WISE+SDSS)33393387321632163216
N obscured74 972 (42.1 per cent)72 587 (41.9 per cent)72 717 (42.1 per cent)69 031 (41.7 per cent)62 715 (41.4 per cent)
N unobscured102 858 (57.9 per cent)101 247 (58.1 per cent)99 894 (57.9 per cent)96 334 (58.3 per cent)88 834 (58.6 per cent)
+CMB lensing mask
Area (WISE)33003189303230323032
N all IR176 243165 044163 356157 036145 238
Area (WISE+SDSS)30783150299429942994
N obscured70 093 (42.1 per cent)67 935 (41.7 per cent)67 839 (42.1 per cent)64 637 (41.7 per cent)58 688 (41.4 per cent)
N unobscured96 401(57.9 per cent)94 885 (58.3 per cent)93 375 (57.9 per cent)90 256 (58.3 per cent)83 187 (58.6 per cent)
ALLSKYALLWISEALLSKY*ALLWISE*Both
Area (WISE)33733422325032503250
N all IR180 860175 911174 597167 364154 862
Area (WISE+SDSS)33393387321632163216
N obscured74 972 (42.1 per cent)72 587 (41.9 per cent)72 717 (42.1 per cent)69 031 (41.7 per cent)62 715 (41.4 per cent)
N unobscured102 858 (57.9 per cent)101 247 (58.1 per cent)99 894 (57.9 per cent)96 334 (58.3 per cent)88 834 (58.6 per cent)
+CMB lensing mask
Area (WISE)33003189303230323032
N all IR176 243165 044163 356157 036145 238
Area (WISE+SDSS)30783150299429942994
N obscured70 093 (42.1 per cent)67 935 (41.7 per cent)67 839 (42.1 per cent)64 637 (41.7 per cent)58 688 (41.4 per cent)
N unobscured96 401(57.9 per cent)94 885 (58.3 per cent)93 375 (57.9 per cent)90 256 (58.3 per cent)83 187 (58.6 per cent)

We explore the impact of our masks on our bias measurements by applying various combinations to the data, including the conservative cases of applying both the AW and AS masks. These samples with both masks applied will be labelled with an asterisk (*) throughout, and the AW* sample is shown in the bottom panel of Fig. 1. The sample labelled ‘both’ contains sources that satisfy our selection criteria in both the AW and AS catalogues, and these make up 92.5 per cent of the full AW* sample. Table 1 summarizes these samples.

2.2 Obscured and unobscured quasars

Hickox et al. (2007) used multiwavelength data in the Böotes field to demonstrate that an optical–IR colour cut at R−[4.5]=6.1 (Vega magnitudes) can robustly separate obscured and unobscured quasar populations. Donoso et al. (2014), D14, and D15 used AS W2 and SDSS r-band fluxes in a similar way to study these populations, with rW2 = 6 separating the samples. We do the same using AW W2 magnitudes here.

The SDSS completely covers our footprint, and reaches 50 per cent completeness at r = 22.6 (York et al. 2000; Abazajian et al. 2009). As in D14/D15, we utilize the SDSS pipeline psfmag values, to isolate as much as possible the contribution from the quasar, as compared to that from the host galaxy, in resolved sources. D14 determined that the use of SDSS modelMags did not affect results significantly, and given the general similarity in the morphologies of AW-selected quasars (see below) there is no reason to believe this will change here. The SDSS r-band magnitudes are corrected for Galactic extinction using extinction values supplied with the SDSS data (based on the dust maps of Schlegel, Finkbeiner & Davis 1998), and converted from AB magnitudes to Vega using mr,AB = mr,Vega + 0.16 (Blanton & Roweis 2007).

For the obscured and unobscured samples (but not the total IR-selected sample) the SDSS bad fields mask6 is applied (e.g. White et al. 2011, 2012). This reduces the AW area to 3387 deg2 and the total number of sources to 173 834. We match the AW-selected quasars to SDSS sources with a radius of 2 arcsec, only keeping objects with 15 < r < 25, and find matches for 83 per cent, the same as was found for AS-selected quasars. Objects without SDSS matches are placed in the obscured sample, resulting in 72 587 (41.9 per cent) obscured and 101 247 (58.1 per cent) unobscured quasars.

Of the objects with an SDSS counterpart, 68 per cent are unresolved (based on the objc_type keyword; Stoughton et al. 2002), marginally higher than the 65 per cent in the AS-selected sample of D14 and significantly higher than the 55 per cent of Donoso et al. (2014). Note however that Donoso et al. (2014) used deeper optical imaging in the COSMOS field to classify source morphologies, which partially explains the lower unresolved fraction. Broken down into subclasses, 36 and 81 per cent of the obscured and unobscured sources are unresolved, respectively, broadly consistent with D14. We note that the additional application of the AS mask to the AW data does not change these ratios.

As in D14/D15, to estimate the redshift distribution (dN/dz) of the AW-selected quasars, which is necessary to understand and interpret the bias measurements, we apply our mask and selection criteria to objects in the Böotes field. This field has been observed in multiple photometric bands, and has extensive follow-up spectroscopy (Brodwin et al. 2006; Hickox et al. 2011; Kochanek et al. 2012). We find 368 with AW data that satisfy our selection in this field, with 145 (39.4 per cent) and 223 (60.6 per cent) obscured and unobscured, respectively. These are consistent with the fractions in our overall sample. All of the unobscured and all but two of the obscured samples have spectroscopic redshifts, and the distributions are shown in Fig. 2. The mean/median/standard deviations of these distributions are 1.02/0.97/0.56 (total), 0.98/0.90/0.54 (obscured), and 1.05/1.04/0.58 (unobscured), consistent with what is found for AS-selected samples.

Redshift distributions of ALLWISE selected quasars, generated by applying our mask and selection criteria to objects with multiwavelength and spectroscopic follow-up in the Böotes field. These distributions are nearly identical as when using the ALLSKY catalogue data (D14, D15).
Figure 2.

Redshift distributions of ALLWISE selected quasars, generated by applying our mask and selection criteria to objects with multiwavelength and spectroscopic follow-up in the Böotes field. These distributions are nearly identical as when using the ALLSKY catalogue data (D14, D15).

We point out that the AW and AS mask do not significantly affect the Böotes region, and so potential differences in the redshift distribution of sources masked in one catalogue and not the other are not accounted for. This will be analysed further in Section 4.4.1. Finally, it is worth noting that the AW- and AS-selected samples have very similar redshift distributions, despite differences in their average photometric properties (see Sections 4.4.1 and 4.4.2).

2.3 Planck CMB lensing maps

The Planck mission (Planck Collaboration I 2011) mapped the CMB at nine frequencies, from 30 to 857 GHz, over the entire sky. The first data release (DR1) in 2013 March (with an update in 2013 December) was based on the nominal mission data with 15 months of observations (Planck Collaboration I 2014) and included a lensing potential (ϕ) map with a lensing signature detected at an overall significance >25σ (Planck Collaboration XVII 2014). This map was cross-correlated with the quasar density by Geach et al. (2013) and D15.

In 2015 a second data set (DR2), using four years of data and improved reduction and calibration pipelines, was released (Planck Collaboration 2015a). This release included an updated map of the lensing potential, with a lensing signal detected at >40σ (Planck Collaboration 2015b). These DR2 data are the highest quality AS CMB maps to date, with sensitivities down to μK and, at the frequencies where most of the lensing information is carried (143 and 217 GHz), resolutions of 7 and 5 arcmin.

Using the healpix routine ISYNFAST, we convert the Planck DR2 ‘alm’ (Alm) file,7 which contains the coefficients of the spherical harmonic transform of the data on the sky, to a healpix map of the lensing convergence (⁠|$\kappa = \frac{1}{2} \Delta ^2 \phi$|⁠) with nside = 2048, which has ∼3 arcmin2 pixels. In Fig. 3 we show the distribution of both the raw and 1° Gaussian-smoothed κ in our region of interest. The updated DR2 distribution is narrower, and has a peak marginally more consistent with zero, which is expected in the absence of systematic effects.

Raw (top) and 1° Gaussian smoothed (bottom) distributions of κ from Planck DR1 (dashed) and DR2 (solid). The DR2 distribution is narrower, and has a peak slightly more consistent with zero.
Figure 3.

Raw (top) and 1° Gaussian smoothed (bottom) distributions of κ from Planck DR1 (dashed) and DR2 (solid). The DR2 distribution is narrower, and has a peak slightly more consistent with zero.

3 METHODS

Our goal is to measure the quasar bias and infer dark mater halo masses using the quasar data alone, as well as cross-correlate the quasar data with CMB lensing maps. We also explore the WISE AS and AW catalogues for possible sources of contamination, to identify which combination of data and masks is the most reliable. In this section we outline the formalisms for these measurements.

3.1 Angular autocorrelations

All of the codes used to calculate the angular autocorrelation functions, including models and bias fitting, are available at https://github.com/mdipompe/angular_clustering.

3.1.1 DM autocorrelation function theory

We utilize the two-point angular (as opposed to real-space, as our sources lack individual redshift measurements) correlation function ω(θ) to analyse how WISE quasars cluster around themselves. The angular correlation function is related to the probability that a given pair of objects (DM haloes hosting quasars, in our case) with mean number density n, separated by a projected angular distance θ, are within a solid angle dΩ (Totsuji & Kihara 1969; Peebles 1980):
(1)
Objects formed in the peaks of a Gaussian random field, such as massive, quasar-hosting galaxies, should cluster more strongly than the underlying typical DM distribution (Kaiser 1984; Bardeen et al. 1986). This excess, or bias, is independent of scale in most models, at least on large scales. The bias of quasars bq is related to the underlying DM autocorrelation by |$\omega _{\rm q}(\theta ) = \omega _{\textrm{dm}}(\theta ) b_{\rm q}^2$|⁠.

To calculate ωdm(θ), we first generate the non-linear matter power spectrum P(k, z) using camb8 (Lewis, Challinor & Lasenby 2000). We note that this procedure is updated from D14, for more consistency with the CMB lensing cross-correlation measurements.

For angular scales θ ≪ 1 radian, as we probe here, we can project the matter power spectrum to an angular autocorrelation in a flat Universe via Limber's approximation (Limber 1953; Peebles 1980; Peacock 1991):
(2)
In this equation, Δ2(k, z) = (k3/2π2)P(k, z) is the dimensionless power spectrum, J0 is the zeroth-order Bessel function of the first kind, χ is the comoving distance along the line of sight, dN/dz is the normalized redshift distribution (taken from a spline fit to the distribution of the appropriate subsample shown in Fig. 2; see Section 5.2), and dz/dχ = Hz/c = (H0/c)[Ωm(1 + z)3 + ΩΛ]1/2. This model of ωdm(θ) can then be rescaled to fit the estimated quasar autocorrelation and measure the effective bias, or the bias integrated over our redshift range. We will explore the effect of other bias models that evolve with z in the discussion section.

3.1.2 Estimating ωqq(θ)

We estimate the quasar autocorrelation ωqq(θ) by comparing the number counts of quasar pairs in annuli of increasing radii with what is expected for a random distribution (Landy & Szalay 1993):
(3)
In this estimator, DD, DR, and RR are the normalized numbers of data–data, data–random, and random–random pairs in each bin of θ [i.e. DD = DD(θ) = Ndata pairs/(NDND)]. The random sample must follow the same angular selection function as the data, which is simple in our case because the WISE selection is uniform over this field with holes described by our mask. We generate a random catalogue that obeys our mask using the mangle function RANSACK. The random catalogue for each measurement is always at least 10 times the size of the data set so that the random counts do not limit the statistical precision. We calculate ωqq using four bins per dex, beginning at ∼12 arcsec and extending to 1| $_{.}^{\circ}$|1.
Errors on the quasar autocorrelations are estimated using inverse-variance-weighted jackknife resampling (e.g. Scranton et al. 2002; Myers et al. 2005, 2007). We divide our full footprint into N = 16 equal-area regions, build N subsamples by iteratively removing a single region, and repeating the autocorrelation measurement using the remaining regions. Denoting each subsample by L, the inverse-variance-weighted covariance matrix |${\boldsymbol{\sf C}}_{ij} = {\boldsymbol{\sf C}}(\theta _i,\theta _j)$| (i and j denote angular size bins) is
(4)
where ω is the angular autocorrelation for all of the quasars and ωL is the angular autocorrelation for subset L. Note here that the RR terms are not normalized by the sizes of the random samples, and account for the different number of counts in each region (though these are very close to the same, given the equal area of each region). The jackknife errors σi are taken from the square root of the diagonal elements of the covariance matrix, though the full matrix is used when performing fits to measure the bias.

3.2 Cross-correlations with systematics

In an ideal data set, the quasar density will not correlate with any observational systematics, but in practice such effects can impact target selection. Below we will cross-correlate WISE AS- and AW-selected quasars with various parameters, to quantify whether one catalogue is superior in this regard, using a pixelization method as in e.g. Scranton et al. (2002). Due to the relatively low density of sources (∼50 deg−2, less for the obscured and unobscured subsamples), we must use somewhat large pixels of 0.66 deg on a side (0.44 deg2 in area), and are thus limited to exploring cross-correlations on scales of this size or larger.

Splitting up our region results in 9360 pixels, with an average of ∼20 sources per pixel. The available area of each pixel is calculated by randomly populating it with 1000 sources, applying our mask, and multiplying the full area by the fraction of random points outside the mask. Only pixels with at least 0.05 deg2 of available area and at least five sources (one in the case of obscured and unobscured subsamples) are used.

For each pixel i, the relative density of quasars is calculated:
(5)
where 〈ρq〉 is the mean quasar density for the whole field. The relative value of systematics |$\delta _i^{\rm s}$|⁠, including W1 and W2 magnitudes, W1 − W2 colour, Galactic reddening Ag, and the WISE Moon level (recall that any region with Moon level >1 is already discarded) is calculated in a similar way for each pixel, using the values of these parameters at the location of each object. Cross-correlations are then calculated as
(6)
where Θij = 1 if the separation between the centres of pixels i and j is within the bin θ, and zero otherwise (the denominator is then a normalization equal to the number of pixel pairs satisfying this criteria). The angular binning here is also four bins per dex, beginning at 0| $_{.}^{\circ}$|6 (the smallest scale possible due to the pixel size), and extending to 10°.

Errors on the angular correlation functions using this method are also estimated via jackknife resampling. In this case, the region is split into N = 25 regions, and for a given iteration in the covariance matrix calculation (using equation 4) any pixel falling within this subregion is excluded. The square root of the diagonal elements of the covariance matrix are adopted as the 1σ errors.

We also use equation (6) to rapidly calculate the quasar autocorrelation on large scales (>1°), where the pair counting method becomes computationally expensive.

3.3 CMB lensing correlations

All of our code to work with the CMB lensing maps and perform the CMB cross-correlations with healpix maps can be found at https://github.com/mdipompe/lensing_xcorr.

3.3.1 CMB lensing-matter cross-correlation theory

The quasar bias bq can also be measured by comparing the cross-correlation of the quasar density with the CMB lensing convergence (⁠|$C_l^{\kappa {\rm q}}$|⁠) with theoretical predictions for a given matter distribution. This formalism is detailed fully elsewhere (Bleem et al. 2012; Sherwin et al. 2012), and we provide a brief summary here.

The lensing convergence (κ) in comoving coordinates (χ) along a line of sight |$\hat{\boldsymbol {n}}$| is the integral over the relative overdensity of matter (⁠|$\delta (\boldsymbol {r},z)$|⁠) multiplied by the lensing kernel Wκ:
(7)
The lensing kernel (Cooray & Hu 2000; Song et al. 2003) is
(8)
where a(χ) = (1 + z(χ))−1 is the scalefactor, and χCMB is the comoving distance to the CMB. Fluctuations in the quasar density are given by
(9)
where Wq(χ) is the quasar host distribution kernel:
(10)
Here, dN/dz is the normalized redshift distribution of the quasar population (again estimated using a spline interpolation; see Section 5.2), which has bias bq, assumed here to be independent of redshift. The cross-power at a Fourier mode l is
(11)
where P(k = l/χ, z) is the matter power spectrum (e.g. Eisenstein & Hu 1999), again generated using camb. Equation (11) gives us the model cross-power spectrum for the underlying distribution of matter when the effective bias is unity (again we will explore other models of the bias in the discussion section). We note that the current defaults of camb have been slightly updated since D15, and these result in a matter power spectrum with slightly higher amplitude that propagates to the final model, resulting in a lower bias measurement. A camb parameter file is included with our supplied code to aid with consistent future measurements.

3.3.2 Measuring CMB lensing auto- and cross-correlations

Leveraging the healpix format of the Planck maps, and the speed with which cross-correlations can be performed using spherical harmonic transforms, we employ routines in the healpix package to measure auto- and cross-correlations with Planck data. The cross-power |$C_l^{\kappa X}$| of a κ map (Mκ) with map X (MX, which could be another κ map or a quasar density map as described below) is measured by taking the Fourier transform of each map and multiplying them:
(12)
where |$\boldsymbol {l}\in l$| describes the binning. We present all results calculated this way with five bins in l per dex, beginning at l = 10.

To cross-correlate the quasar density with the CMB lensing convergence (⁠|$C_l^{\kappa {\rm q}}$|⁠), we first generate a healpix map of the relative quasar density δ at the same resolution as the Planck lensing convergence map (nside = 2048). This requires an estimate of the used area of each pixel, which we find by populating our footprint with 150 million points (for an average of ∼30 per pixel), applying our mask, and using the ratio of points inside and outside the mask per pixel. This area estimation is of course subject to its own errors, and as discussed in D15 can impact the final measurement in several ways. We therefore discard any pixel that overlaps a mask component, which generally removes less than 100 deg2 (see Table 1). The relative quasar density is then calculated with respect to the total number of quasars and area in the remaining complete pixels.

Uncertainties in |$C_l^{\kappa X}$| are estimated by repeating the measurement with several rotations of the κ map. D15 illustrated the consistency of this method with others, such as substituting simulated noise maps. We use 34 rotations, 17 in increments of 20° in Galactic longitude, and another 17 with an additional reflection in latitude about the Galactic equator. We derive covariance matrices (C(li, lj) = Cij):
(13)
where N is the number of rotated cross-correlations and |$C_{l,k}^{\kappa X} - C_{l}^{\kappa X}$| is the cross-correlation from each rotation. We adopt the square root of the diagonal elements of Cij as the 1σ errors on these correlations.

3.4 Fitting procedures

To measure the quasar bias, we fit the models from equations (2) and (11) to the measured autocorrelations and CMB lensing cross-correlations. We use the full covariance matrix to scale our given model fm(x) to the data f(x) with a χ2 minimization:
(14)
Here, f(x) may be ω(θ) or |$C_l^{\kappa {\rm q}}$| and the sums are over bins of θ or l. Both models are only a function of one parameter, the bias bq, and so the errors on the fits are determined by where Δχ2 = 1.

On scales larger than ∼1 Mpc h−1 the clustering of quasars and their parent haloes is sensitive only to the underlying density field and can be fit by the simple models described in the previous sections, while on smaller scales the halo occupation distribution (HOD) and the physics of galaxy formation and evolution become important (e.g. Berlind & Weinberg 2002; Berlind et al. 2003; Richardson et al. 2012; Krumpe, Miyaji & Coil 2014; Eftekharzadeh et al. 2015). At z = 1 (the approximate mean for all of our samples) this linear scale corresponds to ∼0| $_{.}^{\circ}$|04. Therefore, we restrict our fits on the autocorrelation to 0| $_{.}^{\circ}$|04 < θ < 0| $_{.}^{\circ}$|4. For CMB lensing cross-correlations we restrict our fits to 40 < l < 400, which is above the peak of the model cross-correlation and below a possible correlated feature in the Planck DR2 lensing map (which lies in the range 638 < l < 732; Planck Collaboration 2015b). This is also the range with the smallest errors on the cross-correlation (Section 4.2).

3.5 DM halo masses

All of the code used to convert biases to halo masses is provided at https://github.com/mdipompe/halomasses.

Once the quasar bias is determined, it can be converted into a typical DM halo mass (Mh) for a given sample. This is done using the model fits to cosmological simulations of Tinker et al. (2010):9
(15)
where ν = δc/σ(M). The numerator, δc, is the critical density for collapse of a DM halo, and is defined for a Universe containing matter and a cosmological constant Λ by
(16)
where Ωm,z is the density parameter for matter at the redshift under consideration (Navarro, Frenk & White 1997). The denominator, σ(M), is the linear matter variance at the size scale of the halo (⁠|$R_{\rm h} = (3M_{\rm h}/4\pi \bar{\rho }_{\rm m})^{1/3}$|⁠, with |$\bar{\rho }_{\rm m}$| the mean density of matter) calculated by
(17)
We again use camb to obtain the linear P(k, z), and |$\hat{W}(k,R_{\rm h})$| is the Fourier transform of a top-hat window function of radius Rh:
(18)
The constants A, a, B, b, C, and c for equation (15) are taken from table 2 of Tinker et al. (2010) and defined for the overdensity parameter Δ = 200 (y ≡ log10Δ):
(19)
Table 2.

Summary of all bias measurements. Bias measurements for the various sample and mask combinations, using the quasar autocorrelation (top half) and CMB lensing cross-correlations (bottom half). Measurements are made by fitting the model to the data over the range 40 < l < 400 or 0| $_{.}^{\circ}$|04 < θ < 0| $_{.}^{\circ}$|4. Samples with an asterisk have had the masks developed with both ALLSKY and ALLWISE samples applied. Values in parentheses in the bottom half are from measurements with the Planck DR1 data, and all others use the DR2 data. These results are shown in Fig. 6.

ALLSKYALLWISEALLSKY*ALLWISE*Both
Quasar autocorrelation
bq all IR2.20 ± 0.061.93 ± 0.122.01 ± 0.101.87 ± 0.121.92 ± 0.12
bq obscured2.39 ± 0.142.06 ± 0.192.41 ± 0.171.93 ± 0.212.13 ± 0.21
bq unobscured1.87 ± 0.151.86 ± 0.121.91 ± 0.151.83 ± 0.131.88 ± 0.15
CMB lensing cross-correlation
bq all IR(2.18 ± 0.17) 1.94 ± 0.14(1.95 ± 0.15) 1.82 ± 0.14(2.10 ± 0.13)1.86 ± 0.141.87 ± 0.14
bq obscured(2.27 ± 0.25) 2.06 ± 0.21(2.22 ± 0.21) 1.82 ± 0.20(2.25 ± 0.22)2.12 ± 0.222.06 ± 0.22
bq unobscured(1.85 ± 0.24) 1.78 ± 0.17(1.74 ± 0.20) 1.79 ± 0.19(1.93 ± 0.16)1.70 ± 0.191.72 ± 0.18
ALLSKYALLWISEALLSKY*ALLWISE*Both
Quasar autocorrelation
bq all IR2.20 ± 0.061.93 ± 0.122.01 ± 0.101.87 ± 0.121.92 ± 0.12
bq obscured2.39 ± 0.142.06 ± 0.192.41 ± 0.171.93 ± 0.212.13 ± 0.21
bq unobscured1.87 ± 0.151.86 ± 0.121.91 ± 0.151.83 ± 0.131.88 ± 0.15
CMB lensing cross-correlation
bq all IR(2.18 ± 0.17) 1.94 ± 0.14(1.95 ± 0.15) 1.82 ± 0.14(2.10 ± 0.13)1.86 ± 0.141.87 ± 0.14
bq obscured(2.27 ± 0.25) 2.06 ± 0.21(2.22 ± 0.21) 1.82 ± 0.20(2.25 ± 0.22)2.12 ± 0.222.06 ± 0.22
bq unobscured(1.85 ± 0.24) 1.78 ± 0.17(1.74 ± 0.20) 1.79 ± 0.19(1.93 ± 0.16)1.70 ± 0.191.72 ± 0.18
Table 2.

Summary of all bias measurements. Bias measurements for the various sample and mask combinations, using the quasar autocorrelation (top half) and CMB lensing cross-correlations (bottom half). Measurements are made by fitting the model to the data over the range 40 < l < 400 or 0| $_{.}^{\circ}$|04 < θ < 0| $_{.}^{\circ}$|4. Samples with an asterisk have had the masks developed with both ALLSKY and ALLWISE samples applied. Values in parentheses in the bottom half are from measurements with the Planck DR1 data, and all others use the DR2 data. These results are shown in Fig. 6.

ALLSKYALLWISEALLSKY*ALLWISE*Both
Quasar autocorrelation
bq all IR2.20 ± 0.061.93 ± 0.122.01 ± 0.101.87 ± 0.121.92 ± 0.12
bq obscured2.39 ± 0.142.06 ± 0.192.41 ± 0.171.93 ± 0.212.13 ± 0.21
bq unobscured1.87 ± 0.151.86 ± 0.121.91 ± 0.151.83 ± 0.131.88 ± 0.15
CMB lensing cross-correlation
bq all IR(2.18 ± 0.17) 1.94 ± 0.14(1.95 ± 0.15) 1.82 ± 0.14(2.10 ± 0.13)1.86 ± 0.141.87 ± 0.14
bq obscured(2.27 ± 0.25) 2.06 ± 0.21(2.22 ± 0.21) 1.82 ± 0.20(2.25 ± 0.22)2.12 ± 0.222.06 ± 0.22
bq unobscured(1.85 ± 0.24) 1.78 ± 0.17(1.74 ± 0.20) 1.79 ± 0.19(1.93 ± 0.16)1.70 ± 0.191.72 ± 0.18
ALLSKYALLWISEALLSKY*ALLWISE*Both
Quasar autocorrelation
bq all IR2.20 ± 0.061.93 ± 0.122.01 ± 0.101.87 ± 0.121.92 ± 0.12
bq obscured2.39 ± 0.142.06 ± 0.192.41 ± 0.171.93 ± 0.212.13 ± 0.21
bq unobscured1.87 ± 0.151.86 ± 0.121.91 ± 0.151.83 ± 0.131.88 ± 0.15
CMB lensing cross-correlation
bq all IR(2.18 ± 0.17) 1.94 ± 0.14(1.95 ± 0.15) 1.82 ± 0.14(2.10 ± 0.13)1.86 ± 0.141.87 ± 0.14
bq obscured(2.27 ± 0.25) 2.06 ± 0.21(2.22 ± 0.21) 1.82 ± 0.20(2.25 ± 0.22)2.12 ± 0.222.06 ± 0.22
bq unobscured(1.85 ± 0.24) 1.78 ± 0.17(1.74 ± 0.20) 1.79 ± 0.19(1.93 ± 0.16)1.70 ± 0.191.72 ± 0.18

4 RESULTS AND ANALYSIS

4.1 The quasar autocorrelation function

In Fig. 4 we show the quasar angular autocorrelation ω(θ) for various combinations of samples and masks. As a reminder, samples marked with an asterisk have both the AW- and AS-generated masks applied. The final panel of this plot shows our fiducial model, as generated with equation (2) for each of the sample redshift distributions (note that the subtle differences in dN/dz do not strongly affect the models). The final panel also shows the fitting ranges for the bias measurement.

The angular autocorrelation for various data and mask combinations. Note that a sample marked with an asterisk has both the AW and AS masks applied. The final panel shows the fiducial models for each sample, from equation (2), as well as the bias fitting range. The dashed lines in each panel show the model rescaled by measured biases. Small panels under each show the measurements divided by the models, to highlight any scale dependences.
Figure 4.

The angular autocorrelation for various data and mask combinations. Note that a sample marked with an asterisk has both the AW and AS masks applied. The final panel shows the fiducial models for each sample, from equation (2), as well as the bias fitting range. The dashed lines in each panel show the model rescaled by measured biases. Small panels under each show the measurements divided by the models, to highlight any scale dependences.

In each panel, the models rescaled by the measured bias are shown as dashed lines. The panels below each plot show the autocorrelations divided by the models, to highlight any scale dependences of the bias. Bias measurements for these samples are given in the top half of Table 2, and shown in the right-hand side of Fig. 6.

In all cases, over the chosen fitting range the bias is independent of scale. Below ∼0| $_{.}^{\circ}$|04, for the obscured sample particularly, clustering amplitudes anticorrelate with angular scale. This is true for the unobscured sample as well, but generally only below ∼0| $_{.}^{\circ}$|01. The bias of the unobscured sample is very consistent across all samples and mask combinations (Fig. 6). On the other hand, there is a marked decrease in the obscured bias when switching to AW-selected sources, even when the AS mask is applied.

The difference in bias between obscured and unobscured quasars identified by D14 remains present when using the AS data with the additional AW mask components.10 It is also still marginally present in the AW-selected samples, as well as the sample selected from both catalogues – however, in these cases the error bars overlap significantly, reducing the difference to at best 1σ.

4.2 Quasar CMB lensing cross-correlation

Fig. 5 shows the quasar CMB lensing cross-correlations for various combinations of the data, including the AW and AS samples as well as Planck DR1 and DR2 lensing maps. The bottom-centre panel shows the model generated from equation (11) for each dN/dz, and highlights the bias fitting region. Again, the model does not depend strongly on the subtle differences in redshift distribution for the three samples. Dotted lines in each panel show the models rescaled by the bias values, and the small panels below each show the measurement divided by the model to highlight scale dependences. Bias measurements are summarized in the lower half of Table 2, and shown in the left-hand panel of Fig. 6.

The CMB-lensing cross-correlations for various WISE, Planck, and mask combinations. Note that WISE samples marked with an asterisk have both the AW and AS masks applied. The bottom-centre panel shows the fiducial model for each subsample, as well as the fitting range for the bias measurement. The dotted lines in each panel show the model rescaled by the measured bias. Small panels below each measurement show the measurement divided by the model, to highlight any scale dependences.
Figure 5.

The CMB-lensing cross-correlations for various WISE, Planck, and mask combinations. Note that WISE samples marked with an asterisk have both the AW and AS masks applied. The bottom-centre panel shows the fiducial model for each subsample, as well as the fitting range for the bias measurement. The dotted lines in each panel show the model rescaled by the measured bias. Small panels below each measurement show the measurement divided by the model, to highlight any scale dependences.

A comparison of the measured bias via CMB lensing cross-correlations (left) and the quasar autocorrelation (right), using various data sets. In the CMB lensing panel, ‘DR1’ and ‘DR2’ refer to the Planck data release. An asterisk on the WISE catalogue indicates that the masks derived from both catalogues have been applied, and ‘Both’ is the measurement for objects satisfying our criteria in both catalogues. The numbers under each measurement indicate the total area used, in deg2.
Figure 6.

A comparison of the measured bias via CMB lensing cross-correlations (left) and the quasar autocorrelation (right), using various data sets. In the CMB lensing panel, ‘DR1’ and ‘DR2’ refer to the Planck data release. An asterisk on the WISE catalogue indicates that the masks derived from both catalogues have been applied, and ‘Both’ is the measurement for objects satisfying our criteria in both catalogues. The numbers under each measurement indicate the total area used, in deg2.

As seen in the autocorrelation measurements, the lensing cross-correlations are generally scale-independent over the chosen fitting range. However, beyond this range towards smaller scales (larger l), there is a clear increase in cross-correlation power that is particularly strong for the obscured sample, though more prevalent for both samples when Planck DR1 is used. While the measurement is noisy at the largest l, with the Planck DR2 lensing map the unobscured cross-correlation stays consistent with flat out to the smallest scales (∼0| $_{.}^{\circ}$|1, consistent with the clustering measurement). The scale-dependence for the obscured sample is reduced somewhat with DR2, but is still present. The obscured sample cross-correlation with DR2 also shows a prominent bump around l ∼ 600, about the scales for which there may be an unknown correlated feature in the DR2 lensing map (Planck Collaboration 2015b).

The bias values for the unobscured sample are quite consistent across all sample and mask combinations, and with the autocorrelation measurements (Fig. 6; please refer footnote 10). Overall the use of Planck DR2 reduces the measured bias in both samples (compare for example the measurements for DR1+AS and DR2+AS). For the obscured sample, the use of Planck DR2 and the AW-selected sample reduces the bias somewhat, until the AS mask is also applied where the bias is raised again. Overall, the obscured sample tends to have a higher bias, though at a lower significance compared to D15 (at most ∼2σ).

4.3 Planck DR1 versus DR2 lensing convergence

Now that we have seen that the updated Planck CMB lensing map tends to reduce the measured bias, particularly for the obscured sample, we briefly compare the two maps directly to investigate quantitatively whether one map is superior to the other. Of course, the fact that the signature of lensing, at least as averaged over the whole sky, is so much stronger in DR2 is already indicative of the superiority of the newest data. Fig. 3 illustrates the improved behaviour of the DR2 data as well. Finally, DR2 has a lensing potential power spectrum in better agreement with theoretical predictions (see e.g. fig. 6 of Planck Collaboration 2015b).

We also correlate the DR1 and DR2 κ maps with themselves and each other (using the method described in Section 3.3.2), to highlight potential differences. If the maps have similar properties, these auto- and cross-correlations should all be quite similar. The results are shown in Fig. 7. While the DR2–DR2 auto-correlation and DR2–DR1 cross-correlation are nearly identical, the DR1–DR1 autocorrelation has significantly more power (by ∼45 per cent), while having the same shape at all scales. This suggests that the DR2 map contains significantly less correlated noise, while still preserving real features in the data. However, it is not clear how this might affect the obscured sample cross-correlation more than that of the unobscured sample.

Auto- and cross-correlations of the Planck DR1 and DR2 κ maps. The similar shape but higher power at all scales in the DR1 autocorrelation suggests that there is more correlated noise in the DR1 map.
Figure 7.

Auto- and cross-correlations of the Planck DR1 and DR2 κ maps. The similar shape but higher power at all scales in the DR1 autocorrelation suggests that there is more correlated noise in the DR1 map.

4.4 WISE AS- versus AW-selected quasars

In Section 2, many similarities between the AS- and AW-selected quasar samples were discussed. The number densities are similar (three fewer objects per square degree in the AW-selected sample), the obscured and unobscured fractions are indistinguishable (Table 1), the redshift distributions are consistent, and the optical morphological properties are very similar (with a slight increase in unresolved sources in AW). Despite these similarities, the difference in bias measurements when changing from AS- to AW-selected samples (or changing the masks), suggests that there may be some fundamental difference in samples selected from the two catalogues. To help inform the discussion of the bias measurements, we explore the properties of these samples in more detail here.

4.4.1 Photometric properties of AS- and AW-selected quasars

The top row of Fig. 8 shows distributions of W1, W2, and r for the AW- and AS-selected quasars (the full IR-selected samples; these comparisons are very similar for obscured and unobscured subsamples), as well as the difference between the AW and AS W1 and W2 magnitudes for objects that are selected from both catalogues. These comparisons are made after applying both masks to the samples (i.e. AW* and AS*), so any differences are not due to masking (see below). Of the common objects in AS and AW, the vast majority (>99.9 per cent) are matched to the same optical counterpart so differences between r magnitudes for the common sample are not shown as they are generally null.

Top row: comparison of W1, W2, and r distributions for the ALLWISE (solid) and ALLSKY (dashed) selected samples, after applying both masks to both samples (so differences are not due to masking). The top-right panel shows the difference between the W1 (solid) and W2 (dashed) magnitudes for those objects that are common between the samples (r-band magnitudes are excluded since the majority of ALLWISE and ALLSKY sources are matched to the same optical counterpart, and these differences are null). Bottom row: the same, but for colours.
Figure 8.

Top row: comparison of W1, W2, and r distributions for the ALLWISE (solid) and ALLSKY (dashed) selected samples, after applying both masks to both samples (so differences are not due to masking). The top-right panel shows the difference between the W1 (solid) and W2 (dashed) magnitudes for those objects that are common between the samples (r-band magnitudes are excluded since the majority of ALLWISE and ALLSKY sources are matched to the same optical counterpart, and these differences are null). Bottom row: the same, but for colours.

The AW sample shows a subtle shift towards brighter W1 and W2 fluxes, and as noted the r distributions are nearly identical. The reason for this is clear in the top-right panel showing the difference in the common sample. The W1 difference distribution is strongly asymmetric, with an updated AW flux more likely to be brighter, while the W2 difference distribution is much more symmetric about zero and smaller in magnitude. These effects are noted in the AW explanatory supplement, which states that AW photometry is known to be increasingly brighter at magnitudes fainter than W1 ∼ 14 and W2 ∼ 13 (the majority of our sources) due to correction of a faint source underestimation bias in AS. These effects are particularly important for colour-selected samples such as ours.

The second row of Fig. 8 illustrates how these changes propagate through to our colour selection. We see that the larger brightening in W1 compared to W2 leads to an IR-redder sample in AS than in AW. The similarity in the distributions of rW2 illustrates that most of the difference in W1 − W2 colours is driven by the shift in W1.

We also compared the photometric properties of AW-selected quasars that fall within the AS mask (see the middle panel of Fig. 1) to the sources outside of the mask. There are no clear differences that indicate serious problems with the AW objects that fall within the AS mask, and so we do not show them here. Additionally, the obscured and unobscured fractions of these masked sources are consistent with the overall fractions. The same is true for AS-selected quasars that are within the AW mask. These similarities in photometric properties suggest that the fact that the Böotes field is not strongly affected by differences in masking does not have an important impact on our estimates of the redshift distribution.

4.4.2 AS- and AW-only sources

There are 12 544/5776/6632 (total/obscured/unobscured) objects that are only selected by AW and 19 777/9019/10 523 (total/obscured/unobscured) objects that are only selected in AS. These sources are not generally missing from the opposite catalogue, but rather their photometry has been updated such that they no longer meet our selection criteria. Fig. 9 compares the AW and AS photometry and colours for sources that are selected as quasars in only one catalogue (again, only the full samples are shown as the obscured and unobscured trends are similar). In general, objects that do not make the cut in one or the other catalogue are borderline objects, near W2 ∼ 15 or W1 − W2 ∼ 0.8. Most of the difference in the two samples is thus due to the small adjustments to W1 and W2 in AW. The W1 distributions of these objects tend to have a spike around W1 ∼ 15.8, reflective of the sharp W2 cut. On average however, the AW-only sources are somewhat fainter in W1 than the AS-only sources (compare the peaks of either the dashed or solid lines in the top-left and top-right panels). This is also true for W2, as the tail to brighter magnitudes is smaller in the AW-only sources.

A comparison of WISE W1 (top), W2 (middle), and W1 − W2 (bottom) from the ALLWISE and ALLSKY catalogues for objects that satisfy our selection criteria in only one catalogue (right: ALLWISE only, left: ALLSKY only). In the W2 and colour panels the distributions of objects that would not satisfy the other cut (W1 − W2 < 0.08 or W2 < 15.05) are shown in magenta. Clearly, most objects satisfying a cut in one catalogue but not the other are primarily border-line objects in terms of W2, which causes them to also appear generally fainter in W1.
Figure 9.

A comparison of WISE W1 (top), W2 (middle), and W1 − W2 (bottom) from the ALLWISE and ALLSKY catalogues for objects that satisfy our selection criteria in only one catalogue (right: ALLWISE only, left: ALLSKY only). In the W2 and colour panels the distributions of objects that would not satisfy the other cut (W1 − W2 < 0.08 or W2 < 15.05) are shown in magenta. Clearly, most objects satisfying a cut in one catalogue but not the other are primarily border-line objects in terms of W2, which causes them to also appear generally fainter in W1.

Considering that the sources that only meet our selection in AW or AS tend to have borderline photometric properties, can we argue that one catalogue is truly eliminating more contamination, or is this simply just noise near the cuts? In Fig. 10, we plot the relative densities (δ, see Section 3.2) of AS- and AW-only sources as a function of position on the sky, using nside = 2048 healpix pixels smoothed with a 1° Gaussian. We of course expect that intrinsically quasars are uniformly distributed across the field. However, the AS-only selected sources are heavily biased in position, with their density generally increasing towards the Ecliptic plane (the strips masked for Moon contamination are perpendicular to the Ecliptic). AW-only sources are more evenly distributed, though the distribution is still not completely uniform. This suggests that objects selected by AS only may indeed be artefacts, or their selection is biased by a position-dependent factor more-so than in AW (see Section 4.4.3).

The relative densities (δ; see Section 3.2) of ALLSKY-only (right) and ALLWISE-only (left) selected quasars as a function of position. Blue indicates an underdensity relative to the mean, red is overdense. It is clear that the objects selected only by ALLSKY are less evenly distributed, suggesting that many are in fact artefacts or that there is a position-dependent bias affecting their selection. This is greatly reduced, but still present, in ALLWISE.
Figure 10.

The relative densities (δ; see Section 3.2) of ALLSKY-only (right) and ALLWISE-only (left) selected quasars as a function of position. Blue indicates an underdensity relative to the mean, red is overdense. It is clear that the objects selected only by ALLSKY are less evenly distributed, suggesting that many are in fact artefacts or that there is a position-dependent bias affecting their selection. This is greatly reduced, but still present, in ALLWISE.

The clustering and lensing cross-correlation properties of these AS- and AW-only samples may shed some light on their nature, but their position-dependent density complicates this because of the need for a random sample that mimics their distribution. Instead, we measure cross-correlations with the full samples from each catalogue, which are more uniformly distributed and can be normalized with our uniform random sample. This measurement is done via (e.g. Croft, Dalton & Efstathiou 1999)
(20)
where the ‘F’ and ‘O’ subscripts indicate the full and only-in-one catalogue samples, respectively. The results are shown in Fig. 11. For comparison, the autocorrelation measurements for the full AS and AW samples are shown in green in each panel. The AS-only sources clearly show a stronger clustering signal relative to the full sample, increasingly so going from unobscured to obscured objects. There is some indication that the AW-only obscured sources cluster more strongly than the full AW sample, but at much lower significance.
The cross-correlation of objects selected by ALLSKY only (top) and ALLWISE only (bottom) with their respective full samples, after applying both masks to each sample (differences not due to masking). In each panel, the green points show the autocorrelation for the complete sample selected from each catalogue, and the amplitudes of power-law fits to the data are listed in the legend (using a fixed slope of −1). The ALLSKY-only sources clearly cluster more strongly than the full AS sample, and more so for the obscured sample. This is only weakly seen in the AW-only obscured sample. Some of this effect might be explained by the faint objects only selected by one catalogue lying at the higher redshift end of the full-z distribution, but the fact that the behaviour is not the same in both cases suggests additional contamination in the AS sample.
Figure 11.

The cross-correlation of objects selected by ALLSKY only (top) and ALLWISE only (bottom) with their respective full samples, after applying both masks to each sample (differences not due to masking). In each panel, the green points show the autocorrelation for the complete sample selected from each catalogue, and the amplitudes of power-law fits to the data are listed in the legend (using a fixed slope of −1). The ALLSKY-only sources clearly cluster more strongly than the full AS sample, and more so for the obscured sample. This is only weakly seen in the AW-only obscured sample. Some of this effect might be explained by the faint objects only selected by one catalogue lying at the higher redshift end of the full-z distribution, but the fact that the behaviour is not the same in both cases suggests additional contamination in the AS sample.

To quantify this, we fit a simple power-law to the data of the form ωqq(θ) = Aθ−1 (note that we do not fit our model DM autocorrelation to these measurements, as the redshift distributions for the AS- and AW-only samples are not well constrained). The power-law slope of −1 is a typical value for quasar autocorrelations in angular projection (Myers et al. 2006; Shen et al. 2007, 2009; Ross et al. 2009; White et al. 2012, D14), and fits the full sample results here well. The amplitudes of these fits for each sample are given in the figure legend, and highlight the qualitative impression discussed above.

Since objects only selected by one catalogue tend to be faint (see Section 4.4.1), it is likely that at least some of these sources represent the higher redshift end of the distribution, on average, partially explaining their larger clustering signal. However, it is not clear why this would affect one sample more than the other, which suggests that there is some additional signal from contamination present in the AS-only sources that is not present in the AW-only objects. The fact that this is stronger in the obscured subsample is also reflective of why the bias of the obscured sample is affected more significantly by changing samples and masks.

We also cross-correlate the AS- and AW-only samples with the CMB lensing maps, and the results are shown in Fig. 12. The green points in each panel show the result from the full AS and AW samples. To zeroth order, the fact that there is any cross-power here confirms what we found above, that many of these objects are indeed extragalactic, probably at the high-redshift end of our sample. We again fit fixed-slope-power-laws to these results for more direct comparisons, and the amplitudes are given in the legend. Given the amount of noise in these measurements, it is difficult to draw stronger conclusions other than the AS- and AW-only samples have a similar cross-power as the full samples.

A comparison of the CMB lensing cross-correlation for objects selected by ALLSKY only (top) and ALLWISE only (bottom), after applying both masks to each sample (differences not due to masking). In each panel, the green points show the results for the complete sample selected from each catalogue. There is significant noise in all cases due to the small sample sizes, but there is a signal present indicating that these objects are not just artefacts and are indeed extragalactic. However, there is no evidence that the cross-correlation signal is different from the full samples. Without understanding the dndz of these samples, these results are difficult to interpret.
Figure 12.

A comparison of the CMB lensing cross-correlation for objects selected by ALLSKY only (top) and ALLWISE only (bottom), after applying both masks to each sample (differences not due to masking). In each panel, the green points show the results for the complete sample selected from each catalogue. There is significant noise in all cases due to the small sample sizes, but there is a signal present indicating that these objects are not just artefacts and are indeed extragalactic. However, there is no evidence that the cross-correlation signal is different from the full samples. Without understanding the dndz of these samples, these results are difficult to interpret.

4.4.3 Cross-correlations with systematics

To further compare the samples selected from AS and AW, we cross-correlate the quasar density with several systematics, using the method described in Section 3.2. We do this first for the full AW and AS samples, with both masks applied, cross-correlating with W1 and W2 magnitude, W1 − W2, Galactic reddening Ag, and Moon level. The results are shown in Fig. 13.

All panels include a dashed line marking zero correlation as a guide. Top left: the autocorrelation function of the total ALLSKY and ALLWISE samples on large scales, illustrating both the agreement between the two methods used for the calculation (see Sections 3.1 and 3.2) and the fact that the ALLWISE data has less excess signal on large scales. Clockwise from top-centre: cross-correlations of the ALLSKY and ALLWISE data with various systematics – W1 magnitude, W2 magnitude, W1 − W2 colour, extinction in the g band (Ag), and the Moon level in the WISE image tiles. ALLWISE cross-correlations are consistent with zero in most cases, while those with ALLSKY are not, with the exception of the Moon level.
Figure 13.

All panels include a dashed line marking zero correlation as a guide. Top left: the autocorrelation function of the total ALLSKY and ALLWISE samples on large scales, illustrating both the agreement between the two methods used for the calculation (see Sections 3.1 and 3.2) and the fact that the ALLWISE data has less excess signal on large scales. Clockwise from top-centre: cross-correlations of the ALLSKY and ALLWISE data with various systematics – W1 magnitude, W2 magnitude, W1 − W2 colour, extinction in the g band (Ag), and the Moon level in the WISE image tiles. ALLWISE cross-correlations are consistent with zero in most cases, while those with ALLSKY are not, with the exception of the Moon level.

The first panel shows the quasar autocorrelation using the pixelization method (equation 6) as compared to the estimator in equation (3). They agree quite well on scales where they overlap, and the pixelization method allows us to probe larger scales (above ∼2°). On these scales (>50 Mpc h−1 at z = 1), the quasar autocorrelation should approach zero (e.g. Myers et al. 2006; Krumpe et al. 2014). This is true of the AW-selected sample, but the AS sample retains a significant signal even out to ∼7°.

In the absence of systematic observational effects, W1 and W2 should be uncorrelated with the quasars as a function of scale. The next two panels of Fig. 13 show that this is generally true of AW, within errors, but not for AS. This is also true, though to a lesser degree, for W1 − W2.

The next panel shows the cross-correlation of the quasar density with Galactic reddening in the g band, Ag. Again, this is consistent with zero for AW (though there is a slight systematic shift from zero), but not AS. Note that the reddening component of the mask did not change from AS to AW, and so this reflects a change in the data itself. It is also interesting to note the flatness of these cross-correlations, which is most likely a consequence of the fact that the Galactic dust density varies slowly with scale, on average.

Finally, the last panel shows the cross-correlation with the Moon level. Recall that regions with moon_lev >1 in W4 were masked, which does leave some minor Moon contamination possible. In this case, the results are reversed – the AS data do not correlate with the Moon level at all, but AW is anticorrelated. This seems to imply that there may be problems with WISE data calibration due to the Moon present in AW that were not present in AS.

Because the change from AS to AW seems to affect the obscured sample more than the unobscured, we focus in Fig. 14 on cross-correlations of these subsamples of the AW-selected quasars with systematics. The first panel shows agreement with the pair-counting method for calculating the autocorrelation, as well as the fact that both autocorrelations approach zero on larger scales. We omit the cross-correlations with W1, W2, and W1 − W2 as these are null for both subsets of data. However, as seen in the centre panel, obscured sources do correlate with reddening in this sample, while unobscured sources do not. In the final panel, we see that both obscured and unobscured sources correlate in a similar way with the Moon level.

All panels contain a dashed line at zero correlation as a guide. Left: the autocorrelation function of ALLWISE obscured and unobscured quasars, illustrating both the agreement between the two methods (see Sections 3.1 and 3.2) and the lack of signal on large scales for both samples. Center, right: cross-correlation of the samples with Galactic extinction in the g band (Ag; centre), and the Moon level in the WISE imaging tiles (right). The correlation of the data with Ag appears to be primarily driven by the obscured sample, while the correlation with the Moon level is similar for both.
Figure 14.

All panels contain a dashed line at zero correlation as a guide. Left: the autocorrelation function of ALLWISE obscured and unobscured quasars, illustrating both the agreement between the two methods (see Sections 3.1 and 3.2) and the lack of signal on large scales for both samples. Center, right: cross-correlation of the samples with Galactic extinction in the g band (Ag; centre), and the Moon level in the WISE imaging tiles (right). The correlation of the data with Ag appears to be primarily driven by the obscured sample, while the correlation with the Moon level is similar for both.

5 DISCUSSION

5.1 Which measurement is the most reliable?

Based on the analysis of Planck DR1/DR2 and the quasars selected from the WISE AS and AW catalogues, it is clear that the latter in both cases are superior and more reliable data sets. The apparent reduction in correlated noise in Planck DR2, along with the increased lensing signal are obvious improvements. For the WISE selected quasars, the lack of correlations with W1 and W2 magnitudes and Ag, as well as the more uniform distribution of AW-only sources, indicate that this catalogue has certainly improved systematic effects that impact selection and contamination by artefacts. However, Fig. 10 illustrates that this is not completely eliminated with AW, and there is still a slight correlation of obscured quasars with Galactic reddening in AW.

Despite these overall improvements, the correlation of AW quasars with the Moon level is concerning. In order to address this, we perform one final check by repeating our clustering and CMB lensing cross-correlation analysis with any region with moon_lev > 0, in any WISE band, masked out. This removes an additional ∼1000 deg2 of area, and the reduction in sample sizes inflates the uncertainties accordingly. However, we find that the results do not change significantly, though the bias is slightly reduced for both samples (by a few per cent). While it appears that the Moon may have an adverse effect on our selection, we are unable to remove it from our measurement with the current information from WISE.

Given the lack of understanding on how the Moon level is affecting the AW-selected quasars, and the fact that the objects selected only in AW still show a position-dependent relative density, our most conservative approach currently is to use objects selected as quasars by both AW and AS. This naturally includes the mask components of both data sets as well. We adopt the bias and halo masses from this sample as our best current constraints, and these values are summarized in Table 3 for convenience.

Table 3.

Adopted bias and halo mass measurements. The adopted bias and halo masses from the angular autocorrelation (top) and CMB lensing cross-correlations (bottom), using the most conservative sample – those objects satisfying the selection criteria in both ALLSKY and ALLWISE catalogues.

zbqlog (Mh/Mh−1)
Quasar autocorrelation
All IR1.021.92 ± 0.12|$12.79_{-0.11}^{+0.10}$|
Obscured0.982.13 ± 0.21|$13.00_{-0.16}^{+0.14}$|
Unobscured1.051.88 ± 0.15|$12.72_{-0.15}^{+0.13}$|
CMB lensing cross-correlation
All IR1.021.87 ± 0.14|$12.74_{-0.14}^{+0.12}$|
Obscured0.982.06 ± 0.22|$12.94_{-0.18}^{+0.15}$|
Unobscured1.051.72 ± 0.18|$12.56_{-0.21}^{+0.17}$|
zbqlog (Mh/Mh−1)
Quasar autocorrelation
All IR1.021.92 ± 0.12|$12.79_{-0.11}^{+0.10}$|
Obscured0.982.13 ± 0.21|$13.00_{-0.16}^{+0.14}$|
Unobscured1.051.88 ± 0.15|$12.72_{-0.15}^{+0.13}$|
CMB lensing cross-correlation
All IR1.021.87 ± 0.14|$12.74_{-0.14}^{+0.12}$|
Obscured0.982.06 ± 0.22|$12.94_{-0.18}^{+0.15}$|
Unobscured1.051.72 ± 0.18|$12.56_{-0.21}^{+0.17}$|
Table 3.

Adopted bias and halo mass measurements. The adopted bias and halo masses from the angular autocorrelation (top) and CMB lensing cross-correlations (bottom), using the most conservative sample – those objects satisfying the selection criteria in both ALLSKY and ALLWISE catalogues.

zbqlog (Mh/Mh−1)
Quasar autocorrelation
All IR1.021.92 ± 0.12|$12.79_{-0.11}^{+0.10}$|
Obscured0.982.13 ± 0.21|$13.00_{-0.16}^{+0.14}$|
Unobscured1.051.88 ± 0.15|$12.72_{-0.15}^{+0.13}$|
CMB lensing cross-correlation
All IR1.021.87 ± 0.14|$12.74_{-0.14}^{+0.12}$|
Obscured0.982.06 ± 0.22|$12.94_{-0.18}^{+0.15}$|
Unobscured1.051.72 ± 0.18|$12.56_{-0.21}^{+0.17}$|
zbqlog (Mh/Mh−1)
Quasar autocorrelation
All IR1.021.92 ± 0.12|$12.79_{-0.11}^{+0.10}$|
Obscured0.982.13 ± 0.21|$13.00_{-0.16}^{+0.14}$|
Unobscured1.051.88 ± 0.15|$12.72_{-0.15}^{+0.13}$|
CMB lensing cross-correlation
All IR1.021.87 ± 0.14|$12.74_{-0.14}^{+0.12}$|
Obscured0.982.06 ± 0.22|$12.94_{-0.18}^{+0.15}$|
Unobscured1.051.72 ± 0.18|$12.56_{-0.21}^{+0.17}$|

5.2 The redshift distribution and evolution of the bias

There are two approaches we have adopted that can impact the models from equations (2) and (11), and thus the final inferred bias values. The first is in the way we handle the empirical redshift distributions. While the spline interpolation of dN/dz preserves the features seen in Fig. 2, it is possible that some of these are artificial and do not fully reflect the true intrinsic distribution (though the application of our selection criteria to the Böotes field should mitigate some of these issues). Therefore, we also test fit a smoother function to these distributions – a Gaussian with an exponential tail – which results in a smoother overall dN/dz. Propagating these fits through to the models, we find that they shift by at most a few per cent over the scales of interest. The effect is most dramatic for the unobscured sample, due to the spike at low-z, but still only shifts the model by ∼4 per cent, and only at larger scales. Refitting the effective bias, we find that the results shift by ∼1 per cent. Our choice of model for dN/dz does not impact our results significantly.

We have also adopted a simple bias model that does not account for evolution with redshift, and we are really measuring the effective bias, or the bias averaged over our redshift range. Since we lack individual source redshifts, and our error bars are still sufficiently large, we are unable to accurately generate our own empirical b(z) model directly from our data. However, we explore the role of such evolution in our fits, by including in our model calculations the b(z) of Croom et al. (2005): b(z) = 0.53 + 0.289(1 + z)2. We then fit these new models, with an additional scalefactor b0, to our data.

In general, a model of the bias including redshift evolution does not significantly change the quality our fits (based on the values of χ2). In the case of the angular autocorrelation, both obscured and unobscured populations require an additional rescaling, with b0,obsc = 1.34 and b0,unob = 1.26. Including these factors and inserting the mean redshifts of each sample into the model (i.e. b = b0[0.53 + 0.289(1 + 〈z〉)2]), we derive bias values of bobsc = 2.23 ± 0.18 and bunob = 2.21 ± 0.13. These are both increased but roughly consistent with our simpler effective bias measurement, and the values for the two populations are consistent. For the CMB lensing cross-correlations, we measure b0,obsc = 1.28 and b0,unob = 1.04. In this case, the unobscured sample is consistent with the fiducial unobscured Croom et al. (2005) model, while the obscured sample requires additional rescaling. With this model we find bobsc = 2.13 ± 0.14 and bunob = 1.82 ± 0.11, fully consistent with our values assuming a constant bias.

Given that we are unable to currently determine b(z) for these samples directly, the fact that such a model does not improve our fits, and that overall the results are consistent between a constant and evolving bias model, we prefer to adopt the constant bias measurements here as they involve fewer assumptions. Future work will revisit this issue in more detail.

5.3 Comparison with previous results

There have been several measurements of the obscured and unobscured quasar bias and halo masses in the last several years, as it is only recently that samples have grown large enough to do so precisely. There is yet to be convergence on a result.

Our values here are lower than those of D14 and D15 due both to the subtle differences in the samples selected as discussed above, as well as some updates in our procedures. First, the use of camb with its most recent default parameters to produce the model power spectra in a consistent fashion for both the clustering and lensing cross-correlation measurements results in models with slightly larger amplitudes, which naturally decreases the inferred bias. This is seen in the ‘DR1/ALLSKY’ points of Fig. 6, as compared to e.g. fig. 7 of D15. This is a systematic effect that impacts both the obscured and unobscured samples, but not their relative values. However, in Fig. 6 we also illustrate that using the most recent data results in a decrease of the bias as well, for both unobscured and obscured quasars, and more significantly for the clustering measurement. As we argued in the previous section, the new data are quantifiably superior, and these new values should be considered more reliable.

It is difficult however to pinpoint the exact reason for this reduction in the bias. The fact that it is present in both the clustering and CMB lensing cross-correlations, which depend on different systematics, indicates that it is a real effect – AW is picking out objects with lower clustering and lensing cross-correlation amplitudes, especially for the obscured sample. Since the AW W1 − W2 colour distribution is slightly bluer there may be a relationship between IR colour and bias (as opposed to just optical–IR colour, or obscured and unobscured, and bias). However, exploring such a trend will require even larger samples, over extended areas, in order to have large enough subsamples to explore with sufficiently small statistical errors. We will revisit this in future work. The systematic shift in colour in AW may also point to the need for updates to the various WISE colour selection techniques for quasars.

While D14 showed that the factor of 10 larger halo masses for obscured quasars of Donoso et al. (2014) was due to insufficient masking of the WISE data, and found instead a factor of ∼3 difference, that is further reduced here to roughly a factor of 2. However, the significance of this difference is now only ∼2σ.

In Fig. 15 we compare our updated measurements with several from the recent literature. Note that all halo masses in the right-hand panel are calculated using our prescription based on the bias values taken from the respective studies shown on the left (i.e. the halo masses shown may not be the same as those calculated in the original reference). The grey band gives the approximate range of values for optically selected unobscured quasars, primarily from the works referenced in Section 1. We first note that both obscured and unobscured quasars are now roughly consistent with this range, despite the difference in the samples here. While the measurements of Hickox et al. (2011) are at slightly higher redshift, our results remain consistent with theirs. Though our error bars are significantly smaller due to the dramatic increase in area and sample size, the significance of the difference in bias/halo masses is roughly the same once the reduction in the magnitude of the difference is considered.

The adopted bias (left) and halo mass (right) measurements from Table 3 compared to other recent results (comparison with the D14/D15 results can be seen in Fig. 6). Points have been shifted slightly in redshift where necessary for clarity. The grey band represents the range of results typical for optically selected unobscured quasars, largely from the SDSS. The Mendez et al. (2015) results utilize individual redshifts for a projected clustering measurement – note that their halo mass measurements fall outside of the plot range (though the obscured halo mass error bar can be seen), so the actual values are listed. While the Hickox et al. (2011) results agree well with what we find here, the Mendez et al. (2015) estimates are significantly lower. Obscured quasars have halo masses ∼2 times larger than unobscured quasars, with a significance of ∼2σ.
Figure 15.

The adopted bias (left) and halo mass (right) measurements from Table 3 compared to other recent results (comparison with the D14/D15 results can be seen in Fig. 6). Points have been shifted slightly in redshift where necessary for clarity. The grey band represents the range of results typical for optically selected unobscured quasars, largely from the SDSS. The Mendez et al. (2015) results utilize individual redshifts for a projected clustering measurement – note that their halo mass measurements fall outside of the plot range (though the obscured halo mass error bar can be seen), so the actual values are listed. While the Hickox et al. (2011) results agree well with what we find here, the Mendez et al. (2015) estimates are significantly lower. Obscured quasars have halo masses ∼2 times larger than unobscured quasars, with a significance of ∼2σ.

We also show the recent results of Mendez et al. (2015), at slightly lower redshift. Their sample used the IR selection technique of Assef et al. (2013), which extends deeper to W2 < 17.11 with a magnitude-dependent colour criteria. They also used samples with individual redshifts from several fields, and we show their results with the COSMOS field excluded, due to this region containing structures that are particularly overdense as compared to the cosmic average (e.g. Lilly et al. 2007). The deeper W2 limit likely results in a slightly lower average luminosity, though the requirement of spectroscopic redshifts somewhat offsets this effect. Note that the Mendez et al. (2015) halo masses are quite low (the lowest of all the samples they considered) compared to ours and those of other groups – more than two orders of magnitude in the obscured quasar case. Aside from the small halo masses in general, the sense of the difference between obscured and unobscured samples is the opposite of what we find here and was found in Hickox et al. (2011).

The individual redshifts of Mendez et al. (2015) should reduce systematic errors due to differences in the redshift distributions that could be present in our measurements, despite the improved statistical errors we achieve with a larger sample. However, our bias values for the obscured and unobscured samples differ by ∼10 per cent. As shown in fig. 8 of D15, our effective redshift estimates would need to be systematically offset by ∼0.25 to fully account for this difference. It is unlikely that the Böotes field is this poorly representative of the full population. On the other hand, considering that the bias is not weighted equally at each redshift, the presence of a high-z tail in the obscured sample could skew the bias significantly. Such a tail would most likely be present in the sources lacking SDSS counterparts. As a check we repeat our clustering and CMB lensing cross-correlation for only the obscured sources with r-band detections. These objects do not have a significantly different redshift distribution in Böotes, and as in D14/D15, the measurement is completely consistent with that of the entire obscured sample. While incorporating full redshift information into our measurements would certainly be an improvement, it is difficult to see how redshift distribution errors would dramatically impact our findings.

5.4 Interpretation and future work

While the magnitude and significance of the difference in bias and halo mass between obscured and unobscured quasars is reduced in this work, it is still present. The simplest interpretation of this is that obscured quasars reside in higher mass haloes, and are not simple analogues of unobscured quasars seen through a dusty torus. They could instead represent different phases common to the quasar phase of black holes.

It is also possible that our obscured sample is a mix of genuinely obscured objects and lower luminosity unobscured AGN. This ‘contamination’ could imply either that lower luminosity AGN are more biased (which is unexpected, and previous results regarding trends of bias with luminosity have been weak at best; Shen et al. 2009, 2013; Eftekharzadeh et al. 2015; Krolewski & Eisenstein 2015), or that the signal from the true obscured population is being diluted. In addition, it is likely that some fraction of quasars are in fact obscured only due to our particular line of sight, and are intrinsically the same as the unobscured population if seen from a different angle. This would also serve to dilute the bias measurement of sources obscured by other means. Both of these factors imply that the difference we are finding between obscured and unobscured quasars may be a lower limit.

We again use abundance-matching techniques (e.g. Colín et al. 1999; Kravtsov & Klypin 1999; Vale & Ostriker 2004; Shankar et al. 2006; Guo et al. 2010) to estimate the implied lifetimes of obscured and unobscured phases, assuming there is an evolutionary trend. The median Lbol ∼ 1046 erg s−1 (Hickox et al. 2011; Hainline et al. 2014) implies an overall space density at z = 1 of ∼2 × 10−5 Mpc−3 (Hopkins, Richards & Hernquist 2007). This is split by roughly 60 and 40 per cent unobscured and obscured, respectively. For the halo masses measured here (we adopt the results from CMB lensing cross-correlations), and the halo mass function of Tinker et al. (2010), we find host halo space densities of |${\rm d}n/{\rm d}\log _{10}(M) = 8.1_{-5.3}^{+2.8} \times 10^{-4}$| Mpc−3 (unobscured) and |$3.0_{-1.8}^{+1.0} \times 10^{-4}$| Mpc−4 (obscured). Over 0.5 < z < 1.5, or roughly 4 Gyr of cosmic time, this implies lifetimes of |$70_{-37}^{+27}$| Myr for the unobscured phase and |$100_{-52}^{+38}$| Myr for the obscured phase. Each is on the order of 1 per cent of the Hubble time, with the obscured phase lasting about 1.5 times longer than the unobscured phase.

Given the reduced significance of the difference in halo masses found here, it is possible that the large-scale bias of obscured and unobscured quasars is consistent with a pure orientation model. However, there are other indications that even if these sources occupy similar mass haloes overall, there are other differences. For example, the smaller scale signal (below 1 Mpc h−1, or 0| $_{.}^{\circ}$|04 at our redshifts) appears to differ between the samples in both the clustering and CMB-lensing cross-correlations. This could imply a difference in the HODs and in particular the satellite fractions of each population (Zheng et al. 2005; Zheng & Weinberg 2007; White et al. 2012; Chatterjee et al. 2013). However, studying the HOD of obscured quasars in detail is difficult without full redshift information, but such additional information could provide higher-order correlation functions or direct measurement of the mean occupation function of obscured quasars (Chatterjee et al. 2013).

Our next steps are to leverage the AS nature of WISE along with the full footprint of SDSS and other optical surveys (e.g. The Dark Energy Survey) to build the largest IR-selected obscured and unobscured quasar samples possible. This will provide dramatically improved statistical error bars on these bias measurements, to further constrain the potential difference in halo masses. We will also use photometric redshift estimations from multiwavelength photometry and SED fitting (e.g. Hainline et al. 2014; Carroll & Hickox 2015; DiPompeo et al. 2015b) to explore the redshift evolution of obscured quasars.

6 SUMMARY

Using the most recent AW source catalogue from WISE, along with new products from Planck, we present updated measurements of the obscured and unobscured quasar bias via angular autocorrelations and cross-correlation with CMB lensing. We find, as with our previous work, that obscured quasars have a larger bias and therefore halo mass and longer lifetime. However, this difference is reduced with respect to our results using previous data products from WISE and Planck. The inferred typical halo mass of obscured quasars is roughly a factor of 2 larger than that of unobscured quasars, at a significance of ∼1σ–2σ.

In order to explain this change, we have carefully compared the properties of quasars selected from both WISE catalogues using standard colour-cuts. The general properties of these sources – morphology, redshift distribution, unobscured/obscured fractions – are indistinguishable. However, the photometric properties do differ slightly. In particular, the more recent WISE catalogue tends to select brighter and bluer quasars, which could contribute to the change in our results.

We have also carefully explored systematic effects in WISE-selected quasars from both catalogues, and find that the updated version is clearly superior. In particular, the quasar selection in the AW catalogue appears less biased with respect to flux measurements in W1 and W2, as well as with Galactic extinction. However, there may be some complications with respect to Moon contamination in AW, which leads us to rely on the sample that is selected from both catalogues for our current measurement.

Finally, we make available with this paper several sets of codes for making these measurements and generating consistent models.

MAD, RCH, and ADM were partially supported by NASA through ADAP award NNX12AE38G and by the National Science Foundation through grant number 1211096. MAD and ADM were also supported by NSF grant numbers 1211112 and 1515404, and MAD and RHC were also supported by NSF grant number 1515364. RCH also acknowledges support from an Alfred P. Sloan Research Fellowship, and a Dartmouth Class of 1962 Faculty Fellowship.

1

Note that throughout this paper type 1 and type 2 do not refer to spectral classifications based on broad and narrow emission lines, but to optical-to-infrared colours as described in more detail in Section 2.2.

2

See the WISE ALLSKY Release Explanatory Supplement, http://wise2.ipac.caltech.edu/docs/release/allsky/expsup/.

3

See the ALLWISE Release Explanatory Supplement, http://wise2.ipac.caltech.edu/docs/release/allwise/expsup/.

4

This data mask, along with several other files including our samples, is available as a mangle polygon file at http://faraday.uwyo.edu/admyers/wisemask2015/wisemask.html.

8

Code for Anisotropies in the Microwave Background (http://lambda.gsfc.nasa.gov/toolbox/tb_camb_ov.cfm). We use our idl wrapper for camb (camb4idl) which can be found at https://github.com/mdipompe/camb4idl. Our initialization file for camb is available with the other autocorrelation codes linked above, as well as code to parse the camb output and generate the model autocorrelation.

9

Note that here we use the more recent Tinker et al. (2010) model, while in D14/D14 we used the Sheth, Mo & Tormen (2001) model for more direct comparison with previous results. In the context of searching for differences in obscured and unobscured halo masses (the main goal of D14/D15), the choice of model is less critical; however, moving forward we prefer the most precise halo masses for future modelling purposes

10

Note that there is a decrease in the measurements presented here compared to D14/D15 due to our updated model procedure as well as the more restricted fitting range.

REFERENCES

Abazajian
K. N.
et al.
ApJS
2009
182
543

Agarwal
N.
Ho
S.
Shandera
S.
J. Cosmol. Astropart. Phys.
2014
02
038

Alexander
D. M.
Hickox
R. C.
New Astron. Rev.
2012
56
93

Antonucci
R.
ARA&A
1993
31
473

Assef
R. J.
et al.
ApJ
2013
772
26

Assef
R. J.
et al.
ApJ
2015
804
27

Bardeen
J. M.
Bond
J. R.
Kaiser
N.
Szalay
A. S.
ApJ
1986
304
15

Becker
R. H.
White
R. L.
Helfand
D. J.
ApJ
1995
450
559

Berlind
A. A.
Weinberg
D. H.
ApJ
2002
575
587

Berlind
A. A.
et al.
ApJ
2003
593
1

Blanton
M. R.
Roweis
S.
AJ
2007
133
734

Bleem
L. E.
et al.
ApJ
2012
753
L9

Booth
C. M.
Schaye
J.
MNRAS
2010
405
L1

Brodwin
M.
et al.
ApJ
2006
651
791

Carroll
C. M.
Hickox
R. C.
AGN versus Star Formation: The Fate of the Gas in Galaxies
2015
Durham University
A2

Chatterjee
S.
Nguyen
M. L.
Myers
A. D.
Zheng
Z.
ApJ
2013
779
147

Coil
A. L.
Hennawi
J. F.
Newman
J. A.
Cooper
M. C.
Davis
M.
ApJ
2007
654
115

Colín
P.
Klypin
A. A.
Kravtsov
A. V.
Khokhlov
A. M.
ApJ
1999
523
32

Comastri
A.
Setti
G.
Zamorani
G.
Hasinger
G.
A&A
1995
296
1

Cooray
A.
Hu
W.
ApJ
2000
534
533

Croft
R. A. C.
Dalton
G. B.
Efstathiou
G.
MNRAS
1999
305
547

Croom
S. M.
Smith
R. J.
Boyle
B. J.
Shanks
T.
Miller
L.
Outram
P. J.
Loaring
N. S.
MNRAS
2004
349
1397

Croom
S. M.
et al.
MNRAS
2005
356
415

Croton
D. J.
MNRAS
2009
394
1109

da Ângela
J.
et al.
MNRAS
2008
383
565

Das
S.
et al.
Phys. Rev. Lett.
2011
107
21301

DiPompeo
M. A.
Myers
A. D.
Hickox
R. C.
Geach
J. E.
Hainline
K. N.
MNRAS
2014
442
3443
(D14)

DiPompeo
M. A.
Myers
A. D.
Hickox
R. C.
Geach
J. E.
Holder
G.
Hainline
K. N.
Hall
S. W.
MNRAS
2015a
446
3492
(D15)

DiPompeo
M. A.
Bovy
J.
Myers
A. D.
Lang
D.
MNRAS
2015b
452
3124

Donley
J. L.
Rieke
G. H.
Pérez-González
P. G.
Rigby
J. R.
Alonso-Herrero
A.
ApJ
2007
660
167

Donley
J. L.
et al.
ApJ
2012
748
142

Donoso
E.
Yan
L.
Stern
D.
Assef
R. J.
ApJ
2014
789
44

Eftekharzadeh
S.
et al.
MNRAS
2015
453
2779

Eisenstein
D. J.
Hu
W.
ApJ
1999
511
5

Ellison
S. L.
Patton
D. R.
Hickox
R. C.
MNRAS
2015
451
L35

Elvis
M.
et al.
ApJS
1994
95
1

Fan
X.
et al.
AJ
2006
131
1203

Fitzpatrick
E. L.
Massa
D.
ApJ
2009
699
1209

Geach
J. E.
et al.
ApJ
2013
776
L41

Górski
K. M.
Hivon
E.
Banday
A. J.
Wandelt
B. D.
Hansen
F. K.
Reinecke
M.
Bartelmann
M.
ApJ
2005
622
759

Goulding
A. D.
Alexander
D. M.
Bauer
F. E.
Forman
W. R.
Hickox
R. C.
Jones
C.
Mullaney
J. R.
Trichas
M.
ApJ
2012
755
5

Goulding
A. D.
et al.
ApJ
2014
783
40

Guo
Q.
White
S.
Li
C.
Boylan-Kolchin
M.
MNRAS
2010
404
1111

Hainline
K. N.
Hickox
R. C.
Carroll
C. M.
Myers
A. D.
DiPompeo
M. A.
Trouille
L.
ApJ
2014
795
124

Hamilton
A. J. S.
Tegmark
M.
MNRAS
2004
349
115

Helfand
D. J.
White
R. L.
Becker
R. H.
ApJ
2015
801
26

Hewett
P. C.
Foltz
C. B.
Chaffee
F. H.
AJ
1995
109
1498

Hickox
R. C.
et al.
ApJ
2007
671
1365

Hickox
R. C.
et al.
ApJ
2011
731
117

Ho
S.
et al.
ApJ
2012
761
14

Ho
S.
et al.
J. Cosmol. Astropart. Phys.
2015
05
040

Høg
E.
et al.
A&A
2000
355
L27

Hopkins
P. F.
Richards
G. T.
Hernquist
L.
ApJ
2007
654
731

Hopkins
P. F.
Hernquist
L.
Cox
T. J.
Kereš
D.
ApJS
2008
175
356

Hu
W.
ApJ
2001
557
L79

Huterer
D.
Takada
M.
Bernstein
G.
Jain
B.
MNRAS
2006
366
101

Kaiser
N.
ApJ
1984
284
L9

Kochanek
C. S.
et al.
ApJS
2012
200
8

Komatsu
E.
et al.
ApJS
2011
192
18

Kravtsov
A. V.
Klypin
A. A.
ApJ
1999
520
437

Krolewski
A. G.
Eisenstein
D. J.
ApJ
2015
803
4

Krumpe
M.
Miyaji
T.
Coil
A. L.
ApJ
2010
713
558

Krumpe
M.
Miyaji
T.
Coil
A. L.
Acta Polytechnica CTU Proceedings, 1, 71
2014

Lacy
M.
et al.
ApJS
2004
154
166

Lacy
M.
et al.
ApJS
2013
208
24

Lacy
M.
Ridgway
S. E.
Sajina
A.
Petric
A. O.
Gates
E. L.
Urrutia
T.
Storrie-Lombardi
L. J.
ApJ
2015
802
102

Landy
S. D.
Szalay
A. S.
ApJ
1993
412
64

Leistedt
B.
Peiris
H. V.
MNRAS
2014
444
2

Lewis
A.
Challinor
A.
Lasenby
A.
ApJ
2000
538
473

Lilly
S. J.
et al.
ApJS
2007
172
70

Limber
D. N.
ApJ
1953
117
134

Mainzer
A.
et al.
ApJ
2011
743
156

Mateos
S.
et al.
MNRAS
2012
426
3271

Mateos
S.
Alonso-Herrero
A.
Carrera
F. J.
Blain
A.
Severgnini
P.
Caccianiga
A.
Ruiz
A.
MNRAS
2013
434
941

Mendez
A. J.
et al.
2015
preprint (arXiv:e-prints)

Myers
A. D.
Outram
P. J.
Shanks
T.
Boyle
B. J.
Croom
S. M.
Loaring
N. S.
Miller
L.
Smith
R. J.
MNRAS
2005
359
741

Myers
A. D.
et al.
ApJ
2006
638
622

Myers
A. D.
Brunner
R. J.
Nichol
R. C.
Richards
G. T.
Schneider
D. P.
Bahcall
N. A.
ApJ
2007
658
85

Myers
A. D.
et al.
2015
preprint (arXiv:1508.04472)

Navarro
J. F.
Frenk
C. S.
White
S. D. M.
ApJ
1997
490
493

Netzer
H.
ARA&A
2015
53
365

Padmanabhan
N.
White
M.
Norberg
P.
Porciani
C.
MNRAS
2009
397
1862

Peacock
J. A.
MNRAS
1991
253
1

Peebles
P. J. E.
The Large-Scale Structure of the Universe
1980
Princeton, NJ
Princeton Univ. Press

Planck Collaboration
2015a
preprint (arXiv:e-prints)

Planck Collaboration
2015b
preprint (arXiv:e-prints)

Planck Collaboration I
A&A
2011
536
A1

Planck Collaboration I
A&A
2014
571
A1

Planck Collaboration XVII
A&A
2014
571
A17

Porciani
C.
Magliocchetti
M.
Norberg
P.
MNRAS
2004
355
1010

Richards
G. T.
et al.
MNRAS
2005
360
839

Richards
G. T.
et al.
ApJS
2006
166
470

Richardson
J.
Zheng
Z.
Chatterjee
S.
Nagai
D.
Shen
Y.
ApJ
2012
755
30

Ross
N. P.
et al.
ApJ
2009
697
1634

Ross
A. J.
et al.
MNRAS
2011
417
1350

Sabra
B. M.
Saliba
C.
Akl
M. A.
Chahine
G.
ApJ
2015
803
5

Sanders
D. B.
Soifer
B. T.
Elias
J. H.
Madore
B. F.
Matthews
K.
Neugebauer
G.
Scoville
N. Z.
ApJ
1988
325
74

Satyapal
S.
Ellison
S. L.
McAlpine
W.
Hickox
R. C.
Patton
D. R.
Mendel
J. T.
MNRAS
2014
441
1297

Schlegel
D. J.
Finkbeiner
D. P.
Davis
M.
ApJ
1998
500
525

Scranton
R.
et al.
ApJ
2002
579
48

Secrest
N.
Dudik
R.
Dorland
B.
Zacharias
N.
Makarov
V.
Fey
A.
Frouard
J.
Finch
C.
ApJS
2015
221
12

Seljak
U.
Zaldarriaga
M.
Phys. Rev. Lett.
1999
82
2636

Setti
G.
Woltjer
L.
A&A
1989
224
L21

Shankar
F.
et al.
ApJ
2006
643
14

Shen
Y.
et al.
AJ
2007
133
2222

Shen
Y.
et al.
ApJ
2009
697
1656

Shen
Y.
et al.
ApJ
2013
778
98

Sherwin
B. D.
et al.
Phys. Rev. D
2012
86
83006

Sheth
R. K.
Mo
H. J.
Tormen
G.
MNRAS
2001
323
1

Smith
K. L.
Koss
M.
Mushotzky
R. F.
ApJ
2014
794
112

Song
Y.-S.
Song
Y.-S.
Cooray
A.
Cooray
A.
Knox
L.
Knox
L.
Zaldarriaga
M.
Zaldarriaga
M.
ApJ
2003
590
664

Stern
D.
et al.
ApJ
2005
631
163

Stern
D.
et al.
ApJ
2012
753
30

Stoughton
C.
et al.
AJ
2002
123
485

Swanson
M. E. C.
Tegmark
M.
Hamilton
A. J. S.
Hill
J. C.
MNRAS
2008
387
1391

Tinker
J. L.
Robertson
B. E.
Kravtsov
A. V.
Klypin
A.
Warren
M. S.
Yepes
G.
Gottlober
S.
ApJ
2010
724
878

Totsuji
H.
Kihara
T.
PASJ
1969
21
221

Vale
A.
Ostriker
J. P.
MNRAS
2004
353
189

van Engelen
A.
et al.
ApJ
2012
756
142

Volonteri
M.
Natarajan
P.
Gültekin
K.
ApJ
2011
737
50

Warren
W. H. J.
Hoffleit
D.
BAAS
1987
19
733

Werner
M. W.
et al.
ApJS
2004
154
1

White
M.
et al.
ApJ
2011
728
126

White
M.
et al.
MNRAS
2012
424
933

Wright
E. L.
et al.
AJ
2010
140
1868

York
D. G.
et al.
AJ
2000
120
1579

Zheng
Z.
Weinberg
D. H.
ApJ
2007
659
1

Zheng
Z.
et al.
ApJ
2005
633
791