THE SPITZER MID-INFRARED AGN SURVEY. II. THE DEMOGRAPHICS AND COSMIC EVOLUTION OF THE AGN POPULATION

We used the optical through near-infrared photometric data tabulated in Table 2 of Paper I combined with 12 μm flux densities from the Wide field Infrared Survey Explorer (WISE; Wright et al. 2010). The 12 μm flux densities were particularly useful for constraining the mid-infrared spectral energy distributions (SEDs) of the AGNs between the Spitzer 8 μm and 24 μm bands. We matched the objects in Paper I to the WISE all-sky survey using a $3^{\prime\prime}$ match radius, recovering matches to 896 out of 963 (93%) objects. (Some of those not matched had faint detections in WISE 12 μm from which we could estimate flux densities.) These data were used as the basis of the SED fits described in Appendix A.

Figure 1 shows the redshift–luminosity plane for our survey. The rest-frame 5 μm luminosity, ${{L}_{5\,\mu {\rm m}}}$, was calculated by fitting the SED of each object with either a pure type-1 AGN model or a combination of an AGN and a host stellar population, as detailed in Appendix A. The presence of infrared emission from the host galaxy may lead to some subtle selection effects, which are also discussed further in  Appendix A, though we do not believe they have a significant effect on our derived luminosity function, at least at $z\gt 0.4$.

Zoom In Zoom Out Reset image size

The rest-frame 5 μm luminosity of the AGNs (without a stellar host contribution) was then used to calculate ${{L}_{5\,\mu {\rm m}}}$. The 14 samples from which the survey was drawn, with their differing 24 μm flux density limits and areas, account for the relatively wide luminosity range at a given redshift (a factor of $\;\approx 30$) in all but the highest redshift bins. For the analysis presented in this paper, we use the "statistical sample" of Paper I, which comprises the 662 objects in a 90% spectroscopically complete subsample of the entire survey (479 of which are confirmed AGNs). The "statistical sample" was defined by exploiting the close correlation of S₂₄ and emission line flux density. Objects in each of the 14 samples were sorted in order of decreasing S₂₄. The value of S₂₄ at which the survey completeness fell below 90% was then noted, and objects with S₂₄ brighter than this included in the "statistical sample". This procedure prevents the samples from becoming too biassed toward objects which are easy to obtain spectra for (e.g., normal type-1 quasars). As described in Paper I, we classify our AGNs as normal, blue type-1 AGNs, red type-1 AGNs (dust reddened, but still showing broad lines in the rest-frame optical, i.e., with ${{A}_{V}}\sim 1$ toward the AGN), type-2 AGNs (showing narrow lines only, even in the rest-frame optical, and thus with ${{A}_{V}}\;\mathop{_{\sim }}\limits^{\gt }\;5$ toward the AGN), non-AGNs (selected as AGN candidates but showing no signs of AGN activity in their optical spectra), stars, and objects with featureless spectra. Our statistical sample consists of 122 normal type-1 AGNs, 86 red type-1s, and 271 type-2s, 118 non-AGNs, 23 stars, and 42 with featureless spectra.

We performed maximum likelihood fits both to the luminosity functions of the AGN population as a whole, and to the individual luminosity functions for the normal (blue) type-1s, type-2s, and red type-1s. From an initial $1/{{V}_{a}}$ estimate of the luminosity function, we noted that although a power law provided a visually reasonable fit to the data, there was a small but significant flattening at low luminosities, consistent with the luminosity function being a broken power law undergoing luminosity evolution (although our data do not probe below the break at $z\gt 1$, so the faint-end evolution is unconstrained at high redshifts).

We therefore fit a double power-law function of the following form:

$$\begin{eqnarray}&&\phi (L,z)=\frac{{{\phi }^{*}}}{{{(L/{{L}^{*}})}^{\gamma _{1}^{0}}}+{{(L/{{L}^{*}})}^{\gamma _{2}^{0}}}},\end{eqnarray} \tag{ 1 }$$

where ϕ is the comoving space density of the AGN, ${{\phi }^{*}}$ is the characteristic space density, both in units of comoving Mpc⁻³, L is the rest-frame luminosity at 5 μm, and L^* is the break luminosity at 5 μm, both in units of $\;{\rm erg}{{{\rm s}}^{-1}}\;{\rm H}{{{\rm z}}^{-1}}$. The evolution in L^* is the usual cubic expression (e.g., H07):

