THE NuSTAR EXTRAGALACTIC SURVEY: FIRST DIRECT MEASUREMENTS OF THE ≳10 keV X-RAY LUMINOSITY FUNCTION FOR ACTIVE GALACTIC NUCLEI AT z > 0.1

, , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , , and

Published 2015 December 10 © 2015. The American Astronomical Society. All rights reserved.
, , Citation J. Aird et al 2015 ApJ 815 66 DOI 10.1088/0004-637X/815/1/66

0004-637X/815/1/66

ABSTRACT

We present the first direct measurements of the rest-frame 10–40 keV X-ray luminosity function (XLF) of active galactic nuclei (AGNs) based on a sample of 94 sources at 0.1 < z < 3, selected at 8–24 keV energies from sources in the Nuclear Spectroscopic Telescope Array (NuSTAR) extragalactic survey program. Our results are consistent with the strong evolution of the AGN population seen in prior, lower-energy studies of the XLF. However, different models of the intrinsic distribution of absorption, which are used to correct for selection biases, give significantly different predictions for the total number of sources in our sample, leading to small, systematic differences in our binned estimates of the XLF. Adopting a model with a lower intrinsic fraction of Compton-thick sources and a larger population of sources with column densities ${N}_{{\rm{H}}}\sim {10}^{23-24}$ cm−2 or a model with stronger Compton reflection component (with a relative normalization of R ∼ 2 at all luminosities) can bring extrapolations of the XLF from 2–10 keV into agreement with our NuSTAR sample. Ultimately, X-ray spectral analysis of the NuSTAR sources is required to break this degeneracy between the distribution of absorbing column densities and the strength of the Compton reflection component and thus refine our measurements of the XLF. Furthermore, the models that successfully describe the high-redshift population seen by NuSTAR tend to over-predict previous, high-energy measurements of the local XLF, indicating that there is evolution of the AGN population that is not fully captured by the current models.

Export citation and abstract BibTeX RIS

1. INTRODUCTION

Measurements of the luminosity function of active galactic nuclei (AGNs) provide the key observational data available to track the history and distribution of accretion onto supermassive black holes. X-ray surveys have been crucial for performing such measurements as they can efficiently identify AGNs down to low luminosities, over a wide range of redshifts, and where the central regions are obscured by large amounts of gas and dust (see Brandt & Alexander 2015, for a recent review).

A number of studies have presented measurements of the X-ray luminosity function (XLF) of AGNs based on surveys with the Chandra or XMM-Newton X-ray observatories (e.g., Ueda et al. 2003; Barger et al. 2005; Hasinger et al. 2005; Aird et al. 2010; Miyaji et al. 2015). These studies show that the AGN population has evolved substantially over cosmic time. Both the space density of AGNs and the overall accretion density (which traces the total rate of black hole growth) peaked at z ∼ 1−2 and have declined ever since. The evolution also has a strong luminosity dependence, whereby the space density of luminous AGNs peaks at z ≈ 2 whereas lower luminosity AGNs peak at later cosmic times (z ≈ 1). This pattern, reflected in the evolution of the shape of the XLF, provides crucial insights into the underlying distributions of black hole masses and accretion rates (e.g., Aird et al. 2013; Shankar et al. 2013).

A major issue in XLF studies with Chandra and XMM-Newton is the impact of absorption. It is now well established that many AGNs are surrounded by gas and dust that obscures their emission at certain wavelengths (e.g., Martínez-Sansigre et al. 2005; Stern et al. 2005; Tozzi et al. 2006). Soft X-rays (∼0.5–2 keV) will be absorbed by gas with equivalent neutral hydrogen column densities ${N}_{{\rm{H}}}\gtrsim {10}^{22}$ cm−2, whereas higher energy X-rays can penetrate larger column densities. Thus, samples selected at ∼2–10 keV energies are less biased and are typically adopted in XLF studies (e.g., Aird et al. 2010; Miyaji et al. 2015). However, absorption can still suppress the observable flux, especially at higher column densities (${N}_{{\rm{H}}}\gtrsim {10}^{23}$ cm−2) and at lower redshifts (z ≲ 1), where the observed band probes lower rest-frame energies. The emission from Compton-thick sources (${N}_{{\rm{H}}}\gtrsim {10}^{24}$ cm−2) is even more strongly suppressed, although such sources may still be identified at ∼0.5–10 keV energies by their scattered emission and signatures of reflection from the obscuring material (e.g., Georgantopoulos et al. 2013; Brightman et al. 2014). Multiwavelength studies can also identify the signatures of heavy obscuration in X-ray detected AGNs (e.g., Cappi et al. 2006; Alexander et al. 2008; Gilli et al. 2010; Georgantopoulos et al. 2011) or directly identify additional, obscured AGN populations that are not detected in the X-ray band (e.g., Donley et al. 2008; Juneau et al. 2011; Eisenhardt et al. 2012; Del Moro et al. 2015).

Correcting for absorption to accurately recover the XLF is challenging. Luminosity estimates for individual sources can be corrected based on X-ray spectral analysis or hardness ratios (e.g., Ueda et al. 2003; La Franca et al. 2005; Buchner et al. 2015). Alternatively, less-direct statistical approaches can be adopted to constrain the underlying distribution of NH and account for the impact on observed X-ray samples (e.g., Aird et al. 2015, hereafter A15; Miyaji et al. 2015). Despite much progress, the distribution of NH as a function of luminosity and redshift (hereafter, "the NH function") and, most crucially, the fraction of Compton-thick sources remain poorly constrained outside the local Universe.

The cosmic X-ray background (CXB) provides an additional, integral constraint on the fraction of heavily obscured and Compton-thick AGNs. At low energies (≲8 keV), Chandra surveys have resolved the majority (∼70%–90%) of the CXB into discrete point sources, predominantly unabsorbed and moderately absorbed AGNs at z ∼ 0.5–2 (e.g., Worsley et al. 2005; Georgakakis et al. 2008; Lehmer et al. 2012; Xue et al. 2012). However, the peak of the CXB occurs at much higher energies, ∼20–30 keV (e.g., Marshall et al. 1980; Churazov et al. 2007; Ajello et al. 2008). Population synthesis models, based on studies of the XLF and NH function at lower energies, attribute ≳70% of the emission at this peak to absorbed AGNs (${N}_{{\rm{H}}}\gtrsim {10}^{22}$ cm−2 e.g., Gilli et al. 2007; Draper & Ballantyne 2009; Treister et al. 2009), although the required fraction of Compton-thick (${N}_{{\rm{H}}}\gtrsim {10}^{24}$ cm−2) AGNs is still uncertain (e.g., Ballantyne et al. 2011; Akylas et al. 2012, A15).

Until recently, only ∼1%–2% of the CXB at ≳10 keV could be directly resolved into individual objects, due to the limited sensitivity at these energies achieved with non-focusing X-ray missions such as INTEGRAL or Swift (e.g., Krivonos et al. 2007; Tueller et al. 2008; Ajello et al. 2012). Nonetheless, AGN samples identified by these missions—which are less biased against obscured sources—have enabled crucial measurements of the XLF, the NH function, and the fraction of Compton-thick AGNs in the local (z ≲ 0.1) universe (e.g., Beckmann et al. 2009; Burlon et al. 2011; Vasudevan et al. 2013), which are often used to extrapolate to higher redshifts (e.g., Ueda et al. 2014, hereafter U14).

The Nuclear Spectroscopic Telescope Array (hereafter NuSTAR, Harrison et al. 2013) is the first orbiting observatory with >10 keV focusing optics, providing ∼2 orders of magnitude increase in sensitivity compared to previous high-energy observatories. One of the primary objectives of the NuSTAR mission is to identify and characterize the source populations that produce the peak of the CXB. To this end, NuSTAR is executing a multi-layered program of extragalactic surveys. Source catalogs and initial results from the dedicated surveys of the COSMOS and Extended Chandra Deep Field South (ECDFS) regions are presented by Civano et al. (2015, hereafter C15) and Mullaney et al. (2015, hereafter M15), respectively, while Alexander et al. (2013) presented the first results from our ongoing serendipitous survey program. Harrison et al. (2015) present source number counts at 3–8 keV and 8–24 keV energies from the full survey program and show that NuSTAR is directly resolving ∼35% of the CXB emission at 8–24 keV, a factor ∼15–30 times more than previous high-energy X-ray observatories.

