AGNs and Their Host Galaxies in the Local Universe: Two Mass-independent Eddington Ratio Distribution Functions Characterize Black Hole Growth

We assume that active black holes can be observed in two distinct accretion modes: radiative and radio mode. In the radiative mode, AGNs accrete via a geometrically thin, optically thick accretion disc (e.g., Shakura & Sunyaev 1973). The Eddington ratio and the radiative efficiency of these AGNs are of the order of 10%, i.e., the potential energy of the infalling matter is efficiently converted to radiation. In contrast, AGNs that are in the radio mode have a low radiative efficiency and low Eddington ratios. Their accretion flows are geometrically thick and optically thin, which results in a radiative cooling time that is longer than the infall time (Narayan & Yi 1994; Narayan & McClintock 2008). Radio or jet mode AGNs are often found in massive, red elliptical galaxies (e.g., Best et al. 2005), whereas radiative or quasar mode AGNs tend to be hosted by galaxies with bluer colors and higher SFRs (see discussion below and e.g., Nandra et al. 2007; Hickox et al. 2009; Netzer 2009; Schawinski et al. 2009; Treister et al. 2009; Griffith & Stern 2010; Goulding et al. 2014).

The two AGN populations show an overlap in the form of high-excitation radio galaxies. Based on optical emission lines, radio galaxies can be split into high- and low-excitation radio galaxies (HERGs and LERGs, Hardcastle et al. 2006, 2007; Smolčić 2009; Heckman & Best 2014; Ching et al. 2017). The distinction is linked to different AGN accretion modes (Hardcastle et al. 2006; Smolčić 2009): HERGs accrete via radiatively efficient accretion, whereas LERGs accrete via radiatively inefficient, advection-dominated accretion. HERGs, which are similar to Seyferts with low levels of radio emission, can thus also be detected using selection methods based on, for example, X-ray (Wong et al. 2016) or infrared (Hardcastle et al. 2013) emission. Bright HERGs are often referred to as radio-loud quasars.

Heckman & Best (2014) show that below the space density of LERGs is significantly higher than that of HERGs (also see Best & Heckman 2012; Gendre et al. 2013; Heckman & Best 2014; Pracy et al. 2016). While HERGs become dominant at high radio luminosities, their space density is several orders of magnitude lower than that of LERGs at, for example, . The number of radio AGNs that accrete via radiatively efficient accretion is thus negligible.

For our model we neglect this overlap between the radiatively efficient and inefficient AGN populations and consider X-ray and radio AGNs separately. We furthermore assume that they are detected in the hard X-rays and at 1.4 GHz, respectively.

In the local universe, the bimodality that we observe in color–mass and color–magnitude space is well established (Bell et al. 2003; Baldry et al. 2004; Faber et al. 2007; Martin et al. 2007; Schawinski et al. 2014). On one hand, star-forming galaxies, which are predominantly late-type/disc dominated, inhabit the "blue cloud." On the other hand, quiescent galaxies, mostly early type and/or bulge dominated, are found on the "red sequence." Between these two populations lies the transition zone, the "green valley" (Bell et al. 2004; Martin et al. 2007; Fang et al. 2012; Schawinski et al. 2014). Previous work has shown that X-ray-selected AGNs are found in galaxies that show signs of star formation, so either in the blue cloud or the green valley (Hickox et al. 2009; Schawinski et al. 2009; Silverman et al. 2009; Treister et al. 2009; Koss et al. 2011; Rosario et al. 2013; Goulding et al. 2014). For our simple model we assume that most of the X-ray-selected AGNs are hosted by blue and green galaxies and use the stellar mass function of these galaxies as input when constraining the ERDF relative to the XLF.

To be able to compare our predicted to the observed XLF, we need to convert from bolometric to hard X-ray luminosities. Although the spectral energy distributions of AGNs show a high degree of uniformity from the X-rays to the infrared (Elvis et al. 1994; Richards et al. 2006a), previous studies have suggested that the X-ray bolometric correction might be luminosity or Eddington ratio dependent (Marconi et al. 2004; Hopkins et al. 2007; Vasudevan & Fabian 2009; Lusso et al. 2010, 2012). We consider the simplest model and assume a constant bolometric correction. We discuss the effect that a luminosity-dependent bolometric correction would have on our results in Appendix B.2.

Rigby et al. (2009) determine to be 1.34. We use the same value to convert from 15–55 keV to 2–10 keV. Averaging over a wide range of Eddington ratios, we assume based on the results by Vasudevan & Fabian (2009). We conclude:

In our model we assume that the 1.4 GHz RLF is produced by radiatively inefficient AGNs which are hosted by red, quiescent galaxies. The fact that low Eddington ratio radio AGNs are primarily found in massive ellipticals is well established (Matthews et al. 1964; Yee & Green 1987; Best et al. 2005; Hickox et al. 2009). The HERG and LERG distinction has been introduced more recently (Hardcastle et al. 2006, 2007; Smolčić 2009; Heckman & Best 2014). As we discussed above, the RLF is dominated by LERGs which accrete via radiatively inefficient accretion. These LERGs are predominantly hosted by optically quiescent, red galaxies (Smolčić 2009; Smolčić et al. 2009; Janssen et al. 2012; Ching et al. 2017).

To convert our predicted bolometric luminosity functions to 1.4 GHz, we use a constant bolometric correction of

This is a simplified assumption since the connection between the AGN bolometric luminosity and the 1.4 GHz core emission is complex. First of all, radio loudness , i.e., the ratio of radio to optical emission, is a function of the Eddington ratio (Woo & Urry 2002; Ho 2002); the higher λ, the lower the radio loudness. Depending on whether the total or only the core radio luminosity is taken into account, there is furthermore evidence for (Sikora et al. 2007) and against (Broderick & Fender 2011) a radio-loud/radio-quiet dichotomy in space. In addition, radio AGNs release energy in the form of mechanical jet power and so only a fraction of their total energy is emitted in the form of radiation (Falcke & Biermann 1995; Bîrzan et al. 2004, 2008; Körding et al. 2008; Cattaneo & Best 2009; Cavagnolo et al. 2010; Plotkin et al. 2012; Turner & Shabala 2015; Godfrey & Shabala 2016; Mingo et al. 2016).

To constrain the shape of , we need to propose a value. Motivated by the fact that radio AGNs, specifically LERGs, tend to have low Eddington ratios compared to radiatively efficient AGNs or HERGs (Ho 2002; Evans et al. 2006; Hardcastle et al. 2007; Merloni & Heinz 2008; Hickox et al. 2009; Smolčić 2009; Alexander & Hickox 2012; Best & Heckman 2012; Fabian 2012), we choose .

We stress that we do not claim to be precisely −3. This value represents an assumption which allows us to test whether the observed RLF is consistent with a mass-independent ERDF and a constant bolometric correction. While summarizing all complexities that affect into a single constant allows us to constrain the shape of the radio ERDF, we do not claim to be able to determine the absolute value of the ERDF break, as we will discuss in more detail below.

To convert stellar to black hole masses, we use the local scaling relation. We parameterize the relation as

and choose μ = −2.75, β = 1 and σ = 0.3 dex. This μ value lies between the results of Häring & Rix (2004, μ = −2.8), Jahnke et al. (2009, μ = −2.75), Kormendy & Ho (2013b, μ = −2.31), McConnell & Ma (2013, μ = −2.54), Marleau et al. (2013, μ = −3.02) and Reines & Volonteri (2015, μ = −3.55). We do not take into account variations in the relation that might arise due to different levels of star formation, different morphological classifications, or using bulge instead of total stellar mass (Sani et al. 2011; Reines & Volonteri 2015; Savorgnan et al. 2016; Terrazas et al. 2016), and we use stellar mass M for the host mass M_host.

While alternative ERDF shapes have been discussed previously (Kollmeier et al. 2006; Kauffmann & Heckman 2009; Hopkins & Hernquist 2009; Cao 2010; Aird et al. 2012, 2013a; Bongiorno et al. 2012, 2016; Nobuta et al. 2012; Conroy & White 2013; Hickox et al. 2014; Veale et al. 2014; Schulze et al. 2015; Trump et al. 2015; Jones et al. 2016), we assume a broken power-law-shaped ERDF for both radiatively efficient and inefficient AGNs. A broken power-law-shaped ERDF solves the apparent discrepancy between a Schechter function-shaped stellar mass function and an observed broken power-law-shaped luminosity function. Furthermore, a broken power-law-shaped ERDF allows us to change both the high and the low Eddington ratio end slopes. Besides a functional form for the ERDF, we assume that the ERDF is mass independent. Mass-dependent ERDFs have been proposed (Schulze et al. 2015; Bongiorno et al. 2016), yet not introducing a mass dependence constitutes a straightforward first assumption. In Sections 4.3 and 6.1 we discuss a possible mass dependence and alternative ERDF shapes, respectively.

Our first method to predict the bolometric luminosity function is based on random draws and is illustrated in the top panels of Figure 3. We use a stellar mass function and an ERDF as input and perform the following steps.