$$\begin{eqnarray}&&{{{\rm log} }_{10}}{{L}^{*}}(z)={{{\rm log} }_{10}}L_{0}^{*}+{{k}_{1}}\epsilon +{{k}_{2}}{{\epsilon }^{2}}+{{k}_{3}}{{\epsilon }^{3}},\end{eqnarray} \tag{ 2 }$$

where $\epsilon ={{{\rm log} }_{10}}((1+z)/(1+{{z}_{{\rm ref}}}))$, $L_{0}^{*}$ is a free parameter in the fit corresponding to the break luminosity at redshift zero, ${{\gamma }_{1}}$ is the faint end slope, and ${{z}_{{\rm ref}}}$ is fixed at 2.5. We also fit a model with a fixed faint end by setting ${{L}^{*}}=L_{0}^{*}$ in the left-hand term in the denominator of Equation (1) and allowing the bright-end slope, ${{\gamma }_{2}}$, to vary as

$$\begin{eqnarray}&&{{\gamma }_{2}}=\gamma _{2}^{0}+\tau \epsilon .\end{eqnarray} \tag{ 3 }$$

To perform the actual fits, we evaluated the log-likelihood function, ${\rm ln} \ \mathcal{L}$, in the usual way (Marshall et al. 1983), by minimizing $S=-2{\rm ln} \mathcal{L}$ (using the Powell minimization routine from the scipy.optimize library⁸ ):

$$\begin{eqnarray}\begin{array}{rcl} S & = & -2\mathop{\sum }\limits_{i=1,{{n}_{a}}}{\rm ln} [\phi ({{L}_{i}},{{z}_{i}})] \\ {} & {} & +2\mathop{\sum }\limits_{j=1,{{n}_{f}}}\int _{{{L}_{{\rm min} }}}^{{{L}_{{\rm max} }}}\int \phi (L,z){{A}_{j}}\frac{dV}{dz}dzdL \\ \end{array}\end{eqnarray} \tag{ 4 }$$

where the first term sums over the n_a AGN being considered, and the second sums over the n_f = 14 samples, each of area A_j. ${{L}_{{\rm min} }}$ is the minimum luminosity detectable at redshift z in a given sample, and ${{L}_{{\rm max} }}$ the maximum luminosity (for samples where there is an upper flux density limit). For the k corrections in this calculation, mean [5.8]–[8.0] and [8.0]–[24] spectral indices for each AGN type were measured and applied as appropriate to the AGN type and the redshift range of the calculation.

The fitting was performed by optimizing the parameters ${{\gamma }_{1}}$, ${{\gamma }_{2}}$, τ, $\phi _{0}^{*}$, ${{k}_{1}},{{k}_{2}}$, k₃, and ${{{\rm log} }_{10}}(L_{0}^{*})$ for AGNs with $z\gt 0.05$ and $L\gt {{10}^{29.5}}\;{\rm erg}{{{\rm s}}^{-1}}\;{\rm H}{{{\rm z}}^{-1}}$. We fit all three AGN types separately (normal type-1s, type-2s, and red type-1s). We also fit the type-2s and red type-1s together to give a luminosity function for the total obscured population. In order to assess the possible effects of incompleteness, we fit to a set including the objects classified optically as non-AGNs (some fraction of which are, in fact, AGNs; see Paper I) and also containing the featureless spectrum objects which have no redshifts, distributed evenly in comoving volume between $z=0.5-5$, and refer to this as the "Maximal" fit. The best-fit parameters are given in Table 1.

Table 1.  Best-fitting Luminosity Function Parameters

AGN used	na	S	${{{\rm log} }_{10}}({{\phi }^{*}}$	$\gamma _{1}^{0}$	$\gamma _{2}^{0}$	k1	k2	k3	${{{\rm log} }_{10}}(L_{0}^{*}$
			$({\rm Mp}{{{\rm c}}^{-3}}))$						$({\rm erg}{{{\rm s}}^{-1}}{\rm H}{{{\rm z}}^{-1}}))$	τ