In this paper, we present the first measurements of the XLF of AGNs at 0.1 < z < 3 based on direct selection of sources at hard (>8 keV) energies from across the NuSTAR extragalactic survey program. Section 2 describes our data and defines our sample. Section 3 describes our statistical methods to estimate intrinsic luminosities and recover the XLF. In Section 4, we present our measurements of the 10–40 keV XLF and explore the effects of different model assumptions. We discuss our results and future prospects for the NuSTAR survey program in Section 5 and summarize our findings in Section 6. We adopt a flat cosmology with ΩΛ = 0.7 and H0 = 70 km s−1 Mpc−1 throughout this paper.

2. DATA AND SAMPLE SELECTION

The NuSTAR extragalactic survey program (see Harrison et al. 2015, for an overview) consists of three components: (1) a deep (∼400 ks) survey covering both the ECDFS (M15) and Extended Groth Strip (EGS: J. Aird et al. 2016, in preparation) regions; (2) a medium depth (∼100 ks) survey covering the COSMOS field (C15); and (3) a wide-area program searching for serendipitous detections across all NuSTAR observations (Alexander et al. 2013, C. Fuentes et al. 2016, in preparation, G. B. Lansbury et al. 2016, in preparation). In this paper we select sources from across the NuSTAR extragalactic survey program that are directly detected at 8–24 keV energies. NuSTAR provides unprecedented sensitivity at these energies, although the sensitivity of this band is dominated by ≲12 keV energies due to a combination of decreasing effective area and the decrease in source photon flux with increasing energy.

Our overall sample consists of 97 sources. We identify lower-energy X-ray counterparts to the vast majority of these sources (93 out of 97) and identify optical or infrared counterparts (matching to the low-energy X-ray position, if available) for all but one source (NuSTAR J033122-2743.9 in the ECDFS, discussed further below). Reliable (spectroscopic or photometric) redshift estimates are obtained for 91 out of our 97 sources. Additional details are given in Sections 2.1 and 2.2 below, with full details and catalogs provided by M15, C15, or G. B. Lansbury et al. (2016, in prepartation). Figure 1 shows the distribution of the rest-frame 10–40 keV luminosities versus redshift for our sample. Luminosities in this plot are estimated from the 8–24 keV count rates assuming an unabsorbed X-ray spectrum with a photon index Γ = 1.9 and a reflection component with a relative normalization of R = 1, folded through the NuSTAR response. Section 3.1 gives a detailed description of our spectral model and the uncertainties in these luminosity estimates, which are accounted for in our analysis of the XLF. For our estimates of the XLF we use sources in the redshift range 0.1 < z < 3, resulting in a sample of 94 sources (which includes an additional six sources with indeterminate redshifts).

Figure 1.

Figure 1. Rest-frame 10–40 keV X-ray luminosity (not corrected for absorption) vs. redshift for sources in our 8–24 keV selected sample from the various NuSTAR survey components, as indicated. Luminosities are estimated from the 8–24 keV observed fluxes (see Sections 2 and 3.1 for details). The vertical dotted lines show the limits of the redshift bins adopted for our XLF estimates in Section 4 (sources in the shaded regions are excluded). The blue dashed line indicates the characteristic break in the XLF, ${L}_{*},$ based on the 2–10 keV XLF measured by A15.

Standard image High-resolution image

2.1. Dedicated Survey Fields (ECDFS, EGS, and COSMOS)

In the dedicated survey fields (ECDFS, EGS, and COSMOS), we adopt a consistent source detection procedure. We generate background maps using the nuskybgd code (Wik et al. 2014), which we use to calculate the probability that the image counts in a 20'' aperture were produced by a spurious fluctuation of the background (hereafter, the "false probability") at the position of every pixel. We then identify positions where the false probability falls below a set threshold and thus the observed counts can be associated with a real source. See M15 and C15 for full details.

Our final sample of 8–24 keV detections consists of 19, 13, and 32 sources in the ECDFS, EGS, and COSMOS fields, respectively. We identify lower-energy X-ray counterparts to our sources in the deep Chandra or XMM-Newton imaging of these fields (Lehmer et al. 2005; Cappelluti et al. 2009; Puccetti et al. 2009; Xue et al. 2011; Nandra et al. 2015) for all but one of our sources (NuSTAR J033122-2743.9 in the ECDFS, discussed by M15). All of the low-energy counterparts have multiwavelength identifications, as given by M15 and C15 and references therein (for the ECDFS and COSMOS fields) or Nandra et al. (2015) for the EGS field. A very high fraction of these counterparts (86%) have available spectroscopic redshifts; for the remainder we adopt the best photometric redshift estimate given by M15, C15, or Nandra et al. (2015). For two sources in the ECDFS (NuSTAR J033212-2752.3 and NuSTAR J033243-2738.3), which lack photometric estimates in the M15 catalog, we adopt photometric redshifts from Hsu et al. (2014). We retain NuSTAR J033122-2743.9 (which lacks a low-energy or multiwavelength counterpart) in our sample but assume no knowledge of its redshift, effectively adopting a p(z) distribution that is constant in $\mathrm{log}(1+z)$ (see A15 and further discussion in Section 3.1 below). We note that this source could be a spurious detection (and would be consistent with the expected spurious fraction, given our false probability thresholds); this possibility is allowed for by our statistical methodology to determine the XLF.

2.2. Serendipitous Survey

For the serendipitous survey, we adopt the same source detection procedure as in the dedicated fields but use a different method to determine the background maps since in many of the fields a bright target contaminates ∼10%–80% of the NuSTAR field of view. We thus take the original X-ray images and measure the counts in annular apertures of inner radius 30'' and outer radius 90'' centered at each pixel position. We rescale the counts within each annulus to a 20'' radius based on the ratio of the effective exposures. This procedure produces maps giving estimates of the local background level at every pixel based on the observed images. The annular aperture ensures any contribution from a source at that pixel position is excluded from the background estimate. However, any large-scale contribution from a bright target source will be included. We use these background maps, along with the mosaic count images, to generate false-probability maps, and we proceed with source detection as in the dedicated survey fields, adopting a false-probability threshold of <10−6 across all bands and fields. We exclude any detections within $90^{\prime\prime} $ of the target position. We also exclude areas occupied by large, foreground galaxies (based on the optical imaging) or known sources that are associated with the target (but are not at the aimpoint). In addition, we exclude areas where the effective exposure is <30% of the maximal (on-axis) exposure in a given field, which removes unreliable detections close to the edge of the NuSTAR field of view where the background is poorly determined. We consider all fields analyzed as part of the serendipitous program up to 2015 January 1, extending the sample of Alexander et al. (2013). Full details of this extended serendipitous survey program will be given by G. B. Lansbury et al. 2016, (in preparation). For this paper, we apply a number of additional cuts.

  • 1.  
    We exclude fields at Galactic latitudes $| b| \lt 20^\circ ,$ to ensure our sample is dominated by extragalactic sources.
  • 2.  
    We exclude any fields where there are >106 counts within 120'' of the aimpoint; i.e., fields where the target is bright and will substantially contaminate the entire NuSTAR field of view.
  • 3.  
    We only consider fields at declinations > −5°; i.e. accessible from the northern hemisphere.

We do not expect any of these cuts to introduce systematic biases in the sample. The final cut ensures that we have a high spectroscopic redshift completeness,30 thanks to a substantial ongoing follow-up program with Palomar and Keck (PI: Harrison; PI: Stern), in addition to existing redshifts from the literature. Follow-up programs of southern fields are underway using Magellan (PIs: Bauer, Treister) and ESO NTT (PI: Lansbury) but have yet to achieve the high level of spectroscopic completeness required for the present study.

After applying the above cuts, our serendipitious survey spans 106 NuSTAR fields, corresponding to a total area coverage of 4.40 deg2. These NuSTAR observations span a wide range of depths (∼10–1000 ks, although predominantly ≲50 ks), resulting in a wide range of sensitivities (see Section 2.3 below and Figure 2).

Figure 2.

Figure 2. X-ray area curves (area as a function of 8–24 keV flux) for the various NuSTAR extragalactic survey components, as indicated.

Standard image High-resolution image