• 1.  
  Determine N_draw: Before drawing from the stellar mass function, we need to determine the number of values to be drawn, N_draw. This is a simple scaling factor which determines the lowest number densities that we can probe with the random draw method. N_draw can be given a physical meaning by linking it to the stellar mass function. For example N_draw can be defined as:

  

  Here, and denote the minimum and maximum stellar mass values that we consider in the draw, respectively. can be chosen to represent the comoving volume (e.g., Hogg 1999) of a specific survey with solid angle Ω.
• 2.  
  Construct the cumulative distribution function (CDF) of the stellar mass function and draw stellar mass values: After having determined N_draw, we construct the CDF of the stellar mass function by computing:

  

  We randomly draw N_draw numbers between 0 and 1 from a uniform distribution and, by inverting , assign the corresponding stellar mass values.
• 3.  
  Convert all stellar masses to black hole masses by assuming a conversion factor μ with scatter σ: We assume the local relation (see Section 2.6) and draw N_draw conversion factors from a normal distribution with mean and scatter . is then given by the sum of these conversion factors and the previously drawn stellar mass values. Having determined allows us to construct the black hole mass function. We bin in and adjust the normalization so that in bin i is given by

  

  Here, n_i corresponds to the number of simulated systems in bin i, is the volume that we have already used in Equation (8) and represents the bin size.
• 4.  
  Construct the CDF of the ERDF and draw Eddington ratio values : We also assign randomly drawn Eddington ratio values to black hole mass value. In analogy to Equation (9) we construct the CDF for the given ERDF in the range to . The shape of the CDF is not affected by the normalization of the ERDF. However determines how many of the black hole mass values are assigned an Eddington ratio value. If the integral over from to is 1 we draw values from the CDF. However, if the integral over the ERDF is a with , then only black holes are assigned a value.
• 5.  
  Calculate bolometric luminosities to determine : We compute the bolometric luminosities corresponding to the black holes and their corresponding values:

  

  Here, . In analogy to the black hole mass function, we bin in . in bin i is then given by:

  

Zoom In Zoom Out Reset image size

Our second method to predict the shape of the luminosity function is based on convolutions and follows the work of C15. We again use a stellar mass function and an ERDF as input and go through the following steps (see the bottom row of Figure 3):

• 1.  
  Convolve the stellar mass function with a normal distribution to predict the black hole mass function: We convolve the observed stellar mass function with a normal distribution with mean μ and standard deviation σ:

  

  Here, is the normal distribution and is the stellar mass function which is described by either a single or a double Schechter function (see Equations (4) and (5)).
• 2.  
  Convolve the predicted black hole mass function with the given ERDF to determine the shape of the bolometric luminosity function: First, we assume a constant Eddington ratio of and shift from black hole mass to bolometric luminosity space

  

  To take into account the fact that the Eddington ratio is not constant, we convolve with the ERDF, :

  

Besides the ERDF parameters, both the random draw and the convolution method require choices for , , and . We use and since this is the mass range over which the W16 stellar mass functions are constrained. We discuss the effect of this choice on our results in Appendix B.1. We set and to −8 and 1, respectively (also see Table 2). The values of and are degenerate with the normalization of the ERDF. determines what fraction of black holes are assigned an Eddington ratio between and and can therefore be considered as being "on" (see Section 3.5). In the framework of the random draw technique, assigning values to all black hole mass values implies that all black holes in the sample can be considered AGNs. The convolution method produces the same result if the integral over from to is 1. Without a priori knowledge of for the chosen and boundaries only the shape of the bolometric luminosity function, but not its normalization , can be predicted. As we will discuss in more detail below, we determine by rescaling the predicted luminosity function so that the space density over the considered luminosity range matches the observed space density.

Once we have determined the bolometric luminosity function through either the random draw or the convolution method, we use the bolometric corrections and to predict the XLF and the RLF. In both cases, this results in a constant shift toward lower luminosities since our assumed bolometric corrections are constant with luminosity and Eddington ratio. We discuss the effect of a luminosity-dependent bolometric correction in Appendix B.2.

Both methods allow us to predict the bolometric luminosity function once we have assumed the shape of . As we discussed above, a broken power-law shape is an appropriate first assumption for the ERDF since, in contrast to the stellar mass function, the observed AGN luminosity function is also power-law shaped. C15 derive and discuss the properties of the predicted luminosity function if a broken power-law-shaped ERDF is assumed. We summarize these characteristics of and its dependence on and below and in Figure 2.

• 1.  
  The bright end slope : C15 use a simplified ERDF with for to show how the shape affects the predicted luminosity function. They show analytically that the bright end of the luminosity function has the same slope as the high λ end of the ERDF, that is:

  

  As the stellar mass function falls off exponentially at high stellar masses, the shallower slope is necessary to reproduce the observed bright end of the luminosity function. To construct the black hole mass function, we convolve the stellar mass function with a normal distribution. In contrast to the stellar mass function, the black hole mass function therefore no longer falls off exponentially at the high-mass end. In the extreme case, in which is steeper than the exponential cut off of , it is thus the high-mass end of the black hole mass function and not that dominates .
• 2.  
  The faint end slope : The faint end of the luminosity function is determined by either the stellar mass function low-mass end slope α or the ERDF low Eddington ratio end slope, . Using their simplified model C15 conclude that . Here, is the black hole mass function slope. If is steeper than , the ERDF slope determines . The linear relation which we assume between M and M_BH ensures that the stellar and black hole mass functions have the same low-mass end slopes. With we hence conclude⁹

  

  Equation (17) shows that for luminosity functions with we will not be able to constrain well, since it can take any value .
• 3.  
  The break L*: Under the assumption of , the position of the break of the luminosity function is given by:

  

  Here, is a small correction factor which is weakly dependent on the choice of and and varies by less than 0.15 dex.
• 4.  
  The normalization : At L* the normalization of the luminosity function can be predicted using

  

  is the normalization of the ERDF and we have used the fact that in our model . To derive this relation, C15 have again assumed and introduced , a small correction factor which is dependent on and and varies by less than 0.15 dex.

The discussion in the previous section shows that we are likely to find ERDF parameters that allow us to reproduce the observed XLF and RLF. To constrain the best-fitting parameters for and and to quantify the corresponding uncertainties, we introduce an MCMC sampler.

We use the broken power-law ERDF defined in Equation (6) and postulate . This ensures a predicted luminosity function shape similar to the observations and prevents the MCMC sampler from jumping between equivalent solutions during the sampling process. To incorporate this prior, we parameterize the broken power-law ERDF in the following way:

We use the MCMC python package cosmohammer¹⁰ (Akeret et al. 2013) to vary , , and . In each step, the MCMC proposes a new set of ERDF parameters. We use this prediction for the ERDF and the input stellar mass function to estimate the corresponding bolometric luminosity function with the convolution technique. To shift the luminosity function to the hard X-ray or the 1.4 GHz radio regime we apply a constant bolometric correction. We then determine the predicted space densities () in the luminosity bins of the observed luminosity function ().

The normalization of the ERDF, is degenerate with , and . To minimize the number of free parameters and degeneracies, we do not include in the MCMC. Instead we rescale so that the space densities of the predicted and the observed luminosity functions match:

Here and represent the integrals over the predicted and the observed luminosity functions within the observed luminosity function's binning range.

To compute the log-likelihood for each set of new ERDF parameters, we use the observed values and their errors and the predicted values in the bins. We assume that the asymmetric observed errors on follow a log-normal distribution and describe the details of the calculation in Appendix A.1.

Once the MCMC has converged we use the MCMC chain to determine the median , and values and the corresponding 16 and 84 percentiles. Using the definition of in Equation (20), we determine the sum of the and the chains to estimate the median value and its credible intervals.

It is important to acknowledge that the shape of the predicted luminosity function will not be affected by values that lie significantly below . According to Equation (17), will be given by α, the low-mass end of the stellar mass function, if . So, if , we have to interpret as an upper limit. We can only fully constrain , if is significantly steeper than the stellar mass function ().

As we have pointed out in Section 3.3, does not solely rely on , but is also affected by the steepness of the black hole mass function. We illustrate this effect in Figure 4. We vary the scatter in the stellar to black hole mass conversion by convolving the blue+green stellar mass function (panel (a)) with normal distributions with different widths σ. The corresponding black hole mass functions (panel (b)) are then all convolved with the same ERDF (panel (c)) to estimate the corresponding XLFs (panel (d)). Figure 4 illustrates that for low σ values the shape of the black hole mass function does not affect the shape of it is that determines . For high σ values, is however irrelevant and depends on the shallower high-mass end of the black hole mass function.

Zoom In Zoom Out Reset image size

In this section, we discuss how in our model the normalizations of the stellar mass function, black hole mass function, and the ERDF are linked to the black hole occupation fraction, the fraction of active black holes, and the duty cycle.

• 1.  
  Black hole occupation fraction: We define the black hole occupation fraction as the fraction of galaxies that are hosting a black hole. In the local universe most massive galaxies, including our own, host a black hole (Genzel et al. 1996; Ghez et al. 2000, 2008; Magorrian et al. 1998; Schödel et al. 2003). In our model, we thus assume an occupation fraction of 100%. At high redshift the occupation fraction might however be lower than 100% (Treister et al. 2013; Weigel et al. 2015; Trakhtenbrot et al. 2016). Furthermore, depending on the black hole seed formation mechanism, the black hole occupation fraction in dwarf galaxies might be lower than in more massive galaxies (Volonteri et al. 2008; Natarajan 2011; Greene 2012; Reines et al. 2013; Moran et al. 2014; Sartori et al. 2015). If this possible mass dependence is not taken into account, the black hole occupation fraction can be defined in the following way using the stellar mass and black hole mass function:

  