All	446	10997	−4.75 ± 0.02	1.07 ± 0.06	2.48 ± 0.05	1.05 ± 0.03	−4.71 ± 0.13	−0.034 ± 0.19	31.92 ± 0.02	[0]
Normal Type-1	119	3374	−5.18 ± 0.05	0.25 ± 0.18	2.68 ± 0.10	0.537 ± 0.08	−5.48 ± 0.23	0.768 ± 0.76	31.99 ± 0.08	[0]
Type-2	243	6065	−5.04 ± 0.03	1.068 ± 0.12	2.75 ± 0.10	1.01 ± 0.04	−4.26 ± 0.14	−0.241 ± 0.25	31.76 ± 0.02	[0]
Red Type-1	84	2451	−5.39 ± 0.05	0.44 ± 0.20	2.60 ± 0.12	1.36 ± 0.09	−4.96 ± 0.27	0.232 ± 0.56	32.03 ± 0.03	[0]
Type-2+Red Type-1	327	8159	−4.98 ± 0.03	1.09 ± 0.06	2.61 ± 0.08	1.165 ± 0.04	−4.45 ± 0.14	−0.23 ± 0.20	31.91 ± 0.01	[0]
Maximal	588	13830	−4.66 ± 0.02	1.38 ± 0.04	2.47 ± 0.07	1.24 ± 0.05	−4.49 ± 0.07	0.336 ± 0.22	31.94 ± 0.01	[0]
Power law, All	446	11027	−5.28 ± 0.02	⋯	2.16 ± 0.04	1.30 ± 0.03	−4.10 ± 0.10	−0.475 ± 0.035	[32.0]	1.5 ± 0.1
Fixed faint end, All	446	10998	−4.91 ± 0.02	0.57 ± 0.03	2.72 ± 0.05	1.26 ± 0.03	−4.90 ± 0.13	0.539 ± 0.19	32.02 ± 0.02	1.5 ± 0.3

Notes. n_a is the number of AGNs used in the fit. ${{\phi }^{*}},{{\gamma }_{1}},{{\gamma }_{2}}$, and $L_{0}^{*}$ are defined in Equation (1), ${{k}_{1}},{{k}_{2}}$, and k₃ in Equation (2). Parameters fixed in the fit are indicated with square brackets. To convert these to approximate bolometric luminosity functions, add ${{{\rm log} }_{10}}(8.0\nu )=14.68$ to ${{{\rm log} }_{10}}(L_{0}^{*})$.

Download table as:  ASCIITypeset image

For the purposes of comparison with previous work, we also fit a pure power law of the form

$$\begin{eqnarray}&&{{\phi }_{{\rm pl}}}={{\phi }^{*}}{{(L/{{L}^{*}}[z])}^{{{\gamma }_{2}}}}\end{eqnarray} \tag{ 5 }$$

where ${{\gamma }_{2}}$ is allowed to vary according to Equation (3).

To determine whether the luminosity evolution model of Equation (1) was a better fit than the fixed faint-end slope model or the pure power-law model of Equation (5), we examined the likelihood ratio of the fits. For two fits giving values of S₁ and S₂ using models with ${{\nu }_{1}}$ and ${{\nu }_{2}}$ degrees of freedom, respectively, ${\Delta }S={{S}_{1}}-{{S}_{2}}$ is distributed as a ${{\chi }^{2}}$ distribution with ${{\nu }_{1}}-{{\nu }_{2}}$ degrees of freedom (e.g., Freeman et al. 1998). On this basis, the fixed faint-end slope model is slightly, but not significantly worse than the luminosity evolution model (${\Delta }S=1$, with one degree of freedom). The smallness of this difference is perhaps to be expected given the poor constraints on the faint end and its evolution. ${\Delta }S$ for the pure power-law model versus the broken power-law model is 30 (Table 1), and the difference in degrees of freedom is only one, and so the broken power law is formally a much better fit.

To make the binned estimates, we binned the data in five redshift bins ($0.05\lt z\lt 0.4$, $0.4\lt z\lt 0.9$, $0.9\lt z\lt 1.6$, $1.6\lt z\lt 2.5$, and $2.5\lt z\lt 3.8$) assigning a minimum luminosity to each bin at which the bin was approximately complete (unshaded region in Figure 1). The lower-redshift limit was chosen to avoid $z\lt 0.05$ where the Lacy et al. (2004) AGN selection can be confused by the 6.2 μm polycyclic aromatic hydrocarbon (PAH) emission. The bin spacing is in intervals of ${{{\rm log} }_{10}}(1+z)\approx \;0.13$.