Our serendipitous sample contains 33 sources that are detected in the 8–24 keV band. We identify lower-energy counterparts to 30 of these sources, with the majority of counterparts (19/30) identified in the 3XMM source catalog (Watson et al. 2009; Rosen et al. 2015) and the remainder identified manually in archival XMM-Newton, Chandra, or Swift/XRT imaging data. We also identify counterparts in the WISE all-sky survey (Wright et al. 2010) for all but one of our sources (including two of the sources that lack low-energy X-ray counterparts). We identify optical counterparts, matching to the low-energy X-ray position or WISE positions (if available), using imaging from the SDSS (York et al. 2000), USNOB1 (Monet et al. 2003), or our own pre-imaging obtained during the spectroscopic follow-up program. Full details on the cross-matching process for the full serendipitous survey sample will be provided by G. B. Lansbury et al. 2016, (in preparation).

Out of our sample of 33 sources, 28 (85%) have spectroscopic redshifts (all of these sources have an extragalactic origin). We retain the remaining five sources but assume no prior knowledge of their redshift, adopting a p(z) distribution that is constant in $\mathrm{log}(1+z).$ However, after folding these broad p(z) distributions through our models of the XLF a moderate redshift (mean z ∼ 0.7) is generally preferred. We note that ∼3% of sources at $| b| \gt 20^\circ $ with spectroscopic classifications in our full serendipitous sample (including sources detected in the 3–8 keV and 3–24 keV bands) are associated with Galactic sources. Thus, one or more of our five 8–24 keV sources without spectroscopic classifications could have a Galactic origin; given our high overall redshift completeness, this potential contamination will have a negligible impact on our results. A full list of fields and the properties of sources from the serendipitous survey program, marking the fields and sources used in this work, will be given by G. B. Lansbury et al. 2016, (in preparation).

2.3. Sensitivity Analysis

Our source detection procedure (using false probability maps) is essentially identical across all of the different survey components, allowing us to determine the sensitivity in a consistent manner. We construct sensitivity maps following the procedure described by Georgakakis et al. (2008) that accounts for the Poisson nature of the detection. First, for each pixel position in our images we estimate the minimum number of counts, L, in a 20'' extraction aperture that would satisfy our false probability detection threshold given the local background estimate at that position, B. We then estimate the expected counts, s, from a source of flux, F, given the effective exposure and apply an aperture correction31 and fixed flux-to-count rate conversion factor from C15. We then calculate the probability that the combination of the expected source and background counts, s + B, produces a total number of counts that exceeds our threshold for detection, L. Each pixel contributes fractionally to the total area curve—the survey area sensitive to a given flux—in proportion to this probability. We sum over all pixels to determine our overall area curves for a given survey component. Masked areas with low exposure or corresponding to the target (for the serendipitous survey fields) are excluded in this calculation.

We note that the true flux-to-counts conversion factor depends on the spectral shape of the source. In the XLF analysis described below (Section 3), we allow for a range of spectral shapes when converting between luminosities and count rates, which are converted to an equivalent flux to determine the sensitivity from our area curves. Thus, the dependence of the sensitivity on the spectral shape is accounted for in our analysis.

Figure 2 shows the area coverage as a function of the 8–24 keV flux for the various survey components. Our COSMOS and ECDFS area curves are in good agreement with those derived from simulations by C15, verifying our analytic method. The serendipitous survey not only covers the largest overall area but also covers an area comparable to the dedicated deep surveys at fainter fluxes, making it a powerful addition to our study. We note that the dedicated surveys have very deep supporting data at lower X-ray energies and other wavelengths, which enables the high-redshift completeness and will be exploited in future studies of the X-ray and multiwavelength properties of NuSTAR sources.

3. METHODOLOGY

3.1. Luminosity Estimates

For each source in our sample, we must estimate the intrinsic luminosity, corrected for absorption and accounting for any other spectral features (e.g., reflection, scattering), allowing us to trace the accretion power of the AGNs in our sample. We estimate the luminosity in the rest-frame 10–40 keV energy range, ${L}_{10-40\;\mathrm{keV}}.$ This luminosity must be corrected for the effects of absorption along the line of sight and includes the contribution of any reflection component. The rest-frame 10–40 keV energy range roughly corresponds to the observed 8–24 keV NuSTAR band (where we perform detection) for z ≈ 0.3–1 (where the bulk of our sample lies) and is becoming the standard for extragalactic NuSTAR surveys (e.g., Alexander et al. 2013; Del Moro et al. 2014).

To convert between the observed fluxes and L10–40 keV requires knowledge of the X-ray spectrum. X-ray spectral analysis of the NuSTAR sources is underway (see Alexander et al. 2013; Del Moro et al. 2014, A. Del Moro et al. 2016, in preparation and L. Zappacosta et al. 2016, in preparation). For this study, we adopt the statistical approach used by A15. We use a particular X-ray spectral model: an absorbed power law along with a simple modeling of the Compton reflection (pexrav: Magdziarz & Zdziarski 1995) and a soft scattered component. Our model is described in detail in Section 3.1 of A15. We use priors to describe the expected distributions of the spectral parameters. For the photon index, Γ, we apply a normal prior with a mean of 1.9 and standard deviation of 0.2. For the scattered fraction, fscatt, we assume a lognormal prior with mean $\mathrm{log}{f}_{\mathrm{scatt}}=-1.73$ and a standard deviation of 0.8 dex. For the reflection strength, R, we adopt a constant prior in the range 0 < R < 2. These priors allow for our lack of knowledge of the true values for an individual NuSTAR source. This simple spectral model should provide an adequate description of the broad spectral properties of X-ray AGNs in the NuSTAR bandpass. Based on this model, we can derive the joint probability distribution function for the intrinsic luminosity (L10–40 keV), absorption column density (NH), and redshift (z), for an individual source, i, in our sample,

Equation (1)

where di represents any data on the redshift (e.g., optical spectroscopy or a photometric redshift) and Ti and bi represent the NuSTAR data: the total observed 8–24 keV counts in a $20^{\prime\prime} $ aperture, Ti, and the estimated background counts, bi. The final term, $\pi ({L}_{10-40\ \mathrm{keV}},{N}_{{\rm{H}}},z),$ represents any prior knowledge of the expected values of ${L}_{10-40\ \mathrm{keV}}$, NH and z.

In this study we assume $p(z\ | \ {d}_{i})$ is given by a δ-function at the available spectroscopic or photometric redshift value. We neglect the uncertainties in the photo-z for the seven sources where we require such estimates, given their high accuracy (σ≲ 0.04: Luo et al. 2010; Salvato et al. 2011; Hsu et al. 2014; Nandra et al. 2015). For the six sources that lack redshift information completely (one source in ECDFS that lacks a Chandra counterpart and five sources in the serendipitous survey without spectroscopic follow-up) we conservatively adopt a broad redshift distribution spanning 0 < z < 5 that is a constant function of $\mathrm{log}(1+z).$ Our high spectroscopic completeness (86%) ensures that our XLF measurements are not significantly affected by errors in the photometric redshifts or sources with indeterminate redshifts.

We calculate $p({L}_{10-40\ \mathrm{keV}},{N}_{{\rm{H}}}\ | \ z,{T}_{i},{b}_{i})$ by assuming that the observed counts are described by a Poisson process and integrating over the priors on the spectral parameters (Γ, fscatt, and R), as described in Section 3.1 of A15. Our X-ray spectral model is folded through the NuSTAR response function; thus we account for the variations in sensitivity across the 8–24 keV band. The main panel of Figure 3 shows an example of the two-dimensional constraints on ${L}_{10-40\ \mathrm{keV}}$ and NH for a typical source in our sample, assuming no a priori knowledge of these values (adopting log-constant priors for ${L}_{10-40\ \mathrm{keV}}$ and NH). The value of NH is unconstrained, but the observed counts allow us to place constraints on the value of ${L}_{10-40\ \mathrm{keV}}$, depending on the value of NH. If ${N}_{{\rm{H}}}\lesssim {10}^{23}$ cm−2 then, for this example, we estimate that $\mathrm{log}$ (${L}_{10-40\ \mathrm{keV}}$/erg s−1) ≈ 44.35 ± 0.15. The error is a combination of the Poisson uncertainties and the uncertainties on the other spectral parameters (primarily Γ and R) but is not affected by the absorption. For higher values of NH, the same observed counts must correspond to a higher intrinsic luminosity.

Figure 3.