• 2.  
  Active black hole fraction: In the MCMC we initially use an ERDF which is normalized so that the integral over from to results in 1. In terms of the random draw technique this corresponds to assigning a between and value to each black hole (or to each galaxy since the occupation fraction is 1) in the sample. Based on this ERDF we predict the luminosity function. We then rescale the predicted luminosity function so that the space density of the predicted and the observed luminosity functions match. This rescaling factor allows us to determine . Furthermore we can compute the actual fraction of all black holes that have to be assigned a value. We refer to this as the active black hole fraction and define it in the following way:

  

  corresponds to the normalization of the normalized ERDF, and represents the normalization of the ERDF after the rescaling factor has been applied. Note that the active black hole fraction is always linked to the definition of and .

During its lifetime an AGN is expected to change its Eddington ratio and therefore its luminosity (de Vries et al. 2003; MacLeod et al. 2010; Novak et al. 2011; Kasliwal et al. 2015; Schawinski et al. 2015). The ERDF represents the distribution of Eddington ratios for all black hole masses at one moment in time. For instance, the ERDF shape implies that at the time of observation only few galaxies have high Eddington ratios. By postulating that the ERDF does not change significantly over a certain time range, we can further conclude that black holes evolve along the ERDF, spending a small fraction of their lifetime at high Eddington ratios. The ERDF alone does not contain timescale information. We are unable to constrain how quickly black holes change their Eddington ratios and move along the ERDF. Yet by assuming a lifetime model, the ERDF can be constrained from the luminosity function (Hopkins & Hernquist 2009). This leads us to the definition of the AGN duty cycle.

• 1.  
  Duty cycle: The AGN duty cycle is a unitless quantity. It describes the fraction of black holes that, at a given moment in time, have an Eddington ratio above a certain value . Alternatively, the duty cycle corresponds to the fraction of a black hole's lifetime that it is likely to spend at . In the framework of the random draw method this corresponds to the fraction of black holes (or galaxies since the occupation fraction is 1) that are assigned a value and for which . We thus define the duty cycle in the following way:

  

  The duty cycle depends on the definition of , but not the chosen value.

4.1.1. Galaxy Stellar Mass Functions

Our model is based on the stellar mass functions of red, blue and green galaxies in the local universe which we determine using the method and sample of W16. Stellar mass functions are constructed using SDSS DR7 (York et al. 2000; Abazajian et al. 2009) data with apparent magnitudes from the New York Value-Added Galaxy Cataloge (NYU VAGC, Blanton et al. 2005; Padmanabhan et al. 2008) and stellar mass values from the Max Planck Institute for Astrophysics John Hopkins University catalog (MPA JHU, Kauffmann et al. 2003; Brinchmann et al. 2004). The sample is restricted to the redshift range . In the u − r dust (Calzetti et al. 2000; Oh et al. 2011) and k-corrected (Blanton & Roweis 2007) color–mass diagram, galaxies lying above

are referred to as being red and galaxies lying below

are classified as being blue (W16). Galaxies lying between Equations (25) and (26) are part of the green valley and are referred to as green. For more detailed information on the sample see W16.

W16 combine the following three independent methods to generate the stellar mass functions of various subsamples: the classical 1/V_max approach of Schmidt (1968), the non-parametric maximum likelihood method of Efstathiou et al. (1988, SWML), and the parametric maximum likelihood technique of Sandage et al. (1979, STY). To estimate the stellar mass completeness W16 use the method of Pozzetti et al. (2010). In contrast to previous work, W16 do not make any a priori assumptions on which subsamples should be fit with a single or double Schechter function. Instead, a likelihood ratio test is used to determine the better fitting model.

We use the method and the sample presented in W16 to compute the stellar mass functions of the combination of optically blue and green and of optically red galaxies. Due to some randomness in the MCMC, the best-fitting STY parameters for the red stellar mass function are not equivalent to the values reported in W16. They do however lie well within the errors. Both the stellar mass functions of blue+green and of red galaxies are well described by double Schechter functions (see Equation (5)). The best-fitting Schechter function parameters are given in Table 1.

Table 1.  Input Stellar Mass and Luminosity Functions

Stellar Mass Functions	Reference
Blue + green mass function	W16	10.67 ± 0.02	−3.10 ± 0.10	−1.38 ± 0.05	−2.91 ± 0.18	−0.70 ± 0.11
Red mass function	W16	10.77 ± 0.01	−7.12 ± 0.77	−3.08 ± 0.50	−2.67 ± 1.09	−0.46 ± 0.02
Luminosity functions	reference		A/			⋯
Swift/BAT 15–55 keV XLF	A12	43.71 ± 0.12 )		0.79 ± 0.08	2.39 ± 0.12	⋯
Swift/BAT 14–195 keV XLF	Tueller et al. (2008)	43.85 ± 0.26 )				⋯
INTEGRAL17–60 keV XLF	Sazonov et al. (2007)	)				⋯
FIRST/NVSS 1.4 GHz RLF	Pracy et al. (2016)	data points for radio AGN, see Table 2 in Pracy et al. (2016)
NVSS 1.4 GHz RLF	MS07	24.59 ± 0.30 )		0.49 ± 0.04	1.27 ± 0.18	⋯

Note. Overview of the stellar mass functions and luminosity functions that we use as input for our model. To construct the stellar mass function of green and blue galaxies we use the method and sample presented in W16. values are given in units of and not in units of as in W16. For the XLF of A12 we use their non-evolving model fit. For the RLF of MS07 we have converted the normalization A from to . Pracy et al. (2016) do not report functional fits to their radio luminosity functions. In Figure 1, for instance, we thus only show their data points.

Download table as:  ASCIITypeset image

4.1.2. AGN Luminosity Functions

We run the MCMC twice. First, we use the stellar mass function of blue and green galaxies, a logarithmic bolometric correction of and the XLF of A12 to find the best-fitting parameters for , the ERDF of radiatively efficient AGNs. To constrain we consider luminosities between and . Second, for the ERDF of radiatively inefficient AGNs, , we use the red stellar mass function, a logarithmic bolometric correction of and the RLF of MS07. We consider luminosities between and . For both and , we use and . We consider stellar masses between and . This is the stellar mass range that was considered in W16 and we refrain from extrapolating the Schechter function fits to lower stellar masses. In Appendix B.1 we discuss what effect constraining our analysis to this mass range might have on our results.

When running the MCMC we include an initial guess for , and (through ) and constrain the parameters to lie within certain ranges. These initial values and priors are summarized in Table 2. To derive an initial guess for we use Equation (18), neglecting the correction factor . For and we use L* from A12 and MS07 (see Table 1), respectively:

Constraining the allowed parameter ranges when running the MCMC is especially important for the RLF. Equation (18) shows that the break of the luminosity function is dependent on the value. The RLF is shallow and so without a constraint on , the MCMC tries to fit the entire RLF with a single power law: by pushing to low values the entire RLF is fit with the high Eddington ratio end of the ERDF () and the break is ignored. To ensure that instead and are determined by and , respectively, we include a stringent constraint for in the MCMC by restricting it to values between −3.0 and −2.0. In contrast to the RLF, the XLF is steeper. The MCMC thus automatically projects onto . We hence use , testing a wider range of values than for the RLF.

Table 2.  Summary of MCMC Settings

MCMC Parameters	Radiatively Efficient AGNs	Radiatively Inefficient AGNs
Comparison LF	A12	MS07
,	9, 12	9, 12
log. bolometric correction kbol	−1	−3
Fitting ,	41.5, 45.6 )	20.4, 26.4 )
Broken power law
See Section 4.2	Figure 5	Figure 5
Initial	−1.4	−2.5
Initial	0.5	0.5
Initial	2.4	1.3
,	−3.0, 0.0	−3.0, −2.0
,	−1.0, 1.0	0.0, 2.0
,	0.0, 5.0	0.0, 2.0
Schechter function
See Appendix C.1	Figure 14	Figure 15
Initial	−1.3	−2.7
Initial	−1.2	−1.5
,	−2.0, −1.0	−3.0, 0.0
,	−2.0, 1.0	−3.0, 0.0
Log-normal
See Appendix C.1	Figure 14	Figure 15
Initial	−2.8	−2.7
Initial	0.2	0.2
,	−7.0, −2.0	−10.0, −2.0
,	0.1, 1.5	0.0, 2.0

Note. Initial guess and priors for the MCMC runs. In Section 3.4 we discuss the details of the MCMC. We assume a broken power-law ERDF and vary , and . In Appendix C.1 we extend this discussion to log-normal and Schechter function ERDFs.

Download table as:  ASCIITypeset image