The values of the binned estimates were obtained using the formalism of La Franca & Christiani (1997) and Hasinger et al. (2005) where the fitted luminosity function of Equation (1) is used to predict the number of sources in each bin, ${{n}_{{\rm pred}}}$, by integrating the luminosity function over the luminosity and redshift ranges for each bin for each individual field in the survey, properly accounting for the upper and lower flux limits in each field. The value of the luminosity function at the bin center, ${{\psi }_{{\rm mdl}}}$, is then corrected by the ratio of the number of predicted sources in the bin to the actual number observed, ${{n}_{{\rm obs}}}$, i.e.:

$$\begin{eqnarray*}&&\phi _{{\rm b}}^{ij}({{L}_{i}},{{z}_{j}})=\frac{{{\psi }_{{\rm mdl}}}{{n}_{{\rm obs}}}}{{{n}_{{\rm pred}}}}.\end{eqnarray*}$$

Although this has the disadvantage of making the binned estimate somewhat model dependent, the advantage of this technique is that it compensates both for the tendency of objects to be found close to the bottom of a luminosity bin and at the higher-redshift end of redshift bins where the luminosity function is increasing rapidly with redshift (i.e., at $z\lt 2$).

Figure 2 shows the evolution of the luminosity function of all spectroscopically confirmed Spitzer AGNs in our statistical sample with $z\gt 0.05$ and ${{L}_{5\;\mu {\rm m}}}\gt {{10}^{29.5}}\;{\rm erg}{{{\rm s}}^{-1}}\;{\rm H}{{{\rm z}}^{-1}}$. The left-hand panel, Figure 2(a), shows that the luminosity evolution model seems to fit the data well across most of the range of redshifts and luminosities. In the center panel (b), we show the alternative fixed faint-end slope and simple power-law models, both of which seem to represent the data reasonably well (although the single power law is statistically a worse fit than either of the broken power-law models). In the right-hand panel (c), we show the "Maximal" model compared to the luminosity evolution one. Deeper surveys, with AGN identifications at $z\gt 1$, are needed to determine the faint end behavior and evolution, although Juneau et al. (2013) show that the faint end of the luminosity function seems to have evolved little since $z\approx \;1$. The evolution at $z\gt 2$ is relatively poorly constrained, but, as described in Section 2.6, there is evidence for a turnover in the number densities of high-luminosity AGNs at $z\approx \;2.8$. AGNs at the break luminosity contribute most of the luminosity density at a given redshift. Our finding that the break luminosity is a strong function of redshift, decreasing with cosmic time at $z\lt 2.8$, is consistent with the "downsizing" scenario of galaxy evolution (Cowie et al. 1996) in the sense that for a given mean Eddington rate, the characteristic mass of the black holes contributing the bulk of the AGN luminosity density is decreasing with cosmic time.

Zoom In Zoom Out Reset image size

When the optically classified "non-AGNs" and blank spectra with randomly assigned redshifts are added into the survey to produce the Maximal luminosity function, as shown in the right-hand panel of Figure 2 as the dotted–dashed line, we see an increase in the magnitude of the faint-end slope. At high luminosities above the break, however, we see little more than an increase in the normalization (as expected). The larger numbers of low-z AGNs are consistent with the population of optically classified "non-AGNs" being in fact dominated by low-luminosity, heavily obscured AGNs. These conclusions must remain speculative, however, as it is very difficult to rule out the possibility that the some of the "non-AGNs" are interlopers.

In Figure 3, we plot a comparison of our mid-infrared AGN luminosity function with the hard X-ray luminosity functions from La Franca et al. (2005) and Aird et al. (2010; LADE model), and the mid-infrared luminosity function from Matute et al. (2006; based on surveys using the Infrared Space Observatory) at $z\lt 0.8$, and the estimate of the bolometric luminosity function of H07. In Figure 4, we plot the normal (blue) type-1 quasar luminosity function from our survey and compare to both the combined Sloan Digital Sky Survey (SDSS) and 2dF quasar surveys luminosity function (Croom et al. 2009), and the mid-infrared luminosity for type-1 AGNs from Brown et al. (2006) at $z\gt 0.8$ (see also Assef et al. 2011 for a more extensive study of the AGN luminosity function in the field used by Brown et al., including the X-ray AGN).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