Figure 3. Example of the two-dimensional probability distributions of ${L}_{10-40\ \mathrm{keV}}$ and NH for a single NuSTAR 8–24 keV detection (main panel), based purely on the observed 8–24 keV counts and background and assuming no a priori knowledge of ${L}_{10-40\ \mathrm{keV}}$ or NH. Contours indicate the 68.3% (solid) and 95.4% (dashed) confidence intervals (i.e., 1 and 2σ equivalent) in the two-dimensional parameter space. The sub-panels show the marginalized distributions of ${L}_{10-40\ \mathrm{keV}}$ (bottom) and NH (right) for the log-constant prior (black solid line) and applying informative priors based on previous estimates of the XLF and NH function by A15 (blue dashed line) and U14 (red dotted line), extrapolated to the NuSTAR energy band. See Section 3.1 for more details.

Standard image High-resolution image

However, we do not completely lack a priori knowledge of ${L}_{10-40\ \mathrm{keV}}$ or NH. A number of previous studies have presented estimates of the XLF and NH function of AGNs, albeit based on lower energy or lower redshift data (e.g., La Franca et al. 2005; Burlon et al. 2011; Ueda et al. 2014, A15). Such studies tell us that high luminosity sources are significantly rarer than lower luminosity sources (due to the double power-law shape of the XLF) and predict different distributions of NH. We can use these studies to apply informative priors when estimating $p({L}_{10-40\ \mathrm{keV}},{N}_{{\rm{H}}},z\ | \ {d}_{i},{T}_{i},{b}_{i})$ for an individual source.

The sub-panels of Figure 3 illustrate the marginalized distributions of ${L}_{10-40\ \mathrm{keV}}$ and NH (bottom and right panels, respectively). For our non-informative (log-constant) priors on ${L}_{10-40\ \mathrm{keV}}$ and NH (solid black lines), there is a long tail in $p({L}_{10-40\ \mathrm{keV}}),$ corresponding to Compton-thick column densities and correspondingly higher luminosities. We also apply priors based on the XLF and NH functions of A15 (blue dashed lines, see Section 5.1 and Table 9 of A15 for the model specification and parameter values, respectively) and U14 (red dotted lines, see Tables 2 and 4 of U14). These studies give the XLF in terms of the intrinsic (i.e., absorption-corrected) rest-frame 2–10 keV luminosity, which we convert to ${L}_{10-40\ \mathrm{keV}}$ using our X-ray spectral model (marginalizing over the priors on the spectral parameters). With priors based on either the U14 or A15 models, the long tail in $p({L}_{10-40\ \mathrm{keV}})$ is significantly suppressed as high luminosity sources are substantially rarer (according to the XLF). In addition, the peak of $p({L}_{10-40\ \mathrm{keV}})$ is shifted to slightly lower ${L}_{10-40\ \mathrm{keV}}$ than with a constant prior. This shift accounts for the effect of Eddington bias: a detection is more likely to be a lower luminosity (and hence more common) source where the observed flux is a positive fluctuation. Our extensive simulations have shown that this effect is significant in our NuSTAR survey data (see C15).

The distribution of $p({N}_{{\rm{H}}})$ (right panel of Figure 3) is mainly determined by the shape of the NH function at ${N}_{{\rm{H}}}\lesssim {10}^{23}$ cm−2. Thus, in this example, the probability density is slightly higher for ${N}_{{\rm{H}}}={10}^{20-21}$ cm−2 than for ${N}_{{\rm{H}}}={10}^{21-23}$ cm−2. There are only slight differences between the priors based on U14 or A15. At a fixed luminosity, the intrinsic NH function rises at ${N}_{{\rm{H}}}={10}^{23-24}$ cm−2 (for both models), hence $p({N}_{{\rm{H}}})$ increases at ${N}_{{\rm{H}}}={10}^{23}$ cm−2. However, at z = 1.0, absorption by column densities ≳1023 cm−2 starts to suppress the observed 8–24 keV flux; thus, at ${N}_{{\rm{H}}}\gt {10}^{23}$ cm−2 the same observed counts must correspond to a higher ${L}_{10-40\ \mathrm{keV}}$. As higher luminosity sources are rarer, the marginalized $p({N}_{{\rm{H}}})$ decreases as NH increases, indicating that a detected source is less likely to be associated with these higher levels of absorption. An individual detection is very unlikely to be a Compton-thick AGN and the chance of extreme column densities (${N}_{{\rm{H}}}\gtrsim {10}^{25}$ cm−2) is even more strongly suppressed. Nonetheless, the probability of an individual source having a Compton-thick NH is slightly higher for the U14 prior (which has a higher intrinsic Compton-thick fraction32 , ${f}_{\mathrm{CThick}}\approx 50$%) than for the A15 prior (${f}_{\mathrm{CThick}}\approx 25$%).

3.2. Binned Estimates of the XLF

Our NuSTAR sample is relatively small and is limited to ≳L* sources at z > 0.5 (see Figure 1). Thus, we do not attempt to fit an overall parametric model for the shape of the XLF and its evolution. Instead, we produce binned estimates of the XLF in fixed luminosity and redshift bins. We adopt the Nobs/Nmdl method (Miyaji et al. 2001; Aird et al. 2010; Miyaji et al. 2015) and use parametric models from previous, lower-energy studies of the XLF and NH function to correct for underlying selection effects. The binned estimate of the XLF is given by

Equation (2)

where ${\phi }_{\mathrm{mdl}}({L}_{b},{z}_{b})$ is a given model estimate of the XLF, evaluated at Lb and zb33 , and is rescaled by the ratio of the effective observed number of sources in a bin (Nobs) to the predicted number based on the model (Nmdl).

To calculate Nmdl we fold a given model of the XLF and NH function (A15 or U14) through the 8–24 keV sensitivity curves for our surveys (see Section 2.3, Figure 2). We convert between ${L}_{10-40\ \mathrm{keV}}$, NH, and z and the observed 8–24 keV flux based on our X-ray spectral model and the appropriate NuSTAR response function. The Nmdl term accounts for selection biases in our sample, primarily driven by the underlying NH function, and allows for the fact that heavily absorbed sources are expected to be under-represented in our sample.

To calculate Nobs we sum the individual $p({L}_{10-40\ \mathrm{keV}},{N}_{{\rm{H}}},z\ | \ {d}_{i},{T}_{i},{b}_{i})$ distributions of our sources (including priors based on either the U14 or A15 models) and integrate over a given luminosity–redshift bin. As we allow for a range of possible luminosities (and in some cases redshifts), an individual source can make a partial contribution to Nobs for multiple bins.

We estimate errors on our binned XLF based on the approximate Poisson error in the effective observed number of sources in a bin, Nobs. We obtain Poisson errors based on Gehrels (1986), although in many cases Nobs is non-integer and thus these errors are approximations (see also Aird et al. 2010; Miyaji et al. 2015). We only plot bins where ${N}_{\mathrm{obs}}\geqslant 1.$

We note that all three terms in Equation (2) will depend, to varying extents, on the underlying parametric model of the XLF and NH function that is assumed. In particular, changing the assumed ${f}_{\mathrm{CThick}}$ will alter our binned estimates; we expect Compton-thick sources are severely under-represented in our observed sample, which is accounted for in the Nmdl term, and thus our XLF estimates are increased to allow for this missed population. We explore these effects further in Section 4 below.

4. MEASUREMENTS OF THE XLF WITH NUSTAR

In this section we present measurements of the XLF of AGNs based on the NuSTAR 8–24 keV sample. In Figure 4 (top panels) we present binned estimates of the 10–40 keV XLF in three redshift ranges based on the ${N}_{\mathrm{obs}}/{N}_{\mathrm{mdl}}$ method. The green triangles are estimates where we neglect the effects of absorption, providing the most direct estimate of the observed 10–40 keV XLF. We assume all sources have ${N}_{{\rm{H}}}={10}^{20}$ cm−2 and calculate the binned XLF utilising the best-fitting model for the Compton-thin XLF from A15, although these binned estimates are not strongly affected by the assumed XLF model. The red crosses and blue circles are binned estimates where we account for absorption effects using the best-fitting model of the NH function from U14 and A15, respectively, and adopt the total XLF of AGNs (including Compton-thick sources), extraoplated from 2–10 keV energies. These binned estimates are higher than when absorption is neglected (green triangles) due to the corrections for heavily absorbed (${N}_{{\rm{H}}}\gtrsim {10}^{23}$ cm−2) and Compton-thick (${N}_{{\rm{H}}}\gtrsim {10}^{24}$ cm−2) sources that will be under-represented in our observed samples due to selection biases, even at the harder energies probed by NuSTAR. Our binned estimates are listed in Table 1.