We show the MCMC results in Figure 5. The upper and lower rows show the results for the ERDF of radiatively efficient and inefficient AGNs, respectively. The left-hand panels show the three-dimensional probability distributions and the marginalized distributions for , and . In the central panels we illustrate the best-fitting ERDFs. The best-fitting ERDF parameters and the corresponding errors are given within the panels. The right-hand panels show a comparison between the best-fitting predicted and the observed AGN luminosity functions. Also included in the right-hand panels are the residuals of the predicted luminosity functions relative to the observed XLF of A12 (upper panel) and the observed RLF of MS07 (lower panel).

Zoom In Zoom Out Reset image size

The best-fitting broken power-law parameters for the ERDF of radiatively efficient AGNs are:

For radiatively inefficient AGNs we find:

As we mentioned in Section 3.4, the normalization of the ERDF is not a free parameter in the MCMC. Instead we predict the luminosity functions based on a that is normalized to have an integral of 1. We then rescale the predicted so that the integral over the predicted matches the integral over the observed luminosity function. This rescaling factor then allows us to constrain .

Given our and values, and imply active black fractions of ∼17 and ∼1, respectively (see Section 3.5). For radiatively inefficient AGNs this implies that every black hole is assigned an Eddington ratio. The unintuitive value of ∼17 for radiatively efficient AGNs is an artifact of our choice. We chose to use the same values for radiatively efficient and inefficient AGNs to allow for a direct comparison between the resulting ERDFs. However, this shows that to give a physical meaning to the active black hole fraction of X-ray AGNs we would have to choose a higher value. Equivalently we could introduce a cut-off in at the low-λ end.

Figure 5 shows that the shapes of the XLF and the RLF are reflected in the corresponding ERDFs. and have similar values, yet of is significantly steeper than of .

For the XLF δ₁ =  and hence . Here is the slope of the stellar mass function which dominates of blue and green galaxies. The ERDF thus determines the faint end of the XLF (see Equation (17)). The marginalized distribution for has a tail toward lower values. This shows that if , can no longer be constrained since the low-mass end of the stellar mass function determines . Similarly, the marginalized distribution of has a tail toward higher values. As we discussed in Section 3.4, the bright end of the luminosity function not only depends on , i.e., , but also on σ, the assumed scatter in the stellar to black hole mass conversion. For steep values is determined by the high-mass end of the black hole mass function (see Figure 4) and can no longer be constrained.

A similar trend can be seen in the marginalized distribution for . The stellar mass function of red galaxies is dominated by . Since , the is better constrained than the . The marginalized distribution for lacks the tail toward lower values. Due to the stringent constraints, the marginalized distribution is cut off at its maximum.

In our model the error on the input stellar mass function is not taken into account. Table 1 summarizes the errors on M*, , and α. M*, which is degenerate with , is well constrained for both the blue+green and the red population. As we are not constraining in the MCMC, the errors on and are not taken into account. Furthermore, as we discussed, the faint ends of the XLF and the RLF are dominated by and not the slopes of the stellar mass function. The errors on for the blue and green stellar mass function and on for the red stellar mass function thus do not affect our results significantly.

Equations (27) and (28) show that and are degenerate with the bolometric corrections and . As we have discussed in Section 2.5, our choice of is subject to significant uncertainties. Due to it being constant, changing the assumed value causes a shift of radio ERDF, but not a change in the shape. We have thus shown that the observed RLF is consistent with a mass-independent ERDF and a constant bolometric correction. We have determined the shape of , but we have not constrained the absolute value of .

Figure 5 shows that both the observed XLF and RLF are consistent with the simplest form of an Eddington ratio distribution, a mass-independent one. However, mass-dependent ERDF models have been used by, for example, Schulze et al. (2015) and Bongiorno et al. (2016). We now investigate the effect of a mass-dependent ERDF and quantify the allowed ξ variation with . We use the random draw method (see Section 3.1) since it allows us to change the ERDF CDF as a function of stellar mass. We separately vary , , and , while keeping the other two parameters constant. Our tests assume a linear dependence, for instance . Appendix B.3 contains the details of the test. Figures 10–13 show the results.

In summary, our analysis in Appendix B.3 shows that a mild mass dependence of either or is consistent with the observed AGN luminosity functions. For the RLF and the XLF, can be increased by an order of magnitude per order of magnitude in stellar mass.

For , can also be varied by up to ±1 per magnitude in stellar mass. only allows values that decrease by up to 1 per magnitude in stellar mass. For both ERDFs models lead to luminosity functions that are no longer power-law shaped and do not resemble the observed . In Appendix B.1 we show that for our chosen stellar mass range, the RLF and XLF are dominated by galaxies with . This also affects the mass-dependent models which we consider here. A mass dependence of or only leads to agreement with the observations if at M* the mass-dependent ERDF parameter is either equal to its best-fitting mass-independent value. Galaxies with are hence still convolved with the best-fitting mass-independent ERDFs which we determined in Section 4.2.

We conclude that our model predicts XLFs and RLFs that are consistent with the observations. Making the most straightforward assumptions possible and, for instance, ignoring all complexities that might be affecting , we are able to derive these AGN luminosity functions from the galaxy population. The shapes of the XLF and the RLF can be traced back to characteristic, broken power-law-shaped, mass-independent ERDFs. We explored a first-order perturbation to the model and quantified the mild mass dependence of the ERDF that is still consistent with the data. We showed that only certain models are allowed. Specifically, the ERDFs have to resemble our best-fitting mass-independent ERDFs for galaxies. A more extreme dependence of on leads to a deviation from the broken power-law shape of the observed XLF and RLF. Choosing the simplest model possible and not making any assumptions about the mass dependence of , we proceed to discuss the implications of mass-independent ERDFs.

When we introduced our assumption and our model in Section 2 we postulated the need for two ERDFs: one for radiatively efficient and one for radiatively inefficient AGNs. After having discussed our results and how the predicted luminosity function depends on the assumed ERDF (see Section 3.3), we now revisit this fundamental assumption.

The fact that there is no global ERDF which describes the Eddington ratio distribution of both radiatively efficient and inefficient AGNs is primarily due to the significantly different shapes of the XLF and the RLF. With a large difference between and the XLF is steep with a clear break at L*. Compared to the XLF, the RLF is shallow and only has a weak break.

Figures 5 and 6 show that and differ in their values. However, is degenerate with the bolometric correction (see Equation (18)). As we have discussed above, due to the large uncertainties that affect we do not claim to have constrained . We thus do not use the difference in to argue for the need for two ERDFs. Besides the difference in , and have different slopes. determines and is thus steeper for the XLF than for the RLF. Both the XLF and the RLF have a faint end that is steeper than the respective input stellar mass functions. is hence determined by . Unlike , the values that we determine for and are consistent with each other. We conclude that we are unable to reproduce both the steep XLF and the shallow RLF with a single, global ERDF.

Zoom In Zoom Out Reset image size

As both the XLF and the RLF are steeper than their respective input stellar mass functions, the blue+green and the red stellar mass functions do not significantly affect the luminosity function shapes. Their M* values have an effect on and , and their values impact and . Nonetheless and are unaffected by α.

Our results thus show that the different shapes of the XLF and the RLF can be accounted for by using different ERDFs. The fact that and differ in shape is however not due to our using different input stellar mass functions for the X-ray and the radio AGN populations. Using the same stellar mass function, for instance the mass function of the entire galaxy sample, would still result in a steep and a shallow . This is a result of our analysis and not an assumption that we could have made a priori. We hence do not claim to have shown that AGNs in blue+green and red galaxies have different ERDFs. Instead our results imply that X-ray- and radio-selected AGNs must have different ERDF shapes.

A mass-independent ERDF implies that the fraction of galaxies that host AGNs, the AGN fraction, is mass independent. Galaxies which host AGNs are randomly drawn from the galaxy population. Due to a flux or luminosity limit the AGN fraction can however be observed to be stellar mass or black hole mass dependent. At low black hole masses only AGNs with high Eddington ratios will be bright enough to lie above the flux or luminosity limit. The AGN fraction will thus be low. At high black hole masses, we will be able to observe AGNs with a range of Eddington ratios since all of them are bright enough to be detected or included in the sample. Compared to low black hole masses, the AGN fraction at high black hole masses will hence be higher. So, even though the AGN fraction is intrinsically mass independent, a luminosity or flux limit will make it appear mass dependent. Previously, this has for instance been discussed by Aird et al. (2012).

We have shown that to zeroth order the observed XLF and RLF are consistent with mass-independent ERDFs which implies mass-independent AGN fractions. In Section 4.3 we quantified the allowed mass dependence of and . A mild mass dependence of the ERDFs is consistent with the data and would manifest itself in a mass-dependent AGN fraction. Previous studies have reported AGN fractions that vary as a function of stellar mass (Best et al. 2005; Kauffmann et al. 2008; Hickox et al. 2009; Silverman et al. 2009; Haggard et al. 2010; Xue et al. 2010; Tasse et al. 2011; Janssen et al. 2012; Hernán-Caballero et al. 2014; Williams & Röttgering 2015). However, when interpreting these results it is important to ensure that a possible selection effect has been accounted for.