For both simplicity and ease of comparison, a single bolometric correction, independent of redshift and luminosity, was assumed for each population (Table 2). For our mid-infrared AGNs, we assume a uniform bolometric correction of 10. This number was derived as follows. We began with the bolometric correction from 5 μm to the total type-1 quasar emission, $C_{\nu }^{\prime }=8.0$, from Richards et al. (2006). We then corrected $C_{\nu }^{\prime }$ for unobscured type-1s to allow for the fact that the mid-infrared emission is likely to be much more isotropic than the optical/UV, and to have its origin in reprocessed UV light originally emitted away from our line of sight. We thus derive a "true" bolometric correction for type-1 objects, ${{C}_{\nu 1}}\approx \;0.6C_{\nu }^{\prime }$ = 5, where 0.6 is the typical ratio of the UV/optical to the total quasar emission. Considering the type-2s on the other hand, the lower mid-infrared luminosities of type-2 objects at a given radio luminosity (Paper I, Section 6) suggest that type-2 objects have higher mid-infrared extinctions than the type-1s. (This is also manifested as a steeper 8–24 μm spectral index for the type-2s, $\;\approx 1.4$ compared to $\;\approx 1.0$ for the type-1s.) Residual extinction at mid-infrared wavelengths results in an underestimate of the mid-infrared luminosity for the type-2s by a factor of $\;\approx 3$. Thus the effective value of the bolometric correction for the type-2 population will be ${{C}_{\nu 2}}\approx \;3{{C}_{\nu 1}}=15$. For simplicity, we average these two corrections to obtain a value of ${{C}_{\nu }}=10$. However, we will return to this issue in more detail in future work where we will undertake a detailed comparison of type-1 and type-2 AGN SEDs.

Table 2.  Assumed Bolometric Corrections

Wavelength	Correction	Reference
Mid-infrared (5 μm)	10.0	This paper
Hard X-ray (2–10 keV)	40.0	Elvis et al. (1994a)
Optical/UV (1670 $\overset{\circ}{\rm A}$ rest)	4.5	Richards et al. (2006)

Download table as:  ASCIITypeset image

At low redshifts ($z\lt 0.3$), we expect to be incomplete due to contaminating emission from PAHs of up to $\;\approx 50$% (Paper I, Section 2.2), and indeed we find that at $\nu {{L}_{5\,\mu {\rm m}}}\sim {{10}^{44-45}}$ we are below most other AGN luminosity functions in our $0.05\lt z\lt 0.4$ bin, including the mid-infrared one of Matute et al. (2006). At $0.4\lt z\lt 0.9$, however, the luminosity function matches those in other wavebands much better, and we see a small excess over the X-ray luminosity function at low luminosities ($\lesssim {{10}^{44}}\;erg{\rm \ }{{s}^{-1}}$), and a very large one over the optical quasar luminosity function. Considering the blue type-1 objects in our sample alone, however, Figure 4 shows that we agree very well with prior estimates of the normal type-1 luminosity function at all redshifts.

As is already well known (Hasinger 2008; Aird et al. 2010), the hard X-ray luminosity function peaks at much lower redshift than either the mid-infrared or the optical quasar luminosity functions. This difference between the X-ray and other luminosity functions has often been ascribed to a luminosity-dependent bolometric correction (e.g., Steffen et al. 2006; H07, Han et al. 2012). Any redshift-dependent contribution due to absorption effects has been thought to be low, although it has long been known that X-ray absorbers are commonly seen in individual high-redshift quasars (Elvis et al. 1994b). The partial correlation analysis of Steffen et al.(2006) in particular seems to indicate a stronger luminosity than redshift dependence on the ratio between optical and X-ray fluxes, ${{\alpha }_{{\rm OX}}}$, at a given optical or X-ray luminosity (although for X-ray luminosities the result is only significant at 3σ), which would translate into a difference in the bright-end slope of the luminosity function. However, we find no evidence that a strong luminosity-dependent bolometric correction is required based on our comparison of luminosity functions, and we suspect that the apparently high (13.6σ) significance of the anticorrelation ${{\alpha }_{{\rm OX}}}$ on optical luminosity is a by-product of luminosity-dependent obscuration effects, i.e., AGNs with lower luminosities tend to have higher gas columns in the X-ray, increasing the magnitude of ${{\alpha }_{{\rm OX}}}$. In Figure 3, we see that at $0.4\lt z\lt 1.6$, the mid-infrared, optical quasar, and hard X-ray luminosity functions are very similar at the high-luminosity ($\nu {{L}_{5\,\mu {\rm m}}}\gt {{10}^{45}}\;{\rm erg}{{{\rm s}}^{-1}}$) end of the luminosity function. Instead of a difference in the slope of the AGN luminosity functions between the infrared and X-ray, which is what we would expect if there was a strong luminosity dependence to the bolometric correction, what we see is a fall-off in the $z\gt 1.6$ hard X-ray AGN population at all luminosities (particularly when the Aird et al. luminosity function is considered, although this luminosity function predicts somewhat low counts of the most luminous AGNs at $z\gt 3$; see Figure 9 of Kalfountzou et al. 2014). This suggests that the difference in the luminosity functions is predominantly a redshift effect. The lack of a luminosity-dependent bolometric correction in the X-ray is particularly illustrated when the bolometric luminosity function of H07 (which assumes a luminosity-dependent bolometric correction) is plotted in Figure 3, where we see a significant overprediction of the space density of $\nu {{L}_{5\,\mu {\rm m}}}\gt {{10}^{45.5}}\;{\rm erg}{{{\rm s}}^{-1}}$ AGNs at $z\lesssim 2$.