Figure 4.

Figure 4. Top panels: binned estimates of the 10–40 keV XLF of AGNs based on the NuSTAR 8–24 keV selected sample in three redshift ranges. Green triangles are estimates when the effects of absorption are neglected (error bars are omitted for clarity). The red crosses and blue circles are binned estimates of the XLF, where the luminosities and space densities are corrected for absorption effects. We account for absorption-dependent selection effects using the Nobs/Nmdl method, assuming different models of the NH function (including Compton-thick AGNs): U14 and A15 for the red crosses and blue circles, respectively. Errors are based on the Poisson error in the observed number of sources in a bin. There are small discrepancies between the U14 and A15 binned estimates due to differences in the models of the NH function (see Figure 5 and text in Section 4 for details). The black dashed line shows a model of the total XLF of AGNs (from U14, although the A15 XLF model is virtually identical over the range probed by our data), extrapolated from 2–10 keV energies to the 10–40 keV band assuming our baseline, unabsorbed X-ray spectral model (with ${\rm{\Gamma }}=1.9$ and R = 1). Assuming stronger reflection (e.g., R = 2) shifts the model to higher luminosities (gray dotted line) and can thus bring the model into better agreement with the binned estimates based on the U14 NH function. The purple long-dashed lines in all panels show the XLF from Swift/BAT (Ajello et al. 2012) at z ≈ 0, converted to the 10–40 keV band. Bottom panels: ratio between the XLF model and the binned estimates in terms of Nobs/Nmdl.

Standard image High-resolution image

Table 1.  Binned Estimates of the Rest-frame 10–40 keV X-Ray Luminosity Function of AGNs Based on the NuSTAR 8–24 keV Sample

z-range $\langle z\rangle $ a $\mathrm{log}{L}_{10-40\;\mathrm{keV}}$ ϕnoabsb ϕA15c ϕU14d
    (erg s−1) (Mpc−3 dex−1) (Mpc−3 dex−1) (Mpc−3 dex−1)
0.1–0.5 0.28 42.75–43.00 ${1.81}_{-1.08}^{+2.14}\times {10}^{-4}$ ${2.64}_{-1.76}^{+3.86}\times {10}^{-4}$ ${3.42}_{-2.25}^{+4.87}\times {10}^{-4}$
    43.00–43.25 ${9.39}_{-4.67}^{+8.13}\times {10}^{-5}$ ${1.52}_{-0.81}^{+1.46}\times {10}^{-4}$ ${1.74}_{-0.93}^{+1.71}\times {10}^{-4}$
    43.25–43.50 ${5.69}_{-2.27}^{+3.50}\times {10}^{-5}$ ${9.51}_{-3.82}^{+5.91}\times {10}^{-5}$ ${1.02}_{-0.42}^{+0.66}\times {10}^{-4}$
    43.50–43.75 ${2.63}_{-0.98}^{+1.46}\times {10}^{-5}$ ${4.08}_{-1.52}^{+2.28}\times {10}^{-5}$ ${4.60}_{-1.74}^{+2.62}\times {10}^{-5}$
    43.75–44.00 ${8.48}_{-3.85}^{+6.33}\times {10}^{-6}$ ${1.27}_{-0.54}^{+0.86}\times {10}^{-5}$ ${1.59}_{-0.67}^{+1.07}\times {10}^{-5}$
    44.00–44.25 ${3.20}_{-1.66}^{+2.95}\times {10}^{-6}$ ${4.28}_{-2.18}^{+3.84}\times {10}^{-6}$ ${5.73}_{-2.87}^{+5.03}\times {10}^{-6}$
    44.25–44.50 ${1.81}_{-0.95}^{+1.71}\times {10}^{-6}$ ${3.00}_{-1.50}^{+2.63}\times {10}^{-6}$ ${3.83}_{-1.90}^{+3.28}\times {10}^{-6}$
0.5–1.0 0.75 43.50–43.75 ${1.74}_{-1.27}^{+3.09}\times {10}^{-4}$ ${2.12}_{-1.73}^{+4.89}\times {10}^{-4}$ ${2.95}_{-2.44}^{+7.05}\times {10}^{-4}$
    43.75–44.00 ${5.92}_{-3.05}^{+5.42}\times {10}^{-5}$ ${9.49}_{-5.06}^{+9.21}\times {10}^{-5}$ ${1.25}_{-0.66}^{+1.20}\times {10}^{-4}$
    44.00–44.25 ${2.32}_{-0.96}^{+1.50}\times {10}^{-5}$ ${3.21}_{-1.33}^{+2.08}\times {10}^{-5}$ ${3.89}_{-1.59}^{+2.48}\times {10}^{-5}$
    44.25–44.50 ${1.53}_{-0.42}^{+0.57}\times {10}^{-5}$ ${1.99}_{-0.56}^{+0.76}\times {10}^{-5}$ ${2.26}_{-0.64}^{+0.86}\times {10}^{-5}$
    44.50–44.75 ${4.17}_{-1.39}^{+1.98}\times {10}^{-6}$ ${5.89}_{-1.89}^{+2.65}\times {10}^{-6}$ ${6.42}_{-2.08}^{+2.94}\times {10}^{-6}$
    44.75–45.00 ${3.74}_{-2.68}^{+6.38}\times {10}^{-7}$ ${5.50}_{-3.57}^{+7.62}\times {10}^{-7}$ ${5.84}_{-3.88}^{+8.49}\times {10}^{-7}$
1.0–3.0 1.51 44.25–44.50 ${4.21}_{-2.44}^{+4.73}\times {10}^{-5}$ ${5.75}_{-3.58}^{+7.36}\times {10}^{-5}$ ${1.03}_{-0.63}^{+1.28}\times {10}^{-4}$
    44.50–44.75 ${1.61}_{-0.71}^{+1.16}\times {10}^{-5}$ ${2.08}_{-0.95}^{+1.57}\times {10}^{-5}$ ${3.34}_{-1.54}^{+2.56}\times {10}^{-5}$
    44.75–45.00 ${5.69}_{-2.09}^{+3.09}\times {10}^{-6}$ ${6.75}_{-2.49}^{+3.70}\times {10}^{-6}$ ${8.92}_{-3.32}^{+4.95}\times {10}^{-6}$
    45.00–45.25 ${2.00}_{-0.69}^{+1.00}\times {10}^{-6}$ ${2.36}_{-0.81}^{+1.17}\times {10}^{-6}$ ${2.68}_{-0.92}^{+1.33}\times {10}^{-6}$
    45.25–45.50 ${3.76}_{-1.85}^{+3.19}\times {10}^{-7}$ ${4.86}_{-2.27}^{+3.81}\times {10}^{-7}$ ${4.89}_{-2.29}^{+3.84}\times {10}^{-7}$
    45.50–45.75 ${1.56}_{-0.91}^{+1.76}\times {10}^{-7}$ ${1.54}_{-0.87}^{+1.66}\times {10}^{-7}$ ${1.40}_{-0.79}^{+1.51}\times {10}^{-7}$

Notes.

aMean redshift of sources in the redshift bin. bBinned estimates of the XLF, neglecting absorption (i.e., assuming all sources have ${N}_{{\rm{H}}}={10}^{20}$ cm−2). The A15 model for the XLF of Compton-thin AGNs is adopted for the ${N}_{\mathrm{obs}}/{N}_{\mathrm{mdl}}$ method. Errors are based on the Poisson error in the observed number of sources (Nobs). cBinned estimates of the XLF, adopting the A15 model for the NH function and the total XLF (including Compton-thick sources) for the ${N}_{\mathrm{obs}}/{N}_{\mathrm{mdl}}$ method. dBinned estimates of the XLF, adopting the U14 model for the NH function and the total XLF (including Compton-thick sources) for the ${N}_{\mathrm{obs}}/{N}_{\mathrm{mdl}}$ method.

Download table as:  ASCIITypeset image

In Figure 4, we also show the U14 model for the total XLF (including Compton-thick sources), extrapolated from 2–10 keV energies (black dashed line). To extrapolate, we assume our (unabsorbed) X-ray spectral model evaluated at the mean of the priors on the spectral parameters (i.e., Γ = 1.9, fscatt = 1.9%, and R = 1.0), which gives

Equation (3)