AGNs with radiatively efficient and inefficient accretion are thought to play different roles in the quenching of star formation. AGNs with radiatively efficient accretion are considered to be the cause of quenching (Sanders et al. 1988; Di Matteo et al. 2005; Cattaneo et al. 2009; Fabian 2012). AGNs with radiatively inefficient accretion may be keeping their host galaxies quenched (Bower et al. 2006; Croton et al. 2006; Springel et al. 2006; Somerville et al. 2008).

In their phenomenological model, Peng et al. (2010, 2012) showed that the stellar mass function of red and quiescent galaxies can be reproduced by splitting the quenching mechanism into two distinct processes. "Mass quenching" is a mass-dependent, but environment-independent, mechanism. "Environment quenching" summarizes mass-independent, but environment-dependent, processes. While environment quenching could, for example, be associated with ram pressure stripping (Gunn & Gott 1972) or strangulation (Larson et al. 1980; Balogh et al. 2000), AGN feedback could be the physical origin of mass quenching.

One could imagine a simplistic model in which mass quenching is caused by more AGN activity and thus more AGN feedback at high stellar masses. Massive galaxies could for example be more likely to host AGNs with particularly high Eddington ratios. We have shown that the observations are consistent with a mass-independent model. This is inconsistent with this simplest form of mass quenching. We also discussed models which allow for a mild mass dependence of the ERDF. The linear increase of with of up to 1 mag per magnitude in stellar mass is however inconsistent with the exponential increase of the mass quenching probability which the Peng et al. (2010, 2012) model requires. Nonetheless, more sophisticated models could link AGN feedback and mass quenching. For instance, regardless of a mass-independent AGN fraction, AGN feedback in low-mass galaxies could be more effective than in high-mass galaxies due to a shallower potential well. Furthermore, AGN feedback might still be a necessary but not sufficient condition for the quenching of star formation.

For AGNs with radiatively inefficient accretion, the mass-independent ERDF implies that galaxies of all masses are equally likely to host radio AGNs of high and low Eddington ratios. Massive ellipticals contain massive black holes and are thus radio bright, even for low Eddington ratios. Flux-limited radio surveys bias us toward these massive galaxies. Our results imply that radio activity is possible at all stellar and black hole mass scales and not just in massive ellipticals. Furthermore, maintenance mode can occur on all mass scales, if reaching a certain Eddington ratio is a sufficient condition.

After having discussed the need for two ERDFs and the implications of mass-independent ERDFs for the AGN fraction and for the quenching of star formation, we now consider what might cause X-ray- and radio-selected AGNs to have different shapes.

• 1.  
  Host galaxies: As we discussed above, X-ray-selected AGNs tend to be found in optically blue and green galaxies, whereas radio-selected AGNs are often hosted by optically red galaxies (Hickox et al. 2009; Schawinski et al. 2009; Silverman et al. 2009; Smolčić 2009; Smolčić et al. 2009; Treister et al. 2009; Koss et al. 2011; Janssen et al. 2012; Rosario et al. 2013; Goulding et al. 2014; Ching et al. 2017). Star-forming and quiescent galaxies are thought to be accreting in the cold and hot phases (Shabala et al. 2008). A first assumption could thus be that AGNs in red galaxies have to have a different ERDF, since they contain less gas than blue and green galaxies. Black hole growth is however not limited by the amount of available gas, but by the loss of angular momentum (Jogee 2006).¹¹
• 2.  
  Stellar mass: We assumed that both and are mass independent and showed that this choice of ERDF is consistent with the observations. Stellar or black hole mass is hence unlikely to be the physical property that determines the ERDF shape.
• 3.  
  Large-scale properties: Large-scale properties, such as merger events, environment, or halo mass, might be promoting higher or lower Eddington ratios. However, it remains unclear how these physical conditions would be transmitted over multiple orders of magnitude in distance, from the galaxy outskirts to the central black hole accretion disc (see e.g., Alexander & Hickox 2012). The ERDF shape is thus unlikely to be set by large-scale galaxy properties.
• 4.  
  Accretion process: On small scales, it might be the accretion process itself that explains the difference between the two galaxy populations and determines the Eddington ratio distribution. For accretion rates , we would expect accretion via a classical Shakura & Sunyaev (1973) disc. At accretion via advection-dominated accretion flows (Narayan & Yi 1994; Narayan & McClintock 2008) becomes dominant. If the shape of the ERDF, and especially , is purely set by this change between accretion modes, we would expect and to be similar in shape.

None of the aforementioned processes is likely to solely determine the ERDF shape. We hence conclude that the difference between and could be caused by processes on intermediate scales: the shape of the ERDFs could be set by how efficiently the gas is driven from the inner parsecs to the accretion disc, by the properties of the gas itself, or by both.

According to Hardcastle et al. (2007), cold and hot gas lead to radiatively efficient and inefficient accretion, respectively. The host galaxies of X-ray selected AGNs are thus likely to contain cold gas which, due to its temperature, has a clumpy structure. Radio AGN host galaxies are likely to have reservoirs of smooth, hot gas. Radio AGNs could therefore be fuelled by a continuous Bondi-like accretion flow. For X-ray-selected AGNs the clumpy structure of the gas could cause short, episodic bursts of high accretion rate.

Compared to , has a steeper slope . As λ increases we expect the radiative efficiency to decrease once we approach the Eddington limit. This is caused by the increasing significance of the radiation pressure and the transition from a thin to a slim accretion disc (e.g., Abramowicz et al. 1988; Laor & Netzer 1989; Netzer & Trakhtenbrot 2014; Sa̧dowski & Narayan 2016). This change in radiative efficiency could be connected to the steep, high Eddington ratio slope that we observe for the ERDF of blue+green galaxies. A physical limit to the clump size of the cold gas could provide an alternative explanation. In contrast, the shallow slope of the radio AGN ERDF could show that the continuous hot gas accretion flow is not affected by such limits.

Even though we cannot constrain well, it is clear that the two ERDFs overlap at least to some degree. Within the overlapping λ range, a galaxy can thus be part of both the radio and the X-ray population. At high Eddington ratios, such a galaxy would be bright in the radio and the X-rays and we would identify it as a radio-loud quasar. At lower λ values we would classify such a galaxy as a HERG. The fraction of AGNs that are detected both in the X-rays and the radio does not only depend on the overlapping λ range, but also and .

In Section 4.2 we determined the ERDF of all X-ray and all radio AGNs. However, within these two populations there might be groups of AGNs that follow different ERDF shapes. For instance, the ERDF of AGNs in galaxies that are undergoing a galaxy merger might differ from the overall . While such subpopulations are possible, and , which we determined in Section 4.2, have to be recovered when we consider the sum of all X-ray and radio AGNs. For simplicity we focus our discussion on the ERDF of radiatively efficient AGNs. The same arguments can, however, be applied to .

We consider two possibilities for how the ERDFs of such subpopulations might differ from the overall, general ERDF. First, of a subpopulation might still be power-law shaped, but might have , , and values that differ from . Second, could be split into two subpopulations: one group of AGNs that is assigned λ values from the high Eddington ratio end of and one group that covers the low λ range of the best-fitting X-ray ERDF.

In the simplest case could be made up of groups of AGNs that have the same , and values as , but differ in their normalization . Physically this would imply that the active black hole fraction varies among certain groups of galaxies.

Alternatively, could consist of groups of AGNs that have broken power-law-shaped ERDFs with different , and values. Figure 7 shows a generic example: We assume that the population of X-ray-selected AGNs consists of two groups. The first group has an ERDF which has the same , and values as . The ERDF of the second group has the same slope as the ERDF of group 1, but a significantly higher value. We consider two cases. In the top row of Figure 7 we show the stellar mass functions, the ERDFs, and the resulting individual and the total XLF for the assumption that group 1 and group 2 have comparable space densities. The panels in the bottom row show the results for the case in which group 2 has a significantly lower space density than group 1. In the right-hand panels we show the observed XLFs in gray.

Zoom In Zoom Out Reset image size

Figure 7 shows that in the case where group 1 and group 2 have comparable space densities, the significantly higher value of group 2 leads to the combined XLF no longer being consistent with the observations. This is not the case in the bottom row, where the space density of group 2 is significantly lower than the space density of group 1. The very few AGNs in group 2 having a higher value do not alter the shape of the combined XLF.

For Figure 7 we assumed that the AGNs in groups 1 and 2 are drawn from host galaxy populations with similar mass function shapes. Furthermore, we assumed that the active black hole fraction among both groups is the same. If for group 2 we assumed a stellar mass function with a significantly lower M*, this would have accounted for the higher value. Similarly, assuming that group 2 has a lower active black hole fraction than group 1 would have diminished the deviation from the observations which we highlight in the upper panels.

Figure 7 shows that within the population of X-ray AGNs groups of AGNs with ERDFs that differ from can exist. By summing over the ERDFs of all subpopulations and weighing them by the space density of AGNs that they apply to, we have to recover , as it represents the average ERDF of all X-ray AGNs. Introducing subpopulations with different ERDFs still leads to consistency with the observations, if the space density of these AGNs is low. However, if we model the population of all X-ray AGNs as consisting of, for instance, two groups of comparable space densities, it is crucial that the sum of their weighted ERDFs reproduces . For the example shown in the first row of Figure 7 this would require changing the ERDF of group 1.