Figure 5 shows the evolution of the mid-infrared AGN population split up by type. The top two panels, (a) and (b), show the luminosity functions at mean redshifts of 0.65 and 1.25, respectively, for the different types of AGNs with the plotted points showing the binned estimates. The black dotted line indicates the fit for all confirmed AGNs, with the cyan line that for the Maximal model. Compared to the type-1 luminosity functions (shown in blue and red, for normal and red type-1s, respectively), the type-2 luminosity function (green) has a much steeper faint-end slope. In the lower two panels, (c) and (d), the number density of AGNs above luminosities of $\nu {{L}_{5\,\mu {\rm m}}}={{10}^{45}}$ and ${{10}^{46}}\;{\rm erg}{{{\rm s}}^{-1}}$, respectively, are plotted against redshift. Up to $z\approx \;2$, all of the populations shown (all AGNs, blue and red type-1s, type-2s, and the combination of red type-1s and type-2s, representing the obscured population) show a similar evolution. The blue type-1 population shows the same evolution as optically selected quasars, with a peak in number density at $z\approx \;2.8$, but the obscured quasar population (magenta line in panel (d)) shows a peak in number density that is broader and shifted to earlier epochs ($z\gt 3$). This is consistent with some theoretical predictions (Fanidakis et al. 2012), results from studies of the radio galaxy population (Jarvis et al. 2001), and a study of the infrared luminosity function of quasars from SDSS (Vardanyan et al. 2014). Our models are fairly poorly constrained beyond the peak, however, with only 13 objects of any type at $z\gt 2.8$, so the reality of the difference at high redshifts may be affected by small number statistics. Nevertheless, we have performed some tests to investigate whether the difference in the evolution we see might be statistically significant. We take for a null hypothesis that the evolution of obscured AGNs (type-2 and red type-1s) and the unobscured (type-1) population is the same. We then compare the S values for fits which fix the evolution of the obscured population to that of the unobscured one, and vice-versa. Fixing the values of the evolution parameters ${{k}_{1}},{{k}_{2}}$ and k₃ for the obscured population to those fit for the unobscured population, and allowing the remainder of the parameters to be fit results in a difference in the S values of 17, which, with the three degrees of freedom from fixing the ${{k}_{1}}-{{k}_{3}}$ parameter values, allows us to eliminate the null hypothesis that the evolution is the same at $\gt 99.9$% confidence. Forcing the evolution terms for the unobscured population to be the same as those for the obscured one results in a less significant ${\Delta }S=9.5$, $\;\approx 97$% confidence. We also compared to the evolution of all AGNs. As expected, fixing the evolution of the obscured AGNs to match that of all AGNs resulted in no significant difference, but the difference from forcing the type-1 evolution to match that of all AGNs resulted in ruling out the null hypothesis at the 95% confidence level.

Zoom In Zoom Out Reset image size

In summary, we find that the evolution of the obscured and unobscured AGN populations are significantly different, at $\gt 99.9$% confidence using our strongest statistical test. The most obvious difference in the luminosity function evolution is at high redshifts where the number density of high-luminosity obscured quasars appears to peak at a higher redshift than that of the unobscured population. Small number statistics, however, mean that this interpretation should be treated with caution.