where L2–10 keV is the rest-frame 2–10 keV luminosity. Both the U14- and A15-based binned estimates, to first order, are in agreement with this extrapolation of the XLF model, within their errors. We note that the extrapolation of the A15 model for the total XLF (including Compton-thick AGNs) is virtually identical to the U14 XLF model over the range of luminosities and redshifts probed by our data and thus is also in agreement with the binned estimates.

On closer examination, however, there are differences between the binned estimates when the U14 model of the NH function is used to correct for selection biases, rather than A15. The U14 binned estimates are slightly higher than the A15 binned estimates, most noticeably at higher redshifts, and are systematically higher than the extrapolated XLF model. The Nobs/Nmdl ratios, indicating the differences between our binned estimates and the model, are shown in the bottom panels of Figure 4 and also illustrate this pattern: Nobs/Nmdl is systematically higher for the U14-based binned estimates, indicating that the U14 model of the NH function under-predicts the number of sources in our NuSTAR sample (i.e., under-predicts Nmdl). These discrepancies must be due to differences in the underlying model of the NH function, which introduces different corrections for the fraction of heavily absorbed and Compton-thick sources.

To explore this further, in the top panels of Figure 5 we show the intrinsic NH functions from the U14 and A15 models, evaluated at the mean redshift and geometric mean luminosity for NuSTAR sources in each bin. There are clear differences whereby the U14 model has a higher fraction of sources in the Compton-thick regime (${f}_{\mathrm{CThick}}=50$% versus 25% for A15) and a lower relative fraction of sources with ${N}_{{\rm{H}}}={10}^{23-24}$ cm−2 (i.e., heavily obscured, but Compton-thin sources). These differences result in higher binned estimates of the XLF in Figure 4 when the U14 model is used to correct for selection biases—the binned estimates are increased to allow for a larger population of Compton-thick sources that are not represented in our observed samples.

Figure 5.

Figure 5. Model predictions of the intrinsic and observable NH functions. Top panels: intrinsic NH functions evaluated at the mean redshift and geometric mean luminosity of sources in the bin, based on the A15 (blue dashed line) and U14 (red dotted line) models. Bottom panels: predicted numbers of sources in our NuSTAR 8–24 selected sample in each redshift range as a function of NH, based on the A15 and U14 models (folded through our sensitivity and assuming our baseline R = 1 spectral model in both cases). We also give the total predicted numbers (Nmdl) in each redshift range and the observed number of sources (Nobs). The total observed number in a bin (Nobs) is non-integer due to the small number of sources with undetermined redshifts that can thus make a partial contribution to multiple redshift bins.

Standard image High-resolution image

The bottom panels of Figure 5 show the predicted numbers of sources with different NH values in our 8–24 keV selected sample. These are calculated by folding the U14 or A15 models of the NH function (and XLF) through our sensitivity functions (assuming our X-ray spectral model) over the entire redshift bin. The predicted distributions differ substantially from the intrinsic NH functions shown in the top panels; most noticeably the predicted numbers of Compton-thick AGNs in the samples are extremely small. Such small numbers are due to the flux from a Compton-thick AGN being suppressed by a factor of ≳2–10, even in the relatively hard 8–24 keV band that is used for selection, combined with the steep slope of the XLF at the bright luminosities that we probe, which results in strong selection biases against such sources (see Figure 3 and Ballantyne et al. 2011). The predicted numbers of sources based on the U14 model are generally smaller than with the A15 model due to the higher intrinsic ${f}_{\mathrm{CThick}}$ in the U14 model. The total predicted numbers, Nmdl combined over all luminosities, and the effective observed numbers in a bin, Nobs, are given in each panel. Generally, Nmdl based on the A15 model is closer to Nobs in all of the redshift panels, although statistically the Nmdl estimates from both A15 and U14 are consistent with Nobs at the <2.5σ level. Combining the three redshift panels, we find that our total observed number of sources (94) is consistent with the prediction based on A15 (84.5) at the <1σ level but is significantly higher (>3σ) than the U14 prediction (58.9).

Another possible explanation for the U14 model under-predicting the observed number of sources could be our choice of X-ray spectral model. Assuming different prior distributions for Γ, fscatt, and R could change the extrapolation of our models of the XLF to the 10–40 keV band, potentially bringing the model XLF into better agreement with our binned estimates and altering our Nmdl predictions (without requiring any changes to the NH function). Varying fscatt has a minimal impact as the observed 8–24 keV band probes the direct emission and reflection component, rather than the scattered emission that is important at lower energies. Varying Γ has a larger impact, although a very flat, intrinsic photon index (Γ ≈ 1.4 which is then subjected to absorption) is required to bring the U14 model into agreement with our binned estimates. It is now well established that the intrinsic photon index of the power-law emission from AGNs has a mean Γ ≈ 1.9 with scatter ≲0.2 (e.g., Nandra & Pounds 1994; Tozzi et al. 2006; Scott et al. 2011), as assumed in our original spectral model.

Changing our assumed reflection strength, R, has a more significant impact. Our baseline spectral model assumes a flat distribution in the range R = 0–2 and thus an average reflection of $\langle R\rangle =1.$ Assuming stronger reflection increases the luminosity in the 10–40 keV band for the same 2–10 keV luminosity.34 For R = 2, we find

Equation (4)

with our spectral model, which shifts the model XLF to higher ${L}_{10-40\ \mathrm{keV}}$ (gray dotted line in Figure 4) and brings it into good agreement with the U14-based binned estimates. For a fixed R = 2 (for all sources), we predict Nmdl = 99.9 (across all redshift bins) with the U14 model of the NH function, in good agreement with our observed number of sources (94). However, with the A15 model of the NH function we predict Nmdl = 146.2 if we adopt a spectral model with stronger reflection (R = 2), thus over-predicting the observed number of sources.

In conclusion, the U14 model of the NH function predicts fewer sources than observed in our 8–24 keV sample, possibly due to a higher relative fraction of Compton-thick AGNs (compared to A15). Thus, the binned estimates of the XLF are slightly higher. The A15 model, which predicts more sources with ${N}_{{\rm{H}}}\approx {10}^{23-24}$ cm−2, is in better agreement with the observed samples. Alternatively, a stronger reflection component (average R ≈ 2 rather than R ≈ 1) across all luminosities results in better agreement between our measurements of the XLF and the extrapolation of the U14 model. Further study (e.g., X-ray spectral analysis) is required to improve constraints on the NH function, determine the intrinsic distribution of R, and refine our estimates of the XLF.

5. DISCUSSION

5.1. The Evolution of the XLF

We have presented the first measurements of the rest-frame 10–40 keV XLF of AGNs at high redshifts (0.1 < z < 3) based on a sample of sources selected at comparable observed-frame energies (8–24 keV). Selecting at hard X-ray energies—in contrast to prior studies based on lower-energy data—provides estimates of the intrinsic X-ray luminosity that are relatively unaffected by absorption due to Compton-thin column densities (${N}_{{\rm{H}}}\lesssim {10}^{24}$ cm−2). Thus, our work provides the most direct measurements of the XLF currently possible at these redshifts, although our selection remains biased against Compton-thick sources.

Our results are consistent with the strong evolution of the AGN population seen in lower-energy studies of the XLF (e.g., Ueda et al. 2003; Barger et al. 2005) and at longer wavelengths (e.g., Ross et al. 2013; Lacy et al. 2015), characterized by a shift in the luminosity function toward higher luminosities at higher redshifts. This evolution leads to the substantial increase in the space density of luminous (${L}_{10-40\ \mathrm{keV}}\gtrsim {10}^{44}$ erg s−1) AGNs out to z ∼ 2 that is seen in our measurements of the XLF (see Figure 4). No substantial modifications to this overall picture appear to be required based on our NuSTAR data, although the limited range of luminosities probed by NuSTAR means we are unable to measure the faint-end slope of the XLF at z ≳ 0.5 or compare different parameterizations. Furthermore, we find that strong (R ∼ 2) reflection is needed to reconcile the U14 model of the XLF with our NuSTAR sample (see further discussion in Section 5.2 below).