Instead of having subpopulations with power-law-shaped ERDFs, could also be split into different parts at a given λ value. In the most extreme case, one group of X-ray AGNs could have an ERDF which covers . With a cut-off at , the second group could dominate the low-λ end. With similar stellar mass functions, the sum of these two AGN populations would reproduce and the observed XLF.

Such a scenario is possible mathematically, would not break the mass independence of the ERDF and cannot be excluded by our model. A split of the ERDF close to implies that some AGNs are constantly growing, whereas others remain at low accretion rates. If AGNs are unlikely to move between the two groups, we expect their to be a population of non-flickering AGNs (Schawinski et al. 2015) that constantly exhibit high accretion rates. Furthermore, runaway growth for the high-λ group of AGNs and a total stalling of growth among the low-λ group are likely to break local scaling relations. If AGNs move among the two populations, for instance because high accretion rates are caused by temporary events such as starbursts or galaxy mergers, we return to our previous discussion of what shapes the ERDF. While determining if certain conditions promote high accretion rates is an interesting question, there is no need to think of the ERDF as being made up of strictly separate parts. Furthermore, the amount of time that AGNs spend at high and low λ values is directly reflected by the ERDF shape.

Our model does not allow us to constrain the ERDF of such subpopulations. However, by constraining and , we have determined the average ERDF of all X-ray and radio AGNs. We stress that recovering the XLF and the RLF does not require the introduction of subpopulations with different ERDFs and the necessary additional assumptions. Instead, the multifaceted populations of X-ray and radio AGNs can be described by two separate, simple, and global ERDFs. Dissecting these populations into their subpopulations and determining how they contribute to their respective ERDFs, represents the first-order perturbation to our zeroth-order model and requires additional observational constraints.

A broken power-law-shaped ERDF resolves the discrepancy between the exponential decline at the high-mass end of the stellar mass function and the power-law-shaped bright end of the AGN luminosity function. Furthermore, a broken power-law-shaped ERDF allows us to control both and , the high and low Eddington ratio ends of the distribution. Other studies have used alternative functional forms (Kollmeier et al. 2006; Kauffmann & Heckman 2009; Hopkins & Hernquist 2009; Cao 2010; Aird et al. 2012, 2013a; Bongiorno et al. 2012; Nobuta et al. 2012; Conroy & White 2013; Hickox et al. 2014; Veale et al. 2014; Schulze et al. 2015; Trump et al. 2015; Bongiorno et al. 2016; Jones et al. 2016). With an exponential cut-off at high λ, a Schechter function ERDF provides a natural limit to the number of high and super-Eddington ratio sources. However it causes a steep decline at the bright end of the luminosity function. Besides Schechter function-shaped ERDFs, log-normal-shaped distributions have been used in the literature. We adjust the MCMC method and fit the observed XLF and RLF with Schechter and log-normal-shaped ERDFs to investigate the challenges that arise from these λ distributions. The details of how we adjust the MCMC are given in Appendix C.1 and the results are shown in Figures 14 and 15.

While Schechter function and log-normal-shaped ERDFs fail to reproduce the broken power-law-shaped fits to the observed luminosity functions, they lead to broad consistency with the measured values. On one hand, a log-normal-shaped ERDF succeeds in reproducing the values measured by A12. The MS07 RLF values can be reproduced by a Schechter function-shaped ERDF. On the other hand, a Schechter function-shaped and a log-normal-shaped lead to inconsistency between the predictions and the respective observations in the brightest luminosity bins.

Our analysis in Appendix C.1 shows that Schechter and log-normal-shaped ERDFs are likely to lead to significant disagreement between predictions and observations if a wider luminosity range is considered. If AGN luminosity functions are truly broken power-law shaped, Schechter and log-normal-shaped ERDFs are likely to fail at the bright and faint end, respectively. To be able to determine if a broken power-law-, Schechter function- or log-normal-shaped ERDF leads to better agreement between predictions and observations, we thus need luminosity functions that probe a wider luminosity range than those that we considered here.

Our approach and our assumptions make our model subject to caveats.

Observationally, the local black hole mass function is often fit by a Schechter, a modified Schechter, a double power law, or a log-normal function (Aller & Richstone 2002; Marconi et al. 2004; Shankar et al. 2004; Greene & Ho 2007). In our model, the black hole mass function high-mass end depends on the σ of the normal distribution that we assume to convert from stellar to black hole mass. The steepness of this slope affects our MCMC results, especially if σ is large and the black hole mass function is thus shallow at high masses. Our black hole mass function high-mass end might therefore not reflect the observations and might have affected the MCMC results for .

Making the simplest assumptions possible, we chose constant bolometric corrections for the hard X-rays and 1.4 GHz and did not take luminosity or Eddington ratio dependences into account (see Sections 2.3 and 2.5). Due to many underlying complexities affecting especially we were only able to determine shape of . As we discuss in Appendix B.2, a luminosity-dependent hard X-ray bolometric correction results in a change of 's and . The shallower value lies within the uncertainties of the MCMC results where we assumed . For this is not the case. Yet even after adjusting to compensate for the luminosity-dependent bolometric correction, is still steeper than the value of the RLF that we determined for a constant bolometric correction.

We chose not to extrapolate the input stellar mass functions beyond the mass range considered in W16. While stellar mass bins beyond do not affect our results, we are neglecting galaxies with . We quantify this effect in Appendix B.1.

Furthermore, we were unable to constrain the ERDF normalization with our approach. is not a free parameter in the MCMC, since every combination of , and requires its own value. is thus degenerate with the ERDF parameters. Instead of varying , we adjusted the ERDF normalization so that the predicted space density n_pred matched the observed one n_obs (see Equation (21)).

6.3.1. Phenomenological Models

6.3.2. Observational ERDF Constraints

So far we have treated radiatively efficient and inefficient AGNs as separate populations and have not distinguished between AGNs within these two groups. In future work we will investigate these straightforward assumptions further and explore the first-order perturbations to our model. For instance, as we discussed in Section 2.1, HERGs represent the overlap between the populations of X-ray- and radio-selected AGNs. Constraining the HERG ERDF will allow us to test our model further and to expand our discussion of which physical process might be shaping the ERDF. Furthermore, among the population of X-ray-selected AGNs, the ERDF of AGNs in major galaxy mergers is of great interest. In the context of the classical Sanders et al. (1988) model and the results by Treister et al. (2012) it will allow us to investigate the role of AGNs and mergers in the quenching of star formation further.

For AGNs, and especially X-ray-selected AGNs, we observe the effect of cosmic downsizing: the space density of more luminous AGNs seems to peak at higher redshift than the space density of fainter AGNs (Cowie et al. 1996; Ueda et al. 2003; Hasinger et al. 2005) . Within the local volume, the measurement of the XLF is thus limited by the small volume and the low number of luminous AGNs. While the errors on the local XLF are large, the broken power-law shape is supported by luminosity functions of the same shape at higher redshift (Vito et al. 2014; Georgakakis et al. 2015; Miyaji et al. 2015; Fotopoulou et al. 2016). Repeating our analysis for the XLF at higher z is hence possible. However, we need to consider that at higher redshift the stellar mass range within which we can constrain the stellar mass function decreases and the errors on the mass function itself increase (Ilbert et al. 2013; Muzzin et al. 2013).

Our results can be applied to hydrodynamical simulations and semi-analytical models. Without having to resolve or describe the inner few parsecs and the details of the accretion process, the gas temperature could indicate whether the black hole is accreting via efficient or inefficient accretion. Our results thus not only predict the ERDF shape, but also imply whether quasar or radio mode feedback should be used.

Overall, our analysis would benefit from a wider fitting range and smaller errors on the observed luminosity functions for which large surveys are necessary. In the hard X-rays, the upcoming 105 month Swift-BAT data release (flux limit K. Oh et al. 2017, in preparation) will allow for smaller errors on . At lower energies, the Chandra and XMM-Newton coverage of SDSS Stripe 82 (Stripe 82X; flux limit LaMassa et al. 2013a, 2013b) will be probing fainter luminosities than Swift-BAT. In the near future, eROSITA (Merloni et al. 2012), and specifically the eROSITA all-sky survey (eRASS; flux limit ), will be able to significantly increase the number of soft X-ray detected AGNs (Kolodzig et al. 2013). Compared to eRASS, ATHENA will be probing even fainter luminosities over a wider redshift range and area in the far future (flux limit Aird et al. 2013b; Georgakakis et al. 2013).