Figure 6 shows histograms of the number of objects as a function of redshift and luminosity, split up by type. Figure 7 plots the fraction of obscured objects as a function of both luminosity and redshift with the measured points plotted with error bars and the luminosity function models plotted as dashed lines. These plots show that the fraction of obscured objects varies from $\;\approx 90$% for Seyfert galaxies to $\;\approx 50$% at the highest quasar luminosities. In addition, there is evidence that the nature of the obscuration in the survey seems to change as a function of redshift and/or luminosity. At low luminosities and redshifts, most of the obscured objects are highly obscured, type-2 objects, but at high luminosities and redshifts the obscured population is dominated by less obscured red type-1 objects. The right-hand panel suggests that despite the change in the nature of the obscuration, the overall obscured fraction is not a strong function of redshift when comparing bins at $z\lt 0.8$ and $z\gt 0.8$, and, indeed, as discussed in Section 2.6 below, it is not until $z\gt 2$ that we see a possible difference. Of course, selection effects may well be playing a role here as heavily obscured objects become progressively harder to identify at high redshifts even in samples selected in the mid-infrared. As we discuss in Paper I, however, we believe such selection effects are not having a severe effect on our sample, even at $z\gt 3$.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Figure 8 shows the AGN bolometric luminosity density and corresponding accreted black hole mass density since z = 5. Panel (a) shows the luminosity density as a function of redshift, with dusty type-1 objects dominating at high redshifts, normal type-1s at $z\sim 2$, and type-2 objects at $z\lt 1$. As our best-fit luminosity evolution model has a steep faint-end slope, we need to assume a low-luminosity cutoff at high luminosities in order to prevent an excessive amount of AGN luminosity density from being estimated. We assume, as a fiducial model, that the fitted form of the luminosity function we have adopted can be extrapolated a factor of three below the current limit of the survey. This results in luminosity densities in line with those we obtain from our alternative fixed faint-end slope model. With this assumption, we find more AGNs than previous surveys, and integrating to the present day (e.g., Martínez-Sansigre & Taylor 2009) indicates a total energy, U, emitted by all AGNs over the history of the universe of ${{{\rm log} }_{10}}(U/{\rm erg}\ \ {\rm Mp}{{{\rm c}}^{-3}})=37.5\pm 0.2$, where the lower bound comes from cutting off the faint-end of the luminosity function at the luminosity corresponding to the flux limit of the deepest sample in our survey (0.61 mJy at 24 μm), and the upper bound from extrapolating the luminosity function to 10 times fainter than this limit. To obtain the observed energy density today due to AGN, u, we need to divide the integrand by $(1+z)$ (Soltan 1982), resulting in a value of $u=2.8\pm 1.2\times {{10}^{-15}}\;{\rm erg}\ {{{\rm s}}^{-1}}{\rm c}{{{\rm m}}^{-3}}$, or a radiation intensity $b=uc/4\pi =6.7\pm 2.9\;{\rm nW}\ {\rm s}{{{\rm r}}^{-1}}$, within the range of $b=3.5-8\;{\rm nW}\ {\rm s}{{{\rm r}}^{-1}}$ estimated from the X-ray background, and suggesting that $\;\approx 12$% of the total radiation intensity of the universe is contributed by AGNs (Elvis et al. 2002). We therefore need a radiative efficiency of accretion, η at the high end of previous estimates to match the local mass density in black holes ($\eta \approx \;0.15$, Elvis et al. 2002; $\eta \approx \;0.06$, Martínez-Sansigre & Taylor 2009; Figures 8(b) and (c)). We use the compendium of local black hole mass density estimates of Graham & Driver (2007), with a mean of ${{10}^{5.7}}\;{{M}_{\odot }}\;{\rm Mp}{{{\rm c}}^{-3}}$ to estimate a best-fit value of $\eta =0.18_{-0.07}^{+0.12}$.

Zoom In Zoom Out Reset image size

As a function of cosmic time, the evolution of the black hole mass density shows a very rapid rise at $z\gt 3$, then continues to rise but at a smaller rate. Fifty percent of the current mass density has accreted by $z\approx \;1.8$. At $z\gt 3$, we predict that most black hole mass accretion occurs in dusty type-1 AGNs, this switches to normal type-1 quasars at $z\sim 2$, and to type-2 AGN at $z\lt 1$ when powerful quasars become very rare.