Our NuSTAR measurements provide a high-redshift comparison to previous measurements of the XLF at hard (≳10 keV) X-ray energies, which have been restricted to the local (z ≲ 0.1) universe. In Figure 6 we compare the best fit to the XLF of local AGNs in the Swift/BAT 60-month sample (thick purple dashed line, Ajello et al. 2012) with extrapolations of the U14 and A15 models to z = 0.03 (the median redshift of sources in the Swift/BAT sample). Both the A15 model (solid blue line) and the U14 model with a strong (R = 2) reflection component (black dotted line)—either of which is consistent with our higher redshift NuSTAR data—over-predict the local Swift/BAT XLF at all luminosities. The U14 model with R = 1 is in better agreement with the Swift/BAT XLF close to the flux limit (vertical gray dashed line); however, this model significantly under-predicts the number of sources in our high-redshift NuSTAR sample (see Figure 5, discussed in Section 4 above). A similar result is seen in our analysis of the NuSTAR number counts, where the predictions of population synthesis models (based on either the U14 XLF with strong reflection or the A15 model) provide a good agreement with the observed NuSTAR number counts but over-predict the counts at the brighter fluxes probed by Swift/BAT (see Figure 5 of Harrison et al. 2015). Ballantyne (2014) also highlighted the discrepancies between measurements of the local XLF and the extrapolations of high-redshift evolutionary models. These findings indicate that there must be some evolution in the XLF, absorption distribution, or spectral properties of AGNs between z ∼ 0 and the higher redshifts probed by NuSTAR that is not fully accounted for by the current models.

Figure 6.

Figure 6. Comparison of the local XLF of Compton-thin AGNs, based on the best fit to the Swift/BAT 60-month 15–55 keV sample by Ajello et al. (2012; purple long-dashed line, hatched region indicates conversion to 10–40 keV with a photon index between Γ = 1.4 and Γ = 2.1) to extrapolations of the U14 and A15 models (evaluated at z = 0.03, the median redshift of the Swift/BAT sample). The extrapolations of the A15 model with moderate reflection (R = 1) and a relatively high fraction of sources with ${N}_{{\rm{H}}}={10}^{23-24}$ cm−2 (blue solid line) and the U14 model with strong reflection (R = 2, dotted black line), both of which can describe the high-redshift population seen by NuSTAR, over-predict the observed z ∼ 0 XLF at all luminosities. The U14 model with R = 1 is in better agreement with the Swift/BAT measurement close to the flux limit (vertical gray dashed line); however, this model significantly under-predicts the number of sources in our high-redshift NuSTAR sample. These findings indicate that there is evolution in the XLF, absorption distribution, or spectral properties of AGNs between z ∼ 0 and the higher redshifts (z ∼ 0.1–3) probed by NuSTAR that is not accounted for by the current models.

Standard image High-resolution image

5.2. The Absorption Distribution, the Fraction of Compton-thick Sources and the Reflection Strength

Our final binned estimates of the XLF differ slightly depending on the underlying model of the NH function that is used to correct for selection biases, whereby the binned estimates using the U14 model are systematically higher than the binned estimates using the A15 model. Our observed sample appears to favor a higher relative number of heavily absorbed (${N}_{{\rm{H}}}\approx {10}^{23-24}$ cm−2) sources and a lower Compton-thick fraction, as given by the A15 model. We note that the NH function of A15 is based on an indirect method, attempting to reconcile samples of AGNs selected at 0.5–2 keV and 2–7 keV energies from Chandra surveys, spanning out to z ∼ 5. In contrast, U14 determine the NH function in the local universe based on spectral analysis of sources in the Swift/BAT AGN sample (selected at 15–195 keV), which is then extrapolated to higher redshifts. Both of these methods have their limitations. Accurate constraints on the intrinsic NH function in both the local and high-redshift universe are vital to determine the extent of obscured black hole growth and shed light on the physical origin of the obscuring material (e.g., Hopkins et al. 2006; Buchner et al. 2015). While our results provide indirect constraints on the NH function (we do not measure NH for individual sources), our measurements of the 10–40 keV XLF place important constraints on the number densities of absorbed, moderately luminous sources at z ≳ 0.1.

The incidence and strength of Compton reflection, which can substantially boost the observed fluxes at ≳10 keV energies, also has an important impact on our study. Moderate reflection (R ≈ 1) provides a good agreement between our observed NuSTAR sample and the extrapolation of the A15 XLF from lower energies. A stronger (R ≈ 2) reflection component is needed to reconcile our NuSTAR sample with the extrapolated model without requiring changes to the U14 NH function and Compton-thick fraction. Most prior studies have found that R (usually probed indirectly via the strength of the iron K-line) is inversely correlated with luminosity and generally weak (i.e., R ≲ 1) for moderately luminous AGNs (e.g., Iwasawa & Taniguchi 1993; Nandra et al. 1997; Page et al. 2005; Ricci et al. 2011). Spectral analysis of a single source in our NuSTAR sample (NuSTAR J033202-2746.8 in the ECDFS: Del Moro et al. 2014) found a moderate reflection component (R ≈ 0.55) for this high luminosity (${L}_{10-40\ \mathrm{keV}}\approx 6.4\times {10}^{44}$ erg s−1) source. However, Ballantyne (2014) found that strong reflection (R ≈ 1.7) at all luminosities is needed to reconcile different measurements of the local XLF across a wide range of X-ray energies (∼0.5–200 keV). Our measurements of the 10–40 keV XLF indicate that moderate-to-strong reflection (R ∼ 1−2) is required to describe the average spectral characteristics of ${L}_{10-40\ \mathrm{keV}}\sim {10}^{43-46}$ erg s−1 AGNs at z ∼ 0.1–3. The extent, strength, and spectral characteristics of reflection provide insights into the physical nature of the obscuring material and the accretion disk (e.g., García et al. 2013; Falocco et al. 2014; Brightman et al. 2015). Strong reflection could also indicate a substantial population of rapidly spinning black holes in the detected sample; however, a relatively small intrinsic fraction of high-spin sources (∼7%) can potentially dominate the observed number counts at a given flux limit (Brenneman et al. 2011; Vasudevan et al. 2015). Accurately measuring the distribution of reflection is thus an important challenge for future statistical studies of AGN populations.

The Compton-thick fraction and the strength of reflection are also vital parameters for understanding the origin of the CXB, in particular, the peak at ∼20–30 keV (e.g., Gilli et al. 2007; Draper & Ballantyne 2009; Akylas et al. 2012). The degeneracy between these parameters and the consequences for AGN population synthesis models are discussed in detail by Treister et al. (2009). Our measurements of the XLF appear to suffer from a similar degeneracy. Nevertheless, our results constrain the distribution of luminosity and redshift of sources responsible for a large fraction of the CXB emission (∼35% at 8–24 keV, see Harrison et al. 2015), placing an additional constraint on the possible contributions of Compton-thick AGNs or reflection to the CXB.

Recent NuSTAR studies have provided constraints on the distribution of NH of optically selected Type-2 QSOs (Gandhi et al. 2014; Lansbury et al. 2014, 2015) and find that NH is often underestimated for these sources based purely on low-energy (<10 keV) data. Given the small sample sizes and large remaining uncertainties, their recovered NH function is consistent with both the U14 and A15 models adopted in this paper. C15 use a band ratio analysis to identify candidate Compton-thick AGNs among NuSTAR detected sources in the survey of the COSMOS field. Their number of candidates (∼13%–20% of the NuSTAR sources) is significantly higher than our predictions (we expect only ∼0.5 Compton-thick AGNs out of 32 sources detected at 8–24 keV energies in COSMOS), indicating that our models may need updating to a much higher intrinsic ${f}_{\mathrm{CThick}}$. However, X-ray spectral analysis is required to measure NH and confirm the Compton-thick nature of these candidates. Alexander et al. (2013) presented X-ray spectral analysis of the first ten sources from the NuSTAR serendipitous survey, finding a high fraction (≳50%) have ${N}_{{\rm{H}}}\gtrsim {10}^{22}$ cm−2 but none were Compton-thick, consistent with the predictions of our models, albeit for a limited sample size. Del Moro et al. (2014) presented the first spectral analysis from the deep survey of the ECDFS, focusing on a single source.

5.3. Future Prospects