In the radio regime, upcoming surveys and facilities such as, VLASS, MeerKAT, ASKAP and SKA, will not only allow us to constrain the local RLF better, but will also enable us to probe it at much higher redshifts. The VLA Sky Survey¹³ (VLASS) will cover 34,000 at a flux density limit of in the 2–4 GHz range. VLASS will thus reach a similar depth to FIRST with a significantly better resolution of 2.5 arcsec. The Meer Karoo Array Telescope (MeerKAT; Jonas 2009) and the Australian SKA pathfinder (ASKAP; Johnston et al. 2007; DeBoer et al. 2009), are the South African and the Australian SKA pathfinder telescopes, respectively. Planned surveys include the 1.4 GHz MeerKAT international GigaHertz Tiered Extragalactic Exploration survey (MIGHTEE; at rms /beam, Jarvis 2012) and the Evolutionary Map of the universe (EMU; entire southern sky, rms , Norris et al. 2011). Photometric and spectroscopic redshifts for EMU will be provided by SkyMapper (Keller et al. 2007), WALLABY (Koribalski 2012) and TAIPAN (Beutler et al. 2011). Norris et al. (2013) summarize all upcoming radio continuum surveys with the SKA pathfinders. During its first stage of completion, the Square Kilometer Array (SKA; Dewdney et al. 2009) will include two instruments: SKA1-LOW (50–350 MHz) and SKA1-MID (350-14 GHz) (see e.g., Kapinska et al. 2015). Wide and deep radio continuum surveys with the SKA will increase the number of low redshift radio AGNs by orders of magnitude and will allow us to detect radio-quiet AGNs of luminosities as low as out to redshifts of (Smolčić et al. 2015).

When running the MCMC to find the best-fitting ERDF parameters we compute the log-likelihood of each set of tested parameters. Given , and , we determine in the bins of the observed luminosity function () using the convolution method. After rescaling to match the space density of the observed luminosity function, we then use and its corresponding errors to compute .

The errors on are asymmetric. It is thus inappropriate to determine based on a normal distribution and reduced values. Instead we assume that is distributed log-normally, parameterize as x and express the probability density function of x in the following way:

In each bin, we use the observed values and their errors to determine the medians of the log-normal distributions and the values that correspond to the 16 and 84 percentiles:

a represents a shift that ensures that all x values are . It is the same in each MCMC step. and correspond to the lower and the upper errors on , respectively.

To be able to compute we determine μ and σ of the log-normal distributions in each of the bins:

represents the percent point function, the inverse of the cumulative density function. () thus corresponds to the value at which the integral over a normal distribution with mean 0 and standard deviation 1 reaches 16 (84) percent.

if x is truly log-normally distributed. However, the log-normal distribution is only an approximation for the distribution. We thus add and in quadrature to compute σ:

Once μ and σ have been determined in each of the N bins, we compute using :

When applying the MCMC to the observed luminosity functions by A12 and MS07, we only consider stellar masses between and . This is the stellar mass range that was considered for the mass function determination in W16. By not extrapolating the stellar mass functions to lower and higher masses we are introducing a bias in the ERDF determination since the contribution of mass bins beyond is neglected. Here we quantify this effect for the red and the blue+green stellar mass function.

Figure 8 shows the contribution of individual stellar mass bins to the XLF in the top row and to the RLF in the bottom row. The left-hand panels show the input stellar mas functions, the central panels illustrate the assumed ERDF, and the right-hand panels show the predicted luminosity function. For this figure we extrapolated the STY (Sandage et al. 1979) results for the blue+green (top row) and the red (bottom row) stellar mass function to . These additional stellar mass bins were also considered in the convolution with and for which we assumed our results from Section 4.2. In the right-hand panels we show the contribution of individual mass bins with colored lines and the combined luminosity functions with solid black lines. To guide the eye the observed XLF and RLF are shown in gray.

Zoom In Zoom Out Reset image size

The top row of Figure 8 shows that the XLF is dominated by stellar mass bins at . The space density of stellar mass bins beyond is too low to have a significant effect on the XLF. Mass bins below only contribute to the faint end of the XLF. Within our fitting range of masses in the range contribute to the overall space density of the XLF.

While the blue+green stellar mass function is fit with a mild double Schechter function, the red stellar mass function is described by a strong double Schechter function. As the lower panels of Figure 8 show, the steep of the red stellar mass function causes a significant increase in the space density of low-mass objects if we extrapolate to masses below . The faint end of the RLF is hence dominated by the contribution of these low stellar mass bins. Within the range, masses account for of the total space density. Similar to the blue+green stellar mass function, the exponential cut-off of the red stellar mass function results in high-mass bins not affecting the predicted luminosity function.

Extrapolating the red stellar mass function to low stellar masses thus significantly changes the predicted luminosity function and, as a consequence, would affect our ERDF fitting results. The strong double Schechter function shape and the steep value of the red stellar mass function is primarily due to galaxies at the edge of the completeness limit of the W16 sample. An increase in of more than a dex between and is unphysical as previous work by, for instance, Baldry et al. (2012) shows. The effect on the RLF is thus likely to be less than and likely to be comparable to the result that lower stellar mass bins have on the XLF.

We acknowledge that by choosing not to extrapolate the input stellar mass functions we are introducing a bias in our analysis. The cut-off at does not affect our results. However, the lower limit of does result in us missing lower mass galaxies which contribute to the faint end of the predicted luminosity function. Including these lower-mass bins in our MCMC analysis (see Section 4.2) would result in a change in the best-fitting value. These lower-mass bins would, however, not change the fact that, based on mass-independent ERDFs, our simple model is capable of predicting luminosity functions that are consistent with the observations.

For our model we assume a constant bolometric correction for the conversion between bolometric luminosities and hard X-rays. Based on the work by Rigby et al. (2009) and Vasudevan & Fabian (2009) we use . Here we illustrate the effect that a luminosity-dependent bolometric correction has on our results.

To include in our model we use the random draw method (see Section 3.1). We introduce a linear dependence of on :

The slope of 0.72 is motivated by the work of Steffen et al. (2006) who determine 0.72 to be the best-fitting slope for the conversion between the 2500 Å and the 2 keV luminosity. To adapt this relation for our model, we renormalize the relation so that . For we assume the value of A12 (see Table 1).

Figure 9 illustrates our results. The blue dashed and red solid lines show the random draw results including a luminosity-dependent bolometric correction for two different ERDFs. For the blue model corresponds to our best-fitting solution for the X-ray ERDF with (see Section 4.2). The red solid lines show the results of an ERDF with shallower slopes and . The lower panels show the predicted XLFs for these ERDFs and the assumed dependence. Shown in black are the observed XLFs of Sazonov et al. (2007), Tueller et al. (2008), and A12.

Zoom In Zoom Out Reset image size

Instead of a constant shift, a luminosity-dependent bolometric correction causes the XLF to be compressed when we convert from bolometric to hard X-ray luminosities. The lower middle panel of Figure 9 shows that in this case our best-fitting X-ray ERDF from Section 4.2 is no longer consistent with the data. The effect can be compensated for by changing the ERDF slopes and . In the upper central panel the red solid line illustrates a possible modification of the ERDF which leads to consistency with the data given a luminosity-dependent bolometric correction. Compared to our best-fitting ERDF with a constant bolometric correction, we have decreased by 0.15 and by 0.9. For this change lies within the uncertainties that we determined for the best-fitting model. The change in exceeds these uncertainties. However, the shallower value that is necessary to compensate for a luminosity dependent bolometric correction is not consistent with the value that we determined for the radio ERDF. Note that a change in the normalization of Equation (36) would require an additional change of .

We conclude that introducing a luminosity-dependent bolometric correction in our model requires a change of the slopes and to ensure consistency with the data. A luminosity-dependent bolometric correction does however not change our main conclusion that a mass-independent ERDF is consistent with the observed XLF.

For our approach we assume that the broken power-law ERDF is constant and depends neither on stellar nor on black hole mass. In Section 4.2 we showed that, based on this simple assumption, we are able to reproduce the observed XLF and RLF. Aird et al. (2012) argue for a mass-independent ERDF, but Schulze et al. (2015) and Bongiorno et al. (2016) use mass-dependent ERDFs. An ERDF for which, for example, increases with could thus be plausible. We now test whether the observations can also be reproduced with such a and quantify the allowed variation in ξ. In general, we cannot claim that the ERDF is mass independent or that it has to follow a certain mass dependence, simply because we are unable to test all possible ways in which and could depend on . We can only test whether a certain ξ model is consistent with the observations. For simplicity, we limit our analysis to broken power-law ERDFs for which only one of the three parameters is mass dependent.

In the following we will not be examining , a mass-dependent ERDF normalization, which implies a mass-dependent active black hole fraction. We consider the simplest case and assume that is mass independent.

To incorporate a mass-dependent ERDF into our model, we use the random draw technique (see Section 3.1). To make the method computationally less expensive, we pre-determine the ERDF CDF in stellar mass bins (). We then calculate N_draw, draw values and assign . We determine the stellar mass bin for each value and use the corresponding ERDF CDF to draw a value. The values and the corresponding luminosity function are calculated and constructed in the same way as discussed in Section 3.1.

To allow for an easy comparison, we renormalize the luminosity functions that are predicted by mass-dependent ERDFs when plotting our results in the following way:

Here, is the maximum luminosity value that is reached with and for the XLF and for radiatively inefficient AGNs.

B.3.1. A Mass-dependent Break

Figures 10 and 11 show the random draw results for ERDFs with mass-dependent breaks for radiatively efficient and inefficient AGNs, respectively. We show , the resulting for all masses and the predicted luminosity functions in the left-hand, central, and right-hand columns, respectively. In the left-hand columns the vertical dashed lines highlight of the blue+green mass function (Figure 10) and of the red mass function (Figure 11). The horizontal dashed lines mark , which we determined with the MCMC, assuming mass-independent ERDFs (see Section 4.2). For comparison, we illustrate the best-fitting, mass-independent ERDFs as black solid lines in the central columns. In the right-hand column of Figure 10, the black solid lines show the observed XLF of A12. In the right-hand column of Figure 11, the black solid lines illustrate the 1.4 GHz luminosity function of MS07.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We assume and keep and constant, setting them to the values we determined with the MCMC for a mass-independent ERDF. In the top two rows, we vary the y-intercept b, but keep the slope a constant, while in the bottom two rows, we change a, but enforce .

When interpreting these mass-dependent ERDFs, it is important to keep the stellar mass function in mind. Our randomly drawn sample contains more galaxies than galaxies. The values at low stellar masses thus affect the predicted luminosity function more strongly than at higher stellar masses. Furthermore, the simple relations that we derived for the predicted luminosity function in Section 3.3 are no longer valid due to the mass dependence.

For the XLF, Figure 10 shows that for the mass-dependent break the predicted and the observed XLF are in agreement if and . More extreme a values are less likely since they result in L* values that are significantly higher or lower than the observed break.

For the RLF, Figure 11 shows a similar trend. We find an agreement between our predictions and the observations if and .

B.3.2. A Mass-dependent Slope

Our results for mass-dependent ERDF slopes for radiatively efficient and inefficient AGNs are shown in Figures 12 and 13, respectively. In the upper two rows we vary the low Eddington ratio slope . The bottom two rows show . is always constant and set to . Similar to Figures 10 and 11, we show and / as vertical and horizontal dashed lines in the left-hand columns, respectively. The mass-independent ERDFs that we determined with the MCMC are shown as black solid lines in the central columns. In Figure 12, the luminosity function of A12 is shown with black solid lines in the right-hand column. In Figure 13, we show the luminosity function of MS07 in the same way.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In these figures, we additionally highlight the regions in which or in the left-hand column. We previously defined and to be the ERDF's low and high ratio slopes, respectively and postulated . In the gray shaded regions in the left-hand columns, this definition is no longer valid and by varying (), we are changing the high (low) Eddington ratio slope.

Figure 12 shows that a mass-dependent slope can be excluded for the XLF. It leads to XLFs which are no longer broken power-law shaped (first and second row). A leads to agreement between the predicted and the observed XLF if and (third and fourth row).

For the RLF we examine mass-dependent ERDF slopes in Figure 13. The first two rows show that a mass-dependent does not lead to agreement between the predictions and observations. For only with leads to an RLF that resembles the observations (fourth row).

We have not quantified the goodness of fit for these mass-dependent models and have not taken the errors on the observed luminosity functions into account. Furthermore, we fixed two of the three free parameters to the best-fitting values for the mass-independent case which we determined in Section 4.2. Rather than fitting the constant parameters to achieve consistency with the observations for the mass-dependent ERDF model, we thus simply added the mass dependence on top of our best-fitting ERDF results. Nonetheless, our analysis shows that mild mass dependences of the ERDF still result in luminosity functions that resemble the observations. For the XLF and the RLF, mass dependences of order of up to 1 mag in per magnitude in M are consistent with the observations. Changing and as a function of stellar mass leads to more extreme variations of the predicted luminosity function. For both, and , the simple models that we proposed here can be excluded. For can be varied by up to ±1 per magnitude in stellar mass. For the ERDF of radiatively inefficient AGNs a that decreases by up to 1 per magnitude in M still leads to consistency with the observations. For all models which we considered here it is important that at M* the mass-dependent ERDF parameter corresponds to its respective best-fitting mass-independent value. M* galaxies, which within our mass range dominate the XLF and the RLF (see Appendix B.1), are thus still convolved with the best-fitting, mass-independent ERDFs which we determined in Section 4.2.

Our model is based on a broken power-law-shaped ERDF (see Section 3). Yet different ERDF shapes, such as a single power law, a Schechter function, or a log normal function, have been proposed.

A single power-law ERDF has, for instance, been used by Aird et al. (2012). Bongiorno et al. (2012), who also use a simple power-law ERDF, find evidence for a break at high Eddington ratios and Aird et al. (2013a) stress the need for a break in the ERDF to match the observed XLF. In recent work, Bongiorno et al. (2016) thus use a broken power-law fit to the data. Veale et al. (2014) test both a truncated power law and a log-normal-shaped ERDF. Furthermore, a log-normal ERDF is used by Kollmeier et al. (2006) and Conroy & White (2013) ("scattered light bulb"). Kauffmann & Heckman (2009) find a log-normal distribution for galaxies with high Eddington ratios and high star formation rates ("feast mode") and a black hole mass and stellar age-dependent power-law ERDF for galaxies with old stellar populations ("famine mode").

Hopkins & Hernquist (2009), Cao (2010), Hickox et al. (2014), and Trump et al. (2015) use a Schechter function ERDF. Jones et al. (2016) show that an intrinsic Schechter function ERDF can be reduced to the log-normal ERDF observed by Kauffmann & Heckman (2009) if the BPT (Baldwin et al. 1981) selection effects are taken into account. Nobuta et al. (2012) and Schulze et al. (2015) consider both a log-normal and a Schechter function ERDF.

Irrespective of the functional form of the ERDF, we do not distinguish between type-1 and type-2 AGNs. We assume unification (Urry & Padovani 1995; Urry 2004) and thus expect both groups to have the same intrinsic ERDF. In contrast to our model, Tucci & Volonteri (2017), for example, stress the difference between the type-1 (log-normal) and the type-2 AGN (power law) ERDF.

To investigate the differences between a broken power law, a Schechter function, and a log-normal ERDF, we use our MCMC method (see Section 3.4) to fit these additional functional ERDF forms to the observed luminosity functions. For the Schechter function ERDF we vary the low Eddington ratio end slope and the knee . The log-normal ERDF also has two free parameters: the peak and the width .

In analogy to Equation (4), we assume the following functional form for the Schechter function ERDF:

The log-normal ERDF is parameterized as:

In analogy to Section 4.2, we fix the fitting range to for radiatively efficient AGNs and to for radiatively inefficient AGNs.

Our priors are given in Table 2 and our results are shown in Figures 14 and 15. We observe the following.

• 1.  
  A Schechter function ERDF for the XLF: As , where is the low-mass slope of the blue+green stellar mass function, the low-λ end of the ERDF and not determines (see Equation (17)). This is also true for the Schechter function-shaped ERDF which we are considering in Figure 14: is well constrained and steeper than . lies close to our initial guess (see Equation (27)), but is higher than in the broken power-law case (see Section 4.2). Using a Schechter function-shaped ERDF leads to an XLF that is consistent with the majority of observed XLFs. Due to the exponential cut-off, it does, however, fail to reproduce in the brightest luminosity bin.
• 2.  
  A log-normal ERDF for the XLF: When using a log-normal-shaped ERDF, the MCMC fits the entire observed XLF with the declining part of the ERDF. is therefore over one dex lower than in the case of a broken power-law- or Schechter function-shaped ERDF. Due to its broadness, succeeds in reproducing the XLF within the luminosity range that we are considering.
• 3.  
  A Schechter function ERDF for the RLF: As , where is the high-mass end of the red stellar mass function, for the RLF is well constrained and has smaller uncertainties than in the case of the XLF. As in the XLF case, for the Schechter function-shaped ERDF is higher than the best-fitting value for a broken power-law-shaped ERDF. The predicted RLF is consistent with the measured values. However at the bright end it deviates significantly from the best-fitting broken power-law RLF of MS07.
• 4.  
  A log-normal ERDF for the RLF: In analogy to the XLF, reproducing the RLF with a log-normal-shaped ERDF leads to low values. In the case of the RLF, the best-fitting for a log-normal lies over five orders of magnitude below the value for a broken power-law- or Schechter function-shaped ERDF. With the broad log-normal-shaped ERDF we are able to reproduce most of the observed values. However the distribution is too shallow to reproduce in the brightest luminosity bin.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

We conclude that, in comparison to a broken power-law-shaped ERDF, reproducing the observed XLF and RLF with a log-normal or a Schechter function ERDF is more challenging. This is not unexpected since both functional forms have only two free parameters, whereas a broken power-law-shaped ERDF allows us to also vary both and . Using a log-normal-shaped ERDF results in a broad distribution with a low value which reproduces the shape of a power-law- or Schechter function-shaped ERDF. In the case of the XLF, a Schechter function-shaped fails and a log-normal-shaped succeeds in reproducing the observed values. Vice versa, Schechter and log-normal-shaped ERDFs respectively succeed and fail to predict luminosity functions that are consistent with the RLF observations.

To be able to determine whether a broken power-law-, Schechter function-, or log-normal-shaped ERDF leads to better agreement between our predictions and the observations, we need luminosity functions that cover a wider range in . The faint and the bright end are where Schechter function- and log-normal-shaped ERDFs show their weakness. This can already be seen in the analysis that we presented here, but is likely to become more evident at higher and lower luminosities. If the observed luminosity functions do indeed have a broken power-law shape, a Schechter function-shaped ERDF is likely to be too steep to produce the bright end. A log-normal-shaped ERDF is likely to be too shallow to predict a luminosity function that is consistent with the observations at the faint end.