To fully test the conclusions of this paper, compare between the U14 and A15 models, and refine our estimates of the XLF requires accurate measurements of the distribution of NH and R for AGNs across a wide range of luminosities and redshifts. Spectral analysis of sources from across the NuSTAR extragalactic survey program will be the focus of forthcoming papers (A. Del Moro et al. 2016, in preparation, L. Zappacosta et al. 2016, in preparation) and will place crucial constraints on the distribution of NH and R for luminous AGNs out to z ∼ 3. Ultimately though, the ability of the NuSTAR survey program to constrain the XLF, NH function, Compton-thick fraction, or reflection properties of AGNs is limited by both depth and sample size. The ongoing survey program will improve this situation. Our serendipitous sample will roughly double with the inclusion of Southern fields (see G. B. Lansbury et al. 2016, in preparation) and continues to grow as NuSTAR observes targets. The dedicated survey program is also continuing and will initially focus on increasing the area coverage of the deep (∼400 ks) layer via observations of the GOODS-N (Alexander et al. 2003) and UDS (Lawrence et al. 2007; Ueda et al. 2008) fields. These observations will increase the area coverage at the faintest fluxes by $\gtrsim 50$%, improving the number statistics at low luminosities.

Pushing to greater depths, however, is vital to constrain the low-to-moderate luminosity AGN population (≲L*) that corresponds to the bulk of the accretion density. Our current NuSTAR sample is limited to luminous X-ray sources. We do not probe below the break in the XLF (L*) at z ≳ 0.5 and do not place strong constraints on the faint-end slope in any redshift bin. Thus, we are unable to address issues regarding the best parametric description of the evolution of the XLF (i.e., pure luminosity evolution, independent luminosity and density evolution, luminosity-dependent density evolution) or the extent of any evolution in the shape (see Aird et al. 2010; Miyaji et al. 2015, U14, A15). At fainter luminosities the intrinsic Compton-thick fraction is expected to be somewhat higher and the flatter slope of the XLF should reduce the biases against the detection of such sources in the current samples (see Figure 5). Indeed, the observed fraction of Compton-thick AGNs in deep Chandra surveys (∼3%, e.g., Brightman et al. 2014) is similar to the expected fraction in our NuSTAR survey and the absolute numbers (∼100 Compton-thick candidates were identified by Brightman et al. 2014) are much larger, mainly due to the much fainter luminosities that are accessed by Chandra.35 Substantially increasing the nominal depths of the NuSTAR survey is challenging as the observations are background limited for exposures ≳150 ks. An alternative strategy is to consider sources that are detected in the broader 3–24 keV band, providing a sample of sources a factor ≳2 larger and reaching luminosities a factor ∼2 fainter than the 8–24 keV sample considered here. A preliminary analysis of the 3–24 keV sample (using the same, indirect approach to absorption corrections described in this paper) indicates that the XLF is consistent with the U14 and A15 models at ≲L*, although this band is dominated by soft (≲8 keV) photons and larger, uncertain corrections are required for the fraction of absorbed and Compton-thick sources. However, many of the sources detected at 3–8 keV or 3–24 keV in the NuSTAR surveys have statistically significant counts (and thus useful information) at harder energies. Thus, it may be possible to improve constraints on the XLF, NH function, and distribution of reflection via sophisticated analysis of the band ratios or X-ray spectra, or via stacking analyses. Careful consideration of the survey sensitivity and selection biases will be vital in such studies.

6. CONCLUSIONS

  • 1.  
    We have presented the first measurements of the rest-frame 10–40 keV XLF of AGNs at 0.1 < z < 3 based on a sample of 94 sources selected at comparable observed-frame energies (8–24 keV). Our study takes advantage of the unprecedented sensitivity at these energies that is achieved by the NuSTAR survey program.
  • 2.  
    We find that different models of the NH function, used to account for selection biases in our measurements, make significantly different predictions for the total number of sources in our sample, leading to slight differences in our binned estimates of the XLF. The Ueda et al. (2014) model predicts fewer AGNs than observed in our sample, possibly due to a higher ${f}_{\mathrm{CThick}}$, whereas the Aird et al. (2015) model, which predicts more sources with ${N}_{{\rm{H}}}\approx {10}^{23-24}$ cm−2, is in better agreement with the observed samples.
  • 3.  
    Our results are also sensitive to our assumed X-ray spectral model. Stronger reflection (R ≈ 2, compared to our baseline assumption of R ≈ 1) at all luminosities can bring the Ueda et al. (2014) model predictions into better agreement with our NuSTAR sample. However, with R ≈ 2 the Aird et al. (2015) model over-predicts the number of sources in our sample by ≳50%.
  • 4.  
    Our results are consistent with the strong evolution of the AGN population seen in lower-energy studies of the XLF (e.g., Ueda et al. 2003; Barger et al. 2005; Aird et al. 2010), characterized by a shift in the luminosity function toward higher luminosities at higher redshifts. However, the models that successfully describe the high-redshift population detected by NuSTAR tend to over-predict the local (z ≈ 0) XLF measured by Swift/BAT, indicating some evolution of the AGN population that is not fully captured by the current models. Nonetheless, as our sample is limited to luminous (≳L∗) X-ray sources at z ≳ 0.5, we defer an investigation of different parametric descriptions of the evolution of the XLF to future studies.
  • 5.  
    Forthcoming X-ray spectral analysis of the NuSTAR survey should enable us to measure NH and R for the brightest sources, break the degeneracy between the NH function and the average reflection strength, and refine our estimates of the XLF. Including lower-energy NuSTAR detections may enable us to probe a factor ∼2 deeper. The ongoing NuSTAR survey program will also increase our sample size and improve our estimates of the XLF.

This work made use of data from the NuSTAR mission, a project led by the California Institute of Technology, managed by the Jet Propulsion Laboratory, and funded by the National Aeronautics and Space Administration. We thank the NuSTAR Operations, Software, and Calibration teams for support with the execution and analysis of these observations. This research has made use of the NuSTAR Data Analysis Software (NUSTARDAS) jointly developed by the ASI Science Data Center (ASDC, Italy) and the California Institute of Technology (USA). We acknowledge financial support from: ERC Advanced Grant FEEDBACK at the University of Cambridge (J.A., A.C.F.); a COFUND Junior Research Fellowship from the Institute of Advanced Study, Durham University (J.A.); the Science and Technology Facilities Council (STFC) grants ST/I001573/1 (D.M.A. and A.D.M.), ST/K501979/1 (G.B.L.), and ST/J003697/1 (P.G.); the Leverhulme Trust (D.M.A.); the Caltech Kingsley Visitor Program (D.M.A., A.C.); NSF award AST 1008067 (D.R.B.); NASA grants 11-ADAP11-0218 and GO3-14150C (F.C.); an Alfred P. Sloan Research Fellowship and a Dartmouth Class of 1962 Faculty Fellowship (R.C.H.); CONICYT-Chile grants Basal-CATA PFB-06/2007 (F.E.B., E.T.); FONDECYT 1141218 (F.E.B.) and 1120061 (E.T.); "EMBIGGEN" Anillo ACT1101 (F.E.B., E.T.); the Ministry of Economy, Development, and Tourisms Millennium Science Initiative through grant IC120009, awarded to The Millennium Institute of Astrophysics, MAS (F.E.B.); NASA NuSTAR subcontract 44A-1092750 and NASA ADP grant NNX10AC99G (W.N.B., B.L.); ASI/INAF grant I/037/12/0011/13 (A.C., L.Z.); and NASA Earth and Space Science Fellowship Program grant NNX14AQ07H (M.B.).

Footnotes

  • 30 

    The cut on declinations $\gt -5^\circ $ is not applied for the number counts analysis of Harrison et al. (2015), where redshift information is not required, resulting in a larger areal coverage and sample size.

  • 31 

    As discussed by C15, the core of the NuSTAR point-spread function varies by less than a few percent over the field of view. Thus, we can neglect any spatial dependence of the aperture correction.

  • 32 

    We define the intrinsic fraction of Compton-thick AGNs, fCThick, as the ratio of the number of sources with ${N}_{{\rm{H}}}={10}^{24-26}$ cm−2 to all absorbed (${N}_{{\rm{H}}}\gt {10}^{22}$ cm−2) AGNs.

  • 33 

    We fix Lb at the center of the luminosity bin, whereas zb is fixed at the mean redshift of all sources in a given redshift bin.

  • 34 

    The contribution from reflection is negligible at the ∼2–10 keV energies probed by A15 and U14, thus adopting a spectral model with stronger reflection when extrapolating to the NuSTAR band does not contradict the results of these lower-energy studies.

  • 35 

    We note that the higher energies probed by NuSTAR are vital to accurately characterize the spectra of both Compton-thin and Compton-thick AGNs, measure NH, and determine the strength of the reflection component (see Lansbury et al. 2015).

Please wait… references are loading.
10.1088/0004-637X/815/1/66